使用R中的terra返回最接近点的多边形

Question

我在做多边形分析中的一个点

library(terra)
library(rnaturalearth)
  
crdref <- "+proj=longlat +datum=WGS84"

lonlat<- structure(c(-123.115684, -81.391114, -74.026122, -122.629252, 
                      -159.34901, 7.76101, 48.080979, 31.159987, 40.621058, 47.50331, 
                       21.978049, 36.90086), .Dim = c(6L, 2L), 
                      .Dimnames = list(NULL,c("longitude", "latitude")))

pts <- vect(lonlat, crs = crdref)
world_shp <- rnaturalearth::ne_countries()
world_shp <- terra::vect(world_shp, crs = crdref)
world_shp <- terra::project(world_shp, crdref)

plot(world_shp)
points(pts, col = "red", pch = 20)      

所有这些点都位于多边形的边缘，因此当我尝试提取每个点所在的多边形时，我得到一个 NA

e <- terra::extract(world_shp, pts) 
e$sovereignt
NA
  

有什么方法可以使用 terra 包为每个点返回最近的多边形

Answer 1

您可以使用 sf 包中的 st_join:

library(rnaturalearth)
library(tidyverse)
library(sf)
library(rgdal)

devtools::install_github("ropensci/rnaturalearthhires")
library(rnaturalearthhires)


# create points dataframe
points <- tibble(longitude = c(123.115684, -81.391114, -74.026122, -122.629252,
                               -159.34901, 7.76101),
                 latitude = c(48.080979, 31.159987, 40.621058, 47.50331,
                              21.978049, 36.90086)) %>% 
  mutate(id = row_number())

# convert to spatial
points <- points %>% 
  st_as_sf(coords = c("longitude", "latitude"),
           crs = 4326)

# load Natural Earth data
world_shp <- rnaturalearth::ne_countries(returnclass = "sf",
                                         scale = "large")


# Plot points on top of map (just for fun)
ggplot() +
  geom_sf(data = world_shp) +
  geom_sf(data = points, color = "red", size = 2)


# check that the two objects have the same CRS
st_crs(points) == st_crs(world_shp)


# locate for each point the nearest polygon
nearest_polygon <- st_join(points, world_shp,
                           join = st_nearest_feature)

nearest_polygon