R语言NDVI和降水量逐像元相关性分析

给定变量X和Y，它们均有测量数据。如果变量X的高值随变量Y的高值的变化而变化，并且X的低值随Y的低值的变化而变化，那么我们称X和Y是正相关的。如果X的高值随Y的低值的变化而变化，反之亦然，那么我们称变量X和Y是负相关的。

相关系数

相关系数：一种用于测量数值变量间相关程度的度量，取值范围在−1到+1之间。

r值域

对于相关系数的大小，我们一般认为：

•  存在显著性相关
• 高度相关
• 中度相关
• 低度相关
•  不相关

对于相关系数的正负：

• 当 负相关关系
• 当 正相关关系

常用的相关系数计算方法

常用的相关系数计算方法有：

• Pearson
• 一般用于两个呈正态分布的连续变量
• 两个变量观测值是成对的，观测值之间互相独立
• 样本中的极端值影响较大
• Spearman
• 对原始变量分布不作要求
• 适用范围较Pearson更广

P值

• 显著性检验可以用于确定观测到的效果是否落在零假设模型的随机变异范围内。
• 给定一个零假设模型，p 值表示模型所生成的结果与观测到的结果同样极端的概率。
• 一般认为：
• p <0.05显著
• p <0.01非常显著

第一类错误和第二类错误

• 在评估统计显著性时，可能会出现下面两种类型的错误。
• 第一类错误：错误地将仅由随机性导致的效果判定为真。
• 第二类错误：错误地将实际为真的效果判定为假的（即由随机性导致的）。

事实上，第二类错误并不是一种错误，它是由于判断样本的规模过小，而无法检测到效 果。如果p 值不足以表明统计显著性（例如，超过5%），我们应称其为“效果未验证”。增大样本的规模，可能会生成较小的p值。

显著性检验（即假设检验）的基本功能就是防止我们被随机性愚弄。因此，我们通常可以 通过构造显著性检验去最小化第一类错误。

R语言逐像元相关计算

• 方法：terra包栅格计算，Hmisc包计算相关系数和p值
• 读取栅格数据
• 栅格统一分辨率
• 栅格数据裁剪
• 栅格统一范围
• 相关系数逐像元并行计算

library(terra)
gimmslist = list.files(path = "./NDVI", pattern=".tif$")
gimmsdir = paste0("./NDVI/", gimmslist)
gimms = rast(gimmsdir)   #读取中国NDVI数据
pre2015_2017 = rast("pre_2015_2017.nc")    #读取中国降水量
pre2015Mon = pre2015_2017[[1:12]]            #选取2015年的
Pre2015g = resample(pre2015Mon, gimms)           #重采样
gimmsPre = trim(mask(gimms, Pre2015g))            #裁剪
gimmsPre2 = extend(gimmsPre, Pre2015g)

#GIMMS和CDR多年相关
z = c(gimmsPre2, Pre2015g)
#逐像元计算r和P
fun_cor =  function(x) {
  Rs = Hmisc::rcorr(x[1:12], x[13:24], type = "spearman")
  Rx = Rs$r[2]
  Px = Rs$P[2]
  return(c(Rx, Px))
}

#计算GIMMS和CDR的相关和P值
r_NDVI_PRE = app(z, fun_cor, cores=4)
names(r_NDVI_PRE) = c("r", "P-value")
plot(r_NDVI_PRE)
writeRaster(r_NDVI_PRE, filename = "r_NDVI_PRE.tif", names=r_NDVI_PRE@ptr[["names"]])

NDVI和降水量相关性

参考文献

1. https://zhuanlan.zhihu.com/p/27161877
2. https://zhuanlan.zhihu.com/p/350334110
3. https://zhuanlan.zhihu.com/p/102974265
4. 盖磊.面向数据科学家的实用统计学[M]. 人民邮电出版社2018

点击阅读原文即可参加视频课程学习