sp::over()。该点是否属于用 OGRGeoJSON 文件标识的多边形之一？

Question

我正在尝试获取布尔向量，例如，v[i] =1 告诉我第 i 个点（纬度经度对，存在于 train 数据框）属于 OGRGeoJSON 文件标识的地理区域之一。

OGR file的结构大致如下：

1区：24多边形
2区：4个多边形
3区：27多边形
4区：18个多边形
5区：34多边形

这就是我尝试做的。

但是，获得的结果不正确，因为生成的多边形是 OGR 文件中存在的所有不同区域的混合体。

library(rgdal)
library(httr)
library(sp)

r <- GET('https://data.cityofnewyork.us/api/geospatial/tqmj-j8zm?method=export&format=GeoJSON')
nyc_neighborhoods <- readOGR(content(r,'text'), 'OGRGeoJSON', verbose = F)

#New York City polygonal
pol_lat <- c(nyc_neighborhoods_df$lat)
pol_long <- c(nyc_neighborhoods_df$long)
xy <- cbind(pol_lat, pol_long)
p = Polygon(xy)
ps = Polygons(list(p),1)
pol = SpatialPolygons(list(ps))

#Points to analyse (pair of coordinates)
ny_lat <- c(train$pickup_latitude, train$dropoff_latitude)
ny_long <- c(train$pickup_longitude, train$dropoff_longitude)
ny_coord <- cbind(ny_lat, ny_long)
pts <- SpatialPoints(ny_coord)

#Query: Does the point to analyze fall in or out NYC?
over(pts, pol, returnList = TRUE)

我怎样才能解决这个问题以获得正确的结果？

Answer 1

sp 是一个较旧的软件包，正在逐步淘汰以支持较新的 "Simple Features" sf 软件包。让我知道您是否愿意使用 magrittr 包中的管道运算符 %>%，因为它与 sf 包配合得很好（dplyr 和 purrr).

使用 sf，您可以：

library(sf)

# Replace this with the path to the geojson file
geojson_path <- "path/to/file.geojson"

boroughs <- sf::st_read(dsn = geojson_path, stringsAsFactors = FALSE)

现在制作一个非常简单的空间点对象来代替 "trains" 数据。

# Make test data.frame

test_df <- 
  data.frame(
  # Random test point I chose, a couple of blocks from Central Park
      a = "manhattan_point", 
      y = 40.771959, 
      x = -73.964128, 
      stringsAsFactors = FALSE)

# Turn the test_df into a spatial object
test_point <-
  sf::st_as_sf(
    test_df,
    # The coords argument tells the st_as_sf function
    # what columns store the longitude and latitude data
    # which it uses to associate a spatial point to each
    # row in the data.frame
    coords = c("x", "y"), 
    crs = 4326 # WGS84
    )

现在我们可以确定我们的点落在哪个多边形上了：

# Get the sparse binary predicate. This will give a list with as 
# many elements as there are spatial objects in the first argument, 
# in this case, test_point, which has 1 element.
# It also has attributes which detail what the relationship is
# (intersection, in our case) 
sparse_bin_pred <- sf::st_intersects(test_point, boroughs)

# Output the boro_name that matched. I think the package purrr
# offers some more intuitive ways to do this, but
lapply(
  sparse_bin_pred, 
  function(x) boroughs$boro_name[x]
  )

最后一部分输出：

[[1]]
[1] "Manhattan"

sp::over()。该点是否属于用 OGRGeoJSON 文件标识的多边形之一？

sp::over(). Does the dot belong to one of the polygons identified with an OGRGeoJSON file?

r

geojson

spatial

sp

ogr