当前位置:文档之家› 北太平洋SST(1980-1992)的EOF分解及初步分析

北太平洋SST(1980-1992)的EOF分解及初步分析

北太平洋SST(1980-1992)的EOF分解及初步分析数据来源:本文利用来自英国气象局Hadley气候预测和研究中心的HadISST海温资料,对北太平洋(20N-80N 130E-90W)1980至1992年156个月的SST进行研究,。

本文选取下载了1870年1月至2011年3月的全球SST数据,数据格式为*.netcdf,分辨率为1°。

以下给出了数据下载的途径和链接:①/hadobs/hadisst/index.html点击页面中心处蓝色“main data page”:图1②/hadobs/hadisst/data/download.html选择SST数据下载:图2数据处理:本文选取下载了1870年1月至2011年3月的全球SST数据,数据格式为*.netcdf,数据名为HadISST_sst.nc,截取北太平洋(20N-80N 130E-90W)1980年1月至1992年12月156个月份的SST数据进行经验正交函数(Empirical Othorgnal Function)分解,简记为EOF分解,得到该区域该时段的海温时空特征。

在编写Matlab程序过程中,应特别注意:⑴剔除与其它站点相关系数小的站点的数据~简单的认为剔除陆地和冬季结冰点的数据;⑵求距平值的协方差矩阵时,要进行逐月平均求距平,而不能是156个月的平均值,否则会导致第1模态的方差贡献率很大;⑶当变量数m远大于观测样本数n时,导致协方差矩阵mRm=(nXm’)*(nXm)的阶数较大,可先求(nXm)*(nXm’)矩阵的特征值和特征向量,再求(nXm’)*(nXm)的特征值和特征向量,这叫做时空转换;⑷M文件编写时要尽量减少循环量,提高运算速度;⑸EOF分析能够有效地体现物理场主要信息,保留次要信息,并排除外来的随机干扰。

数据分析:用主成分分析(Principal Component Analysis)的方法,即PCA对结果进行分析:主成分分析是多元统计分析中一个非常重要的内容,它是一种从多个变量化为少数变量的统计方法。

由于多个变量之间是相互影响的,它们之间的关系是非常复杂的,为简化分析又不损失信息,并提取它们之间相互关系的主要特征,主成分分析利用多个变量之间的相互关系构造一些新变量,这些新变量不仅能综合反映原来多个变量的信息,而且彼此之间是相互独立的,同时是按方差贡献大小排列的。

方差贡献率小的变量通常规律性很差,其实际物理意义也不清晰,因此在实际分析过程中常常视为误差量或噪声而忽略,只取方差贡献率大的变量来研究,从而达到降维分析的目的。

通过对相应数据处理分析,前13个主成分的累积方差贡献率占总方差的0.810567965759345,对前13个主模态的方差贡献率和累积方差贡献率列表格:图3现仅列出北太平洋前5个主模态的空间分布填色图及时间序列,并对第1和第3模态进行分析:北太平洋第1模态填色图及时间序列图4 北太平洋第2模态填色图及时间序列图5 北太平洋第3模态填色图及时间序列图6 北太平洋第4模态填色图及时间序列图7北太平洋第5模态填色图及时间序列图8 对第1主模态进行分析:北太平洋洋流图9图10(a)图10(b)图10(a)是第1模态空间分布型,它解释海温场总方差的22.9%,此型在北太平洋西、中部被一片强负值控制,负中心约在170°E,40°N和150°W,40°N附近,而北太平洋东部和北美沿岸为较弱的正值区,说明北太平洋西、中部海温与东部海温是反相关关系,负区与北太平洋西风漂流区(如图9)吻合。

由美国海洋学家斯蒂文•黑尔于1996年发现的太平洋年代际振荡(PDO)被科学研究的初步结果表明其与厄尔尼诺(El Niño)和拉尼娜(La Nina)现象有着极其密切的关系。

该型可以反映和PDO有关的大尺度分布特征,因此这种分布型是全球海洋与大气相互作用的一个重要组成部分,它是北太平洋海温非季节变化的最重要的型式。

Monthly values for the PDO index: 1900—January2008 图11图10(c)观察发现图11(Monthly values for the PDO index: 1900—January2008)1980年至1992年时间段的指数和第1模态的时间序列图10(c)有很好的对应关系,可以验证北太平洋海表面温度第1模态空间分布型确实与PDO有很强的相关性。

