R - Find "n" closest points to each point in SpatialPointsDataFrame

Question

I am working with a dataset that is similar to this one:

library(raster)
library(sf)

# Define CRSs
crs.latlon <- "+init=epsg:4326"
crs.sirgas <- "+init=epsg:5880"

# define edges in somewhere within brazil
(spbb <- matrix(c(-53,-44,-25,-20), nrow=2, byrow=F))

# create spatial points in lat/lon
spbb <- SpatialPoints(spbb, CRS(crs.latlon))

# convert to sirgas 2000
(spbb <- spTransform(spbb, CRS(crs.sirgas)))

# make a raster extent object rounded to km:
ext <- extent(5200000, 6046000, 7234000, 7756000)

# create raster
r <- raster(ext, ncol=252, nrow=156, crs=crs.sirgas) # ~3 km resolution
r[] <- runif(ncell(r)) * 10

# create station points
pts <- sampleRandom(r, 500, na.rm=TRUE, sp=TRUE, cells=T) #%>%
#plot(st_geometry(sf::st_as_sfc(pts)), col=sf.colors(4), axes = TRUE, cex=2, pch=19)

#check it
plot(r)
plot(pts, add=T, col="red", pch=19)

For each point in pts, I would like to find the n (i.e. 5) closest surrounding ones. The resulting data frame would be something roughly similar to:

id      nearest     distance
pt1       pt3         14607
pt1       pt204       8540
pt1       pt301       6306
pt1       pt455       4956
pt1       pt337       2145

And so on for each point in the dataset.

What is the fastest way to achieve that?

Answer 1

The geos functions provided in sf, combined with some matrix operations, will do the trick. First, convert and transform your points.

library(sf)

spbb                             # define your points here

spbb <- st_as_sf(spbb)           # convert to simple features
spbb <- st_transform(spbb, 5880) # transform to your desired proj (unit = m)

Then, use st_distance to create a matrix of distances between each point in the data frame. Replace the diagonal with NAs (because distance of a point to itself is 0).

dist_matrix   <- st_distance(spbb, spbb)           # creates a matrix of distances
diag(dist_matrix) <- NA                            # replaces 0s with NA

Lastly, use some matrix operations to get the smallest value of each row.

spbb$distance <- rowMins(dist_matrix)              # get the dist of nearest element
spbb$nearest  <- rowMins(dist_matrix, value = T)   # get the index of the nearest element

Answer 2

As I am not completely understanding what the exact output is you want (you want new "features" with every time the point and the five nearest (=500 new features), or do you want just to eliminate the points that are never the nearest, or just a dataframe with indicating the index of the 5 nearest point for every point, or....?), I will tell you how it is possible to calculate the distance between every point:

You can convert your pts to an sf object using st_as_sf(pts) and then calculate the distance between the points using:

dist = as.data.frame(st_distance(pts,pts))

Then you can extract for every point the 5 nearest ones. If you specify what the exact output is you want, I can perhaps give a more precise answer.

Answer 3

I ended up writing a code that works as expected. Starting from the example provided in the original post, the following should work:

# convert points to sf objects
pts <- st_as_sf(pts)

# creates a matrix of distances between the points
dist_matrix <- st_distance(pts)

# replaces 0s with NA
diag(dist_matrix) <- NA

# convert matrix to data frame and set column and row names
dist_matrix <- data.frame(dist_matrix)
names(dist_matrix) <- pts$cell
rownames(dist_matrix) <- pts$cell

# find the 5 nearest stations and create new data frame
library(tidyr)
library(dplyr)

near <- dist_matrix %>% 
  mutate(ID=rownames(.)) %>% 
  gather('closest','dist',-ID) %>% 
  filter(!is.na(dist)) %>% 
  group_by(ID) %>% 
  arrange(dist) %>% 
  slice(1:5) %>% 
  mutate(dist_rank=1:5)

Hope this is useful to someone else.