`poly.counts`的等价物用来计算落在多边形内的sf包中的纬度/长度对

Question

sf包提供了一种处理地理要素的好方法，但是我找不到与poly.counts功能从GISTools包需要sp对象。

poly.counts从SpatialPointsDataFrame秋天SpatialPolygonsDataFrame的多边形内计算的点的数量，并且可以如下使用：

数据

每个多边形

我生成两个

## Libraries 
library("GISTools") 
library("tidyverse") 
library("sf") 
library("sp") 
library("rgdal") 
## Obtain shapefiles 
download.file(url = "https://www2.census.gov/geo/tiger/TIGER2016/STATE/tl_2016_us_state.zip", destfile = "data-raw/states.zip") 
unzip(zipfile = "data-raw/states.zip", exdir = "data-raw/states") 
sf_us_states <- read_sf("data-raw/states") 

## Our observations: 
observations_tibble <- tribble(
    ~lat, ~long, 
31.968599, -99.901813, 
35.263266, -80.854385, 
35.149534, -90.04898, 
41.897547, -84.037166, 
34.596759, -86.965563, 
42.652579, -73.756232, 
43.670406, -93.575858 
) 

计算点我的sp物品：

sp_us_states <- as(sf_us_states, "Spatial") 

observations_spdf <- observations_tibble %>% 
    select(long, lat) %>% # SPDF want long, lat pairs 
    SpatialPointsDataFrame(coords = ., 
         data = ., 
         proj4string = sp_us_states@proj4string) 

现在我可以用poly.counts

points_in_states <- 
    poly.counts(pts = observations_spdf, polys = sp_us_states) 

一下添加到sp对象：

sp_us_states$points.in.state <- points_in_states 

现在我已经完成了我会转换回sf对象，并可以将其想象为：

library("leaflet") 
updated_sf <- st_as_sf(sp_us_states) 
updated_sf %>% 
    filter(points.in.state > 0) %>% 
    leaflet() %>% 
    addPolygons() %>% 
    addCircleMarkers(
    data = observations_tibble 
) 

问题

我可以在sf和sp对象之间进行繁琐的转换吗？

Answer 1

尝试以下操作：

sf_obs = st_as_sf(observations_tibble, coords = c("long", "lat"), 
    crs = st_crs(sf_us_states)) 
lengths(st_covers(sf_us_states, sf_obs)) 
# check: 
summary(points_in_states - lengths(st_covers(sf_us_states, sf_obs))) 

st_covers返回与由每个国家覆盖点的索引中的列表; lengths返回这些向量长度的向量或点数。您将看到的警告表明，尽管您拥有地理坐标，但底层软件假设它们是笛卡尔式的（在这种情况下，这很可能不会产生问题;如果要以正确的方式摆脱它，则转到投影坐标） ）