使用底格里斯从纬度/经度获取人口普查区域

Question

我有相对大量的坐标，我想获取人口普查区域（除了 FIPS 代码）。我知道我可以使用call_geolocator_latlon（如这里所做的那样）查找单个纬度/经度对，但这对我的目的来说似乎不切实际，因为该函数对人口普查局的 API 发出一次调用，我想这需要很长时间在我的约 200,000 双上运行。

有没有更快的方法来做到这一点，也许是通过使用block_groups函数下载每个州的 shapefile并从那里映射到经纬度到人口普查区？

Answer 1

这不使用tigris，而是sf::st_within()用于检查重叠区域的点数据框。

我在tidycensus这里使用将加利福尼亚大片的地图导入 R。

library(sf)

ca <- tidycensus::get_acs(state = "CA", geography = "tract",
              variables = "B19013_001", geometry = TRUE)

现在来模拟一些数据：

bbox <- st_bbox(ca)

my_points <- data.frame(
  x = runif(100, bbox[1], bbox[3]),
  y = runif(100, bbox[2], bbox[4])
  ) %>%
  # convert the points to same CRS
  st_as_sf(coords = c("x", "y"),
           crs = st_crs(ca))

我在这里做了 100 分以ggplot()获得结果，但是 1e6 的重叠计算速度很快，在我的笔记本电脑上只有几秒钟。

my_points$tract <- as.numeric(st_within(my_points, ca)) # this is fast for 1e6 points

结果：

head(my_points) # tract is the row-index for overlapping census tract record in 'ca'

# but part would take forever with 1e6 points
library(ggplot2)

ggplot(ca) +
  geom_sf() +
  geom_sf(data = my_points, aes(color = is.na(tract)))