如何在将大型空间多边形数据框转换为光栅堆栈时提高速度

Question

我正在处理空间多边形数据框（SPDF）数据集。该数据集包含一个月度气候数据的时间序列。我的目标是将SPDF转换为1000米分辨率的栅格堆叠，以进行进一步的统计分析。我已经编写了R代码，但它非常缓慢，转换一个列需要很长时间。如果您有任何建议可以使其更快，我将不胜感激。

hru<-readOGR("E:\\Tade\\HRU\\ubn_merge.shp",layer="ubn_merge") # spatial polygon
spdf<-merge(hru,spdf.2000,by.x="HRU",by.y="HRU",all.x=T,sort=F) # spdf nrow=565 ncol=375
# convert sp to raster
hru.ras<-raster(hru,resolution=1000) # raster hru shape to 1km
for (i in 1:length(spdf){
  et.ras<-rasterize(spdf,hru.ras,field=paste("m",1,sep="")) # rasterize
  et.stack<-stack(et.stack,et.ras)
}

谢谢

- Tade

---

根据我的经验，如果使用R和rasterize()函数来完成此任务，你将会面对许多问题。通过我的测试，直接使用gdal_rasterize工具将矢量转换为栅格更易管理。我使用系统调用来调用它，例如system(paste0("gdal_rasterize -burn 1 -l census_buffer ", census_buffer_path, " /tmp/census_mask.tif"))等，这当然假定了您的路径中包含gdal实用程序。 - Forrest R. Stevens

@ForrestR.Stevens 谢谢。 - Tade

Answer 1

正如Forrest所说（而你已经经历过的那样），rasterize有点慢，但你可以做得比现在更好。你根本不需要使用循环：

r <- raster(spdf, resolution=1000)
et.ras <-rasterize(spdf, r, field=paste0("m",1:ncol(spdf)))

这将创建一个带有栅格属性表的单个 RasterLayer。要创建 RasterStack，请执行以下操作：

s <- deratify(et.ras)

如果您要使用循环，请仅使用rasterize一次，获取多边形ID，然后使用subs获取实际感兴趣的变量。

当然，这还可以避免外部依赖带来的痛苦。

- Robert Hijmans

---

谢谢。速度显著提高（几个小时到几分钟）。 - Tade

一如既往，Robert，你提供了很多值得思考的内容（以非常好的方式）。我之前没有想过你所概述的方法，实际上我认为它会在某些用例中提高性能，而不是发出多个调用到gdal_rasterize。谢谢！ - Forrest R. Stevens

Answer 2

截至2020年3月，您现在可以考虑使用这个新的包。请参考以下链接以了解如何使用 fasterize 包：https://cran.r-project.org/web/packages/fasterize/vignettes/using-fasterize.html。

fasterize()

这个功能完全符合

rasterize()

但速度快100-1000倍。我认为这可以帮助很多人。

- Antoine C.