使用 sf,选择包含至少一个选择的点特征的区域特征?
With sf, selecting area features that contain at least one of a selection of point features?
例如,假设我有一个 sf
对象,其中包含美国本土的 4 个城市及其坐标。然后我有一个具有 48 个特征的 sf
对象(每个可能的状态一个)。有没有办法 select 包含指定城市的州子集?类似于:
cities_sf
state_sf %>%
filter(states s.t. there exists x in cities_sf s.t. x in states_sf) +
ggplot() +
...
编辑: st_within(my_cities, my_states)
给了我
structure(list(290L, 378L, 51L, integer(0), 283L, 478L, 415L,
380L, 489L, 64L, 189L, 184L, 311L, 488L, 66L, 73L, 49L, 1L,
359L, 111L, 502L, 489L, 272L, 115L, 352L, 241L), predicate = "within",
region.id = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11",
"12", "13", "14", "15", "16", "17", "18", "19", "20", "21", "22", "23",
"24", "25", "26"), ncol = 544L, class = "sgbp")
我可以看出这 26 个索引对应于 my_states
中包含城市的多面体,但我不确定如何使用该 SGBD("sparse geometry binary predicate",根据文档) ggplot
/geom_sf
条款
中的对象
编辑 2:我最终使用了 slice(states_sf, unlist(st_within(cities_sf, states_sf)))
,它给出了一个 sf
对象,它是我需要的子集
没有样本数据,这是我能想到的最好的。
library( sf )
#find intersecting points/polygons
intersect <- st_intersection(x = polygons, y = points)
#and go further from there
更新
使用@Spacedman 在他的回答中提供的样本数据。
library(dplyr)
library(sf)
states %>%
#create ID's for the states (if they don't have one already)
#state ID should be equal to rownumber (fot the filter later on)
mutate( id = row_number() ) %>%
#filter out states that do not have any intersetcions with the points/cities
filter( id %in% unlist( st_intersects(cities, states) ) ) %>%
#plot
mapview::mapview()
使用 USAboundaries
包中的 us_states 函数,让我们制作一小组状态:
> states <- us_states(map_date = "2000-01-01", resolution = "high", states = c("CA", "OR", "WA","NV","NM","UT","CO","ID","AZ"))
这是我创建的一些点:
> pts
Simple feature collection with 4 features and 0 fields
geometry type: POINT
dimension: XY
bbox: xmin: -121.7663 ymin: 34.86508 xmax: -110.7263 ymax: 46.65593
epsg (SRID): 4326
proj4string: +proj=longlat +datum=WGS84 +no_defs
geometry
1 POINT (-110.7263 34.86508)
2 POINT (-111.7345 38.64123)
3 POINT (-120.1531 46.65593)
4 POINT (-121.7663 39.37335)
要测试交集:
> st_intersects(states, pts)
although coordinates are longitude/latitude, st_intersects assumes that they are planar
Sparse geometry binary predicate list of length 9, where the predicate was `intersects'
1: 1
2: 4
3: (empty)
4: (empty)
5: (empty)
6: (empty)
7: (empty)
8: 2
9: 3
该对象是一个列表,因此您可以获得元素的长度并找到大于零的元素 - 即其中包含一些东西:
> lengths(st_intersects(states, pts)) > 0
although coordinates are longitude/latitude, st_intersects assumes that they are planar
[1] TRUE TRUE FALSE FALSE FALSE FALSE FALSE TRUE TRUE
然后以常规方式对空间多边形进行子集化:
> plot(st_geometry(states[lengths(st_intersects(states, pts)) > 0,]))
绘制了具有四个点的四个状态。
创建子集并将其提供给 ggplot
如果这就是您绘制地图的方式。
例如,假设我有一个 sf
对象,其中包含美国本土的 4 个城市及其坐标。然后我有一个具有 48 个特征的 sf
对象(每个可能的状态一个)。有没有办法 select 包含指定城市的州子集?类似于:
cities_sf
state_sf %>%
filter(states s.t. there exists x in cities_sf s.t. x in states_sf) +
ggplot() +
...
编辑: st_within(my_cities, my_states)
给了我
structure(list(290L, 378L, 51L, integer(0), 283L, 478L, 415L,
380L, 489L, 64L, 189L, 184L, 311L, 488L, 66L, 73L, 49L, 1L,
359L, 111L, 502L, 489L, 272L, 115L, 352L, 241L), predicate = "within",
region.id = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11",
"12", "13", "14", "15", "16", "17", "18", "19", "20", "21", "22", "23",
"24", "25", "26"), ncol = 544L, class = "sgbp")
我可以看出这 26 个索引对应于 my_states
中包含城市的多面体,但我不确定如何使用该 SGBD("sparse geometry binary predicate",根据文档) ggplot
/geom_sf
条款
编辑 2:我最终使用了 slice(states_sf, unlist(st_within(cities_sf, states_sf)))
,它给出了一个 sf
对象,它是我需要的子集
没有样本数据,这是我能想到的最好的。
library( sf )
#find intersecting points/polygons
intersect <- st_intersection(x = polygons, y = points)
#and go further from there
更新
使用@Spacedman 在他的回答中提供的样本数据。
library(dplyr)
library(sf)
states %>%
#create ID's for the states (if they don't have one already)
#state ID should be equal to rownumber (fot the filter later on)
mutate( id = row_number() ) %>%
#filter out states that do not have any intersetcions with the points/cities
filter( id %in% unlist( st_intersects(cities, states) ) ) %>%
#plot
mapview::mapview()
使用 USAboundaries
包中的 us_states 函数,让我们制作一小组状态:
> states <- us_states(map_date = "2000-01-01", resolution = "high", states = c("CA", "OR", "WA","NV","NM","UT","CO","ID","AZ"))
这是我创建的一些点:
> pts
Simple feature collection with 4 features and 0 fields
geometry type: POINT
dimension: XY
bbox: xmin: -121.7663 ymin: 34.86508 xmax: -110.7263 ymax: 46.65593
epsg (SRID): 4326
proj4string: +proj=longlat +datum=WGS84 +no_defs
geometry
1 POINT (-110.7263 34.86508)
2 POINT (-111.7345 38.64123)
3 POINT (-120.1531 46.65593)
4 POINT (-121.7663 39.37335)
要测试交集:
> st_intersects(states, pts)
although coordinates are longitude/latitude, st_intersects assumes that they are planar
Sparse geometry binary predicate list of length 9, where the predicate was `intersects'
1: 1
2: 4
3: (empty)
4: (empty)
5: (empty)
6: (empty)
7: (empty)
8: 2
9: 3
该对象是一个列表,因此您可以获得元素的长度并找到大于零的元素 - 即其中包含一些东西:
> lengths(st_intersects(states, pts)) > 0
although coordinates are longitude/latitude, st_intersects assumes that they are planar
[1] TRUE TRUE FALSE FALSE FALSE FALSE FALSE TRUE TRUE
然后以常规方式对空间多边形进行子集化:
> plot(st_geometry(states[lengths(st_intersects(states, pts)) > 0,]))
绘制了具有四个点的四个状态。
创建子集并将其提供给 ggplot
如果这就是您绘制地图的方式。