问题描述
我是地理数据的新手，所以也许这是一个愚蠢的问题。 我有两个多边形，我想要做的是测试每个“病房”多边形包含在哪个“选区”多边形中。基本上，查看每个病房在哪个选区。
但是，我遇到了 st_within() 函数 sf 的问题，因为它仅在第一个几何体完全在第二个几何体中时才返回 TRUE。当较小的单元接触较大单元的边界时，这是一个问题。这是一个例子
我有一个叫做“选区”的多边形和一个叫做“病房”的多边形。您可以在此处找到它们：https://www.dropbox.com/sh/xu08wm79rym00zz/AADFDpyPe0EuDSDSY-6SvqpYa?dl=0
选区中有 2 个对象：
plot(constituencies$geometry)
wards <- st_transform(wards,crs = st_crs(constituencies))
test98 <- st_join(wards,constituencies,join=st_within)
head(test98)
Simple feature collection with 6 features and 4 fields
geometry type: MULTIpolyGON
dimension: XY
bBox: xmin: 531936.6 ymin: 180569.5 xmax: 533595.1 ymax: 182082.5
projected CRS: OSGB 1936 / British National Grid
WD98CD WD98NM name label geometry
1 00AAFA Aldersgate <NA> <NA> MULTIpolyGON (((532104.9 18...
2 00AAFB Aldgate <NA> <NA> MULTIpolyGON (((533319.2 18...
3 00AAFC Bassishaw Cities of London and Westminster 107 MULTIpolyGON (((532587.3 18...
4 00AAFD Billingsgate Cities of London and Westminster 107 MULTIpolyGON (((533167.9 18...
5 00AAFE Bishopsgate <NA> <NA> MULTIpolyGON (((533410.7 18...
6 00AAFF Bread Street Cities of London and Westminster 107 MULTIpolyGON (((532300.3 18...
当我绘制它们以查看哪些名称中包含 NA（红色）时，您可以看到它是触及选区边界的那些：
plot(constituencies_short$geometry)
plot(test98$geometry[!is.na(test98$name),],col = "green",add = T)
plot(test98$geometry[is.na(test98$name),col = "red",add = T)
我假设发生了以下两种情况之一：1) 接触边界不算作完全在几何图形内，或者 2) 选区和选区的边界可能不完全匹配，因此选区并不完全在几何图形内选区。
我的问题是：是否有一种方法可以完美匹配边界，然后测试哪些选区在哪些选区范围内，或者是否有一种方法可以测试哪些选区不是整个选区而是选区边界的大部分？
非常感谢！
解决方法
您的边界没有完全对齐。让我们看看当我们将第一个选区与两个选区相交时会发生什么：
> st_intersection(wards$geom[1],constituencies$geom)
Geometry set for 2 features
geometry type: GEOMETRY
dimension: XY
bbox: xmin: 531936.6 ymin: 181262.6 xmax: 532309 ymax: 182012
projected CRS: OSGB 1936 / British National Grid
POLYGON ((532104.9 182011.9,532126.3 181948,5...
MULTIPOLYGON (((532022.4 181893.5,532022.4 181...
我们得到两个输出特征，因为病房在两个选区都有一些。多少？让我们看看：
> st_area(st_intersection(wards$geom[1],constituencies$geom))
Units: [m^2]
[1] 1.298829e+05 9.542473e-01
所以，一个是 129882m^2，另一个是 0.95m^2。您的边界未对齐，因此“错误”选区中有 1 平方米的病房！
解决方法是根据病房面积的值或比例计算完整交叉点集的面积和阈值。
您可以使用 st_intersects 获取所有重叠的列表：
> st_intersects(wards,constituencies)
Sparse geometry binary predicate list of length 25,where the predicate was `intersects'
first 10 elements:
1: 1,2
2: 1
3: 1
4: 1
5: 1,2
6: 1
该列表中具有一个元素的任何元素都是与一个且仅一个选区重叠的选区，否则您需要计算该选区与这些选区重叠的面积比例，以确定该选区主要覆盖哪些选区。
,由于我们只处理细微的差异，一种方法是在加入之前在病房几何形状周围创建一个轻微的负/插入缓冲区。
准备
library(sf)
library(tmap)
wards_transformed <- st_transform(wards,st_crs(constituencies))
无缓冲：
wards_transformed %>%
st_join(constituencies,st_within)
带负缓冲：
wards_transformed %>%
st_buffer(-10) %>%
st_join(constituencies,st_within)
