大气统计作业Word文档下载推荐.doc
《大气统计作业Word文档下载推荐.doc》由会员分享，可在线阅读，更多相关《大气统计作业Word文档下载推荐.doc（19页珍藏版）》请在冰豆网上搜索。
%存储Y中各元素的排行
%计算Xrank中的各个值
fori=1:
N
cont1=1;
%记录大于特定元素的元素个数
cont2=-1;
%记录与特定元素相同的元素个数
forj=1:
ifX（i）<
X（j）
cont1=cont1+1;
elseifX（i）==X（j）
cont2=cont2+1;
end
end
Xrank（i）=cont1+mean（[0:
cont2]）;
%计算Yrank中的各个值
forj=1:
ifY（i）<
Y（j）
cont1=cont1+1;
elseifY（i）==Y（j）
end
Yrank（i）=cont1+mean（[0:
%利用差分等级（或排行）序列计算斯皮尔曼等级相关系数
fenzi=6*sum（（Xrank-Yrank）.^2）;
fenmu=N*（N^2-1）;
coeff=1-fenzi/fenmu;
end%函数mySpearman结束}
计算相关系数：
A=load（'
RBJsum.txt'
B=load（'
TBJsum.txt'
A1=A（:
2）;
B1=B（:
%相关系数
R=corrcoef（A1,B1）;
R=
1.0000-0.4332
-0.43321.0000
%秩相关系数
X=A1;
Y=B1;
coeff=mySpearman（X,Y）;
coeff
coeff=
-0.4693
（2）显著性检验：
相关系数自由度n=52-2
查表可知，rc=0.273<
0.4332
计算得到的相关系数大于0.273，则在显著水平0.05是显著的。
秩相关系数：
同理，对秩相关系数检验可得：
rc=0.273<
0.4693，因此秩相关系数在显著水平0.05上是显著的，通过检验。
（2）将降水量作为预报因子，气温作为预报量，试给出回归方程，并说明降水每增加100mm，气温大致下降多少º
C？
legend（'
降水量-气温散点图'
scatter（x,y）
xlabel（'
降水量'
ylabel（'
气温'
title（'
%建立回归方程：
stepwise（x,y）；
b=25.995
-0.0019776
x1=[ones（52,1）,x];
y1=x1*b;
plot（x,y,'
*'
x1,y1）;
拟合曲线及散点图'
那么降水每增加100mm，气温大致下降0.19776º
C
二、EOF分析
（1）进行EOF分析，给出前十个模态的方差解释率，并根据North准则（绘出ErrorBar图），说明需要选择前几个模态进行分析？
fid=fopen（'
ssta.dat'
'
r'
a=fread（fid,’float’）;
%a=reshape（a,56*11,336）;
b=a’;
b=zeros（336,616）;
c=zeros（1,616）;
j=0;
forn=1:
336;
c=a（j+1:
j+616）;
b（n,:
）=reshape（c,1,616）;
j=j+616;
%b1=fread（fid,[616336]）;
b=b’;
end
b=b'
;
C=b*b'
/616;
[EOF,E]=eig（C）;
PC=EOF'
*b;
E=fliplr（flipud（E））;
lambda=diag（E）;
EOF=fliplr（EOF）;
PC=flipud（PC）;
sum_d=sum（lambda）;
count=0;
fori=1:
count=count+lambda（i）;
G1（i）=count/sum_d;
%累计方差贡献率
G2（i）=lambda（i）/sum_d;
%方差贡献率
方差贡献率
累计方差贡献率
第一模态
0.716516015743107
第二模态
0.147473083981738
0.863989099724845
第三模态
0.050392162180530
0.914381261905375
第四模态
0.024518677922468
0.938899939827844
第五模态
0.014502739245485
0.953402679073329
第六模态
0.012702830971315
0.966105510044644
第七模态
0.007148842058117
0.973254352102761
第八模态
0.006019078625006
0.979273430727767
第九模态
0.004340892184153
0.983614322911921
第十模态
0.003048921846347
0.986663244758267
%North准则
f=lambda（1:
10）;
errorf=f*sqrt（2/（336-2））;
number=1:
1:
10;
errorbar（number,f,errorf）;
number'
eigenvalue'
northerrorbar'
只有对第一、二、三、四、五模态进行分析。
因为第五模态后误差重合了，无法区分开来。
（2）给出需要分析的前几个模态的空间分布和时间序列。
x=1:
plot（x,PC（1,:
）,'
1979-2006年336个月'
INDEX'
赤道中东太平洋第一模态1979-2006年的SST时间序列'
gridon;
lon=[160:
2:
270];
lat=[-10:
10];
m_proj（'
EquidistantCylindrical'
lat'
[-1010],'
lon'
[160270]）;
m_contourf（lon,lat,rot90（reshape（EOF（:
1）,56,11）'
2）,'
linestyle'
none'
colorbar;
m_grid（'
box'
fancy'
tickdir'
in'
longitude'
fontsize'
15,'
fontname'
comicsansms'
latitude'
北太平洋第1模态填色图'
幼圆'
%MATLAB安装m_map绘图工具才能调用以上程序，如果没有，直接调用eof1=reshape（eof（:
n）,56,11）;
eof1=eof1'
contour（eof1）;
n表示不同的模态
三、聚类分析
x=
1.0e+003*
0.01621.42902.0000-0.00820.0062
0.01570.97002.2090-0.02060.0019
0.01631.26002.0850-0.01730.0028
0.01721.44201.7260-0.00950.0046
0.01881.87401.7090-0.00490.0080
0.01791.69801.8480-0.00450.0075
0.01630.97601.2390-0.00460.0056
x1=zscore（x）;
y1=pdist（x2）;
y1=pdist（x1）;
y1=pdist（x1,'
mahalanobis'
z1=linkage（y1）;
c1=cophenet（z1,y1）;
T=cluster（z1,6）;
H=dendrogram（z1）;
c1
c1=
0.5540
T
T=
3
1
2
4
5
6
四、功率谱分析
（1）太阳黑子数的离散功率谱，绘出谱图。
sunspot=importdata（'
sunspot.dat'
year=sunspot（:
1）;
wolfer=sunspot（:
figure
plot（year,wolfer）
xlabel（'
Years'
SunspotData'
）
plot（year（1:
50）,wolfer（1:
50）,'
b.-'
Atthefirst50years'
figure
n=length（Y）;
power=abs（Y（1:
n/2））.^2;
nyquist=1/2;
%取采样频率为0.5；
freq=（1:
n/2）/
升级会员 