Returns the indices of layer y
which are nearest neighbors of each feature of layer x
. The number of nearest neighbors k
and the search radius maxdist
can be modified.
The function has three modes of operation:
lon-lat points—Calculation using C code from GeographicLib
, similar to sf::st_distance
projected points—Calculation using nabor::knn
, a fast search method based on the libnabo C++ library
lines or polygons, either lon-lat or projected—Calculation based on sf::st_distance
st_nn( x, y, sparse = TRUE, k = 1, maxdist = Inf, returnDist = FALSE, progress = TRUE, parallel = 1 )
x | Object of class |
---|---|
y | Object of class |
sparse |
|
k | The maximum number of nearest neighbors to compute. Default is |
maxdist | Search radius (in meters). Points farther than search radius are not considered. Default is |
returnDist |
|
progress | Display progress bar? The default is |
parallel | Number of parallel processes. The default |
If sparse=TRUE
(the default), a sparse list
with list element i
being a numeric vector with the indices j
of neighboring features from y
for the feature x[i,]
, or an empty vector (integer(0)
) in case there are no neighbors.
If sparse=FALSE
, a logical
matrix with element [i,j]
being TRUE
when y[j,]
is a neighbor of x[i]
.
If returnDists=TRUE
the function returns a list
, with the first element as specified above, and the second element a sparse list
with the distances (as units
vectors, in meters) between each pair of detected neighbors corresponding to the sparse list
of indices.
C. F. F. Karney, GeographicLib, Version 1.49 (2017-mm-dd), https://geographiclib.sourceforge.io/1.49/
data(cities) data(towns) cities = st_transform(cities, 32636) towns = st_transform(towns, 32636) water = st_transform(water, 32636) # Nearest town st_nn(cities, towns, progress = FALSE)#>#> [[1]] #> [1] 70 #> #> [[2]] #> [1] 145 #> #> [[3]] #> [1] 59 #># Using 'sfc' objects st_nn(st_geometry(cities), st_geometry(towns), progress = FALSE)#>#> [[1]] #> [1] 70 #> #> [[2]] #> [1] 145 #> #> [[3]] #> [1] 59 #>st_nn(cities, st_geometry(towns), progress = FALSE)#>#> [[1]] #> [1] 70 #> #> [[2]] #> [1] 145 #> #> [[3]] #> [1] 59 #>st_nn(st_geometry(cities), towns, progress = FALSE)#>#> [[1]] #> [1] 70 #> #> [[2]] #> [1] 145 #> #> [[3]] #> [1] 59 #># With distances st_nn(cities, towns, returnDist = TRUE, progress = FALSE)#>#> $nn #> $nn[[1]] #> [1] 70 #> #> $nn[[2]] #> [1] 145 #> #> $nn[[3]] #> [1] 59 #> #> #> $dist #> $dist[[1]] #> [1] 1425.896 #> #> $dist[[2]] #> [1] 1818.766 #> #> $dist[[3]] #> [1] 2878.653 #> #>if (FALSE) { # Distance limit st_nn(cities, towns, maxdist = 7200) st_nn(cities, towns, k = 3, maxdist = 12000) st_nn(cities, towns, k = 3, maxdist = 12000, returnDist = TRUE) # 3 nearest towns st_nn(cities, towns, k = 3) # Spatial join st_join(cities, towns, st_nn, k = 1) st_join(cities, towns, st_nn, k = 2) st_join(cities, towns, st_nn, k = 1, maxdist = 7200) st_join(towns, cities, st_nn, k = 1) # Polygons to polygons st_nn(water, towns, k = 4) # Large example - Geo points n = 1000 x = data.frame( lon = (runif(n) * 2 - 1) * 70, lat = (runif(n) * 2 - 1) * 70 ) x = st_as_sf(x, coords = c("lon", "lat"), crs = 4326) start = Sys.time() result1 = st_nn(x, x, k = 3) end = Sys.time() end - start # Large example - Geo points - Parallel processing start = Sys.time() result2 = st_nn(x, x, k = 3, parallel = 4) end = Sys.time() end - start all.equal(result1, result2) # Large example - Proj points n = 1000 x = data.frame( x = (runif(n) * 2 - 1) * 70, y = (runif(n) * 2 - 1) * 70 ) x = st_as_sf(x, coords = c("x", "y"), crs = 4326) x = st_transform(x, 32630) start = Sys.time() result = st_nn(x, x, k = 3) end = Sys.time() end - start # Large example - Polygons set.seed(1) n = 150 x = data.frame( lon = (runif(n) * 2 - 1) * 70, lat = (runif(n) * 2 - 1) * 70 ) x = st_as_sf(x, coords = c("lon", "lat"), crs = 4326) x = st_transform(x, 32630) x = st_buffer(x, 1000000) start = Sys.time() result1 = st_nn(x, x, k = 3) end = Sys.time() end - start # Large example - Polygons - Parallel processing start = Sys.time() result2 = st_nn(x, x, k = 3, parallel = 4) end = Sys.time() end - start all.equal(result1, result2) }