eof文档
的eof函数给出了固有模式的地图变化和相应的主成分时空数据的时间序列分析。这个函数是专门为3 d matricies海洋表面温度等数据的维度1和2维空间(例如,纬度和经度;经度和纬度;x和y,等等),第三维度代表不同的片或快照的数据。
参见:reof。
内容
语法
eof_maps = eof (A) eof_maps = eof (A, n) eof_maps = eof(…,“面具”,面具)[eof_maps、pc、expvar] = eof (…)
描述
eof_maps = eof (A)计算所有模式的变化一个,在那里一个是一个三维矩阵的第一个二维空间,第三维是时间,数据被认为是等距的。输出eof_maps有相同的尺寸吗一个沿着第三维度,每个地图代表一个可变性重要性顺序模式。
eof_maps = eof (n),只计算第一n模式的变化。对于大型数据集,它的计算速度只计算模式你需要的数量。如果n没有指定,所有计算eof(每个时间片)。
eof_maps = eof(…,“面具”,面具)只执行EOF分析网格细胞为代表的逻辑面具的维度对应维度1和2的一个。此选项提供了防止解决不需要解决的事情,或者让你做分析独立于另一个地区。默认情况下,任何网格细胞一个含有任何nan掩盖了。
[eof_maps、pc、expvar] = eof (…)返回主分量时间序列个人电脑的每一行代表一个不同的模式从1到n和列对应于时间的步骤。例如,电脑(1:)时间序列的第一个varibility(主导)模式。第三输出expvar是方差的百分比用每种模式来解释。请参见下面的说明解释expvar。
一个简单的例子
这里的如何使用一个简单的例子eof函数。适当的EOF分析需要消除趋势和deseasoning数据计算EOF之前,这些步骤中描述下面的教程,但是现在让我们假装这个示例数据集已经准备好进行分析。加载示例数据,然后计算第一个EOF。
这里我策划第一EOF地图使用cmoceanδcolormap (Thyng et al ., 2016)“主”参数,以确保为中心是零。
%样本数据加载:负载pacific_sst.mat%计算第一个EOF海面温度及其%主分量时间序列:(eofmap, pc) = eof (sst, 1);%画出第一个EOF地图:imagescn(经度、纬度、eofmap);轴xy图像从%可选:使用cmocean colormap:cmocean (“δ”,“主”,0)
这是第一个EOF的SST数据集,但由于我们没有季节性周期,第一个EOF主要是季节性变化。作为证据,上面的模式与季节性周期有关,看看对应的主成分时间序列。
图异常(t, pc)轴紧xlim ([datenum (1990年1月1日的)datenum (1995年1月1日的)))datetick (“x”,“keeplimits”)
我看起来很季节性。
教程:从原始气候再分析数据来ENSO, PDO等等。
CDT提供了样本数据集pacific_sst.mat哈德利中心的,这是一个downsampled子集的HadISST海面温度数据集。本教程的最后有一个部分描述如何我原始NetCDF数据导入到Matlab和过程我用来子集。如果你跟随本教程从上到下你应该能够应用EOF分析,任何类似的数据集。
如果你还没有加载示例数据集,现在装入并了解其内容通过检查变量的名称和大小:
负载pacific_sst.mat谁
类属性名称大小字节eofmap 60 x55 26400双纬度60 x1 480双经度55 x1 440海温60双x802 pc 6416双x55x802 21172800双t 802 x1 6416双
我们有一个3 d风场矩阵维度对应纬度x经度的x。什么时间范围,什么是时间步,你问?让我们看看第一个和最后一个日期,和平均时间步骤:
结束datestr (t([1])的意思是(diff (t))
ans = 2×20字符数组的15 - 1月- 1950 12:00:00的15 - 10月- 2016 12:00:00 ans = 30.44
海洋表面平均温度
好吧,这是月度数据,集中在每个月的15日,从1950年到2016年。得到一个数据集的样子,显示平均温度。我使用imagescn自动使NaN值透明,但您可以使用显示亮度图像,pcolor,或任何等价映射工具箱函数如果你喜欢。我也在使用cmocean热colormap (Thyng et al ., 2016):
图imagescn(经度,纬度,意味着(sst, 3));轴xy从cb = colorbar;ylabel (cb、“平均气温C{\保监会}”)cmocean热
全球变暖
全球变暖是真实的吗?的趋势函数让我们容易温度的线性趋势从1950年到2016年。一定要这一趋势乘以10 * 365.25将从每十年度每天度:
imagescn(经度、纬度、10 * 365.25 *趋势(海温,t, 3))轴xy从cb = colorbar;ylabel (cb、每十年的温度趋势{\保监会}C”)cmocean (“平衡”,“主”)
消除全球变暖的信号
全球变暖趋势是有趣,但EOF分析都是瑞士,没有长期趋势,所以我们必须删除的趋势detrend3:
海温= detrend3 (sst、t);
消除季节性周期
如果你再绘制温度趋势,你会发现这都是减少到零,或许几eps数值噪音。现在这是一个SST数据集,即使安东尼·瓦会批准的。
我们现在去趋势对海温数据集(也删除了的意思),但它仍然包含相当多的季节性变化,应该删除之前EOF分析季节性信号。因为我们不感兴趣
海温= deseason (sst、t);
所以现在我们风场数据集已经去趋势,意味着移除和季节性周期移除。剩下的风场异常——事情改变,但不是长期趋势或短期的年度周期。这是其余的方差风场异常数据集:
图imagescn(经度、纬度、var (sst, [], 3));轴xy从colorbar标题(温度的差异)colormap(飞机)%飞机是不可原谅的,除了当重建旧的情节caxis ([0 1])
上面的地图线路与图2的一个很好Messie和查韦斯(2011)告诉我们,我们在正确的轨道上。
计算eof
EOF分析不仅可以让我们了解事情各有不同,但多长时间,什么地区倾向于彼此在一起或的阶段不同。与我们的去趋势,deseasoned风场数据集,这个简单的EOF分析eof功能:
[eof_maps, pc, expv] = eof (sst);%画第一个模式:图显示亮度图像(经度、纬度、eof_maps(:,: 1)轴xy图像cmocean (“卷”,“主”)标题“第一个EOF模式!”
特征向量分析有一个有趣的行为可以产生EOF地图是积极的还是消极的,和解决方案可以不同的每次使用相同的输入。金宝搏官方网站积极的和消极的解决方案同样有效,认为鼓面的振金宝搏官方网站动模式,一些地区的鼓头上升下降,而其他区域,然后开关,同样海温变化的特征值的解决方案可能是积极的还是消极的。唯一重要的是,当我们重构时间序列的EOF的解决方案,我们每一个EOF映射乘以相应的主成分(个人电脑)。
的eof函数写入产生一致的结果每次运行相同的数据,但别担心如果解决不匹配的符号的符号别人的结果,如果他们使用了不同的程序来计算eof——这就意味着他们的程序反号的解决方案,这是非常好。
正如EOF地图可以有积极或消极的解决方案,都是同样有效的,有一些灵活性EOF地图的大小是如何显金宝搏官方网站示的。你可以把一个EOF地图的大小乘以任何你想要的价值,只要你相应的主成分时间序列除以相同的值。让我们看看前三的时间序列模式的可变性,策划subsubplot:
图subsubplot(3、1、1)的阴谋(t, pc (1:))从轴紧ylabel“pc1”标题前三个主成分的subsubplot(3, 1, 2)情节(t, pc(2:))盒子从轴紧集(gca),“yaxislocation”,“对”)ylabel“pc2”subsubplot(3、1, 3)情节(t, pc (3:))从轴紧ylabel“生物”datetick (“x”,“keeplimits”)
可选的主成分和EOF地图缩放
这些主成分时间序列是好方法,但是有些人更喜欢每个时间序列跨度规模所需的范围。查看图5的Messie和查韦斯(2011)看来,他们选择规模每个主分量时间序列,使其涵盖范围1到1。让我们做同样的事情,每个主分量时间序列除以它的最大价值和不要忘记乘以相应的EOF映射相同的值:
为k = 1:尺寸(pc, 1)%的时间序列中的最大值每个主要组件:maxval = max (abs (pc (k,:)));%将时间序列的最大值:电脑(k,:) = pc (k,:) / maxval;%乘以相应的EOF地图:eof_maps (:,:, k) = eof_maps (:,:, k) * maxval;结束
厄尔尼诺南方涛动(ENSO)时间序列
第一个模式去趋势,deseasoned太平洋与ENSO assoiciated。我们可以再次画出时间序列作为一个简单的情节,但异常的情节往往填充。让我们使用异常情节第一个模式,乘以1与图5的迹象Messie和查韦斯(2011)。
图(“pos”,100 100 600 250)异常(t - pc (1:)“topcolor”rgb (“泡泡糖”),…“bottomcolor”rgb (“玉黍螺蓝”))%第一主成分是enso盒子从轴紧datetick (“x”,“keeplimits”)文本((724316 729713 736290)。95 .99结果),“厄尔尼诺”,“水平的”,“中心”)
果然,一些最强的厄尔尼诺事件记录发生在1982 - 1983,1997 - 1998,2014 - 2016。更多的方法来调查ENSO,查看enso函数。
ENSO在频域
有时我们听到厄尔尼诺现象的特征频率每隔五年,或五到七年,有时你听到每两到七年。很难看到在时间序列,所以我们把第一个主成分在频域plotpsd,指定每年12个样品的采样频率,标注在日志x轴,在单位的x值λ(年)而不是频率:
图plotpsd (pc(1:), 12日“计算lnx”,“λ”)包含“周期(年)”集(gca),“xtick”33,[1:7])
如您所见,ENSO信号没有大幅共振频率定义,但整个两到七年的能量范围。我也贴上33年周期性因为奈奎斯特的这个特定的数据集——任何能量比尼奎斯特更长一段(或接近)应该被认为是垃圾。
可变性的地图
eof不仅仅是对时间序列——他们对空间的变化模式。每个模式都有一个特征变化模式就像一个鼓的不同振动模式。在任何给定的时间,不同的模式可以总结为创建一个总温度异常的照片。的正交经验正交函数的一部分意味着每个模式倾向于做自己的事情,独立于其他模式。让我们看看第一个6模式的再现图4Messie和查韦斯(2011)。我的一些模式乘以- 1因为我想匹配他们的迹象,记住,我们可以这样做。
s = [1 1 1 1 1 1];%(签署乘数匹配Messie和查韦斯2011)图(“pos”,100 100 500 700)为k = 1:6次要情节(3 2 k) imagescn(经度、纬度、eof_maps (:,:, k) * (k));轴xy从标题([“模式”num2str (k),“(”num2str (expv (k),' % 0.1 f '),“%)”2)caxis ([2])结束colormap飞机
百分比差异解释为每种模式不匹配Messie查韦斯,因为我们用更短的时间序列相比,我们也使用一个空间世界数据的子集。尽管如此,普遍认同模式。
的飞机colormap不是完全相同的一个使用Messie查韦斯,这解释了为什么一些模式高于Messie和查韦斯看起来略有不同。但由于我们谈论colormaps,彩虹是非常不擅长代表数值数据(Thyng et al, 2016年)。科学也是一种耻辱,我们不能完全复制情节不知道什么颜色是用于发布的版本,但是我跑题了…
鉴于这些地图代表异常,他们应该由发散colormap表示出同等重量的0。设的colormaps次要情节在这个图更平衡:
colormap (gcf cmocean (“平衡”))
制作电影的海温变化eof
在任何给定的时间,海面温度异常与ENSO相关的快照可以通过绘制的地图模式1所示,乘以相应的主成分(向量电脑(1:))。同样,你可以得到一幅全球海面温度异常在给定的时间总结所有的EOF地图,每个乘以相应的主成分。这样我们可以构建一个越来越多的海温异常在我们完整的电影包括越来越多的更多的模式变化。
例如,海温异常与前三的地图模式的可变性,指定的时间,可以通过求和eof每个模式的地图,乘以相应的pc值。你可以手动做加法,说,1990年是这样的:
%指数电影的开始和结束日期:startind =找到(t > = datenum(1990年1月1日),1,'第一次');endind =找到(t < = datenum(1999年12月31日),1,'最后');
%的地图海温异常在开始前三个模式:地图= eof_maps (:,: 1) * pc (1, startind) +…% 1模式,1990年1月eof_maps (:,: 2) * pc (2, startind) +…%模式2,1990年1月eof_maps (:,:, 3) * pc (3 startind);%模式3,1990年1月
然而,这是比必要的更复杂,因为reof函数是用来为我们做加法。让前三的海温异常时间序列模式是这样的:
sst_f = reof (eof_maps、pc、1:3);
现在1990年1月绘制地图作为第一帧电影:
ind_1990s = 481:3:600;%(每三值减少gif大小)
图h = imagescn (sst_f (:,:, ind_1990s (1)));caxis (2 [2]) cmocean平衡cb = colorbar;ylabel (cb、温度异常(模式1 - 3))标题(datestr (t (ind_1990s (1)),“yyyy”))
gif (“SSTs_1990s.gif”、“框架”,gcf) %第一帧写道
k = 2:长度(ind_1990s) h。CData = sst_f (:,:, ind_1990s (k));标题(datestr (t (ind_1990s (k),“yyyy”)) gif %将这个框架添加到gif
你可能注意到的第一件事是,1990年代的海温异常时间序列是由ENSO,并检查出1997 - 1998信号!难怪这是这样的一个热门话题新闻那一年。尽管如此,重要的是要记住,上面的电影并不是一个完整的重建对海温异常,而是只有前三种模式,占
sum (expv (1:3))
ans = 51.88
…超过一半的总方差的SST数据集。重建的绝对温度场而不是异常的前三个模式,需要包括所有的EOF地图,你还需要添加在平均SST地图,趋势和季节性周期。
我得到了样本数据如何
上面所示的示例数据集来自HadISST哈德利中心找到在这里(雷纳et al ., 2003),在超过200 MB。如果你想执行同样的在世界不同地区的分析,您可以下载HadISST_sst。数控数据集,将其导入Matlab。将采样数据集或构造子集是你:
%加载完整的SST数据集:lat =双(ncread (“HadISST_sst.nc”、“纬度”));朗=双(ncread (“HadISST_sst.nc”、“经度”));双(t = ncread (“HadISST_sst.nc”,“时间”)+ datenum(1870年,1,0));海温= ncread (“HadISST_sst.nc”、“场”);
%季度样本数据集的大小,我对每一个网格点大致downsample:海温(sst(1:2:结束,1:2:最终,);lat =纬度(1:2:结束);朗=经度(1:2:结束);
%,以进一步降低尺寸,我剪的背阔肌和经度,保持只有1950后数据:行=朗< -70;朗=朗(行);关口= lat > = -60 & lat < = 60;lat =纬度(关口);* = t > = datenum(1950年1月1日);t = t(次);海温(sst(行、关口倍);海温(sst < -50) =南;
%我发现它更容易重新排列纬度x经度x:海温=排列(风场(2 1 3));
%保存样本数据:保存(“PacOcean.mat”、“纬度”,“朗”、“t”,“场”)
引用
Monique Messie,旧金山查韦斯。“全球海洋表面温度变化模式与区域气候指标。”Journal of Climate 24.16 (2011): 4314-4331.jcli3941.1 doi: 10.1175/2011。
雷纳:A。,Parker, D. E., Horton, E. B., Folland, C. K., Alexander, L. V., Rowell, D. P., Kent, E. C., Kaplan, A. (2003). Global analyses of sea surface temperature, sea ice, and night marine air temperature since the late nineteenth century. J. Geophys. Res.Vol. 108, No. D14, 4407jd002670 doi: 10.1029/2002。
Thyng K.M.,C.A. Greene, R.D. Hetland, H.M. Zimmerle, and S.F. DiMarco. 2016. True colors of oceanography: Guidelines for effective and accurate colormap selection. Oceanography 29(3):9-13,doi: 10.5670 / oceanog.2016.66。
作者信息
的eof函数是由乍得a。格林德克萨斯大学的地球物理研究所(UTIG) 2017年1月,但靠Guillame迷宫的caleof函数从他PCATool的贡献。本教程是乍得格林的帮助下写的Kaustubh Thirumalai。