右 | sf：我有点和每个点周围的 2 个缓冲区。如果较大的缓冲区重叠（但较小的缓冲区不重叠），如何将这些点组合成单个单元？
R | sf: I have points and 2 buffers around each point. How to combine the points into single units if larger buffer overlaps (but smaller does not)?
我目前在地图（学校）上有一些点。每个点（学校）周围都有两个缓冲区。一种是 450m，另一种是 250m。如果点重叠，我想将它们视为一个单元（因为否则事情会变得复杂），但我希望它们保持它们覆盖的 geometries/area。
所以在此处给出的示例地图上，我希望前三个 schools/points 合并为一个“单元”。我希望他们保留他们覆盖的区域，但只将 R 算作一个单元。如果我使用“st_union”函数，我必须对较小和较大的缓冲区都执行此操作。但这会导致单元总数的差异，因为更多的较大缓冲区重叠。所以基本上，如果较大的缓冲区重叠并因此合并，我希望较小的缓冲区与其他较小的缓冲区合并（成为一个单元）。
最终目标是我想计算一家商店在距学校 250 米以内与 450 米以内的次数。我认为最简单的方法是如上所述合并，因为学校之间存在大量重叠。
我的代码如下，但请注意，由于上述解释的原因，它会导致 schools/units 的差异数。
county.sf <- get_acs(state = "MO",
county = c( "St. Louis City"),
geography = "tract",
variables = "B03002_001",
output="wide",
geometry = TRUE) %>%
sf::st_transform(crs = "ESRI:102003")
class(county.sf)
# School data
school <- read.csv("C:\myfile1.csv")
school.sf <- st_as_sf(school, coords = c("long", "lat"), crs = "epsg:4326")
school.sf.utm <- st_transform(school.sf, crs = "ESRI:102003")
# Store data
store <- import("C:\myfile2.csv")
store.sf <- st_as_sf(store, coords = c("XCoord", "YCoord"), crs = "ESRI:102696")
store.sf.utm <- st_transform(store.sf, crs = "ESRI:102003")
elem.buff <-st_buffer(school.sf.utm, 250)
elem.buff2 <-st_buffer(school.sf.utm, 450)
pts_com<-st_union(elem.buff)
pts_pol<-st_cast(pts_com, "POLYGON")
pts_com2<-st_union(elem.buff2)
pts_pol2<-st_cast(pts_com2, "POLYGON")
#unmerged zone map
ex.map<- tm_shape(county.sf) +
tm_polygons() +
tm_shape(elem.buff) +
tm_borders(col="red") +
tm_shape(school.sf.utm) +
tm_dots(col = "red") +
tm_shape(elem.buff2) +
tm_borders(col="blue") +
tm_shape(pts_pol) +
tm_borders(col="black") +
tm_shape(store.sf.utm) +
tm_dots()
ex.map
#merged zones map
ex.map<- tm_shape(county.sf) +
tm_polygons() +
#(elem.buff) +
#tm_borders(col="red") +
tm_shape(school.sf.utm) +
tm_dots(col = "red") +
#tm_shape(elem.buff2) +
#tm_borders(col="blue") +
tm_shape(pts_pol) +
tm_borders(col="red") +
tm_shape(store.sf.utm) +
tm_dots() +
tm_shape(pts_pol2) +
tm_borders(col="blue")
ex.map
lengths(st_intersects(elem.buff, store.sf.utm)) #gives per school but ignores overlapping
lengths(st_intersects(pts_pol, store.sf.utm))
lengths(st_intersects(elem.buff2, store.sf.utm)))
lengths(st_intersects(pts_pol2, store.sf.utm)))
最初我建议使用 st_intersection()，但开始变得复杂，所以我改用了这种取决于 st_join() 的方法。
首先，加载库和数据。
library(readxl)
library(tidyverse)
library(sf)
library(tmap)
tmap_mode('view')
read_xlsx('Schools and Stores_all.xlsx', sheet = 1) %>%
st_as_sf(., coords = c("long", "lat"), crs = "epsg:4326") %>%
st_transform(crs = "ESRI:102003") %>%
{. ->> schools}
schools
# Simple feature collection with 156 features and 1 field
# Geometry type: POINT
# Dimension: XY
# Bounding box: xmin: 208163.7 ymin: -81868.36 xmax: 565039.3 ymax: 151922.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 2
# School geometry
# * <chr> <POINT [m]>
# 1 AcademyOf Envt Sci/math Elementary School (497569.1 143116)
# 2 AcademyOf Envt Sci/math Middle School (498610.1 140067.7)
# 3 Adams Elementary School (495000.6 141023.2)
# 4 Ames Visual/perf. Arts (500120.2 144483.3)
# 5 Ashland Elementary And Br. (496514.9 146392.7)
# 6 Aspire Academy (495458 149014.6)
# 7 Beaumont Cte High School (497748.6 145417.7)
# 8 Bryan Hill Elementary School (495926.2 144709.6)
# 9 Buder Elementary School (492994.6 136888.8)
# 10 Busch Ms Character Athletics (491906.9 135569.1)
# # ... with 146 more rows
然后，在学校周围做缓冲区。
schools %>%
st_buffer(250) %>%
{. ->> schools_250}
schools_250
# Simple feature collection with 156 features and 1 field
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 207913.7 ymin: -82118.36 xmax: 565289.3 ymax: 152172.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 2
# School geometry
# * <chr> <POLYGON [m]>
# 1 AcademyOf Envt Sci/math Elementary School ((497819.1 143116, 497818.8 143103, 497817.7 143089.9, 497816 ~
# 2 AcademyOf Envt Sci/math Middle School ((498860.1 140067.7, 498859.7 140054.7, 498858.7 140041.6, 498~
# 3 Adams Elementary School ((495250.6 141023.2, 495250.3 141010.1, 495249.3 140997, 49524~
# 4 Ames Visual/perf. Arts ((500370.2 144483.3, 500369.9 144470.2, 500368.9 144457.2, 500~
# 5 Ashland Elementary And Br. ((496764.9 146392.7, 496764.6 146379.6, 496763.6 146366.5, 496~
# 6 Aspire Academy ((495708 149014.6, 495707.6 149001.5, 495706.6 148988.4, 49570~
# 7 Beaumont Cte High School ((497998.6 145417.7, 497998.3 145404.7, 497997.3 145391.6, 497~
# 8 Bryan Hill Elementary School ((496176.2 144709.6, 496175.9 144696.5, 496174.9 144683.5, 496~
# 9 Buder Elementary School ((493244.6 136888.8, 493244.2 136875.7, 493243.2 136862.7, 493~
# 10 Busch Ms Character Athletics ((492156.9 135569.1, 492156.6 135556, 492155.5 135543, 492153.~
# # ... with 146 more rows
schools %>%
st_buffer(450) %>%
{. ->> schools_450}
schools_450
# Simple feature collection with 156 features and 1 field
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 207713.7 ymin: -82318.36 xmax: 565489.3 ymax: 152372.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 2
# School geometry
# * <chr> <POLYGON [m]>
# 1 AcademyOf Envt Sci/math Elementary School ((498019.1 143116, 498018.5 143092.5, 498016.7 143069, 498013.~
# 2 AcademyOf Envt Sci/math Middle School ((499060.1 140067.7, 499059.5 140044.2, 499057.6 140020.7, 499~
# 3 Adams Elementary School ((495450.6 141023.2, 495450 140999.6, 495448.2 140976.1, 49544~
# 4 Ames Visual/perf. Arts ((500570.2 144483.3, 500569.6 144459.8, 500567.8 144436.3, 500~
# 5 Ashland Elementary And Br. ((496964.9 146392.7, 496964.3 146369.1, 496962.5 146345.6, 496~
# 6 Aspire Academy ((495908 149014.6, 495907.4 148991, 495905.5 148967.5, 495902.~
# 7 Beaumont Cte High School ((498198.6 145417.7, 498198 145394.2, 498196.2 145370.7, 49819~
# 8 Bryan Hill Elementary School ((496376.2 144709.6, 496375.6 144686.1, 496373.8 144662.6, 496~
# 9 Buder Elementary School ((493444.6 136888.8, 493444 136865.2, 493442.1 136841.8, 49343~
# 10 Busch Ms Character Athletics ((492356.9 135569.1, 492356.3 135545.5, 492354.4 135522, 49235~
# # ... with 146 more rows
使用st_union()加入所有450米的缓冲区，以确定“单位”的身份。
schools_450 %>%
st_union %>%
st_cast('POLYGON') %>%
st_sf %>%
mutate(
unit = row_number()
) %>%
{. ->> schools_450_units}
schools_450_units
# Simple feature collection with 39 features and 1 field
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 207713.7 ymin: -82318.36 xmax: 565489.3 ymax: 152372.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# First 10 features:
# geometry unit
# 1 POLYGON ((565450.4 -82051.3... 1
# 2 POLYGON ((499821.5 138524.2... 2
# 3 POLYGON ((498299.1 138974.4... 3
# 4 POLYGON ((500735.6 142090.4... 4
# 5 POLYGON ((499872.4 142927.1... 5
# 6 POLYGON ((496870.2 135641.5... 6
# 7 POLYGON ((495561.7 135210, ... 7
# 8 POLYGON ((498291.9 141786.4... 8
# 9 POLYGON ((496337.3 144526.6... 9
# 10 POLYGON ((208574.8 -53382.3... 10
请注意，在 156 所学校中，我们只有 39 个单元。当然，当多个 450 米缓冲区链接在一起时，这些单元可以达到相当远的范围 - 请参阅最后的地图。
然后，我们使用st_join()来确定每个单元内有哪些学校（或它们的缓冲区）。指定 join = st_within 是这里的关键。
schools %>%
st_join(schools_450_units, join = st_within) %>%
{. ->> schools_units}
schools_units
# Simple feature collection with 156 features and 2 fields
# Geometry type: POINT
# Dimension: XY
# Bounding box: xmin: 208163.7 ymin: -81868.36 xmax: 565039.3 ymax: 151922.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 3
# School geometry unit
# * <chr> <POINT [m]> <int>
# 1 AcademyOf Envt Sci/math Elementary School (497569.1 143116) 5
# 2 AcademyOf Envt Sci/math Middle School (498610.1 140067.7) 3
# 3 Adams Elementary School (495000.6 141023.2) 3
# 4 Ames Visual/perf. Arts (500120.2 144483.3) 5
# 5 Ashland Elementary And Br. (496514.9 146392.7) 32
# 6 Aspire Academy (495458 149014.6) 37
# 7 Beaumont Cte High School (497748.6 145417.7) 33
# 8 Bryan Hill Elementary School (495926.2 144709.6) 9
# 9 Buder Elementary School (492994.6 136888.8) 18
# 10 Busch Ms Character Athletics (491906.9 135569.1) 13
# # ... with 146 more rows
然后将同一单元的 250 m 缓冲区（即使它们不接触）合并为一个 MULTIPOLYGON，我们使用 group_by(unit) 并使用 summarise(do_union = TRUE)（默认）。这类似于 st_union()，但尊重 group_by()（实际上可能使用 st_union() 获得完全相同的结果，但这也有效）。
schools_units %>%
`st_geometry<-`(schools_250 %>% st_geometry) %>%
group_by(unit) %>%
summarise(do_union = TRUE) %>%
{. ->> schools_250_units}
schools_250_units
# Simple feature collection with 39 features and 1 field
# Geometry type: GEOMETRY
# Dimension: XY
# Bounding box: xmin: 207913.7 ymin: -82118.36 xmax: 565289.3 ymax: 152172.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 39 x 2
# unit geometry
# <int> <GEOMETRY [m]>
# 1 1 POLYGON ((565289.3 -81868.36, 565289 -81881.45, 565288 -81894.49, 565286.2 -81907.47, 565283.9 -8...
# 2 2 MULTIPOLYGON (((499659 138681.1, 499657.3 138668.1, 499654.9 138655.3, 499651.9 138642.5, 499648....
# 3 3 MULTIPOLYGON (((496446.7 137764.8, 496442.9 137752.3, 496438.6 137739.9, 496433.6 137727.9, 49642...
# 4 4 MULTIPOLYGON (((500574.2 142260.4, 500573.1 142247.3, 500571.4 142234.4, 500569 142221.5, 500566 ...
# 5 5 MULTIPOLYGON (((497408.5 142703.5, 497404.7 142690.9, 497400.4 142678.6, 497395.4 142666.5, 49738...
# 6 6 MULTIPOLYGON (((496993.6 135888.6, 496982.8 135881.2, 496971.6 135874.4, 496960.1 135868.1, 49694...
# 7 7 MULTIPOLYGON (((496252.1 134426.9, 496250.4 134413.9, 496248 134401.1, 496244.9 134388.3, 496241....
# 8 8 POLYGON ((498130.8 141969.4, 498130.4 141956.3, 498129.4 141943.3, 498127.7 141930.3, 498125.3 14...
# 9 9 MULTIPOLYGON (((496175.9 144696.5, 496174.9 144683.5, 496173.1 144670.5, 496170.8 144657.6, 49616...
# 10 10 POLYGON ((208413.7 -53199.34, 208413.3 -53212.42, 208412.3 -53225.47, 208410.6 -53238.45, 208408....
# # ... with 29 more rows
希望这已经回答了你的问题，因为现在我们有：
- 每个单位属于哪个学校
schools_units
- 联合 250 米缓冲区
schools_250_units
- 联合 450 米缓冲区
schools_450_units
就查找靠近 units/schools 的商店而言，例如在您的 中，您可能会考虑这些单位的分布范围 - 当 450 米的缓冲区链接在一起时。以下图为例。
tm_shape(schools_450_units, bbox = schools_450_units %>% filter(unit == 19) %>% st_buffer(2000) %>% st_bbox)+
tm_polygons(col = 'unit', border.col = 'blue', legend.show = FALSE)+
tm_shape(schools_250_units)+
tm_polygons(col = 'unit', border.col = 'red', legend.show = FALSE)+
tm_shape(schools_units)+
tm_text(text = 'unit')
看看第 3 单元 - 它非常大，几乎吞没了第 22 单元。由你决定你想做什么以及这是否合适，但这是我基于快速地图的最初想法之一。
我目前在地图（学校）上有一些点。每个点（学校）周围都有两个缓冲区。一种是 450m，另一种是 250m。如果点重叠，我想将它们视为一个单元（因为否则事情会变得复杂），但我希望它们保持它们覆盖的 geometries/area。
所以在此处给出的示例地图上，我希望前三个 schools/points 合并为一个“单元”。我希望他们保留他们覆盖的区域，但只将 R 算作一个单元。如果我使用“st_union”函数，我必须对较小和较大的缓冲区都执行此操作。但这会导致单元总数的差异，因为更多的较大缓冲区重叠。所以基本上，如果较大的缓冲区重叠并因此合并，我希望较小的缓冲区与其他较小的缓冲区合并（成为一个单元）。
最终目标是我想计算一家商店在距学校 250 米以内与 450 米以内的次数。我认为最简单的方法是如上所述合并，因为学校之间存在大量重叠。
我的代码如下，但请注意，由于上述解释的原因，它会导致 schools/units 的差异数。
county.sf <- get_acs(state = "MO",
county = c( "St. Louis City"),
geography = "tract",
variables = "B03002_001",
output="wide",
geometry = TRUE) %>%
sf::st_transform(crs = "ESRI:102003")
class(county.sf)
# School data
school <- read.csv("C:\myfile1.csv")
school.sf <- st_as_sf(school, coords = c("long", "lat"), crs = "epsg:4326")
school.sf.utm <- st_transform(school.sf, crs = "ESRI:102003")
# Store data
store <- import("C:\myfile2.csv")
store.sf <- st_as_sf(store, coords = c("XCoord", "YCoord"), crs = "ESRI:102696")
store.sf.utm <- st_transform(store.sf, crs = "ESRI:102003")
elem.buff <-st_buffer(school.sf.utm, 250)
elem.buff2 <-st_buffer(school.sf.utm, 450)
pts_com<-st_union(elem.buff)
pts_pol<-st_cast(pts_com, "POLYGON")
pts_com2<-st_union(elem.buff2)
pts_pol2<-st_cast(pts_com2, "POLYGON")
#unmerged zone map
ex.map<- tm_shape(county.sf) +
tm_polygons() +
tm_shape(elem.buff) +
tm_borders(col="red") +
tm_shape(school.sf.utm) +
tm_dots(col = "red") +
tm_shape(elem.buff2) +
tm_borders(col="blue") +
tm_shape(pts_pol) +
tm_borders(col="black") +
tm_shape(store.sf.utm) +
tm_dots()
ex.map
#merged zones map
ex.map<- tm_shape(county.sf) +
tm_polygons() +
#(elem.buff) +
#tm_borders(col="red") +
tm_shape(school.sf.utm) +
tm_dots(col = "red") +
#tm_shape(elem.buff2) +
#tm_borders(col="blue") +
tm_shape(pts_pol) +
tm_borders(col="red") +
tm_shape(store.sf.utm) +
tm_dots() +
tm_shape(pts_pol2) +
tm_borders(col="blue")
ex.map
lengths(st_intersects(elem.buff, store.sf.utm)) #gives per school but ignores overlapping
lengths(st_intersects(pts_pol, store.sf.utm))
lengths(st_intersects(elem.buff2, store.sf.utm)))
lengths(st_intersects(pts_pol2, store.sf.utm)))
最初我建议使用 st_intersection()，但开始变得复杂，所以我改用了这种取决于 st_join() 的方法。
首先，加载库和数据。
library(readxl)
library(tidyverse)
library(sf)
library(tmap)
tmap_mode('view')
read_xlsx('Schools and Stores_all.xlsx', sheet = 1) %>%
st_as_sf(., coords = c("long", "lat"), crs = "epsg:4326") %>%
st_transform(crs = "ESRI:102003") %>%
{. ->> schools}
schools
# Simple feature collection with 156 features and 1 field
# Geometry type: POINT
# Dimension: XY
# Bounding box: xmin: 208163.7 ymin: -81868.36 xmax: 565039.3 ymax: 151922.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 2
# School geometry
# * <chr> <POINT [m]>
# 1 AcademyOf Envt Sci/math Elementary School (497569.1 143116)
# 2 AcademyOf Envt Sci/math Middle School (498610.1 140067.7)
# 3 Adams Elementary School (495000.6 141023.2)
# 4 Ames Visual/perf. Arts (500120.2 144483.3)
# 5 Ashland Elementary And Br. (496514.9 146392.7)
# 6 Aspire Academy (495458 149014.6)
# 7 Beaumont Cte High School (497748.6 145417.7)
# 8 Bryan Hill Elementary School (495926.2 144709.6)
# 9 Buder Elementary School (492994.6 136888.8)
# 10 Busch Ms Character Athletics (491906.9 135569.1)
# # ... with 146 more rows
然后，在学校周围做缓冲区。
schools %>%
st_buffer(250) %>%
{. ->> schools_250}
schools_250
# Simple feature collection with 156 features and 1 field
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 207913.7 ymin: -82118.36 xmax: 565289.3 ymax: 152172.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 2
# School geometry
# * <chr> <POLYGON [m]>
# 1 AcademyOf Envt Sci/math Elementary School ((497819.1 143116, 497818.8 143103, 497817.7 143089.9, 497816 ~
# 2 AcademyOf Envt Sci/math Middle School ((498860.1 140067.7, 498859.7 140054.7, 498858.7 140041.6, 498~
# 3 Adams Elementary School ((495250.6 141023.2, 495250.3 141010.1, 495249.3 140997, 49524~
# 4 Ames Visual/perf. Arts ((500370.2 144483.3, 500369.9 144470.2, 500368.9 144457.2, 500~
# 5 Ashland Elementary And Br. ((496764.9 146392.7, 496764.6 146379.6, 496763.6 146366.5, 496~
# 6 Aspire Academy ((495708 149014.6, 495707.6 149001.5, 495706.6 148988.4, 49570~
# 7 Beaumont Cte High School ((497998.6 145417.7, 497998.3 145404.7, 497997.3 145391.6, 497~
# 8 Bryan Hill Elementary School ((496176.2 144709.6, 496175.9 144696.5, 496174.9 144683.5, 496~
# 9 Buder Elementary School ((493244.6 136888.8, 493244.2 136875.7, 493243.2 136862.7, 493~
# 10 Busch Ms Character Athletics ((492156.9 135569.1, 492156.6 135556, 492155.5 135543, 492153.~
# # ... with 146 more rows
schools %>%
st_buffer(450) %>%
{. ->> schools_450}
schools_450
# Simple feature collection with 156 features and 1 field
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 207713.7 ymin: -82318.36 xmax: 565489.3 ymax: 152372.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 2
# School geometry
# * <chr> <POLYGON [m]>
# 1 AcademyOf Envt Sci/math Elementary School ((498019.1 143116, 498018.5 143092.5, 498016.7 143069, 498013.~
# 2 AcademyOf Envt Sci/math Middle School ((499060.1 140067.7, 499059.5 140044.2, 499057.6 140020.7, 499~
# 3 Adams Elementary School ((495450.6 141023.2, 495450 140999.6, 495448.2 140976.1, 49544~
# 4 Ames Visual/perf. Arts ((500570.2 144483.3, 500569.6 144459.8, 500567.8 144436.3, 500~
# 5 Ashland Elementary And Br. ((496964.9 146392.7, 496964.3 146369.1, 496962.5 146345.6, 496~
# 6 Aspire Academy ((495908 149014.6, 495907.4 148991, 495905.5 148967.5, 495902.~
# 7 Beaumont Cte High School ((498198.6 145417.7, 498198 145394.2, 498196.2 145370.7, 49819~
# 8 Bryan Hill Elementary School ((496376.2 144709.6, 496375.6 144686.1, 496373.8 144662.6, 496~
# 9 Buder Elementary School ((493444.6 136888.8, 493444 136865.2, 493442.1 136841.8, 49343~
# 10 Busch Ms Character Athletics ((492356.9 135569.1, 492356.3 135545.5, 492354.4 135522, 49235~
# # ... with 146 more rows
使用st_union()加入所有450米的缓冲区，以确定“单位”的身份。
schools_450 %>%
st_union %>%
st_cast('POLYGON') %>%
st_sf %>%
mutate(
unit = row_number()
) %>%
{. ->> schools_450_units}
schools_450_units
# Simple feature collection with 39 features and 1 field
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 207713.7 ymin: -82318.36 xmax: 565489.3 ymax: 152372.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# First 10 features:
# geometry unit
# 1 POLYGON ((565450.4 -82051.3... 1
# 2 POLYGON ((499821.5 138524.2... 2
# 3 POLYGON ((498299.1 138974.4... 3
# 4 POLYGON ((500735.6 142090.4... 4
# 5 POLYGON ((499872.4 142927.1... 5
# 6 POLYGON ((496870.2 135641.5... 6
# 7 POLYGON ((495561.7 135210, ... 7
# 8 POLYGON ((498291.9 141786.4... 8
# 9 POLYGON ((496337.3 144526.6... 9
# 10 POLYGON ((208574.8 -53382.3... 10
请注意，在 156 所学校中，我们只有 39 个单元。当然，当多个 450 米缓冲区链接在一起时，这些单元可以达到相当远的范围 - 请参阅最后的地图。
然后，我们使用st_join()来确定每个单元内有哪些学校（或它们的缓冲区）。指定 join = st_within 是这里的关键。
schools %>%
st_join(schools_450_units, join = st_within) %>%
{. ->> schools_units}
schools_units
# Simple feature collection with 156 features and 2 fields
# Geometry type: POINT
# Dimension: XY
# Bounding box: xmin: 208163.7 ymin: -81868.36 xmax: 565039.3 ymax: 151922.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 156 x 3
# School geometry unit
# * <chr> <POINT [m]> <int>
# 1 AcademyOf Envt Sci/math Elementary School (497569.1 143116) 5
# 2 AcademyOf Envt Sci/math Middle School (498610.1 140067.7) 3
# 3 Adams Elementary School (495000.6 141023.2) 3
# 4 Ames Visual/perf. Arts (500120.2 144483.3) 5
# 5 Ashland Elementary And Br. (496514.9 146392.7) 32
# 6 Aspire Academy (495458 149014.6) 37
# 7 Beaumont Cte High School (497748.6 145417.7) 33
# 8 Bryan Hill Elementary School (495926.2 144709.6) 9
# 9 Buder Elementary School (492994.6 136888.8) 18
# 10 Busch Ms Character Athletics (491906.9 135569.1) 13
# # ... with 146 more rows
然后将同一单元的 250 m 缓冲区（即使它们不接触）合并为一个 MULTIPOLYGON，我们使用 group_by(unit) 并使用 summarise(do_union = TRUE)（默认）。这类似于 st_union()，但尊重 group_by()（实际上可能使用 st_union() 获得完全相同的结果，但这也有效）。
schools_units %>%
`st_geometry<-`(schools_250 %>% st_geometry) %>%
group_by(unit) %>%
summarise(do_union = TRUE) %>%
{. ->> schools_250_units}
schools_250_units
# Simple feature collection with 39 features and 1 field
# Geometry type: GEOMETRY
# Dimension: XY
# Bounding box: xmin: 207913.7 ymin: -82118.36 xmax: 565289.3 ymax: 152172.7
# Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
# # A tibble: 39 x 2
# unit geometry
# <int> <GEOMETRY [m]>
# 1 1 POLYGON ((565289.3 -81868.36, 565289 -81881.45, 565288 -81894.49, 565286.2 -81907.47, 565283.9 -8...
# 2 2 MULTIPOLYGON (((499659 138681.1, 499657.3 138668.1, 499654.9 138655.3, 499651.9 138642.5, 499648....
# 3 3 MULTIPOLYGON (((496446.7 137764.8, 496442.9 137752.3, 496438.6 137739.9, 496433.6 137727.9, 49642...
# 4 4 MULTIPOLYGON (((500574.2 142260.4, 500573.1 142247.3, 500571.4 142234.4, 500569 142221.5, 500566 ...
# 5 5 MULTIPOLYGON (((497408.5 142703.5, 497404.7 142690.9, 497400.4 142678.6, 497395.4 142666.5, 49738...
# 6 6 MULTIPOLYGON (((496993.6 135888.6, 496982.8 135881.2, 496971.6 135874.4, 496960.1 135868.1, 49694...
# 7 7 MULTIPOLYGON (((496252.1 134426.9, 496250.4 134413.9, 496248 134401.1, 496244.9 134388.3, 496241....
# 8 8 POLYGON ((498130.8 141969.4, 498130.4 141956.3, 498129.4 141943.3, 498127.7 141930.3, 498125.3 14...
# 9 9 MULTIPOLYGON (((496175.9 144696.5, 496174.9 144683.5, 496173.1 144670.5, 496170.8 144657.6, 49616...
# 10 10 POLYGON ((208413.7 -53199.34, 208413.3 -53212.42, 208412.3 -53225.47, 208410.6 -53238.45, 208408....
# # ... with 29 more rows
希望这已经回答了你的问题，因为现在我们有：
- 每个单位属于哪个学校
schools_units - 联合 250 米缓冲区
schools_250_units - 联合 450 米缓冲区
schools_450_units
就查找靠近 units/schools 的商店而言，例如在您的
tm_shape(schools_450_units, bbox = schools_450_units %>% filter(unit == 19) %>% st_buffer(2000) %>% st_bbox)+
tm_polygons(col = 'unit', border.col = 'blue', legend.show = FALSE)+
tm_shape(schools_250_units)+
tm_polygons(col = 'unit', border.col = 'red', legend.show = FALSE)+
tm_shape(schools_units)+
tm_text(text = 'unit')
看看第 3 单元 - 它非常大，几乎吞没了第 22 单元。由你决定你想做什么以及这是否合适，但这是我基于快速地图的最初想法之一。