用Matlab进行MK检验Word文档下载推荐.doc
《用Matlab进行MK检验Word文档下载推荐.doc》由会员分享,可在线阅读,更多相关《用Matlab进行MK检验Word文档下载推荐.doc(5页珍藏版)》请在冰豆网上搜索。
1<
j<
i<
n,β所代表的是斜率,正的则是上升,负的代表下降,值的大大小代表趋势明显与否。
Mann-Kendall检验如下所示:
零假设H0:
β=0
当,拒绝H0假设。
为标准正态方差,α为显著性检验水平。
2、Mann-kendall突变检验
气候系统变化是一个不稳定且不连续的变化过程,而检验其变化的常用方法之一就是Mann-kendall突变检验方法[13],该方法对于变化要素从一个相对稳定状态变化到另一个状态的变化检验非常有效。
且广泛应用于水文,气候,化学,矿物成分检验等各个方面。
Mann-kendall突变检验方法如下:
对于具有n个样本量的时间序列x,构造一个秩序列:
(1.6)
其中
()(1.7)
可见,秩序列是第时刻数值大于j时刻数值个数的累计数。
在时间序列随机独立的假定下,定义统计量:
(1.8)
其中UF1=0,,是累计数的均值和方差,在相互独立,且有相同连续分布时,它们可由下式算出:
(1.9)
(1.10)
为标准正态分布,它是按时间序列x顺序计算出的统计量序列,给定显著性水平a,查正态分布表,若,则表明序列存在明显的趋势变化。
把此方法引用到时间序列的逆序序列中,按xn,xn-1,…x1,再重复上述过程,同时使UFk=-UBk,k=n,n-1,…1,UB=0。
给定显著性水平α,将UFk和UBk两个统计量曲线和显著性水平线绘在同一个图上,若UFk和UBk的值大于0,则表明序列呈上升趋势,小于0则呈下降趋势。
当超过临界直线时,表明上升或下降趋势显著,超过临界线的范围确定为突变的时间区域。
如果UFk和UBk两条曲线出现交点,且交点在临界线之间,那么交点对应的时刻便是突变开始的时间。
二、M-K程序介绍
1、M-K趋势检验
function[slope,zc,za,sign]=MannKendall(x)
%计算S%
s=0;
len=size(x,2);
form=1:
len-1
forn=m+1:
len
ifx(n)>
x(m)
s=s+1;
elseifx(n)==x(m)
s=s+0;
else
s=s-1;
end
end
end
%计算vars和e%
vars=len*(len-1)*(2*len+5)/18;
%计算zc%
ifs>
zc=(s-1)/sqrt(vars);
else
zc=(s+1)/sqrt(vars);
end
%计算za
za=var(x);
sign=0;
zc1=abs(zc);
ifzc1>
=za
sign=1;
else
sign=0;
end
%计算倾斜度
ndash=len*(len-1)/2;
slope1=zeros(ndash,1);
m=1;
fork=1:
len-1,
forj=k+1:
len,
slope1(m)=(x(j)-x(k))/(j-k);
%
m=m+1;
end;
slope=median(slope1);
2、M-K突变检验
clc;
y=data(1:
50,6);
%输入数据
n=length(y);
Sk=zeros(size(y));
UFk=zeros(size(y));
s1=0;
fori=2:
n%
forj=1:
i
ify(i)>
y(j)
s1=s1+1;
s1=s1+0;
%()
end;
end;
Sk(i)=s1;
%
E=i*(i-1)/4;
Var=i*(i-1)*(2*i+5)/72;
UFk(i)=(Sk(i)-E)/sqrt(Var);
%
y2=zeros(size(y));
Sk2=zeros(size(y));
UBk=zeros(size(y));
s2=0;
fori=1:
n
y2(i)=y(n-i+1);
forj=1:
ify2(i)>
y2(j)
s2=s2+1;
s2=s2+0;
end;
Sk2(i)=s2;
E=i*(i-1)/4;
%Sk2(i)的均值
Var=i*(i-1)*(2*i+5)/72;
%Sk2(i)的方差%
UBk(i)=0-(Sk2(i)-E)/sqrt(Var);
UBk2=zeros(size(y));
UBk2(i)=UBk(n-i+1);
year=1961:
2010;
year=year'
;
xlswrite('
D:
\test.xls'
year,'
Sheet1'
'
A1'
);
UFk,'
B1'
UBk2,'
C1'