资料显示,近100多年来,PDO已出现了两个完整的周期:第一周期的“冷位相”发生于1890年至1924年,而1925年至1946年为“暖位相”;第二周期的“冷位相”出现于1947年至1976年,1977年至90年代后期为“暖位相”。

当PDO现象以“暖位相”形式出现时,北美大陆附近海面的水温就会异常升高,而北太平洋洋面温度却异常下降。

并且,在20-30年的冷、暖位相中,会存在短期的反向指数。

由时间序列可知:1980年至1988年底,时间序列指数基本为正值,说明图10(a)中蓝色区域海表面温度低于红色区域海表面温度,即北太平洋西、中部海温低于东部海温;1989年初至1992年初,时间序列指数为负值,说明10(a)中蓝色区域海表面温度高于红色区域海表面温度,即北太平洋西、中部海温高于东部海温。

此分析和历史资料相吻合。

对第3主模态进行分析:图12(a)图11(a)是第3模态空间分布型,它解释海温场总方差的10.0%,此型在北太平洋阿留申群岛南部被一片强正值控制,正值中心约在150°W,40°N附近,而其西南部和日本海海域为一片较强的负区与之相互补偿,其东南部北美洲沿岸为较弱的负值区,说明北太平洋中部海表面温度与东、西两侧海表面温度呈反相关。

正、负中心基本上与北太平洋大气活动中心(阿留申低压与西太平洋高压)对应。

表明正、负区域是海气相互作用最活跃的区域。

查阅资料可知,这种模态的变化,对我国的天气与气候有明显的影响。

图12(b)结合第3模态的时间序列可知:1980年至1983年6月、1987年6月至1988年底及1990年,时间序列指数基本为负值,说明图10(a)中蓝色区域海表面温度高于红色区域海表面温度,即北太平洋西、中部海温高于东部海温;1983年7月至1987年初、1989年及1991年初至1992年底,时间序列指数为正值,说明10(a)中蓝色区域海表面温度低于红色区域海表面温度,即北太平洋西、中部海温高于东部海温。

总结:本文通过对北太平洋1980年至1992年SST数据处理、EOF分解和初步分析,掌握了主成分分析、时空转换的原理和方法,提高了Matlab的编程和绘图技巧并对主模态和时间序列的分析进行了联系。

在此过程中克服了很多困难,受益匪浅。

但和很多同学相比,仍然有很大差距,将会更加虚心请教,刻苦钻研,以取得不断进步。

