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

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

``````library(raster)
library(sf)

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)
``````

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?

• The `knn` function in the spatialEco package will do exactly this, with the added ability to expand the problem into multivariate space (ie., add covariates in addition to distance). Feb 28, 2020 at 18:53

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
``````
• Cool, I look forward to testing your suggestion! Which package is the function `rowMins` from? Feb 27, 2020 at 18:43
• matrixStats::rowMins Aug 11, 2021 at 5:54
• How does this give "n" closest points? Aug 16, 2021 at 13:51

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.

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.