Matlab完成MK检验_matlab做mk突变分析

MK突变检测需要，查阅资料，发现口径不一，传说中的DPS软件不能用（正版的要钱）。大多数代码只能自己看图说话，找UB、UF交点作为突变点，但好多没有格网，放大都看不了，还有，作图上没有加For循环，一次画一个图，倒不是缺点，只是处理的数据多，懒，所以，根据网上说明自己修改了代码，感谢相关博主（见附录）。另，这个检测还可以用一个神奇额R包：Trend，进行相应检验。这里不详细论述。若有其他好方法请指点。

代码：

clear;clc;
 
Data=xlsread('AnnualIntensity_16.csv');
 
Year=Data(:,1);
 
[a,b]=size(Year);
 
Sig1=1.96*ones(a,b);
 
Sig2=-1.96*ones(a,b);
 
for i =1:17
 
Runoff=Data(:,i);
 
UF=MannKendall(Runoff);
 
Runoff2=flipud(Runoff);
 
UF2=MannKendall(Runoff2);
 
UB=-flipud(UF2);
 
figure(i)
 
plot(Year,UF,'k');
 
hold on
 
plot(Year,UB,'r');
 
plot(Year,Sig1,'--');
 
plot(Year,Sig2,'--');
 
xlabel('Year');
 
ylabel('统计量');
 
legend('UF','UB','α<0.05');
 
hold off
 
title={'Year','UF','UB','α','α'};
 
shuju=[Year UF UB Sig1 Sig2];
 
hold on
 
%hold on
 
grid(gca,'minor')% 画出次格网
 
Runoff=[] %清空数组
 
end;
 
%xlswrite('test01.xlsx',title,'Sheet1','A1');
 
%xlswrite('test01.xlsx',shuju,'Sheet1','A2');

            ------------------------------分-----------------------割----------------------线---------------------------------

增加保存图像功能：

clear;clc
 
Data=xlsread('numofExtreme_16.csv');
 
StationALL=xlsread('StationID.xlsx');
 
Year=Data(:,1);
 
[a,b]=size(Year);
 
Sig1=1.96*ones(a,b);
 
Sig2=-1.96*ones(a,b);
 
for i =2:17
 
Runoff=Data(:,i);
 
StationID=StationALL(i-1);
 
UF=MannKendall(Runoff);
 
Runoff2=flipud(Runoff);
 
UF2=MannKendall(Runoff2);
 
UB=-flipud(UF2);
 
figure(i)
 
plot(Year,UF,'k');
 
hold on
 
plot(Year,UB,'r');
 
plot(Year,Sig1,'--');
 
plot(Year,Sig2,'--');
 
xlabel('Year');
 
ylabel('统计量');
 
legend('UF','UB','α<0.05');
 
hold off
 
title={'Year','UF','UB','α','α'};
 
%shuju=[Year UF UB Sig1 Sig2];
 
 
 
hold on;%hold on
 
grid(gca,'minor');% 画出次格网
 
Runoff=[]; %清空数组
 
savepath=['D:\Figures\EPF',num2str(StationID),'.png'];%保存路径
 
saveas(i,savepath);%保存为png文件
 
end;
 
%xlswrite('test01.xlsx',title,'Sheet1','A1');
 
%xlswrite('test01.xlsx',shuju,'Sheet1','A2');

---------------------------------------------------------------------------------------------------------------------------------------------------------------------

应邀加个R的trend包趋势分析代码，需要先安装trend包。为防止乱码，在这里写一下，仅供参考

#EPI
rm(list=ls());
library(trend);
data <- ts(data=c(64.9434375,68.8955625,75.9810625,59.627625,53.6349375,58.2556875,53.319375,53.4730625,59.842,61.718125,56.6345,55.6553125,52.0409375,58.372875,56.81325,60.070625,58.4079375,52.10657143,54.649125,53.1128125,63.21625,57.8413125,55.9164375,56.5718,57.909375,54.8268125,55.968625,55.69225,60.05925,61.7584375,63.6616875,59.565375,57.103625,61.4910625,60.4355,58.321875,63.4930625,54.5916,63.313,56.3408,61.45193333,56.2614,56.25506667,58.99106667,71.18806667,54.0644,72.23966667,58.05253333,66.1172,54.56493333,62.973,67.73566667),start=1961);
res <- pettitt.test(data);res

参考：

Mann-Kendall突变检测(mk突变检测) - CSDN博客  http://blog.csdn.net/liyanzhong/article/details/41867859

科学网—matlab版本的MK程序，包含突变和趋势检验 - 张乐乐的博文  http://blog.sciencenet.cn/blog-1103122-851049.html

Mann-Kendall突变检测(mk突变检测) – MATLAB中文论坛  http://www.ilovematlab.cn/thread-246993-1-1.html