从R中的多边形计算人口普查区域的面积

Question

我正在尝试计算美国每个人口普查区的面积。我使用 tidycensus 为每个人口普查区域及其坐标提取多面体。

library(tidycensus)
library(purrr)

census_api_key("KEY")

us <- unique(fips_codes$state)[1:51]
geometry <- map_df(us, function(x) {
  get_acs(geography = "tract", variables = "B19013_001", geometry = TRUE, 
          state = x)
})

我找到的用于计算多面体面积的方法都不能用于列表，这是 tidycensus 存储坐标的格式。还有其他更容易使用的格式吗？有没有办法使用这种格式的坐标计算面积？

我尝试过 rgeos 包中的 gArea()、raster 包中的 area()、sf 包中的 st_area() 和几何包中的 polyarea()。

Answer 1

使用sf::st_area() 计算多面体的面积对我来说效果很好。这是从tidycensus 下载数据并计算佛蒙特州所有人口普查区面积的示例。

library(tidycensus)
library(sf)
library(dplyr)

vt_tracts <- get_acs(
  geography = "tract",
  variables = "B19013_001",
  geometry = TRUE,
  state = "VT"
)

vt_tracts %>% 
  mutate(area = st_area(.))
#> Simple feature collection with 184 features and 6 fields
#> geometry type:  POLYGON
#> dimension:      XY
#> bbox:           xmin: -73.43774 ymin: 42.72685 xmax: -71.46455 ymax: 45.01666
#> CRS:            4269
#> First 10 features:
#>          GEOID                                       NAME   variable estimate
#> 1  50001960100 Census Tract 9601, Addison County, Vermont B19013_001    78750
#> 2  50001960200 Census Tract 9602, Addison County, Vermont B19013_001    74958
#> 3  50001960300 Census Tract 9603, Addison County, Vermont B19013_001    59464
#> 4  50001960400 Census Tract 9604, Addison County, Vermont B19013_001    78672
#> 5  50001960500 Census Tract 9605, Addison County, Vermont B19013_001    56066
#> 6  50001960600 Census Tract 9606, Addison County, Vermont B19013_001    58343
#> 7  50001960700 Census Tract 9607, Addison County, Vermont B19013_001    53910
#> 8  50001960800 Census Tract 9608, Addison County, Vermont B19013_001    66713
#> 9  50001960900 Census Tract 9609, Addison County, Vermont B19013_001    65515
#> 10 50001961000 Census Tract 9610, Addison County, Vermont B19013_001    62763
#>      moe                       geometry            area
#> 1   8764 POLYGON ((-73.1968 44.26663... 212702063 [m^2]
#> 2   8731 POLYGON ((-73.39963 44.1553... 158329687 [m^2]
#> 3  10068 POLYGON ((-73.27275 44.1768...   6614856 [m^2]
#> 4   4139 POLYGON ((-73.43774 44.0450... 359804506 [m^2]
#> 5  10897 POLYGON ((-73.1495 44.17669... 109378523 [m^2]
#> 6   3601 POLYGON ((-73.06466 44.0408... 530756434 [m^2]
#> 7   7404 POLYGON ((-73.16763 44.0289...  70741373 [m^2]
#> 8  10440 POLYGON ((-73.19061 44.0189...  30961721 [m^2]
#> 9   4949 POLYGON ((-73.41258 43.9827... 478027734 [m^2]
#> 10  7192 POLYGON ((-73.18767 43.8964... 134156214 [m^2]

虽然要回答有关查找美国人口普查区域面积的更狭窄的问题，但您可以使用 tigris pacakge（这是 tidycensus 在引擎盖下使用来拉几何），除了几何之外，它还返回直接来自人口普查文件的每个区域（或任何其他请求的地理区域）的土地面积 (ALAND) 和水域 (AWATER)。

tigris::tracts(state = "VT", class = "sf")

#> Simple feature collection with 184 features and 12 fields
#> geometry type:  MULTIPOLYGON
#> dimension:      XY
#> bbox:           xmin: -73.4379 ymin: 42.72693 xmax: -71.46505 ymax: 45.01666
#> CRS:            4269
#> First 10 features:
#>    STATEFP COUNTYFP TRACTCE       GEOID NAME          NAMELSAD MTFCC FUNCSTAT
#> 1       50      007  000600 50007000600    6    Census Tract 6 G5020        S
#> 2       50      007  000800 50007000800    8    Census Tract 8 G5020        S
#> 3       50      007  000900 50007000900    9    Census Tract 9 G5020        S
#> 4       50      007  001100 50007001100   11   Census Tract 11 G5020        S
#> 5       50      001  960400 50001960400 9604 Census Tract 9604 G5020        S
#> 6       50      001  960900 50001960900 9609 Census Tract 9609 G5020        S
#> 7       50      001  960500 50001960500 9605 Census Tract 9605 G5020        S
#> 8       50      001  961000 50001961000 9610 Census Tract 9610 G5020        S
#> 9       50      001  960600 50001960600 9606 Census Tract 9606 G5020        S
#> 10      50      001  960300 50001960300 9603 Census Tract 9603 G5020        S
#>        ALAND   AWATER    INTPTLAT     INTPTLON                       geometry
#> 1    2240206   147061 +44.4833940 -073.1936217 MULTIPOLYGON (((-73.20609 4...
#> 2    1349586        0 +44.4609464 -073.2086195 MULTIPOLYGON (((-73.21524 4...
#> 3     602429        0 +44.4714041 -073.2088310 MULTIPOLYGON (((-73.21393 4...
#> 4    1736310        0 +44.4509355 -073.2197727 MULTIPOLYGON (((-73.23253 4...
#> 5  321773486 38567561 +44.0870615 -073.2671642 MULTIPOLYGON (((-73.43771 4...
#> 6  456817599 21511403 +43.9055842 -073.2904483 MULTIPOLYGON (((-73.41251 4...
#> 7  107525416  1714148 +44.0978095 -073.0485314 MULTIPOLYGON (((-73.1487 44...
#> 8  128778705  5371192 +43.9008582 -073.1035817 MULTIPOLYGON (((-73.18875 4...
#> 9  528821667  1880822 +43.9987500 -072.9398701 MULTIPOLYGON (((-73.06466 4...
#> 10   6413118   203405 +44.1632873 -073.2485362 MULTIPOLYGON (((-73.27275 4...