基于Matlab的逐像元Hurst指数（R/S）分析以全球19822015年NDVI变化为例。

自相似性和长期依赖性是自然界普遍存在的现象，并在水文、气候、地质和地震等领域广泛运用，Hurst 指数是描述该现象的有效方法 。目前，Hurst 指数的估算方法有多种，如绝对值法、聚合方差法、周期图法、小波分析法、残差分析法和 Ｒ/S 分析法等，有关研究表明:Ｒ/S 分析法和小波分析法估算的 Hurst 指数要比其他方法估算的结果更可靠 。基于重标极差(Ｒ/S)分析方法的 Hurst 指数最早是由英国水文学家 Hurst 在研究尼罗河水库流量和储存能力的关系时提出 。本文采用常用的Ｒ/S 分析法，以中国2000-2015年NDVI变化状况为例，进行逐像元分析，并形成一副栅格图层。
基于R/S的技术流程可自行参考文末的文献，本文利用Matlab来进行实现。主要代码如下：

clear
[aa,R]=geotiffread('H:\Global\NDVI3g\GIMMSraster\raster\最大合成\GIMMS_NDVI2015.tif');%先投影信息
info=geotiffinfo('H:\Global\NDVI3g\GIMMSraster\raster\最大合成\GIMMS_NDVI2015.tif');
ndvisum=zeros(size(aa,1)*size(aa,2),34);%34期数据
for year=1982:2015
filename=strcat('H:\Global\NDVI3g\GIMMSraster\raster\最大合成\GIMMS_NDVI',int2str(year),'.tif');
ndvi=importdata(filename);
ndvi=reshape(ndvi,size(ndvi,1)*size(ndvi,2),1);
ndvisum(:,year-1981)=ndvi;
end
hsum=zeros(size(aa,1),size(aa,2))+NaN;
for kk=1:size(ndvisum,1);
ndvi=ndvisum(kk,:);
if min(ndvi)>0
ndvi_cf=[];
for i=1:length(ndvi)-1
ndvi_cf1=ndvi(i+1)-ndvi(i);
ndvi_cf=[ndvi_cf,ndvi_cf1];
end
M=[];
for i=1:size(ndvi_cf,2)
M1=mean(ndvi_cf(1:i));
M=[M,M1];
end
S=[];
for i=1:size(ndvi_cf,2)
S1=std(ndvi_cf(1:i))*sqrt((i-1)/i);
S=[S,S1];
end
for i=1:size(ndvi_cf,2)
for j=1:i
der(j)=ndvi_cf(1,j)-M(1,i);
cum=cumsum(der);
RR(i)=max(cum)-min(cum);
end
end
RS=S(2:size(ndvi_cf,2)).\RR(2:size(ndvi_cf,2));
T=[];
for i=1:size(ndvi_cf,2)
T1=i;
T=[T,T1];
end
lag=T(2:size(ndvi_cf,2));
g=polyfit(log(lag/2),log(RS),1);
H=g(1);
hsum(kk)=H;
clear der
end
end
geotiffwrite('1982-2015年全球NDVI_Hurst指数.tif',hsum,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);


结果如下所示：

Hurst指数.jpg

参考文章：
严恩萍,林辉,党永峰,夏朝宗. 2000—2012年京津风沙源治理区植被覆盖时空演变特征[J]. 生态学报,2014,34(17):5007-5020.