参考文献:[1]左军成.海洋水文环境要素的分析方法和预报[2]胡基福.气象统计原理与方法[3]黄嘉佑.气象统计分析与预报方法[4]杜凌.海洋要素计算(2011)PPT[5]姜霞.气象统计原理与方法(2011)PPTM文件:clear;clc;close alladdress='E:\oceanelement\HadISST_sst.nc';fid=netcdf.open(address,'NC_NOWRITE');sstid = netcdf.inqVarid(fid,'sst');sst=netcdf.getVar(fid,sstid); % 读取nc格式数据%**************************************************************************sst1=sst(1:90,11:70,1320:1475); % 选取所需要区域的数据sst2=sst(311:360,11:70,1320:1475);sst3=zeros(140,60,156);sst3(90:-1:1,1:60,1:156)=sst1;sst3(140:-1:91,1:60,1:156)=sst2;sst=sst3;%**************************************************************************sst_area1=zeros(156,8400); % zeros全零数组for i=1:156;squ=squeeze(sst(:,:,i)); % 执行该指令后sst数据转换为二维数组sst_area1(i,:)=reshape(squ,1,8400); % 将数据转变为二维end%**************************************************************************% 剔除与其它站点相关系数小的站点的数据~简单的认为剔除陆地和冬季结冰点的数据sst_area1(sst_area1<-10000)=NaN; % 陆地和冰点的填充值为-1.00000001504747e+30~将此值定义为NaN% i=1;% for j=1:8400% if sst_area1(i,j)==-1.79999995231628% sst_area1(i,j)=NaN; % 冰点的填充值为-1.79999995231628 % i=i+1;% end% endsst_nan=isnan(sst_area1);i=0;for j=1:8400if sum(sst_nan(:,j))==0;i=i+1;sst_region(:,i)=sst_area1(:,j);endend%**************************************************************************% 求距平~注意季节的变换X=zeros(size(sst_region)); % 学者给的程序for i=1:12X(i:12:144+i,:)=sst_region(i:12:end,:) - repmat( mean(sst_region(i:12:end,:),1) , size(sst_region(i:12:end,:),1), 1);end% s=0; % 自己的~很多循环计算慢% for k=1:12;% for j=1:4344;% for i=k:12:144+k% s=s+sst_region(i,j);% end% s_ave=s/13;% for i=k:12:144+k% ave_region(i,j)=s_ave;% end% s=0;% end% end% X=sst_region-ave_region; % X1是156*8400的矩阵%**************************************************************************% 由于变量数(8400)>观测样本数~协方差矩阵的阶数较大~因此可用时空转换的方法提高计算速度R=X*X'; % 协方差矩阵R=X*X'是8400*8400的方阵~现定义矩阵R=X'*X是156*156的矩阵[v,d]=eig(R); % 进行EOF分解~因为X'*X 与X*X'的秩相同所以特征值相同~d为x的特征值组成的对角阵~v为X*X'的特征向量~v=fliplr(v); % 矩阵作左右翻转d=rot90(d,2); % 矩阵上下翻转后再左右翻转(查看生成的对角阵是由小到大排列的~此指令可使其由大到小排列~fliplr(flipud(d))=rot90(d,2))diagonal=diag(d);spacef=X'*v;for i=1:156;spacef(:,i)=spacef(:,i)/sqrt(diagonal(i)); % 空间本征函数endtimef=X*spacef; % 时间本征函数sum_d=sum(diagonal);count=0;for i=1:156;count=count+diagonal(i);G1(i)=count/sum_d; % G1(i)是累积方差贡献率endfor i=1:156;G2(i)=diagonal(i)/sum_d; % G2(i)是方差贡献率end%**************************************************************************% 将删去的陆地与冰点的填充值补回sst_area2=zeros(156,8400);sst_area2(:,:)=NaN;i=1;spacef2=spacef';for j=1:8400if sum(sst_nan(:,j))==0;sst_area2(:,j)=spacef2(:,i);i=i+1;endendsst_area3=sst_area2';%************************************************************************** % 画图% subplot(2,1,2) % 将绘图窗口划分为2*1个子窗口,并在第2个子窗口中绘图x=1:156;plot(x,timef(:,6),'g');ylim([ -80 80 ]);% xlabel('1980-1992年156个月','fontsize',15,'fontname','隶书')ylabel('INDEX','fontsize',12,'fontname','黑体')set(gca,'xtick',(1:6:162))set(gca,'xticklabel',{'1980','','1981','','1982','','1983','','1984','','1985','','1986','','1987','','1988','','1989','','1990','','1991','','1992','','1993' })title('北太平洋第6模态1980至1992年SST时间序列', 'color', 'k','fontsize',15,'fontname','幼圆')grid onhold off% %subplot(2,1,1)% lon=[130.5:269.5];% lat=[20.5:79.5];% m_proj('Equidistant Cylindrical','lat',[20.5 79.5],'lon',[130.5 269.5]);%m_contourf(lon,lat,rot90(reshape(sst_area3(:,6),140,60)',2),30,'linestyle','non e');% colorbar% m_coast('patch',[.95 .95 .95]);% m_coast('color','k');% m_grid('linestyle','none','tickdir','out','linewidth',1.5);% % xlabel('longitude','fontsize',15,'fontname','comic sans ms');% % ylabel('latitude','fontsize',15,'fontname','comic sans ms');% title(['北太平洋第6模态填色图'],'fontsize',15,'fontname','幼圆');lon=[130.5:269.5];lat=[20.5:79.5];m_proj('lambert','long',[130.5 269.5],'lat',[20.5 79.5]);m_contourf(lon,lat,rot90(reshape(sst_area3(:,1),140,60)',2));% colorbar;m_coast('patch',[1 .85 .7]);m_elev('contourf',[500:500:6000]);m_grid('box','fancy','tickdir','in');colormap(flipud(copper));% xlabel('longitude','fontsize',15,'fontname','幼圆');% ylabel('latitude','fontsize',15,'fontname','幼圆');title(['北太平洋第1模态填色图'],'fontsize',15,'fontname','幼圆');。

相关主题