deldir | 生成泰森多边形的一段R语言代码

看到一段生成在R语言中生成泰森多边形的代码，感觉很有用。这段代码的最初来源是Carson Farmer的个人网站，网址是：https://carsonfarmer.com/2009/09/voronoi-polygons-with-r/。

voronoipolygons = function(layer) {
    require(deldir)
    require(sp)
    crds = layer@coords
    z = deldir(crds[,1], crds[,2])
    w = tile.list(z)
    polys = vector(mode = 'list', length = length(w))
    for (i in seq(along = polys)) {
        pcrds = cbind(w[[i]]$x, w[[i]]$y)
        pcrds = rbind(pcrds, pcrds[1,])
        polys[[i]] = Polygons(list(Polygon(pcrds)),
                              ID = as.character(i))
    }
    SP = SpatialPolygons(polys)
    voronoi = SpatialPolygonsDataFrame(SP, data = data.frame(x = crds[,1], 
        y = crds[,2], row.names = sapply(slot(SP, 'polygons'), 
        function(x) slot(x, 'ID'))))
}

这个函数只有一个参数layer，用起来很方便。不过它是基于sp结构的对象运行的，如果使用的矢量数据是sf结构的，需要进行相关转换。但是直接使用as函数进行转换会报错：

library(sf)
library(sp)
usa <- albersusa::counties_sf(proj = "laea")
usa.sp <- as(usa, "Spatial")
v <- voronoipolygons(usa.sp) # error

经过摸索，发现下面的代码是可行的：

library(sf)
library(sp)
usa <- albersusa::counties_sf(proj = "laea")
# to sp
coord <- st_coordinates(st_centroid(usa))
usa.data <- st_drop_geometry(usa)
usa.sp <- SpatialPointsDataFrame(coords = coord, data = usa.data)
# 生成泰森多边形
v <- voronoipolygons(usa.sp)
# 可视化
plot(v)

网址https://www.it1352.com/790552.html也提供了两个使用示例。

---

往期推荐阅读：

• 《数据处理通识》专辑-base | 使用apply族函数进行向量化运算

• 《制表与可视化》专辑-ggplot2 | ggplot2作图语法入门

• 《数学模型》专辑-car | 线性回归（三）——残差分析和异常点检验

• 《地理计算与分析》专辑-spdep | 如何在R语言中计算空间自相关指数