如何在sf对象中使用"anti-join" - 子集找出不在多边形内的点

4

我有两个空间对象,一个是点nc_point(n=50),另一个是多边形sample_polygons(n=20)。其中一些点位于多边形中,可以使用nc_point[sample_polygons , ]函数进行识别。

我的问题是如何获取所有未落在任何多边形中的点的子集,即其余的点?

# points layer contain 50 points
nc_point <- st_read(system.file("shape/nc.shp", package="sf")) %>% 
  st_centroid()

# the number of polygons is 20
set.seed(1)
sample_polygons <- st_read(system.file("shape/nc.shp", package="sf")) %>% 
  sample_n(20) %>% 
  select(geometry) # to mimic a situation where points are only identified using spatial correlation. 



# points that fall in polygons can be identified using: 
points_in <- nc_point[sample_polygons , ]

# how to find out points that are not fallen in any polygons? 

感谢您,Phil。

1
没有一个可重现的例子,很难给出具体细节。但是,sf::st_join 的文档向您展示了用于连接 sf 对象的默认函数是 st_intersects,但您可以提供其他函数来进行连接。文档提供了一些建议列表,供您使用其他函数。 - camille
谢谢 @camille,我已经包含了一个可重现的例子。 - Phil Wu
2个回答

3
您可以返回一个TRUE/FALSE输出来选择您的要点。例如,在列表中查找长度为0的元素:
outside <- sapply(st_intersects(nc_point, sample_polygons),function(x){length(x)==0})

那会给你一个逻辑向量,你可以用它来子集:
points_out <- nc_point[outside, ]
points_out

Simple feature collection with 80 features and 14 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: -84.05976 ymin: 34.07663 xmax: -75.80982 ymax: 36.49101
epsg (SRID):    4267
proj4string:    +proj=longlat +datum=NAD27 +no_defs
First 10 features:
    AREA PERIMETER CNTY_ CNTY_ID        NAME  FIPS FIPSNO CRESS_ID BIR74 SID74 NWBIR74 BIR79 SID79 NWBIR79                   geometry
1  0.114     1.442  1825    1825        Ashe 37009  37009        5  1091     1      10  1364     0      19  POINT (-81.49826 36.4314)
2  0.061     1.231  1827    1827   Alleghany 37005  37005        3   487     0      10   542     3      12 POINT (-81.12515 36.49101)
3  0.143     1.630  1828    1828       Surry 37171  37171       86  3188     5     208  3616     6     260 POINT (-80.68575 36.41252)
4  0.070     2.968  1831    1831   Currituck 37053  37053       27   508     1     123   830     2     145  POINT (-76.0275 36.40728)
5  0.153     2.206  1832    1832 Northampton 37131  37131       66  1421     9    1066  1606     3    1197 POINT (-77.41056 36.42228)
7  0.062     1.547  1834    1834      Camden 37029  37029       15   286     0     115   350     2     139  POINT (-76.23435 36.4012)
8  0.091     1.284  1835    1835       Gates 37073  37073       37   420     0     254   594     2     371 POINT (-76.70448 36.44423)
9  0.118     1.421  1836    1836      Warren 37185  37185       93   968     4     748  1190     2     844 POINT (-78.11043 36.39697)
10 0.124     1.428  1837    1837      Stokes 37169  37169       85  1612     1     160  2038     5     176 POINT (-80.23428 36.40034)
11 0.114     1.352  1838    1838     Caswell 37033  37033       17  1035     2     550  1253     2     597 POINT (-79.33477 36.39347)

2

类似于Carlos的答案,但更简单一些,利用函数lengths来获取相交列表中每个元素的长度:

outside <- nc_point[!lengths(st_intersects(nc_point, sample_polygons)), ]
# same as `outside <- nc_point[lengths(st_intersects(nc_point, sample_polygons)) == 0, ]`
outside

Simple feature collection with 80 features and 14 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: -84.05976 ymin: 34.07663 xmax: -75.80982 ymax: 36.49101
epsg (SRID):    4267
proj4string:    +proj=longlat +datum=NAD27 +no_defs
First 10 features:
    AREA PERIMETER CNTY_ CNTY_ID        NAME  FIPS FIPSNO CRESS_ID BIR74 SID74 NWBIR74 BIR79 SID79 NWBIR79
1  0.114     1.442  1825    1825        Ashe 37009  37009        5  1091     1      10  1364     0      19
2  0.061     1.231  1827    1827   Alleghany 37005  37005        3   487     0      10   542     3      12
3  0.143     1.630  1828    1828       Surry 37171  37171       86  3188     5     208  3616     6     260
4  0.070     2.968  1831    1831   Currituck 37053  37053       27   508     1     123   830     2     145
5  0.153     2.206  1832    1832 Northampton 37131  37131       66  1421     9    1066  1606     3    1197
7  0.062     1.547  1834    1834      Camden 37029  37029       15   286     0     115   350     2     139
8  0.091     1.284  1835    1835       Gates 37073  37073       37   420     0     254   594     2     371
9  0.118     1.421  1836    1836      Warren 37185  37185       93   968     4     748  1190     2     844
10 0.124     1.428  1837    1837      Stokes 37169  37169       85  1612     1     160  2038     5     176
11 0.114     1.352  1838    1838     Caswell 37033  37033       17  1035     2     550  1253     2     597
                     geometry
1   POINT (-81.49826 36.4314)
2  POINT (-81.12515 36.49101)
3  POINT (-80.68575 36.41252)
4   POINT (-76.0275 36.40728)
5  POINT (-77.41056 36.42228)
7   POINT (-76.23435 36.4012)
8  POINT (-76.70448 36.44423)
9  POINT (-78.11043 36.39697)
10 POINT (-80.23428 36.40034)
11 POINT (-79.33477 36.39347)


使用lengths()非常整洁!效果很好。 - Phil Wu

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接