有没有一种有效的方法来计算 R 中 sf 中所有成对的多边形交集（不使用 for 循环）？

Question

我想根据数据框中的变量计算多边形的成对交集。

这就是我要找的东西：

set.seed(131)
library(sf)
m = rbind(c(0,0), c(1,0), c(1,1), c(0,1), c(0,0))
p = st_polygon(list(m))
n = 5
l = vector("list", n)
for (i in 1:n)
  l[[i]] = p + 2 * runif(2)
s = st_sfc(l)
s.f = st_sf(s)
s.f$id = c(1,1,2,2,3)
plot(s.f.2, col = s.f.2$id) 
    
s.f.2 = s.f %>% group_by(id) %>%   summarise(geometry = sf::st_union(s)) # %>% summarise(area = st_area(s))

s.f.2$area = st_area(s.f.2)


all.ids.pol = unique(s.f.2$id)
df = data.frame(NULL)
for (i in 1:length(all.ids.pol)) {
  for (j in 1:length(all.ids.pol)) {
    ol1 = st_intersection(s.f.2[c(all.ids.pol[i],all.ids.pol[j]),])
    if (dim(ol1)[1]>2 | i == j) {
      #add in an area count column (area of each arable poly in intersect layer)
      ol1$areaover <- st_area(ol1$geometry)
    } else {ol1$areaover = 0}
    #run the intersect function, converting the output to a tibble in the process
    if (i == j) {
      ol1.1 <- as_tibble(ol1)[1,]  
    } else {ol1.1 <- as_tibble(ol1)[2,]  }
    
    id.names = c(all.ids.pol[i],all.ids.pol[j])
    my.df =data.frame(area = ol1.1$areaover,id1 = id.names[1],id2 =id.names[2])
    df = rbind(df,my.df)
  }    
}
intersected.areas = df %>% 
  arrange(id1) %>%
  mutate(area.ha = units::set_units(area, ha)) %>% 
  select(-area) %>% 
  pivot_wider(names_from = id1, values_from = area.ha) %>%  
  arrange(id2) %>%
  rename(ID = id2)

最后，我想要一个矩阵，在行和列中显示“ID”，矩阵中的每个单元格显示多边形之间的交叉区域。

# A tibble: 3 × 4
     ID   `1`   `2`   `3`
  <dbl>  [ha]  [ha]  [ha]
1     1 1.59  0.501     0
2     2 0.501 1.86      0
3     3 0     0         1

所以代码可以工作，但我想知道 sf 中是否存在一个函数，它可以比使用 for 循环更有效地做到这一点（这对于大型数据集变得非常不切实际）

Answer 1

library(sf)
library(tidyverse)
s.f.2 %>%
  st_intersection(s.f.2) %>%
  mutate(intersect_area = st_area(.)) %>%
  st_drop_geometry() %>%
  pivot_wider(id_cols = id, names_from = id.1, 
              values_from = intersect_area,
              values_fill = 0)

# # A tibble: 3 x 4
#      id   `1`   `2`   `3`
#   <dbl> <dbl> <dbl> <dbl>
# 1     1 1.59  0.501     0
# 2     2 0.501 1.86      0
# 3     3 0     0         1