st_intersection有没有办法返回0值而不是切割不相交的特征？

Question

我有一个缓冲区形状文件和一个道路形状文件，我试图计算出与每个缓冲区相交的道路的总长度。有一些缓冲区不与任何道路相交，因此它们从结果中被删除。我仍然需要这些缓冲区的值，即使它们等于零。任何见解将不胜感激，谢谢。我的代码如下：

###  Read in the buffer shapefile
buffer_only <- st_read("1kmbuff.shp")

### read in roads shapefile: 
  ## Road shapefile accessed from state of michigan GIS opendata
  ## https://gis-michigan.opendata.arcgis.com/datasets/Michigan::all-roads-v17a/about

roads <- read_sf(dsn = "All_Roads_(v17a).shp")
  # Distance (length) unit of measure is meters.
roads
roads.prep <- st_transform(roads, crs(buffer_only))

roads.buff  <- sf::st_intersection(roads.prep, buffer_only) %>%
  mutate(length_m = st_length(geometry))

我尝试使用 st_intersection 但它剪掉了不相交的行。

Answer 1

我建议使用

sf::st_intersects()

来计算与缓冲区相交的要素（道路）的数量。

让我们先准备一些数据（您将用您的形状文件替换它）：

bb <- osmdata::getbb("Oborniki Śląskie", format_out = "sf_polygon") |>
  sf::st_centroid() |>
  sf::st_buffer(dist = 200) |>
  sf::st_bbox() #|>
  as.array()
#> Error in as.array.default(): argument "x" is missing, with no default
  osmdata::getbb()
#> Error in get_nominatim_query(place_name, featuretype, is_polygon, display_name_contains, : argument "place_name" is missing, with no default

h <- osmdata::opq(bb) |>
  osmdata::add_osm_feature(key = "highway") |>
  osmdata::osmdata_sf()

h <- h$osm_lines

buffer_only <- h[1:10, ] |>
  sf::st_buffer(dist = 8)

tmap::tmap_mode("plot")
#> tmap mode set to 'plot'
tmap::tm_shape(buffer_only) +
  tmap::tm_polygons() +
  tmap::tm_shape(h[30:45, ]) +
  tmap::tm_lines(lwd = 2)

我们有 10 个多边形（缓冲区）和 16 条道路，使用

st_intersects()

我们可以检查哪些（以及其中有多少）与我们的多边形相交。

st_intersects()

创建一个稀疏几何二元谓词列表，例如：

sf::st_intersects(buffer_only, h[30:45,])
#> Sparse geometry binary predicate list of length 10, where the predicate
#> was `intersects'
#>  1: 5, 15
#>  2: (empty)
#>  3: 15, 16
#>  4: (empty)
#>  5: (empty)
#>  6: 15
#>  7: 15
#>  8: 10, 15
#>  9: 1
#>  10: (empty)

我们可以看到第一个缓冲区与 2 条路相交，第二个缓冲区与 0 条相交，等等。使用

sparse

参数，我们可以创建

FALSE/TRUE

矩阵，可用于计算相交道路的数量：

sf::st_intersects(buffer_only, h[30:45,], sparse = FALSE)
#>        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12]
#>  [1,] FALSE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
#>  [2,] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
#>  [3,] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
[...]

让我们计算每个缓冲区的相交道路数量：

mat <- sf::st_intersects(buffer_only, h[30:45,], sparse = FALSE)

x <- unlist(lapply(1:10, function(i) length(which(mat[i, ] == TRUE))))

x
#>  [1] 2 0 2 0 0 1 1 2 1 0

您将使用

st_length()

单独计算长度，但是您可以使用

as.data.frame(st_intersects())

的结果将缓冲区数量和相交道路数量作为数据框中的单独行，这可能很有用。

^{创建于 2024 年 10 月 30 日，使用 reprex v2.1.0}