R|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)?

查看:40
本文介绍了R|SF:我有点,每个点周围有2个缓冲区。如果较大的缓冲区重叠(但较小的缓冲区不重叠),如何将这些点合并为单个单元?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我目前在地图上有点(学校)。每个人在点(学校)周围都有两个缓冲区。一个是450米,一个是250米。如果点重叠,我希望将它们视为单个单元(因为否则事情会变得复杂),但我希望它们保持其覆盖的几何图形/面积。

因此,在这里给出的示例地图上,我希望将排名前三的学校/点合并为一个单位。&我希望它们保留其覆盖的面积,但仅将R作为一个单位计算。如果我使用";st_Union";函数,则必须同时对较小和较大的缓冲区执行此操作。但这会导致单元总数的不同,因为更多较大的缓冲区会重叠。因此,基本上,如果较大的缓冲区重叠从而合并,我希望较小的缓冲区与其他较小的缓冲区合并(合并为一个单元)。

最终目标是我想要计算一下,与450米相比,一家商店落在离学校250米的范围内有多少次。我认为,如上所述的合并是最容易的,因为学校之间有很大的重叠。

我有sample data here.

我的代码如下,但请注意,由于上述原因,它会导致学校/单元数量不同。

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米缓冲区,以确定&Quot;单位&Quot;的身份。

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个单元。当然,当多个450m缓冲器链接在一起时,这些单元可能会影响相当远-请参阅末尾的地图。

然后,我们使用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
然后,要将同一单元的250M缓冲区(即使它们没有接触)合并为一个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
  • 统一的250M缓冲区schools_250_units
  • 统一的450M缓冲区schools_450_units

在查找靠近单元/学校的商店方面,就像在您的other question中一样,您可能会考虑单元的分布范围有多广-当450M缓冲区链接在一起时。以下面的地图为例。

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单元。您可以决定要做什么以及这是否合适,但这只是我基于快速地图的最初想法之一。

这篇关于R|SF:我有点,每个点周围有2个缓冲区。如果较大的缓冲区重叠(但较小的缓冲区不重叠),如何将这些点合并为单个单元?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