高效批量创建多边形：R语言sp与terra包性能优化探究

我需要在apply函数内基于坐标矩阵快速批量创建多边形，发现这是代码中耗时较长的环节。我测试了两种实现方式：

方法1：先创建sp多边形再转Terra格式 #

library(sp)
library(terra)    
t0 <- Sys.time()

poly_list <- apply(matrix(1:10000), 1, function(idx){      
  # 生成坐标
  coords <- cbind(rnorm(100), rnorm(100))      
  # 创建Polygons对象
  Polygons(list(Polygon(coords)), idx)
})
  
# 转换为Terra多边形
poly_terra <- vect(SpatialPolygons(poly_list))

# 输出耗时
print(Sys.time() - t0)   
# 耗时：2.082166 秒

方法2：直接用Terra创建多边形 #

t0 <- Sys.time()    
poly_list <- apply(matrix(1:10000), 1, function(idx){      
  # 生成坐标
  coords <- cbind(rnorm(100), rnorm(100))      
  # 创建Terra多边形
  vect(coords, type = "polygon")
})

# 合并为单个Terra多边形对象
poly_terra <- vect(poly_list)

print(Sys.time() - t0)
# 耗时：16.38044 秒

疑问 #

为什么先创建sp多边形再转换为terra格式，比直接用terra创建多边形更快？而且vect(SpatialPolygons(Polygons(list(Polygon(coords)), idx)))这类代码写法繁琐，是否存在更快或更简洁的实现方式？

编辑补充：目前我找到的最快方案，尽管感觉不够规范：

t0 <- Sys.time()    
dummy <- Polygons(list(Polygon(cbind(rep(0,4), rep(0,4)))), "0")   
poly_list <- apply(matrix(1:10000), 1, function(idx){
  
  # 生成坐标
  coords <- cbind(rnorm(100), rnorm(100))     
  # 复用dummy对象修改属性
  new <- dummy
  new@ID <- as.character(idx)
  new@Polygons[[1]]@coords <- coords
  return(new)
})

# 转换为Terra多边形
poly_terra <- vect(SpatialPolygons(poly_list))
print(Sys.time() - t0)
# 耗时：0.7147191 秒

解答 #

性能差异的原因 #

1. 底层实现逻辑不同：sp包的Polygon/Polygons/SpatialPolygons是纯R层面的S4对象构造，逻辑简单直接，仅做基础的对象属性赋值，没有额外的几何校验、空间参考处理或跨语言接口调用开销；而terra::vect创建多边形时，会调用底层GDAL库完成几何对象初始化、拓扑验证（比如检查多边形是否闭合、是否存在自相交），还会默认处理空间参考等元数据，这些额外操作会拉高单次创建的耗时，当循环执行10000次时，累计开销会被放大数倍。
2. 批量处理的效率差异：方法2中每次调用vect都会生成一个独立的SpatVector对象，后续合并时又要重复处理每个对象的元数据；而方法1是先批量构造轻量的Polygons对象列表，一次性转换为SpatialPolygons再转SpatVector，减少了多次初始化和元数据处理的重复开销。

更高效/简洁的实现方式 #

1. 优化“复用对象”方案（保留高性能，代码更整洁）

直接修改S4对象槽位的方式虽然看起来“不规范”，但只要保证坐标格式合法（比如多边形闭合），就是性能最优的方案之一。可以封装成函数让代码更易维护：

library(sp)
library(terra)

fast_build_polygon <- function(coords, id) {
  # 创建一个dummy模板对象
  dummy <- Polygons(list(Polygon(cbind(0,0,0,0,0,0,0,0))), "0")
  # 修改ID和坐标
  dummy@ID <- as.character(id)
  dummy@Polygons[[1]]@coords <- coords
  dummy
}

t0 <- Sys.time()
# 用lapply替代apply（性能略优）
poly_list <- lapply(1:10000, function(idx) {
  coords <- cbind(rnorm(100), rnorm(100))
  # 确保多边形闭合（如果原始坐标未闭合）
  if (!all(coords[1,] == coords[nrow(coords),])) {
    coords <- rbind(coords, coords[1,])
  }
  fast_build_polygon(coords, idx)
})
poly_terra <- vect(SpatialPolygons(poly_list))
print(Sys.time() - t0)

2. 批量构造Terra多边形（无需依赖sp）

如果能把所有多边形的坐标整理成统一格式，可以直接用terra::vect批量创建，彻底避免循环开销——因为只需要做一次GDAL初始化和几何验证，性能远高于循环调用：

library(terra)

t0 <- Sys.time()
# 批量生成所有多边形的坐标（确保每个多边形闭合）
all_coords <- do.call(rbind, lapply(1:10000, function(idx) {
  coords <- cbind(rnorm(100), rnorm(100))
  if (!all(coords[1,] == coords[nrow(coords),])) {
    coords <- rbind(coords, coords[1,])
  }
  coords
}))
# 指定每个多边形的顶点数量（这里每个多边形是101个点，含闭合点）
part_counts <- rep(nrow(all_coords)/10000, 10000)
# 批量创建多边形
poly_terra <- vect(all_coords, type = "polygon", part = part_counts)
print(Sys.time() - t0)

内容的提问来源于stack exchange，提问作者Zoe