分享一下滑动t检验检测突变的代码，根据魏凤英书P58修改

orchiq

%滑动t检验
%程序来源-http://bbs.06climate.com/forum.php?mod=viewthread&tid=41568
%根据魏凤英书P58修改，输入时间序列。
clc;
clear;
load('test.mat');
%以上数据中有data变量
%data变量共二维，结构为2*2000
%data(1,:)储存了1,2,3...2000
%data(2,:)储存了每年对应的数值

years=length(data);%数据的时长
data=data';
x=data(:,1);%year
y=data(:,2);%data
t=zeros(years,1);
checkline99=zeros(years,1);
checkline999=zeros(years,1);

%MMT方法检验%
%先设定子序列的长度
sublen=20;
n1=sublen;n2=sublen;
m=sqrt((1.0/n1)+(1.0/n2));

%开始计算每个基准点的统计量序列t
for i=sublen:years-sublen
    x1avg=mean(y(i-sublen+1:i));
    x2avg=mean(y(i+1:i+sublen));
    x1var=var(y(i-sublen+1:i));  
    x2var=var(y(i+1:i+sublen));
    s=sqrt((n1*x1var+n2*x2var)/(n1+n2-2));
    t(i-sublen+1)=(x2avg-x1avg)/(s*m);  
end

plot(x,t,'r-');
checkline99(:)=2.62;                   %n=sublen*2-2时，显著性检验标准,0.01

hold on;
plot(x,checkline99,'k:');
hold on;
plot(x,-checkline99,'k:');
axis([sublen,years-sublen,round(min(t))-1,round(max(t))+1]);
xlabel('year');ylabel('t');title('滑动t检验');

orchiq

此楼附件为程序中用到的数据。

小温馨岁月

谢谢楼主分享

青椒酱

感谢楼主的无私分享

zihan6688

谢谢楼主分享，要是写成函数会更好

垫底辣妹

请问子序列的长度怎么设定？

lzwbeifeng

初学者，进步中

诗遇上歌

正好做这个

My93131

楼主你好，倒数第二行中，round是什么意思呢

shower107

t(i-sublen+1)=(x2avg-x1avg)/(s*m);  要改为：t(i+1)=(x2avg-x1avg)/(s*m); 否则断点对不上