terra 最近（）仅对所有 `to_id` 列返回 NA

Question

我有一些海上的地理点和一张瑞典海域的地图。我想知道每个点位于哪个海域，对于海域之外的点，到最近海域的距离，包括最近海域的ID。

library(terra)
#> terra 1.7.71

map <- vect("/vsizip//vsicurl/https://www.smhi.se/polopoly_fs/1.140307!/Havsomr_SVAR_2016_3b.zip")

pts <- data.frame(id = 11:25,
                  lon = c(15.29416656, 15.10333347, 15.27083302, 15.27250004, 15.1291666, 
                          15.31333351, 15.24166679, 15.30716705, 15.32333374, 15.38833332, 
                          15.41416645, 15.40666676, 15.27166653, 15.26083374, 15.1916666),
                  lat = c(56.1566658, 56.15000153, 56.15250015, 56.13999939, 56.16383362, 
                          56.15333176, 56.1558342, 56.17416763, 56.12366486, 56.14250183, 
                          56.15916824, 56.17666626, 56.17666626, 56.17300034, 56.1570015)) |> 
  vect(geom = c("lon", "lat"), crs = "EPSG:4326") |> 
  project("EPSG:3006")

正如您所看到的，当您运行这个简单的代码时，

to_id

始终是

NA

并且属于逻辑类型。我想要将

"HID"

中的

map

列放入

to_id

列（或地图中的其他 ID）。我可以通过某种方式得到它吗？

npol <- nearest(pts, map[, "HID"], centroids = FALSE)
npol
#>  class       : SpatVector 
#>  geometry    : points 
#>  dimensions  : 15, 5  (geometries, attributes)
#>  extent      : 506419.7, 525716.6, 6219890, 6225817  (xmin, xmax, ymin, ymax)
#>  coord. ref. : SWEREF99 TM (EPSG:3006) 
#>  names       : from_id    from_x    from_y     to_id distance
#>  type        :   <int>     <num>     <num> <logical>    <num>
#>  values      :       1 5.183e+05 6.224e+06      <NA>        0
#>                      2 5.064e+05 6.223e+06      <NA>        0
#>                      3 5.168e+05 6.223e+06      <NA>        0

Answer 1

不幸的是，

?terra::nearest

的文档没有提供

nearest()

的任何示例，尤其是

centroids = FALSE

。就我个人而言，如果该点位于多边形内，我预计距离为零，否则对应于到多边形最近节点的距离 - 但我无法确定这意味着什么。显然这似乎不符合我的期望。

但是，对于这些数据，我们可以使用

terra::distance()

简单地计算完整的距离矩阵，然后选择我们感兴趣的信息：相关多边形的索引（您可以使用它来索引来自

"HID"

的

map

列） ）和距离本身：

distm <- distance(pts, map[, "HID"])

pts[, "to_id"] <- apply(distm, MARGIN = 1, FUN = which.min)
pts[, "to_HID"] <- map[pts$to_id, ]$HID
pts[, "distance"] <- apply(distm, MARGIN = 1, FUN = min)

pts
#>  class       : SpatVector 
#>  geometry    : points 
#>  dimensions  : 15, 4  (geometries, attributes)
#>  extent      : 506419.7, 525724.3, 6219890, 6225817  (xmin, xmax, ymin, ymax)
#>  coord. ref. : SWEREF99 TM (EPSG:3006) 
#>  names       :    id to_id        to_HID distance
#>  type        : <int> <int>         <chr>    <num>
#>  values      :    11   279 560940-151740        0
#>                   12   288 560850-150580        0
#>                   13   279 560940-151740        0

根据此输出，您的 15 个点中有 3 个在岸上。

编辑：

只是出于好奇而使用

{sf}

来重现这一点：

library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUE

map_sf <- st_as_sf(map)

pts_sf <- st_as_sf(pts)

st_nearest_feature(pts_sf, map_sf)
#>  [1] 279 288 279 371 289 279 327 279 371 244 280 280 279 279 327

st_distance(pts_sf, map_sf) |> apply(MARGIN = 1, FUN = min)
#>  [1]  0.000000  0.000000  0.000000  0.000000 58.765524  0.000000  0.000000
#>  [8] 27.636855  0.000000  0.000000  7.755143  0.000000  0.000000  0.000000
#> [15]  0.000000

最后一次使用

{nngeo}

：

library(nngeo)

pts_nn <- st_nn(pts_sf, map_sf, returnDist = TRUE)
lapply(pts_nn, unlist)
#> $nn
#>  [1] 279 288 279 371 289 279 327 279 371 244 280 280 279 279 327
#> 
#> $dist
#>  [1]  0.000000  0.000000  0.000000  0.000000 58.765524  0.000000  0.000000
#>  [8] 27.636855  0.000000  0.000000  7.755143  0.000000  0.000000  0.000000
#> [15]  0.000000

Answer 2

我的数据库中有 15000 个点，所以速度是一个问题。对不同方法进行基准测试显示最接近的优势。我只能希望 terra::nearest 能够更新，因为它具有如此巨大的速度优势（例如比 nngeo 方法快 515 倍）。