R gis：基于邻居值添加多边形属性？

Question

我有一组林分（SpatialPolygonDataFrame），它们随机分布在整个景观上，即分散并聚集在一起。对于每个多边形，我想确定它是否具有开放边缘。 多边形具有开放边缘如果：

• 没有邻居
• 至少在一侧没有邻居；
• 有邻居，但是林分和树高之间的差异邻居超过5

我徘徊如何将属性open_edge = TRUE/FALSE添加到各个多边形？

我虽然（（[[pseudocode）：

每个子集一一站立（在for循环中]
• 创建周围的缓冲区
• 通过buffer与支架的重叠将支架子集包围
• [如果有邻居->比较高度。如果差异> 5，则open_edge = TRUE

但是，这种方法没有考虑摊位所具有的条件，比如说仅在3面的邻居，例如说是新邻居。 poly2nb工具似乎很有前途，但是如何为单个看台添加属性？

这是我的虚拟方法，但我想知道您是否有更有效的解决方案？

创建虚拟数据：

library(ggplot2)  # for choropleth map plot
library(broom) # to convert spatial data to dataframe
library(mapproj)
library(spdep)    # neighbours
library(rgdal)
library(rgeos)
library(sf)
library(raster)
library(dplyr)
library(spData)
library(sf)

r <- raster(nrow=6, ncol=6, crs = "+init=epsg:2957")
values(r) <- matrix(data = c(9,  NA, NA, NA, NA, NA,
                             NA, NA, NA, NA, NA, NA, 
                             NA, NA, NA, NA, NA, NA, 
                             NA, NA, NA, 1, 1, 1,
                             NA, NA, NA, 1, 9, 1,
                             NA, NA, NA, 1, 1, 1),
                    nrow = 6,
                    ncol = 6, 
                    byrow = TRUE)



# Convert raster to polygon
polys <- rasterToPolygons(r)

[以一个架子为例，确定架子是否具有开放边缘：
# Subset first row in SpatialPolygonDataFrame
i = 10
one = polys[i, ]

# Keep the remaining polygons
left = polys[-i,]

# Create buffer within distance
buff = buffer(one, width = 100)



# subset set of neighbours by spatial overlap
nbrs <- left[which(gContains(sp::geometry(buff),
                                    sp::geometry(left), byid = T)),]


# Compare if the values are different
height.one  = rep(one$layer[1], nrow(nbrs))
height.nbrs = nbrs$layer

# Get the differences between the neighbouring stands
difference = height.one - height.nbrs

# If the difference in at least one stand is 
# in more than 5, set open_edge = TRUE 
# or if no neighbours find
one$open_edge <- any(difference > 5)

Answer 1

让您开始使用spdep::poly2nb

library(raster)

r <- raster(nrow=6, ncol=6, crs = "+init=epsg:2957")
values(r) <- matrix(data = c(9,  NA, NA, NA, NA, NA,
                             NA, NA, NA, NA, NA, NA, 
                             NA, NA, NA, NA, NA, NA, 
                             NA, NA, NA, 1, 1, 1,
                             NA, NA, NA, 1, 9, 1,
                             NA, NA, NA, 1, 1, 1),
                    nrow = 6,
                    ncol = 6, 
                    byrow = TRUE)

# Convert raster to polygon
polys <- rasterToPolygons(r)

library(spdep)
nb <- poly2nb(polys)

plot(polys)
text(polys, 1:10)

str(nb)
#List of 10
# $ : int 0
# $ : int [1:3] 3 5 6
# $ : int [1:5] 2 4 5 6 7
# $ : int [1:3] 3 6 7
# ...

所以poly 1没有邻居，poly 2有邻居3、5、6，等等。>
现在您可以使用sapply计算事物。例如邻居数
sapply(nb, function(i) length(i[i>0]))
#[1] 0 3 5 3 5 8 5 3 5 3