是在 R 中使用 sf 的 MULTIPOLYGON 中的一个点

Question

我有一个数据框 (addresses)，其中包含 6 个位置的纬度/经度信息。我还添加了（出于说明目的）这个位置实际属于的国家（我的真实数据中没有这个）。然后，我使用 st_as_sf 函数将此数据帧转换为 sf 对象 (addresses_sf)。

然后我通过 rnaturalearth 包获取加拿大边界的 sf 对象并将其保存为 can。

我想确定 address_sf 中的哪些条目在加拿大，所以我使用 sf 包中的 st_intersection 函数。

不幸的是，当它运行时，它仅将第 2 个条目（45.41126，-75.70591）标识为在加拿大边界内，而无法在结果中标识其余 3 个条目。

我在下面附上了一些可重现的代码，希望它可以识别我做错了什么...

## Load necessary libraries
library(tidyverse)
library(rnaturalearth)
library(sf)

## Create data
addresses <- tibble::tribble(
                   ~lat,      ~lon, ~country,
               43.77264,  -88.4462,    "usa",
               45.41126, -75.70591, "canada",
                 42.235, -83.06181, "canada",
               45.38441, -70.84368, "canada",
               42.11125,  -83.1107, "canada",
               44.74441, -74.86597,    "usa"
               )

## Convert to sf
addresses_sf <- addresses %>% 
  rowid_to_column() %>% 
  st_as_sf(coords = c("lon", "lat"), crs = 4326)
  
## Get Canadian borders
can <- rnaturalearth::ne_countries(country = "Canada", returnclass = "sf") %>% 
  select(name) %>% 
  st_transform(crs = 4326)

## Determine which entries are in Canada
sf::st_intersection(addresses_sf, can)

更新：

感谢@Ray 将上述脚本中的错误突出显示为没有足够分辨率的错误之一。我用下面的行更新了对 rnaturalearth 的调用，代码现在按预期运行。

## Get Canadian borders
can <- rnaturalearth::ne_countries(country = "Canada", scale = "large", returnclass = "sf") %>% 
   select(name) %>% 
   st_transform(crs = 4326)

Answer 1

@Tony Machado，您运行遇到了解决问题。来自 naturalearth 的多边形（即加拿大）的边界线的分辨率可能会削减边界线。如果你绘制手头的问题，你会看到一些城市靠近边界。

library(ggplot2)
ggplot() + 
  geom_sf(data = can, fill = "lightblue") + 
  geom_sf(data = addresses_sf, color = "red", size = 2)

你可以

获得更高分辨率的边界（例如，在您的 ne_countries() 调用中添加 scale = "large"）；或工作
有缓冲区。

为了缓冲，将 WGS84 重新投影到 UTM 并使用更精细的距离测量可能是值得的。在 WGS84 中，dist 以弧度为单位。那可能很痛苦。但是从下面，您可以看到当我们缓冲加拿大边界线时，您增加了 includes。不幸的是，这也抓住了“底特律”:(

希望这能让你继续前进。

sf::st_intersection(addresses_sf, st_buffer(can, dist = 0.001))
Simple feature collection with 5 features and 3 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: -83.1107 ymin: 42.11125 xmax: -70.84368 ymax: 45.41126
Geodetic CRS:  WGS 84
# A tibble: 5 x 4
  rowid country name               geometry
* <int> <chr>   <chr>           <POINT [°]>
1     2 canada  Canada (-75.70591 45.41126)
2     3 canada  Canada   (-83.06181 42.235)
3     4 canada  Canada (-70.84368 45.38441)
4     5 canada  Canada  (-83.1107 42.11125)
5     6 usa     Canada (-74.86597 44.74441)

是在 R 中使用 sf 的 MULTIPOLYGON 中的一个点

Is a POINT within a MULTIPOLYGON using sf in R

r

sf