优化代码：将大型全球矢量数据集拆分至栅格单元独立文件

针对将大型全球矢量数据按栅格单元拆分的性能瓶颈，以下是几个高效的优化方向，结合代码调整说明：

1. 移除不必要的栅格文件IO操作 #

当前代码中先将每个栅格单元写入GPKG再读取获取边界框，会产生大量磁盘IO开销。直接在内存中处理栅格的边界信息即可：

library("sf")
library("terra")
library("glue")
library("rnaturalearth")
library("tidyverse")

dir.create("tmp", showWarnings = FALSE)

# 生成全球栅格（直接在内存中处理，无需写入文件）
world_grid <- rast(nrows=1, ncols=1, xmin=-180, xmax=180, ymin=-90, ymax=90, crs="epsg:4326") %>%
  st_bbox() %>%
  st_as_sfc() %>%
  st_make_grid(cellsize = 10) %>%
  st_as_sf() %>%
  mutate(grid_id = row_number(),
         bbox = map(geometry, st_bbox)) # 直接在内存中存储每个栅格的边界框

# 生成每个栅格对应的输出文件名，无需写入栅格文件
grid_info <- world_grid %>%
  transmute(output_name = glue("tmp/worldwide_poly_{grid_id}.gpkg"),
            bbox = bbox)

2. 用原生空间操作替代系统调用ogr2ogr #

反复调用system()启动ogr2ogr进程会产生大量进程创建销毁开销，改用sf/terra的原生空间裁剪函数，直接在内存中处理后写入文件：

# 加载源数据到内存（无需写入临时文件）
world_data <- ne_countries(type = "countries", scale = "large", returnclass = "sf") %>% 
  select(iso_a2)

# 拆分函数：直接用st_crop裁剪并写入文件
split_world_to_tiles <- function(output_name, bbox, source_data) {
  # 将bbox转换为sfc对象用于裁剪
  crop_extent <- st_as_sfc(bbox)
  # 裁剪数据
  cropped_data <- st_crop(source_data, crop_extent)
  # 写入文件（如果裁剪结果非空）
  if(nrow(cropped_data) > 0) {
    st_write(cropped_data, output_name, append = FALSE, quiet = TRUE)
  }
  return(output_name)
}

# 应用拆分（后续可替换为并行版本）
pwalk(grid_info, split_world_to_tiles, source_data = world_data)

3. 为源数据建立空间索引 #

给大型矢量数据建立空间索引能大幅加速裁剪查询的速度，内存中处理时可手动创建，或写入文件时自动生成：

# 为内存中的源数据创建空间索引
st_geometry(world_data) <- st_geometry(world_data) %>% st_sfc(crs = st_crs(world_data))
st_set_geometry(world_data, st_geometry(world_data))

# 超大型数据建议写入文件时创建索引再读取
st_write(world_data, "tmp/world_polygons.gpkg", layer_options = "SPATIAL_INDEX=YES")
world_data <- st_read("tmp/world_polygons.gpkg", quiet = TRUE)

4. 优化并行处理策略 #

使用并行时避免过度并行导致磁盘IO竞争，用更高效的并行框架（如furrr）：

library(furrr)
# 设置并行工作进程数（建议为CPU核心数的70%-80%）
plan(multisession, workers = parallel::detectCores() - 2)

# 并行执行拆分
future_pwalk(grid_info, split_world_to_tiles, source_data = world_data)

5. 调整ogr2ogr参数（若仍需使用系统调用） #

如果必须用ogr2ogr，优化参数减少不必要操作：

• 移除-nlt GEOMETRYCOLLECTION，保持原几何类型更高效
• 使用-sql在OGR层面直接过滤，减少数据传输
• 批量生成命令减少进程启动次数

ogr2ogr_path <- "C://Program Files//QGIS 3.26.1//bin//ogr2ogr.exe"
# 优化后的ogr2ogr调用命令
command <- glue('"{ogr2ogr_path}" -spat {bbox$xmin} {bbox$ymin} {bbox$xmax} {bbox$ymax} -f GPKG {output_name} {worldwide_data}')

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