使用R语言并行提取栅格数据

在处理栅格数据时，经常需要将栅格数据中的数据按照shapefile边界或者中心点进行提取，当数据量较大时耗时很久，且容易造成内存溢出。因此考虑通过并行的方式对栅格数据进行提取，从而保证处理速度。

NOTE：目前对栅格数据的处理以raster包为主，同时terra包也可以进行，但实测表明terra包的处理速度在非并行条件下没有优势，因此仅对并行下两包的运算速度进行对比。

软件：R

电脑配置：Intel(R) Core(TM) i7-10700 CPU @ 2.90GHz ； 16核 32GB

计算背景：法国本土共有96个departments，现在需要提取每个departments的气温时间序列数据。气温数据来自ERA5，为nc格式，涵盖2019-01-01至2021-09-01这一时间段内每天的逐小时数据。提取方式是简单提取每个departments所覆盖的格点，计算均值作为该department的观测值，从而形成每个department的逐小时气温时间序列数据，并进行后续分析。

导入包及数据

library(terra);library(tidyverse);library(sf);library(ncdf4);library(rgdal)
library(ggplot2);library(lubridate);library(humidity);library(tsModel);library(raster)## nc 数据
era5_t2m  <- "ERA5_FRANCE_2019_21.nc"## read France shapefile
fr_shp <- vect("FRA_adm/FRA_adm2.shp")
fr_shp <- fr_shp[,c("ID_2")]

非并行使用raster

## LOAD ERA5 DATA
era5_t2m <- brick("ERA5_FRANCE_2019_21.nc",varname="t2m")## RUN AND CALCULATE TIME
start.time <- Sys.time()
dt.t2m = raster::extract(era5_t2m,fr_shp,fun=mean,na.rm=T)
dt.t2m <- t(dt.t2m[,-1])
Sys.time()-start.time

并行使用raster

library(foreach);library(doParallel)## SET PARAMETERS FOR PARALLEL
cores <- detectCores()
cl <- makeCluster(cores-2)
registerDoParallel(cl)start.time <- Sys.time()
dt.raster <- foreach(i = 1:96,.combine=cbind, ## i is the number of points or polygons.packages = c("raster","sf")) %dopar%{## READ NC FILEera5_t2m <- brick("ERA5_FRANCE_2019_21.nc",varname="t2m")## READ SHAPEFILEfr_shp <- st_read("FRA_adm/FRA_adm2.shp") fr_shp <- fr_shp[,c("ID_2")]    ## EXTRACT FILEtmp = raster::extract(era5_t2m,fr_shp[i,],fun=mean,na.rm=T)tmp <- t(tmp)
}
Sys.time()-start.timestopCluster(cl)

并行使用terra

library(foreach);library(doParallel)## SET PARAMETERS FOR PARALLEL
cores <- detectCores()
cl <- makeCluster(cores-2)
registerDoParallel(cl)start.time <- Sys.time()
dt.parallel <- foreach(i = 1:96,.combine=cbind,.packages = c("terra")) %dopar%{era5_t2m <- rast('NETCDF:"ERA5_FRANCE_2019_21.nc":t2m')fr_shp <- vect("FRA_adm/FRA_adm2.shp") fr_shp <- fr_shp[,c("ID_2")]                  tmp = terra::extract(era5_t2m,fr_shp[i,],fun=mean,na.rm=T)tmp <- t(tmp)
}Sys.time()-start.timestopCluster(cl)

运行结果

计算方式	运行时间
非并行使用raster	25 分钟
并行使用raster	7.5 分钟
并行使用terra	19.9分钟

讨论

上述结果表明，使用并行运算可以大幅度提升raster的数据提取时间，而对terra包的提升空间不大。

使用R语言并行提取栅格数据相关推荐

R语言批量提取excel当中的数据
使用R语言批量提取excel文件前言步骤总结前言在工作当中,为了方便统计单项数据,经常会遇上同一格式的几百份数据项需要提取到同一个表格的情况,博主非计算机专业,R语言代码不一定规范,但是经验 ...
R语言入门——读取栅格数据参数解读
读取栅格数据引言 1.数据读入 2.参数解读 2.1 class (类) 2.2 dimensions(维度) 2.3 resolution (像素) 2.4 extent(范围) 2.5 name ...
R语言中提取文件名和文件名
0. 背景需要对文件夹进行操作,需要提取文件的路径,以及文件的名称,搜索了半天,找到了答案,总结如下,希望对后来者有所帮助. 文件夹+文件名例子: nn [1] "f1_blue/y1/P ...
R语言学习-提取igraph的节点和边
网络分析的时候,可能需要提取出网络中的节点或者边,igraph包中其实提供了很多可用的函数. #创建网络方法之一:data.frame data<-data.frame(id1=c(1,1,2, ...
R语言——字符串提取函数
1.substr函数:提取指定位置的字符根据身份证号计算年龄: # 读入数据 id_number <- readxl::read_excel("D:/身份证.xlsx",s ...
raster | R语言中的空间栅格对象及其基本处理方法（Ⅰ）
前面的系列推文已经完成了对R语言中的两个管理空间矢量数据的工具包(sf和sp)的介绍,以及空间自相关.空间插值等空间分析方法. 这里小编再推出一个系列来介绍R语言中管理空间栅格数据的工具包:raste ...
R-GIS: 如何用R语言实现GIS地理空间分析及模型预测
前言:随着地理信息系统(GIS)和大尺度研究的发展,空间数据的管理.统计与制图变得越来越重要.R语言在数据分析.挖掘和可视化中发挥着重要的作用,其中在空间分析方面扮演着重要角色,与空间相关的包的数量也 ...
R语言笔记3：提取R对象的子集
R语言基础系列前情提要: 1数据类型(向量.数组.矩阵. 列表和数据框) 2读写数据所需的主要函数.与外部环境交互 Subsetting R Objects 取子集的三种基本方法 [ :"单 ...
R语言使用psych包的fa函数对指定数据集进行因子分析（输入数据为相关性矩阵）、使用rotate参数指定进行斜交旋转提取因子、使用fa.diagram函数可视化斜交旋转因子分析、并解读可视化图形
R语言使用psych包的fa函数对指定数据集进行因子分析(输入数据为相关性矩阵).使用rotate参数指定进行斜交旋转提取因子.使用fa.diagram函数可视化斜交旋转因子分析.并解读可视化图形目 ...