第五章 主成分分析.

1、第五章 主成分分析主成分的概念主成分的导出与计算经验正交函数分解第五章1 主成分概念两个变量的主成分的概念两个变量的主成分的概念设研究对象为两两个变量个变量x1和x2，两个变量存在相关性两个变量存在相关性，共观测了n次， 样本数据排列为X：1112121222.=.nnxxxxxxX分别以x1和x2为横纵坐标轴，把n个点绘成点聚图(黄色圆点)：这n个圆点表征了x1和x2两个变量两个变量的综合信息。010.900.91NN 这里，设两个变量都来自服从二维正态分布的总体 ：（，），总体协方差为0.9，自方差为1.-202-202x1x2由二维正态分布的概率密度函数可知，这些散点的密集程度可用一组同

2、心椭圆来表示（同心椭圆即二维概率密度函数值的等值线）如果以这条红线红线为横坐标轴，那么各点在红线上的横坐标值（记为y1 ）就可以在很大程度上在很大程度上反映x1和x2的变化信息。-202-202x1x2即： x1和x2两个变量两个变量的信息，可以近似地用1个新变量个新变量 y1来表达，这样就降低了变量的维度，使问题得到简化。观察散点图可以发现，n个点基本都散布在椭圆的长长轴轴（即右图的红色直线）红色直线）附近，y1可以反映x1和x2随时间变化的大部分大部分信息，若要想表达散点的全部全部信息，需再考虑以椭圆的短轴短轴作为纵坐标，获得纵坐标值的信息，记作y2.qy2的离散程度（方差）要远远小于y1

3、的离散程度（方差），qy2与y1相互独立相互独立(相关系数为0)。-202-202x1x2y1y2问题： y1与y2与原变量x1、x2有何关系？该问题即“将坐标轴旋转一定角度后，新坐标与原坐标的关系”。已知A点在原坐标系下，坐标为(x1, x2)；那么，在新坐标系下，坐标(y1，y2)该如何用x1和x2表示？1x2x2x1x2211( ,),A():xyxy1122cossinsincosyxyx写成矩阵形式，即：211212cos+sinsin+cos-=yyxxxx根据右图可知，新坐新坐标标与原坐标原坐标之间的关系为： 111212122122cos-sin=,=,=sincosllllL

