对时间序列栅格数据进行趋势分析，最后得到趋势图，可以用来分析栅格数据在每个栅格位置的变化趋势。
数据准备：时间序列栅格数据
ps: 最好按照时间序列规则命名
代码：
tic;clc;clear;s1 = 'sm\2000.tif';% 栅格数据存储在sm文件夹下[A,R] = geotiffread(s1);info = geotiffinfo(s1);%获取坐标系信息x = (1:21); %自变量[row,col] = size(A); %行列数 M = zeros(21,row,col)+NaN; %建立空数据立方体%将时间序列数据读入空数据立方体for year = 2000:s2 = ['sm\',int2str(year),'.tif'];d = geotiffread(s2);d(d<-200) = NaN; % 本数据的填充值为负数，所以要去掉填充值，不然会参与运算因此在计算之前要先查看 数据内容，如果有填充值就要根据实际情况去掉M(year-1999,:,:) = d;% 年份放在第一个位置endSlope = zeros(row,col);% 预先建立斜率空数组 Rsq = zeros(row,col); % R2Pvalue = zeros(row,col)+1; % P值B = zeros(row,col); % 截距%一元线性回归for i = 1:rowfor j = 1:colY = M(:,i,j);%因变量X = [ones(length(Y),1),x'];%自变量[b,bint,r,rint,stats] = regress(Y,X);Pvalue(i,j) = stats(3); %P值Slope(i,j) = b(2);%斜率B(i,j) = b(1);%截距Rsq(i,j) = stats(1); %R2enddisp(row); endSlope(isnan(Slope)) = 0; Rsq(isnan(Rsq)) = 0; Pvalue(isnan(Pvalue)) = 1; %替换空值%设置存储路径f1 = 'Trend\Slope.tif';%f2 = 'Trend\Pvalue.tif';f3 = 'Trend\Rsq.tif';f4 = 'Trend\B.tif';%导出数据geotiffwrite(f1,Slope,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);geotiffwrite(f2,Pvalue,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);geotiffwrite(f3,Rsq,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);geotiffwrite(f4,B,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);toc;
PS：代码和时间序列数据文件要放在同一个文件夹下，保证路径正确（如下图：sm中存储的是时间序列数据，右边为matlab代码，中间是趋势分析结果存放的位置，在运行程序前确保三个文件都处于同一文件夹下），其次，Trend文件夹需要在运行代码前先建立。
