0

I have a csv file that contains the GPS coordinates of about 100 villages and their unique ID. I have an other csv file with the GPS coordinates of some potential destinations and their unique ID. I found a way to measure the shortest distance between each village and the nearest destination using get.knnx. However the get.knnx function only provides me with the distance, it does not provide me with the ID of the village nor the ID of the destination. For now I tryed this:

g1 = get.knnx(coordinates(points_period1$X), coordinates(villages$hid), k=1)

Here is my code :

library(sp)
villages <- read.csv("/path/villages.csv", header=TRUE)

head(villages,15)
     hid exact_lat exact_lon
1      1  16.97587  7.992230
2      2  16.97371  7.998500
3      3  16.96942  7.984560
4      4  16.95704  7.965890
5      5  18.74667  7.400229
6      6  18.73254  7.401260
7      7  18.71110  7.336340
8   7020  18.71110  7.336340
9   7060  18.71110  7.336340
10     8  16.31698  8.013590
11     9  15.62224  7.892040
12 10040  15.52151  7.034985
13 10050  15.52151  7.034985
14 10060  15.52151  7.034985
15 10070  15.52151  7.034985

coordinates(villages) <- c('exact_lon', 'exact_lat')
proj4string(villages) <- CRS("+init=epsg:32633")
villages = spTransform(villages, CRS("+init=epsg:32632"))

head(villages,15)
     hid
1      1
2      2
3      3
4      4
5      5
6      6
7      7
8   7020
9   7060
10     8
11     9
12 10040
13 10050
14 10060
15 10070

names(villages)
[1] "hid"

points_period1 <- read.csv("/path/points_p1.csv", header=TRUE)

head(points_period1,15)
     X layer        x        y optional
1    1     1 6.625000 17.56250     TRUE
2    2     1 6.633929 17.56250     TRUE
3   79     1 6.616071 15.59821     TRUE
4   80     1 6.625000 15.59821     TRUE
5   81     1 6.633929 15.59821     TRUE
6   84     1 6.633929 15.58929     TRUE
7  117     1 5.294643 15.50893     TRUE
8  128     1 5.116071 15.50000     TRUE
9  129     1 5.125000 15.50000     TRUE
10 130     1 5.133929 15.50000     TRUE
11 133     1 5.285714 15.50000     TRUE
12 136     1 6.160714 15.50000     TRUE
13 144     1 5.107143 15.49107     TRUE
14 145     1 5.116071 15.49107     TRUE
15 146     1 5.125000 15.49107     TRUE

coordinates(points_period1) <- c('x', 'y')
proj4string(points_period1) <- CRS("+init=epsg:32633")
points_period1 = spTransform(points_period1, CRS("+init=epsg:32632"))

head(points_period1,15)
     X layer optional
1    1     1     TRUE
2    2     1     TRUE
3   79     1     TRUE
4   80     1     TRUE
5   81     1     TRUE
6   84     1     TRUE
7  117     1     TRUE
8  128     1     TRUE
9  129     1     TRUE
10 130     1     TRUE
11 133     1     TRUE
12 136     1     TRUE
13 144     1     TRUE
14 145     1     TRUE
15 146     1     TRUE

names(points_period1)
[1] "X" "layer" "optional"


require(FNN)
g1 = get.knnx(coordinates(points_period1), coordinates(villages), k=1)
[1] "nn.index" "nn.dist"

write.csv(g1, file = "/path/dist_2011.csv")
0

Okay, here's a fully reproducible (cut and paste into a new R session) to illustrate:

Need these packages:

library(sp)
library(FNN)

Create a simple villages data set of 10 villages with a hid column:

villages = data.frame(
   hid=c(1,2,3,4,5,7007,9009,1001,123,4453),
   lat=runif(10), lon=runif(10))
coordinates(villages)=~lon+lat
proj4string(villages)=CRS("+init=epsg:27700")

Create a simple points_period1 data set with 5 points and some columns a bit like yours:

points_period1 = data.frame(X=letters[1:5], layer=1, optional=runif(5)>0.5, lat=runif(5), lon=runif(5))
coordinates(points_period1)=~lon+lat
proj4string(points_period1)=CRS("+init=epsg:27700")

That's all set up. Now find the nearest neighbour in points_period1 to each point in villages:

nn = get.knnx(coordinates(points_period1), coordinates(villages),k=1)

The $nn.index component is a matrix with a column for each nearest neighbour. In this case one nearest neighbour so one column:

> nn$nn.index
      [,1]
 [1,]    5
 [2,]    1
 [3,]    4
 [4,]    5
 [5,]    3
 [6,]    2
 [7,]    3
 [8,]    3
 [9,]    2
[10,]    5

This is a number from 1 to 5 for each of the ten villages showing which of the points_period1 is closest. To get the ten records, extract using the first column to extract rows:

points_period1[nn$nn.index[,1],]

which gives:

                coordinates X layer optional
5     (0.900748, 0.3987913) e     1    FALSE
1    (0.4582524, 0.5365765) a     1     TRUE
4   (0.02994982, 0.4127038) d     1    FALSE
5.1   (0.900748, 0.3987913) e     1    FALSE
3    (0.5235299, 0.7821527) c     1     TRUE
2    (0.3874531, 0.7162511) b     1    FALSE
3.1  (0.5235299, 0.7821527) c     1     TRUE
3.2  (0.5235299, 0.7821527) c     1     TRUE
2.1  (0.3874531, 0.7162511) b     1    FALSE
5.2   (0.900748, 0.3987913) e     1    FALSE

If you want to put that data onto the village records, use cbind:

> cbind(villages, points_period1[nn$nn.index[,1],])
                coordinates  hid X layer optional
5   (0.6365308, 0.02949442)    1 e     1    FALSE
1    (0.4260667, 0.5436705)    2 a     1     TRUE
4   (0.01214433, 0.1087677)    3 d     1    FALSE
5.1  (0.9500509, 0.3324811)    4 e     1    FALSE
3    (0.7907734, 0.8053311)    5 c     1     TRUE
2    (0.1087975, 0.7607013) 7007 b     1    FALSE
3.1  (0.5175766, 0.6968651) 9009 c     1     TRUE
3.2  (0.5229172, 0.9995429) 1001 c     1     TRUE
2.1   (0.2870176, 0.738454)  123 b     1    FALSE
5.2    (0.89698, 0.3847053) 4453 e     1    FALSE

Note the coordinates are of the village (the points_period1 coordinates are gone but you can always look them up again). So village with hid=7007 is nearest points_period1 with X="b" and so on.

  • Oh yeah absolutely. However, the issue is that instead of knowing the nearest point to p2[1,] I want to know the nearest point to p2[101,] with 101 being the ID of the village and not the line number as it appears on the data. I'd like the CSV file to include these columns "village ID","nn.index","nn.dist" – Marcel Campion Apr 25 at 10:41
  • Extract them from the relevant spatial data frame using that integer as an index. If that's not clear then please create a reproducible example data set (like mine - small but with an ID column) and I'll walk through it. – Spacedman Apr 25 at 11:02
  • I edited the question a little (I added a short description of the dataset). How would you extract the relevant spatial data frame? – Marcel Campion Apr 26 at 9:51
  • Adding a description of the dataset is not as useful as adding some code to create a sample dataset that we can use - I've done this in my question now and re-edited it from scratch to show how to merge two datasets by nearest-neighbour. – Spacedman Apr 26 at 12:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.