4、ll因此，L= l1, l2就是新坐标基 1, 2 1112122211-1121222,=xyxyyxyx LL坐标转换定理：某点在以为基的坐标系下的坐标为：，在以为基的坐标系下的坐标为：如果两个坐标基满足关系式：，那么，有：12121 0,=0 1cos-sin,=sincos对于本例，原基为：，新坐标基(在 坐标系下)的坐标为11-1T1T211222cossin=-sincosxxxxyyxxLLLLL，易知， 即 为正交矩阵， 因此有：，从另一种角度得出y1与y2与原变量x1、x2的关系 12121 0,=0 1 LL L同时还有：1x2x1(1,0)1(cos ,sin )2( s

5、in , cos )2(0,1)222212121122221212221212222122212121121221212121211112(=(+)(+)()()(nnnnttttnntttnnttttnntttttttttttynyyyyyynyxnxxyxxxnyxxn yynxxxx222222122212212212(cossin)(cosin)+sxxxxyyyy由和的 表 达 式 易 知 ：2212221121211121=()()()nntttnntttttyyyxxyxx下面考察n个点在新旧坐标系中的离散程度（总离差平方和）：121212222212121222211212=

6、cos+sin= -sin()+cos(+)(+)()=nnttttttyxxyxyyn yynxxxxx对于第一项，有：对于第二项：把和代入，易得：所以，坐标系转换后，散点总的离差平方和（或方差）不变：211212cos+sinsin+cos-=yyxxxx观察右图易知，长轴方向的离散程度要远大于短轴方向的离散程度，即y1的方差要远大于y2的方差，如果对本例计算计算变换前后各变量占总离差平方和的比重比重为：2212121122121211()()52%()()9548=5%nnttttnnttttyyyxyxxx2212221121211121=()()()nntttnntttttyyyxx

7、yxx-202-202x1x2y1y2于是，原变量x1和x2的方差在新变量y1和y2中得到了重新分配重新分配，其中，y1占绝大部分, 解释方差解释方差为95%，称为原变量x1和x2的“第一主成分”，y2为“第二主成分”；由于y1的变化占了95%的方差比重，因此我们可以只分析只分析y1来替代原变量x1和x2 ，从而达到降维分析降维分析的目的；这就是 主成分分析主成分分析 (Principal Component Analysis; PCA) 。即：q各主成分的系数向量都是单位向量；q不同主成分的系数向量相互正交。-202-202x1x2y1y2主成分y1和y2的一个重要性质是：y1和y2相互相互

8、独立独立（如右图所示相关系数为0）。T1,0,ijijij12lll l，是一组新的正交坐标基，因此有如下性质：，22221121122211 1221 2210lllll ll l1T111 1212112112=+=xyl xl xllxl x1T2121222122222=+=xyl xl xllxl x111211221222cossin=,sincosllxllxllx其中，系数向量112T212cossin=-sincosyyxxxxL，把主成分的各分量用系数向量l来表示：以上分析是当x有2个变量（二维）时的情形。下面给出当x具有m个变量时，主成分的一般定义一般定义：设x= x1,

9、 x2, , xm T是一个由m个随机变量组成的随机向量，设x的数学期望为零向量0，x的第i个主成分的定义为：TT(1,1,2,.,)iiiiyiml xl l且满足以下条件：TT111,2,.,=iiyimyl xl x(1)在一切中()，方差最大者称为第1主成分。只要有一组观测（例如n次），x成为m行n列的矩阵，yi就是1行n列的行向量，可以对yi计算方差。T11,2,.,.iiyimyl x(2)第2主成分是指在一切()中，与第1主成分无关（相互独立），并且方差最大者。T1211,2,.,-1,iikkyimkyyyl x(k)第 个主成分是指在一切()中，与前个主成分,.,都无关（相互

10、独立），并且方差最大者。T111112121T221212222T1122=.=.=.mmmmmmmmmmmyl xl xlxyl xlxlxylxlxlxl xl xl x把m个主成分的表达式展开：写成矩阵形式：11211112122222T12.mmmmmmmmlllyxylllxyxlllMMMMMyL x： 由主成分和系数向量的定义或性质可知，主成分与正负符号无关。即：如果l1是第一主成分的系数向量，y1=l1Tx是第一主成分，那么，-l1也是第一主成分的系数向量，-y1= -l1Tx也可称为x的第一主成分。(1)易知 -l1和l1具有相同的性质：(2) 由于方差与符号无关，所以y1和

11、-y1都符合主成分的定义，具有相同的性质。T1,0,ijijijl l1112112212221112., , .,.mmmmmmmmlllyylllylllMMMM其中，yLlll如果有n次观测，那么x有n列：111212122212.=.nnmmmnxxxxxxxxxMMMMX11121111212122221222TT1212.=.nnnnmmmnmmmnnxxxyyyxxxyyyxxxyyyMMMMMMMM于是， 也有 列：yYL XL第1主成分PC1第2主成分PC2第m主成分TyL x第五章2 主成分的导出与计算问题：如何求出一组(m个)变量的主成分？主成分是原变量的线性函数，yi=

12、liTx, i=1,2,m, 因此，确定主成分就是要确定系数向量li依据什么准则求解系数向量依据什么准则求解系数向量li ？根据主成分的定义，各系数向量li具有以下性质：(1). 系数向量为单位向量(模等于1)，即liTli=1(2). 不同主成分的系数向量应该相互正交。即liTlj=0在以上约束条件下，把使得主成分的方差达到最大的系数向量作为第一主成分的系数向量l1。第一主成分第一主成分y1的导出的导出111212122212.=.nnmmmnxxxxxxxxxMMMMX10niiittxxx即：对于矩阵的每一行 ，有这时的y1=l1TX为一行向量，根据主成分的定义，第1主成分y1应该在l1

13、Tl1=1 的条件下方差最大方差最大, 即： 1121T111 12121111.,mmmllyl xl xl xlM其中，l xl对x进行n次观测，然后计算出m行n列的资料阵X，12TTTTTTT11111111T21111()()1ymiisly yl Xl Xl XX ll Sll l达到最大约束条件：严格意义的协方差阵定义应为 S=XXT/(n-1), 此处为简便，暂且将XXT也称为协方差阵，他们只差一个倍数。T=mS XX其中，为 个变量的协方差阵m维随机变量x=x1, x2, , xmT，第一主成分为：因为X是距平，所以y1=l1TX的均值也为0，所以，y1Ty1就是y1的方差(整

14、数倍)。111111111()0nnmmnmkijikkjikkjikijjkkjkyyl xlxl xnnn该问题为多元函数的条件极值问题多元函数的条件极值问题（m元）， 可以应用拉格朗日乘数法拉格朗日乘数法。什么是拉格朗日乘数法？什么是拉格朗日乘数法？(,)(,)= 0(,) =(,) +(,)0(,)(,)(,)0(,)(,)(,)0(,)0fx yx yQx yfx yx yQxyQx yfx yx yxxxQx yfx yx yyyyx y以 二 元 函 数 为 例 ， 寻 找在的 约 束 条 件 下 的 极 值 ，先 构 造 拉 格 朗 日 函 数， 其 中 参 数为 拉 格 朗

15、日 乘 数 。依 次 令对、和的 偏 导 数 等 于， 得 ：三 个 方 程 三 个 未 知 数 ， 可 以 该 方 程 组 求 解 。对于本例来说，要求解m元函数l1TSl1(自变量为l11, l21, , lm1)在约束条件l1Tl1-1=0下的极值, 因此，构造拉格朗日函数拉格朗日函数：TT1111=- (1)Ql Sll lQ是关于l11, l21, , lm1 和的m+1元函数，先依次对l11, l21, , lm1求偏导，写成矩阵形式112111mQlQQlQllTTTTT()()()()+()ddddddddddddddddABxABABBAxxxvu,vxuvvuuvvxxxu

16、uxx矩阵乘积对列向量求导的一般规则：矩阵与的乘积对列向量 求导：为两个列向量，它们的乘积对列向量 求导：TT111111TT111111TTTT111111111111TT11111TTTT1111111TTT11111111- (1)()()()()22222Ql Sll llll Sll llllSlllSlllllllll SISllIlllSSlSllllllSlS llSlIS lllSl112lT1111()nmmnlIl的行向量对的列向量求导办法：把行向量中的各元素分别对列向量求导，得到的矩阵。因此有：单位阵111111=22=QQlSlllSll根据极值原理，令等于零向量，

17、即：，因此，0 所以，就是S的特征值特征值，l1是对应于的特征向量特征向量。第一主成分y1的方差： l1TSl1 =l1Tl1= ，因此，要使y1的方差达最大， 应取S的特征值里面最大的那一个: 1 ，系数向量l1就是对应于最大特征值1的特征向量。分析S的特征值：协方差阵S (m阶)一一定定是半正定矩阵，所以有m个不小于0的特征值，一般情况下，m个变量(X阵中的各行)是线性无关的(S非奇异)， 这时的S成为正定矩阵，正定矩阵S有m个大于0的实特征值；因为对于任意非零列向量l，有:lTSl = lTXXTl= (XTl)T(XTl) 0所以， S=XXT为半正定矩阵，这时n行m列(设nm)的矩阵

18、XT的秩等于m，因此不存在一个非零向量l作为齐次方程XTl=0的解，所以，恒有lTSl = lTXXTl= (XTl)T(XTl) 0,即S为正定矩阵。12.0m就是S的m个特征值之一，l1是对应于的特征向量最后，根据约束条件，将1的特征向量l1单位化，就得到了第一主成分的系数向量l1 的最终结果。特征向量与符号无关(l1与l1都是正确结果。)但是， S有m个特征值， 应取哪一个？第二主成分第二主成分 y2=l2Tx 的确定的确定l2的约束条件： 在以上两个约束条件下，寻找l2使y2的方差： 2TTTTTTT222222222()()ysy yl Xl Xl XX ll Sl达到最大依据两个约

19、束条件，构造拉格朗日乘数法函数Q：TTT22222 1(1)()Ql Sll ll l1 2TTTTTTTT212121211 2 1()() =0y ysy yl Xl Xl XX ll Sll lT112 1=00=0l l无意义，因此，所以需要有：（1）l2为单位向量： l2Tl2=1（2）第二主成分y2与第一主成分y1无关，即：T22T2 1=1=0l ll l两个约束条件：TTT2 1212222122222221()2222220Ql lllSllSllllllllSlllT1TTT1212110ll Sll ll l用 左乘上式，得：TTTTTT21211 2 112121 12

20、T22=0=0=1y yl Sll ly yl Sll ll l约束条件：，0利用约束条件，整理上式得：22222220,=QSlllSll因此，即：所以，l2与l1类似，也是S的特征向量，此处的也是S的一个特征值。2TT2222221T222 1=0yysyl Sll lSll l要使 的方差达到除 方差之外的最大，因此 应取 的第二大特征值 ， 就是相应的特征向量且线性代数定理线性代数定理：对于实对称阵，属于不同特征值的特征向量是正交的最后，根据约束条件(2)，将2的特征向量单位化，即得到第二主成分的系数向量l2的最终结果。那么，应是S的哪一个特征值？与第一、二主成分类似， 我们可以推导出

21、所有的主成分，m个变量的主成分的系数向量就是其协方差阵S=XXT的特征值1 2 m所对应的单位化的特征向量 l1, l2, , lm。第k个特征值k ，就是第k主成分yk的方差 ykykT。因此，求解主成分问题，就是求解协方差阵的特征值与特征向量问题，m个主成分的重要性是按其方差（即协方差阵的特征值）的大小来排列的。严格来说， ykykT是yk方差的整数倍，如果令S=XXT/(n-1), 并计算特征值，那么，这时S的特征值k就是yk严格意义的方差：2TTTTTTTT111()()-1-1-1kykkkkkkkkkk kksnnny yl X l XlXXll Slll（1） 不同主成分是相互正

22、交的(协方差为0)，各个主成分的方差就是S的特征值，所以，m个主成分y =y1, y2, , ymT的协方差阵为对角阵，即：11 212 122122122T2.0.0.0.0=00.mmmmmyy yy yy yyy yyymy yy yysssssssssMMMMMMSYY（2）m个主成分的方差总和方差总和与原m个变量的方差总和方差总和相等，即 11mmkiikisiisXS其中， 是 的协方差 主对角线元素。因为线性代数定理：n阶方阵S的迹(迹的定义为主对角线各元素之和)等于S的所有特征值的和。T=TTTYYL XX LL SL或写成：（3）第一、二、 主成分的方差大小是由大到小顺序排列

23、的, 也即按S的特征值由大到小排列。定义两个常用的指标，用以衡量主成分的方差贡献率：第k个主成分的方差贡献率为：前p个主成分的累积方差贡献率为：1kkmiiU11( )pkkmkkG p1222212myyymsss ，即：因此，原m个变量的方差总和，通过主成分得到了重新分配，其中，前几个主成分贡献了大部分的方差。主成分的各种形式前文所介绍的主成分在计算时都是使用距平变量，X中的m个分量都已经扣除了各自的平均值，即导出第2主成分时，要求第2主成分y2和第1主成分y1有如下关系：1 1T221111111=()() =nnjjy yjjyyysy y达最大110niijjxxn于是，对于任一主成

24、分yi=liTX, 其均值也为 0：TTTTTTT11111111()()yy yl Xl Xl XX ll Sl1在导出第一主成分 时， 总是要求达到最大111111111()0nnmmnmkijikkjikkjikijjkkjkyyl xlxl xnnnyi yiT就是yi的方差方差（整数倍），达最大TTT TTT2121211 2 1=0y yl XX ll Sll l1 2T1212121211=()()=0nnjjjjy yjjy yyyyysy yyi yjT就是yi和yj的协方差协方差（整数倍），等于0距平变量距平变量原始变量原始变量如果X为原始变量，使用原始变量计算主成分，得到

25、S=XXT的特征值与特征向量，注意此时的S不是X的协方差阵TTTTTTT11111111()()yy yl X l Xl XX ll Sl1第一主成分 ， 依然会保证：达最大这时，对于任一主成分yi=liTX, 不能保证其均值也为 0：1 1T221111111=()() =nnjjy yjjyyysy y但此时的y1 y1T达最大，但这并非达最大，但这并非y1的方方差差。第2主成分y2和第1主成分y1依然会有如下关系：TTT TTT2121211 2 1=0y yl XX ll Sll l1 2T1212121211=()()=nnjjjjy yjjy yyyyysy y但，0yi yjT=

26、0，但这不是不是yi和yj的协方协方差差（整数倍）。 yi和yj的协方差协方差不等于不等于0标准化变量标准化变量由于标准化变量同时也是距平变量，因此，对于距平变量主成分的性质，标准化变量的主成分也都具备。但是， 二者是否相等？由于标准化变量的协方差阵(即相关矩阵R)与距平变量的协方差阵S是不同的矩阵，于是，R与S具有不同不同的特征值与特征向量l，所以，标准化变量的主成分（及系数向量）与距平变量的主成分(及相关系数)是不同不同的。由于标准化数据的协方差阵R的主对角线各元素都等于1， 因此，根据性质2可知：R的m个特征值之和等于m1mkkmT12322.922.95xxNN 设两个变量 都来自服从

27、二维正态分布的总体 ：（，），协方差为2.9，方差分别为2和5.12300300 xnxxX对进行次随机采样()，得到原始资料阵 （2行300列）：1112121222.=.nnxxxxxxXXXTXdXdTX*X*T记X的距平资料阵为：Xd； 记X的标准化资料阵为：X*-10-50510-50510 x1x2-10-50510-10-50510 x1x2-10-50510-10-50510 x1x2通过上图可以看出，如果对原始变量计算XXT的特征向量，将无法达到主成分分析的目的，因此，如果要分析原始变量的主成分，应该对协方差阵S=XdXdT或（S=（XdXdT)/n）计算的特征值与特征向量：

28、原始变量和距平变量具有相同的协方差阵， 所以，这时的特征向量li也是距平变量距平变量主成分的系数向量。将系数向量li (i=1,2)与原始变量原始变量相乘， yi= liTX， 即为下图结果(散点在新坐标轴(红色)中的坐标就是y)。11 22T11112T1212()()()()0iyy ysyysyyyyyyy这时， 满足主成分的性质：第一主成分方差：最大不同主成分的协方差为：总结：总结：在计算主成分时，要使用XdXdT或协方差阵S=(XdXdT)/n (或相关系数阵R= (X*X*T)/n)来计算特征值与特征向量(这时原始变量与距平变量的主成分具有相同的系数向量)，而不要计算XXT的特征值

29、与特征向量。XdXdT-10-50510-10-50510 x1x2总结 主成分的计算步骤：1. m个变量n次观测，根据样本资料阵计算协方差矩阵S如果想对标准化变量计算主成分，则需要计算相关系数阵R2. 计算S的m个特征值1, 2, , m及其对应的单位化的特征向量l1, l2, , lm3. 主成分的表达式： yi= liTx, 可把x的n次观测代入即得到主成分的时间序列；还可计算各主成分的方差贡献，以及前几个主成分的累积方差贡献。注意！ 教材显示：本例中，利用距平资料距平资料和标准化资料标准化资料计算的前两主成分的方差贡献方差贡献率相同率相同( (都是0.60和0.85)，但严格来说二者并

30、不相等，只是接近。利用距平资料距平资料和标标准化资料准化资料算得的方差贡献率并不具有可比性方差贡献率并不具有可比性！本例中之所以二者方差贡献率几乎相同是因为该例中的三个变量（3、4、5月气温）的标准差非常接近，导致协方差阵与相关阵几乎成比例。主成分的计算 图例设有30个变量，进行了50次观测(m=30, n=50)，计算距平，得到30行50列的距平资料阵X：把X绘成填色和等值线图如下所示：n=50m=30 102030405051015202530-4-2024下面将绘出前两个主成分的系数向量 (m行1列)以及时间序列 (1行n列)111212122212.=.nnmmmnxxxxxxxxxM

31、MMMXn=50m=30 102030405051015202530-4-2024-0.500.5102030l2m=301020304050-20-1001020n=50y2第2主成分的系数向量l2第2主成分y2的时间序列（PC2）PC2方差贡献：19%把各主成分与X进行比较，观察它们分别反映了X的怎样的特征？第1主成分的系数向量l1第1主成分y1的时间序列(PC1)-0.500.5102030l1m=301020304050-20-1001020n=50y1PC1方差贡献：73%注意：系数向量和主成分都无关符号,即：二者可以同时同时乘以负号。第五章3 经验正交函数分解12T,kkmxxmk

32、yxM个变量：， 第 主成分为：xl x111212122212.=.nnmmmnxxxxxxxxxMMMMXT1111211T221222T211T12., , .,.mmmmmmmmmlllyylllmylllMMMMM个主成分构成一列向量(m行1列)：其中，l xl xyL xLllll x如果对这m个变量进行了n次观测，那么x变成m行n列的资料阵X这时，m个主成分为m行n列的矩阵Y：112111112111121122222122221222T121212.mnnmnnmmmmmmmnmmmnlllxxxyyylllxxxyyylllxxxyyyMMMMMMMMMMYL X12=kkk

33、mklllM其中lTYL X因为L是正交矩阵（即LTL=LLT=I），所以，将上式左乘L，可得：于是，多变量资料阵X被分解分解为两个矩阵乘积的形式：X = LY设X的m个变量是同一个气象要素在m个空间点(或台站)上的取值，根据分解公式，第i个变量的第j次观测，可写为：11221.mijijijimmjikkjkxl ylylyl y可以发现， lik 仅依赖于空间点i，不随时间j变化，称为“空间函数”； ykj 只与时间j有关，与空间点i无关，称为“时间函数(系数)”111211112112212222122212()() ()1212.= , .,.mnmnmm nm mm nmmmmmmm

34、mnlllyyylllyyylllyyyyyXLYllly12121, .,mmkkkkimmnyyX = LYllll yl yy， 其中，每个都是 行 列的矩阵。将资料阵X分解为空间函数和时间系数相乘的形式，称为：X的“经验正交函数分解” (Empirical Orthogonal Function (EOF) decomposition)，用于考察气象要素的时空分布时空分布特征。“Empirical(经验)”的含义：“以实验或观察为依据的”， 因为EOF分解要依赖于对研究对象的n次观测。“Orthogonal(正交)”的含义：m个空间函数相互之间是正交的，m个时间函数也是相互正交的。主成

35、分分析(PCA)与经验正交函数分解(EOF)拥有相同的计算方法，只是物理解释的角度略有不同，对lk 和yk 的称呼也有所不同：主成分分析（主成分分析（PCA）经验正交函数（经验正交函数（EOF）分解）分解l k (m行1列)第k主成分的系数向量第k模态模态的空间函数(空间向量)y k (1行n列) 第k主成分的时间序列第k模态模态的时间函数(时间系数)()mnmX整个行 列 资料阵可分解为 个矩阵相加的形式：1112111121111211212222212222122212()() ()121212.= , .,.nmnnmnmm nm mm nmmmmnmmmmmmxxxlllyyyxxx

36、lllyyyxxxlllyyyyXLYllly.mny空间函数值的性质空间函数值的性质12TT121 12 2000011= , ., .,000mm mmnnXYLYYllllllT1jjjnjXyl 所以，第 模态的空间向量：T21ijjijx yji jjysnilsx yl 的第 个元素TTTTT1211, .,mnnXYXYXyXyXy把的各列展开：T111()()ijniijikikx yikxxyysnnx y2T1(=)jjyjddsnSX X这里的是的特征值2ijjx yijiyjslisjxy因此，就是以“第 个变量的观测值 （1行向量）”为预报量，以“第 主成分 （1行向

37、量）”为因子的一元回归方程的回归系数。TT1jjjjnXYXyl取的第 列得： T1nX = LYY对两边右乘，得：气象要素场的估计气象要素场：某气象要素，m个空间点，n次时间观测，气象要素场的资料可表示为X：111212122212.=.nnmmmnxxxxxxxxxXMMMMX = LY由EOF分解可知，X可表示为两个矩阵的乘积：111211212222121 122()() ()12.= , .,+.+.nnmmmm nm mm nmmmmnxxxxxxxxxyyXLYllll y l yl yyMMMM M由于li是m行1列的向量，yi是一行m列的向量，因此，每一个liyi都是一个与X

38、形状相同的矩阵：1111111211111212111211221111121111111211.nniinmmmmnl yl yl ylll yl yl yyyyllylylyl yMMMM将其记为X(i)11121111212122221222121212121122()1112221112221112= , .,+.+.nnnnmmmnmmmmmm nmxxxxxxxxxxxxxxxxxyyXLYllll yl yl yy（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）111212122212()()()()()()22()()()(1)(

39、2)().=nnmnmmmnmmmmmmmmmmxxxxxxxxxxX+ X+ .+ X（ ）（ ）每个li都是单位化的，但时间序列yi却是按方差从大到小排列的，即越前面的X (i) 就越重要（总方差越大），如果取前p个EOF模态(主成分)，即L前p列, Y取前p行(p=m), 这时就得到对原资料阵X的一种估计估计(重构)：记为( )11=ppiiiiiXXl ymnXm原资料阵 就是 个与之形状相同( 行 列)的矩阵相加所构成的：对X估计的精确程度，可采用实测值与估计值之间的残差平方和来度量，即：T2111=() ()()nnmjjijjjijjjiQxxxxxx$jjjxxXX其中，和分别

40、是 和 的第 列。QXX即为 和 两矩阵作差然后求各元素的平方和。P152-153-Qm p可以证明()， 可以用剩下的个特征值之和来表示：1mkkm由于原要素场 个变量的总方差可表示为：，所以可用“1减去两者的比值”来代表估计场对原始场拟合程度的百分率：1111= 1 -=mpkkkpkmmkkkkG p （ ）111pmmkkkkpkkQ 此即“累积方差贡献率”。-0.500.5102030l1m=301020304050-20-1001020n=50y1-0.500.5102030l2m=301020304050-20-1001020n=50y230个变量，50次观测(m=30, n=5

41、0)，得距平资料阵X：对其进行EOF分解，取前三个模态，绘出各有关结果：空间向量和时间系数可以同时乘以负号！前三个模态的方差贡献率分别为：67.6%，19.2%，6.7%, 累积方差贡献率分别为：67.6%，86.8%， 93.5%。l1y1-5-4-3-2-1012345n=50m=30102030405051015202530-5-4-3-2-1012345n=50m=30102030405051015202530X(1)=l1 y1X(1)l2y2-5-4-3-2-1012345n=50m=30102030405051015202530-5-4-3-2-1012345n=50m=3010

42、2030405051015202530X(2)=l2 y2X(1)+ X(2)-5-4-3-2-1012345n=50m=30102030405051015202530-5-4-3-2-1012345n=50m=30102030405051015202530-0.500.5102030l3m=301020304050-20-1001020n=50y3l3y3X(3)=l3 y3X(1)+ X(2)+X(3)-5-4-3-2-1012345n=50m=301020304050510152025300.010.020.030.040.050200400600800100012001400m=140

43、0l1010020030040050060070080090010002000300040005000600070008000n=1000y1对风景图的明度进行EOF分解明度即明亮程度，白色最大，黑色最小。高度 宽度：1400 1000像素资料阵X：m=1400, n=1000.第一特征向量l1表现为在300-400行左右有极大值，对应于瀑布落地处的大片白色。做EOF之前，肉眼观察该图的从左到右最显著的特征是： “暗明暗”，推测第一主成分应该会反映此特征。由瀑布从上而下的形态推测特征向量不会有明显起伏，可能会在瀑布落地处所在的行呈现峰值。在特征向量l1的坐标系下，第一主成分y1反映出从左到右“

44、暗明暗”的特征，与肉眼观察推测一致。First 1 EOFs (55%)First 2 EOFs (76%)First 10 EOFs (93%)First 50 EOFs (98%)利用前p个EOF模态( p=1, 2, 10, 50) 对原图进行重构(估计):计算说明：根据协方差阵S=XdXdT计算特征向量，然后利用原始资料阵yi = liTX计算主成分；最后重构的是原始资料阵X(非距平阵)，因为图片的色值都是正值。高圆圆、赵又廷 11月28日在台北举行婚礼照片从左到右，明度最大的变化是由新娘的亮色转为新郎的暗色，因此第一主成分y1体现了该特征。第一特征向量表现为在下部100-200行是高

45、值区，350-450行为低值区，基本反映下半部婚纱与鲜花的亮色，而上面为皮肤与头发所代表的相对暗色。00.020.040.060.0850100150200250300350400450500550m=550l1501001502002503003504004505000200040006000n=500y1利用前p ( p=1, 2, 10, 50)个EOF模态对原图进行重构(估计):实际应用 热带印度洋海表面温度年际异常的EOF分解1960-2013年月平均海表面温度数据，选取热带印度洋区域(20S-20N，30E-120E) ，通过EOF分析考察海表面温度年际异常的时空分布特征时空分布特

46、征。空间：分辨率（纬度经度）为11，共有6292个格点，陆地格点数为1805，海洋格点数为3899；时间：54年共648个月。热带印度洋的海温异常第一模态（方差贡献33% ）为海盆范围的一致增暖（或变冷）的模态，常称之为：“海盆（一致）模”（Indian Ocean Basin-wide mode； IOB）。通过分析时间序列(PC-1)，可以看到该模态在哪些年份或季节比较显著，也可将其作为一个变量的序列进行各种分析(相关或回归分析、周期分析、判别等)。EOF-1 30oE 60oE 90oE 120oE 20oS 10oS 0o 10oN 20oN -0.03-0.02-0.0100.010

47、.020.03196019701980199020002010-50050PC-1(33%)因此，X资料阵的形状：m=3899, n = 648注意：分析EOF结果时，空间向量与时间系数总是要成对地同时呈现成对地同时呈现。做EOF前准备资料阵X时，要把空间场从二维(经纬度)按一定规则顺序转换成一维列向量，EOF做完后再按同一样的规则转化为二维，以便对空间向量绘图。 30oE 60oE 90oE 120oE 20oS 10oS 0o 10oN 20oN EOF-2 (13%)-0.03-0.02-0.0100.010.020.03196019701980199020002010-50050PC-

48、2热带印度洋海温EOF第二模态热带印度洋海温异常EOF2（方差贡献率为13% ）呈现为东西符号相反的异常分布形态，称为：“偶极子模态”（Indian Ocean Dipole mode；IOD）。EOF的空间向量和时间系数的数数值值分别有什么物理意义？EOF第i模态的空间向量li 是单位化的，各元素的平方和等于1：22212.1iimilll于是，各元素的数值与格点的总数量格点的总数量有关，因此，空间向量的数值没有物理意义。 所以，li乘以资料阵X所得的时间序列时间序列yi = liT X 的数值数值也与格点数有关与格点数有关，没有物理没有物理意义意义（尽管主成分的标准差存在相对大小：y1 的

49、标准差大于y2.）没有物理意义的空间向量和时间系数数值限制了我们利用EOF分析问题的能力，是否可以使它们具备一定物理意义？是否可使空间向量或时间系数数值具备物理意义？第j模态对气象要素场X的估计：111212122212( ).nnmmmnjjjjjjjjjjjjxxxxxxxxxXl y（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）（ ）而导致lj和yj的数值没有物理意义的根源是： 单位化的lj的数值受格点总数的影响。X(j)的每个元素和原气象要素具有相同的单位，因此, X(j)是具有物理意义的，,jjjyy把时间系数 进行标准化，即：将 除以( ):()(=)jjjjjjjjjjXlXl

50、 yly（ ）为保持不变，需再将空间向量 乘以可进行调整：不再保持空间向量的单位化，而是使时间序列时间序列yj标准化标准化。时间序列yj的方差是？TTTTTT2TT2TTTTTTTT=()=11=1=111()()jjjjjjjjjjjjjyjjjjyjjjjjjjjjjjjsnnnsnnnSXXy ylXXll Slllyy ySXXyy yl X l XlXXll Slll如果对计算特征值与特征向量，那么，因此，的方差：如果对计算特征值与特征向量（），那么，的方：推荐差T1jjnS =XXy为了方便，建议计算的特征值与特征向量,这时 的标准差为：)/(jjjjjXyl（）如此一来，不同模态

51、的时间系数的标准差都统一为1，因此，在量值上将主要依赖于空间向量()。/jjjjly在绘图时，分别把和分别作为空间向量与时间系数绘出，T21=ijjijx yji jijjysnilxysx yl前面的学习已知： 的第 个元素是以 为预报量，以为因子的回归系数,j ijiixxl因此，可认为与 带有相同的单位，具备了物理意义：表示当标准化的时间序列数值等于1时，该模态所反映的 的异常值。jjl 的本质是？TTT11()1()1jij i jjjjjjjijjijij i jj iijyynxxlnnllyxyyxx所以，就是以 为预报量，以标准化的时间序列（方差等于1）为因子的回归系数,同时，

52、也是 和的协方差。1iixxs X特别地，如果原变量 是标准化的资料阵： 即 的标准差，ijijjjxly这时， 和都是标准化的，所以它们的协方差就是二者的相关系数。EOF-1 30oE 60oE 90oE 120oE 20oS 10oS 0o 10oN 20oN -0.3-0.2-0.100.10.20.3196019701980199020002010-4-2024PC-1调整前后的EOF-1(印度洋海温异常)EOF-1 30oE 60oE 90oE 120oE 20oS 10oS 0o 10oN 20oN -0.03-0.02-0.0100.010.020.031960197019801

53、99020002010-50050PC-1调整前与调整后的空间和时间函数在形态上相同，只是在量值上相差一个常数。调整后的空间向量的数值可认为具有单位：“摄氏度”，表示该模态在标准状态下（即当PC值为1的时）所反映的海温异常的幅度。 30oE 60oE 90oE 120oE 20oS 10oS 0o 10oN 20oN EOF-1-1-0.500.51196019701980199020002010-4-2024PC-1调整前调整后调整后标准化标准化变量的EOF调整后，如果采用标准化变量做EOF，这时的空间分布图中各格点数值的含义为：原变量与主成分的相关系数。利用EOF对我国汛期降水分型P156，图5.2第1特