一、概念和原理

经验正交函数分析方法(Empirical Orthogonal Function,缩写为EOF)，也称特征向量分析(Eigenvector Analysis)，或者主成分分析(Principal Component Analysis,缩写PCA)，是一种分析矩阵数据中的结构特征，提取主要数据特征量的一种方法。Lorenz在1950年代首次将其引入气象和气候研究，在地学及水声学等其他学科中得到了非常广泛的应用。地学数据分析中通常特征向量对应的是空间样本，所以也称空间特征向量或者空间模态；主成分对应的是时间变化，也称时间系数。因此地学中也将EOF分析称为时空分解。

E0F分析方法能够把随时间变化的变量场分解为不随时间变化的空间函数部分以及只依赖时间变化的时间函数部分。空间函数部分概括场的地域分布特点，而时间函数部分则是由场的空间点的变量线性组合所构成，称为主要分量。这些分量的头几个占有原场内空间点所有变量的总方差的很大部分，这就相当于把原来场的主要信息浓缩在几个主要分量上，因而研究主要分量随时间变化的规律就可以代替场的时间变化研究，且可以通过这一分析得出的结果来解释场的物理变化特征。

二、计算方法

具体的计算步骤通常包括：（1）对数据进行距平处理，形成一个距平数据矩阵；（2）计算距平数据的协方差矩阵；（3）计算协方差矩阵的特征值和特征向量，特征向量即为EOF；（4）将EOF投影到原始数据矩阵上，得到所有空间特征向量对应的时间系数（即主成分）。EOF分析的核心是计算协方差矩阵的特征值和特征向量，有以下两种计算方法：

1、直接计算特征值和特征向量

利用MeteoInfoLab的numeric包中的线性代数子包linalg中的eig函数直接计算矩阵的特征值和特征向量。

生成一个简单的数据矩阵X（二维数组）：

对数据进行距平处理：

生成协方差矩阵C：

计算特征值（E）和特征向量（EOF）：

计算主成分（PC）：

对结果矩阵进行顺序调整：

得到的结果：

检验EOF和PC的正交性：

EOF和PC的矩阵乘积可以恢复X的距平值：

2、SVD分解

利用linalg包中的svd函数进行SVD分解，程序如下：

计算结果：

SVD分解对内存要求较少，计算速度较快。

为了方便EOF计算，MeteoInfoLab的meteolib包提供了一个eof子包来直接进行EOF分析，相关代码主要从Python的eofs库（https://ajdawson.github.io/eofs/latest/）移植。

三、气象场EOF分析

下面通过两个例子演示如何对实际气象场进行EOF分析。

1、太平洋海温距平数据

冬半年太平洋海温的变化和厄尔尼诺-南方振荡现象（ENSO）有关，示例文件包含了多年11月至次年3月太平洋海表温度距平数据。首先从文件中读出海表温度（sst）和经纬度（lon, lat）数据，sst是三维数组（时间、纬度、经度），lon和lat是一维数组。

由于已经是距平数据，省去数据距平的步骤。对于等经纬度网格数据，低纬度网格面积更大，对计算结果影响也应该更大，因此对数据进行纬度加权操作，用纬度余弦平方根计算权重，然后创建一个meteolib包中的EOF类的对象solver，参数是sst数组和权重：

然后从solver中计算表示为相关系数的EOF第一模态，与空间第一模态对应的时间系数，以及第一模态的方差贡献率。

还有一个需要注意的问题，气象格点数据可能包含缺测值，这个示例数据里sst就有缺测值（NaN，数据区域内的陆地范围），在SVD等线性代数计算中不能出现NaN，因此必须在EOF计算之前去掉NaN值，分析完成后再重新添加NaN值生成规则格点数据。为了方便使用，MeteoInfoLab meteolib包中的eof函数包含了对缺测值的自动处理功能。下面的代码计算EOF和PC，其中e1是第一个空间模态的方差贡献率：

然后可以绘图显示第一个模态的空间分布和时间序列，完整代码如下：

可以看出上图空间模态是典型的厄尔尼诺，相关的时间序列显示了厄尔尼诺和拉尼娜事件的峰值和低谷。第二、第三等模态数据的处理和绘图也类似。

2、冬季北大西洋位势高度数据

北大西洋涛动（North Atlantic oscillation，NAO）于1920年由Sir Gilbert Walker发现，三大涛动之一，指亚速尔高压和冰岛低压之间气压的反向变化关系，即当亚速尔地区气压偏高时，冰岛地区气压偏低，反之亦然。NAO是北大西洋地区大气最显著的模态。其气候影响最突出的主要是北美及欧洲，但也可能对其它地区如亚洲的气候有一定影响。

数据处理过程多了一个距平，其他都类似：

四、REOF分析

旋转正交经验分解（REOF）也即旋转主因子分析，REOF克服了EOF分离出的空间场不能清晰表示不同地理区域特征及存在较大取样误差（取整体或部分会造成不同的分析结果）的缺陷。旋转后的典型场空间结构清晰，不但可以较好地反映不同区域的变化特征，还可以反映不同区域的相关分布。Varimax是最常用的旋转算法，计算出EOF后用varimax方法对其进行旋转即可获得REOF结果。