产品
最近在做工作的时候，想到NDVI数据受到降水量、温度的双重影响，所以在进行分析时，有些up主的工作只是对降水量和NDVI或者是温度和NDVI的相关分析。但是其实在做的时候应该考虑使用控制变量法来进行相关的研究，一般的研究只是假设NDVI值受到三个因素的影响，即与温度、降水量和人类活动因素有关，NDVI值、温度、降水量都可以获得相应的数据，人类活动因素则是使用以上的三个数据做残差分析后得到的。具体的过程大家可以在小破站上找到，此处不再赘述。
于是呢，我就在网上也学习了一下，加了一点自己的理解写了一个程序。我自己的水平也不是很高，公开出来呢只是为了方便后面的人在做到和我类似的工作的时候能够减少一点在程序上面花费的时间，更好的去分析解决自己所遇到的问题，如果帮到你了就给我一个点赞吧，有问题的话大家也可以在评论区里留言，大家一起交流进步，我也会不定时的回复大家的。废话不多说了，程序就放在下面了：
clc
clear
%导入投影信息
[a,R]=geotiffread('G:\毕业论文材料\裁剪数据\年度数据\年NDVI值\2000.tif');
info=geotiffinfo('G:\毕业论文材料\裁剪数据\年度数据\年NDVI值\2000.tif');
[m,n]=size(a);
%导入NDVI数据
NDVI=zeros(m*n,21);
for year=1:21
filename=strcat('G:\毕业论文材料\裁剪数据\2000-2020年归一化数据按月份、季度排序\二月\','NDVI',int2str(year),'.tif');%此处要修改,将年份修改为总共的年份。
data=importdata(filename);
data=reshape(data,m*n,1);
NDVI(:,year-0)=data;
end
%导入降水数据
PRE=zeros(m*n,21);
for year=1:21
filename=strcat('G:\毕业论文材料\裁剪数据\2000-2020年降水量裁剪数据按月份、季节排序\2月\','CNPre',int2str(year),'.tif');%此处要修改,将年份修改为总共的年份。
data=importdata(filename);
data=reshape(data,m*n,1);
PRE(:,year-0)=data;
end
%导入温度数据
TEMP=zeros(m*n,21);
for year=2000:2020
filename=strcat('G:\毕业论文材料\裁剪数据\2000-2020年温度裁剪数据按月份、季节排序\2月\','CNtmp_',int2str(year),'_M2_1x1.tif');%此处要修改,将年份修改为总共的年份。
data=importdata(filename);
data=reshape(data,m*n,1);
TEMP(:,year-1999)=data;
end
%温度的相关性和显著性
NDVI_TEMP_xgx=zeros(m*n,1);
NDVI_TEMP_p=zeros(m*n,1);
for i=1:length(NDVI)
ndvi=NDVI(i,:);
ndvi=ndvi';
temp=TEMP(i,:)';
pre=PRE(i,:)';
[r,p]=partialcorr(ndvi,pre,temp);
NDVI_TEMP_xgx(i)=r;
NDVI_TEMP_p(i)=p;
end
NDVI_TEMP_xgx=reshape(NDVI_TEMP_xgx,m,n);
NDVI_TEMP_p=reshape(NDVI_TEMP_p,m,n);
filename1=('G:\毕业论文材料\裁剪数据\温度降水量偏相关分析\降水量_NDVI\降水量_NDVI_2月数据偏相关分析.tif');
filename2=('G:\毕业论文材料\裁剪数据\温度降水量偏相关分析\降水量_NDVI\降水量_NDVI_2月数据偏相关分析显著性.tif');
geotiffwrite(filename1,NDVI_TEMP_xgx,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);
geotiffwrite(filename2,NDVI_TEMP_p,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);
%显著性检验
NDVI_TEMP_xgx=reshape(NDVI_TEMP_xgx,m*n,1);
NDVI_TEMP_p=reshape(NDVI_TEMP_p,m*n,1);
trend=zeros(m*n,1);
for i=1:length(NDVI)
if NDVI_TEMP_xgx(i)>0
if NDVI_TEMP_p(i)<0.05
trend(i)=1;
else NDVI_TEMP_p(i)>0.05;
trend(i)=2;
end
elseif NDVI_TEMP_xgx(i)<0
if NDVI_TEMP_p(i)<0.05
trend(i)=4;
else NDVI_TEMP_p(i)>0.05;
trend(i)=5;
end
else NDVI_TEMP_xgx(i)=0;
trend(i)=3;
end
end
trend=reshape(trend,m,n);
filename3=('G:\毕业论文材料\裁剪数据\温度降水量偏相关分析\降水量_NDVI\降水量_NDVI_2月数据通过显著性检验的偏相关分析.tif');
geotiffwrite(filename3,trend,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);
这个程序会输出三个图，第一个图是相关系数的图像，第二张图是显著性的图像，第三张图是通过显著性检查的出结论的图像，这里重点说说第三张图。
有些人在做完显著性检验之后，会发现自己的图一片黑或者是偶尔有几个点，所以我就把显著性检验的图改成了5个数，分别为1、2、3、4、5，分别代表着显著正向关、不显著正相关、非线性相关、显著负相关和不显著负相关，这些需要在整个图出来以后在ArcGIS里进行重分类，重分类之后就是你想要的结果了。重分类的步骤大家可以在网上找相关的文章，这里也不进行赘述。
最后说几个比较重要的问题：
155-2731-8020
