如何从人口普查 Shapefile（邮政编码级别）中删除所有小岛？答案

【问题标题】：how to remove all the small islands from the Census Shapefile (zip code level)?如何从人口普查 Shapefile（邮政编码级别）中删除所有小岛？
【发布时间】：2018-10-26 18:54:25
【问题描述】：

我已经下载了邮政编码级别的人口普查 shapefile，cb_2017_us_zcta510_500k.shp (https://www.census.gov/geo/maps-data/data/cbf/cbf_zcta.html)

我还下载了映射文件，允许我添加相应的STATE 变量（https://www.census.gov/geo/maps-data/data/zcta_rel_download.html）

我将两者合并得到：

library(sf)
library(dplyr)

big_df 

Simple feature collection with 44434 features and 2 fields
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: -176.6847 ymin: -14.37374 xmax: 145.8304 ymax: 71.34122
epsg (SRID):    4269
proj4string:    +proj=longlat +datum=NAD83 +no_defs
First 10 features:
   ZCTA5CE10 STATE                       geometry
1      35442     1 MULTIPOLYGON (((-88.25262 3...
2      35442     1 MULTIPOLYGON (((-88.25262 3...
3      35442     1 MULTIPOLYGON (((-88.25262 3...

现在，我尝试过滤所有小岛和阿拉斯加：

remove_list <-  c("02", "15", "72", "66", "78", "60", "69",
"64", "68", "70", "74", "81", "84", "86", "87", "89", "71", "76",
"95", "79")



big_df %>% filter(!STATE %in% map(remove_list, as.integer)) %>% 
  tm_shape(.) + tm_polygons('pt_count',palette = "Reds", 
                            style = "quantile", n = 10, 
                            title = "counts")

但我仍然有一些小岛。

我在这里缺少什么？谢谢！

【问题讨论】：

您可以使用mapedit::editMap(x = mapview::mapview(st_bbox(big_df))) 在美国大陆周围绘制一个框或多边形，然后使用st_crop 或st_intersection 我想这样做。
非常有趣...您介意发布一个可行的解决方案吗？很乐意接受这一点。谢谢！！
它不会是可重现的，因为它涉及手工绘制特征......因此，我认为它并不适合答案。
也许有办法使用sf 说，现在放下所有被水包围的东西？
你怎么知道什么被水包围，什么不是？那么湖上的岛屿呢？

标签： r shapefile sf census

【解决方案1】：

这是获取美国大陆几何（国家轮廓）的一种方法：

library(raster)
library(sf)
library(dplyr)

us = getData('GADM', country='USA', level=0) %>%
  st_as_sf() %>%
  st_cast("POLYGON") %>%
  mutate(area = st_area(.)) %>%
  arrange(desc(area)) %>%
  slice(1) # mainland US should be the largest

然后您可以使用它来运行st_intersection(big_df, us) 以仅提取big_df 在us 中的部分。请注意，首先在us 周围创建st_buffer 或st_convex_hull 以确保您的big_df 不会在其边界的某处被剪裁，这可能会有所回报。

【讨论】：

太棒了！为了完整起见，您可以尝试使其与示例数据一起使用吗？这些链接是官方人口普查的链接，非常简单......这将非常有用！！！
缺少什么使其与示例数据一起工作？ st_intersection(big_df, us)还不够吗？
这取决于你我的朋友。但是发布输出（最终情节）会很棒！我认为人们正在寻找这样的规范答案。如果那太难了，请告诉我。再次感谢 fkr 的帮助！
我很乐意提供一个完整的解决方案，如果您更新您的问题以使其可重现，因为数据是通过例如自动检索的。 tidycensus 或通过download.file 到tempfile() 或类似的。此外，数据的合并是不可重现的。
是的，你是对的。让我尽快更新。感谢您指出这一点