气候、水文、气象等时间序列突变检测方法

hykamw123

谢谢楼主分享，{:eb502:}{:eb502:}{:eb502:}

江湖小鱼儿

感谢大佬分享~学习中

Jude

楼主大大，我下载了你的Pettitt程序，发现和我自己编的不一样，计算出的突变点有所不同啊

kingjames

Jude 发表于 2018-9-22 17:19
楼主大大，我下载了你的Pettitt程序，发现和我自己编的不一样，计算出的突变点有所不同啊

？？？这样啊 方便发你的代码上来吗？我对比验证下结果

Jude

kingjames 发表于 2018-9-23 09:19
？？？这样啊 方便发你的代码上来吗？我对比验证下结果

楼主大大，我是做水资源方向的，我用自己编的程序，结合66年数据资料，检测出突变点在44年附近，而用楼主的程序计算出的突变点在66年，这显然有点问题吧。

kingjames

Jude 发表于 2018-9-23 09:32
楼主大大，我是做水资源方向的，我用自己编的程序，结合66年数据资料，检测出突变点在44年附近，而用楼主 ...

我忘记以前的代码了，写论文就用过一次，感觉挺好的呀。你要不试一试这个

%% function [Ut,chg_point]=pettitt_change(x,data)
% Pettitt test for detection of change point
A=xlsread('shuju2008.xlsx');
x=A(:,1);   % 时间
data=A(:,2);% 年平均气温
n=length(data);
s=0;
Ut=zeros(n-1,1);
for i=1:n-1;
  for j=i+1:n;
      s=s+sign(data(i)-data(j));
  end
  Ut(i)=s;
end
Kmax=max(abs(Ut));
chg_point=find(abs(Ut)==Kmax);%找出最大值对应的序列号
coef=2*exp(-6*Kmax^2/(n^3+n^2));
figure(1)
hold on
plot(x(1:n-1),Ut,'K-','linewidth',1.5);%画出整个曲线的变化趋势
plot([x(chg_point),x(chg_point)],[-Kmax,Kmax],'r--','linewidth',1.5);%画出该竖线即为突变点
Klim=Kmax-Kmax*0.1;
plot([x(1),x(n)],[-Klim,-Klim],'b-');
plot([x(1),x(n)],[Klim,Klim],'b-');

王当当

非常感谢楼主的分享

Jude

kingjames 发表于 2018-9-26 10:51
我忘记以前的代码了，写论文就用过一次，感觉挺好的呀。你要不试一试这个

%% function =pettitt_chang ...

这个应该和我编的，计算结果是差不多的。

qinjiancheng01

谢谢楼主分享

zihan6688

很实用的方法，但没法看啊