气象编程|利用Python对夏季降水和同期大西洋海温进行SVD分析

作者 | 摸鱼咯

链接 | https://www.jianshu.com/p/354b7cd158ec

犹豫了很久要不要讲SVD的python实现，今天还是写了吧，纠结的原因在于我对SVD也是新手理解，很多东西也是一知半解，怕误导大家。所以这篇文章也只是我个人的理解，不对之处还请大家评论区帮我纠正。

SVD

气象中所说的SVD方法是干吗用的？

首先还是先讲讲什么是SVD，在气象中的作用是什么，不想看我废话的直接拉下去看代码。我们通常计算两个变量之间的相关关系时往往使用相关系数计算，然而相关系数只能用于两个序列(两个一维变量)或者一个序列和一个场(一个一维变量和一个三维变量)之间的相关关系，那如果我们想找到两个变量场之前的相关性，这时就要用到SVD了。前边的文章有讲到EOF经验正交分解，核心思想是提取变量场的几个主要模态，并且计算各个模态的时间系数。那么SVD则可以理解为，对两个场分别提取模态，看两个变量的模态之间的协同变化关系。

SVD方法中几个特殊名词的概念？

SVD分析中有一些特殊名词，我用通俗的话来解释。比如说：左场，右场，分别对应我们的要求的两个变量场，谁左谁右并不重要。左场提取的模态称为左奇异向量，右场提取的模态为右奇异向量。需要注意的是，各个场的奇异向量均为相互正交的。第一左奇异向量和第一右奇异向量及其各自的时间系数共同构成了SVD的第一模态，也可以叫第一对模态。还有三个重要的名词要掌握。第一个是总体相关系数，指的是一对奇异向量对应的左右时间系数的相关系数，用来看左场第一模态和右场第一模态的相关性（总体相关系数是一个数）。第二个名词是同性相关系数，表示原场和原场某一模态的时间序列的相关系数（为一个场），在一定程度上可以反应该变量的一个遥相关型。第三个名词是异性相关系数，代表原场（比如左场）和对立场（比如右场）某个模态的时间序列的相关系数（为一个场），表是一个场对另一个场的影响关键区。这段话说的十分绕口，等下结合例子再认真想想几个名词的含义应该会清楚一些。那么就开始吧。

核心代码

首先是讲下核心代码， 官方文档参考[1]

git clone https://github.com/Yefee/xMCA.git
cd xMCA
python setup.py install

这是官方提供的安装方式，能否用CONDA安装大家可以自己试一下。

from xMCA import xMCA
#a,b 为左场和右场
svd = xMCA(a,b)
svd.solver()
#lp,rp分别为左场和右场的模态，n=2，即为求取前两模态
lp, rp = svd.patterns(n=2)
#le,re为对应的时间系数序列
le, re = svd.expansionCoefs(n=2)
#方差贡献
frac = svd.covFracs(n=2)

以上代码即为求解SVD的核心代码，接下来我以一个例子来说明。

举个栗子

我准备的数据是1961-2016年夏季的中国平均降水和同期大西洋海温，分别来自中国CN0.51格点资料和英国哈德来中心的海温资料。以下是经过处理后的资料格式（直接用于求解SVD）：

输入数据维度

svd = xMCA(sst_summer, pre_summer)
svd.solver()
lp, rp = svd.patterns(n=2)
le, re = svd.expansionCoefs(n=2)
frac = svd.covFracs(n=2)
print(frac)
#<xarray.DataArray 'frac' (n: 2)>
#array([0.4099809 , 0.16947176], dtype=float32)
#Coordinates:
#  * n        (n) int64 0 1
#Attributes:
#    long_name:  Fractions explained of the covariance matrix between left and...

得出第一模态方差贡献约41%，然后是画图部分，画图细节可参考先前几篇文章，详细讲过了画图部分。

fig, (ax1, ax2) = plt.subplots(2, 2, figsize=(12, 5))
lp[0].plot(ax=ax1[0])
le[0].plot(ax=ax1[1])
rp[0].plot(ax=ax2[0])
re[0].plot(ax=ax2[1])

使用该代码即可快速出图，然而会发现子图之间排列密集，且没有地图地图，在确认计算结果无误后可根据需要调整图形美观度。下图是我修改后的结果

SVD第一模态

两条时间序列相关系数(就是前文提到的总体相关系数为0.67)说明这对空间分布型有密切的关系。我国夏季降水第一分布型在内蒙为正相关区，长江流域附近为负相关区，而大西洋中部海温为一致的负相关型分布，说明当大西洋海温偏低时，我国内蒙地区降水偏多，长江流域附近降水偏少。当然，分析时还要结合时间序列进行分析，分析方法同EOF。我们前边还提到了同性相关和异性相关，十分方便的是，这个库的作者提供了直接计算两者的方法。并且是自带双尾t检验的，据说作者正在开发蒙特卡洛检验，先期待下。

#homogeneous 同性
lho, rho, lphot, rphot = svd.homogeneousPatterns(n=1, statistical_test=True)
#heterogeneous 异性
lhe, rhe, lphet, rphet = svd.heterogeneousPatterns(n=1, statistical_test=True)

这样求出来的8个变量均是二维场，

le, re, lphet, rphet = svd.heterogeneousPatterns(n=1, statistical_test=True)