其他
R语言并行计算Sen+MK趋势分析
并行计算Sen+MK趋势分析
前面的推文中讲过Sen+MK分析
上面的方法有个很大的缺点就是R语言计算过程太慢,对此我改进了前面的代码,使用terra
包重写,开启并行计算,提高计算速度。
并行计算Sen+MK
rast
读取栅格数据function
中长度、年份需要根据实际情况修改输出结果包含三个波段,Z值、slope和p值(详细看参考文献推文)
library(terra)
library(trend)
#输入一个文件夹内的单波段TIFF数据,在这里是历年的NDVI年最大值
flnames <- list.files(path = './ChinaYearMean/', pattern = '.tif$')
fl <- paste0("./ChinaYearMean/", flnames)
firs <- rast(fl)
#Sen+MK计算
fun_sen <- function(x){
if(length(na.omit(x)) <34) return(c(NA, NA, NA)) #删除数据不连续含有NA的像元
MK_estimate <- trend::sens.slope(ts(na.omit(x), start = 1982, end = 2015, frequency = 1), conf.level = 0.95) #Sen斜率估计
slope <- MK_estimate$estimate
MK_test <- MK_estimate$p.value
Zs <- MK_estimate$statistic
return(c(Zs, slope, MK_test))
}
firs_sen = app(firs, fun_sen, cores=4 )
names(firs_sen) = c("Z", "slope", "p-value")
#显示输出
plot(firs_sen)
writeRaster(firs_sen, filename = "./firs_sen.tif", names=firs_sen@ptr[["names"]])