CN112016588B - 一种面向遥相关模式的空间自相关聚类方法 - Google Patents

Abstract

本发明提供了一种面向遥相关模式的空间自相关聚类方法，其通过考虑每一个空间网格单元与相邻单元遥相关的相关程度，在局部Moran指数定义的基础上，采用相关系数原始数值而不作中心化处理，从而改进局部Moran指数计算公式，得到新的空间自相关局部指标LISAAC，实现显著正或负遥相关集聚范围的探测，同时实现异常值(即显著正值区域中出现不显著或负值网格，显著负值区域中出现不显著或正值网格)的识别。实验结果表明，本发明能根据遥相关系数本身的标准化性质，实现对不同类型遥相关的空间聚类，结果便于对不同季节的遥相关程度进行横向对比。

Description

一种面向遥相关模式的空间自相关聚类方法

技术领域

本发明涉及气象与气候学领域，更具体的，涉及一种面向遥相关模式的空间自相关聚类方法。

背景技术

大尺度的气象因子与全球各地降水异常的相互关系一直受到人们的重视，利用相关的气象因子指标和各地月或季节尺度降水的相关系数，可以量化这种遥相关作用的强度和模式。然而在实际计算中，相邻区域的降水往往并非独立，而是存在较强的相关关系。为了探测这种空间上的自相关性，往往采用空间自相关局部指标(LISA)进行空间自相关分析，探测变量是否存在空间格局的有关信息。然而对于以相关系数衡量的遥相关模式而言，系数本身已经具有标准化的性质，而常用的LISA指标如Moran指数，根据《计量地理学》第四章中所描述的定义，其计算过程中将对变量进行标准化处理，因此所得结果只能反映在同一时间及指标下，相对的高值和低值分布，而难以对不同时间或指标的计算结果进行横向对比。而大尺度气象因子往往具有较强的季节间或年际周期性，不同季节下的遥相关强度和模式通常具有较大的差别。因此，有必要对传统的LISA指标进行改进，以便于对不同季节、不同指标下遥相关模式的空间聚类模式。

发明内容

本发明为解决传统的LISA指标得到的结果只能反映在同一时间及指标下，相对的高值和低值分布，而难以对不同时间或指标的计算结果进行横向对比的技术缺陷，提供一种面向遥相关模式的空间自相关聚类方法。

为实现以上发明目的，而采用的技术手段是：

一种面向遥相关模式的空间自相关聚类方法，包括以下步骤：

S1：获取研究区域空间网格坐标信息，依据坐标信息计算空间权重矩阵；

S2：获取研究区域网格尺度降雨数据及相同时间范围的大尺度气象因子指标，得到降雨-气象指标时间序列；

S3：依据所获得的降雨-气象指标时间序列，逐网格计算降雨-气象指标相关系数r；

S4：根据相关系数r及空间权重矩阵，计算逐网格的空间自相关局部指标LISAAC；

S5：重排步骤S2中的降雨及气象指标时间序列，获得新的降雨-气象指标时间序列，重复执行步骤S3-S4，直到达到预设的迭代次数n；

S6：根据步骤S4得到的n组随机LISAAC，构建参考经验分布H；

S7：根据观测LISAAC、观测相关系数r以及经验分布H，得到在指定显著性水平下各网格的分类结果。

其中，在所述步骤S1的空间权重矩阵计算中，距离的表示方法为欧几里得距离。

其中，在所述步骤S1的空间权重矩阵计算中，权重的表示方法为距离平方的倒数。

其中，在所述步骤S2中，记观测的气象指标为：

X＝[x_t]

其中，x_t为第t年观测到的指标数值，X为x_t组成的集合；记网格总数为N，以i为网格索引，类似地，记观测降雨为：

Y＝[y_t,i]

其中，y_t,i为第t年在网格i观测到的指标数值，Y为y_t,i组成的集合。

其中，在所述步骤S3中，对于给定的网格i，由X和Y序列得到气象因子和降雨之间的相关系数：

式中，rx_t(ry_t)代表第t年气象因子，即观测降雨在原始序列中的排序，而则代表序列rx_t(ry_t)的均值；r的值介于-1与1之间，当两个变量完全单调相关时，r的值为1或-1；对于网格i，所计算的相关系数可记为r_i，由此可得到相关系数组成的集合R＝[r_i]。

其中，在所述步骤S4中，所述空间自相关局部指标LISAAC的计算公式具体为：

其中，w_i,j是用于将r_i与邻近网格r_j进行联系的权重系数；为网格i邻近网格的空间加权相关系数；C_i描述r_i和其邻近网格r_j之间的关系的强度；式中的权重系数是影响C_i的重要参数。

上述方案中，所述步骤S4的LISAAC计算采用相关系数本身不作中心化处理。

其中，为了使权重系数随距离增大而衰减，采用反距离权重法，以欧几里得法计算距离，并以该距离的平方的倒数作为两网格之间的距离权重系数：

式中d(i,j)为第i和j各网格之间的欧式距离；当d(i,j)增大，邻近网格相关系数r_j对C_i的影响越小，通常根据数据尺度预设距离阈值以减少计算量。

其中，在所述步骤S5中，所述重排计算为直接随机重排每个网格的观测降雨-气象指标时间序列。

其中，所述步骤S6具体为：

判断遥相关的空间聚类显著性，其对应原假设为气象指标与观测降雨不存在显著的相关性，即遥相关不显著；因此，通过随机排列第i(j)个网格历史气象指标和观测降雨序列，重新计算r_i(r_j)，及相应的C_i；通过多次重复计算，构建C_i的参考经验分布H。

其中，所述步骤S7具体为：

依据参考经验分布H，获得观测C_i对应的p_i值，以表示原假设为真时的强度,根据观测r_i和相应的p_i值，计算得到网格i的分类：

式中PP表示某一网格点和其周围网格点的r都显著为正；PN表示某一网格的r为正，但其周围网格的r值偏低或为负，ns代表网格r_i及其周围网格r_j均不显著；NP代表负r_i被正的r_j围绕；最后一种情况NN表示的是负r_i被负r_j围绕；PP和NN表示r值具有较高的空间正相关性，提示存在区域集聚；而PN和NN则反映了r的空间分布存在异质性。

上述方案中，本发明最后的网格分类结果具体为PP(显著正遥相关网格被相同符号遥相关的网格围绕)、PN(异常值)、ns(不显著网格)、NP(异常值)、PP(显著负遥相关网格被相同符号网格围绕)。

与现有技术相比，本发明技术方案的有益效果是：

本发明提供的一种面向遥相关模式的空间自相关聚类方法，通过考虑每一个空间网格单元与相邻单元遥相关的相关程度，在局部Moran指数定义的基础上，改进局部Moran指数计算公式，得到新的空间自相关局部指标LISAAC，实现显著正或负遥相关集聚范围的探测，同时实现异常值的识别，能根据遥相关系数本身的标准化性质，实现对不同类型遥相关的空间聚类，结果便于对不同季节的遥相关程度进行横向对比。

附图说明

图1为本发明所述方法的流程示意图；

图2为三种方法下四个季度各类型网格占比堆叠柱状图。

具体实施方式

附图仅用于示例性说明，不能理解为对本专利的限制；

为了更好说明本实施例，附图某些部件会有省略、放大或缩小，并不代表实际产品的尺寸；

对于本领域技术人员来说，附图中某些公知结构及其说明可能省略是可以理解的。

下面结合附图和实施例对本发明的技术方案做进一步的说明。

实施例1

如图1所示，一种面向遥相关模式的空间自相关聚类方法，包括以下步骤：

S1：获取研究区域空间网格坐标信息，依据坐标信息计算空间权重矩阵；

S2：获取研究区域网格尺度降雨数据及相同时间范围的大尺度气象因子指标，得到降雨-气象指标时间序列；

S3：依据所获得的降雨-气象指标时间序列，逐网格计算降雨-气象指标相关系数r；

S4：根据相关系数r及空间权重矩阵，计算逐网格的空间自相关局部指标LISAAC；

S5：重排步骤S2中的降雨及气象指标时间序列，获得新的降雨-气象指标时间序列，重复执行步骤S3-S4，直到达到预设的迭代次数n；

S6：根据步骤S4得到的n组随机LISAAC，构建参考经验分布H；

S7：根据观测LISAAC、观测相关系数r以及经验分布H，得到在指定显著性水平下各网格的分类结果。

更具体的，在所述步骤S1的空间权重矩阵计算中，距离的表示方法为欧几里得距离。

更具体的，在所述步骤S1的空间权重矩阵计算中，权重的表示方法为距离平方的倒数。

更具体的，在所述步骤S2中，记观测的气象指标为：

X＝[x_t]

其中，x_t为第t年观测到的指标数值，X为x_t组成的集合；记网格总数为N，以i为网格索引，类似地，记观测降雨为：

Y＝[y_t,i]

其中，y_t,i为第t年在网格i观测到的指标数值，Y为y_t,i组成的集合。

更具体的，在所述步骤S3中，对于给定的网格i，由X和Y序列得到气象因子和降雨之间的相关系数：

式中，rx_t(ry_t)代表第t年气象因子，即观测降雨在原始序列中的排序，而则代表序列rx_t(ry_t)的均值；r的值介于-1与1之间，当两个变量完全单调相关时，r的值为1或-1；对于网格i，所计算的相关系数可记为r_i，由此可得到相关系数组成的集合R＝[r_i]。

更具体的，在所述步骤S4中，所述空间自相关局部指标LISAAC的计算公式具体为：

其中，w_i,j是用于将r_i与邻近网格r_j进行联系的权重系数；为网格i邻近网格的空间加权相关系数；C_i描述r_i和其邻近网格r_j之间的关系的强度；式中的权重系数是影响C_i的重要参数。

在具体实施过程中，所述步骤S4的LISAAC计算采用相关系数本身不作中心化处理。

更具体的，为了使权重系数随距离增大而衰减，采用反距离权重法，以欧几里得法计算距离，并以该距离的平方的倒数作为两网格之间的距离权重系数：

式中d(i,j)为第i和j各网格之间的欧式距离；当d(i,j)增大，邻近网格相关系数r_j对C_i的影响越小，通常根据数据尺度预设距离阈值以减少计算量。

更具体的，在所述步骤S5中，所述重排计算为直接随机重排每个网格的观测降雨-气象指标时间序列。

更具体的，所述步骤S6具体为：

判断遥相关的空间聚类显著性，其对应原假设为气象指标与观测降雨不存在显著的相关性，即遥相关不显著；因此，通过随机排列第i(j)个网格历史气象指标和观测降雨序列，重新计算r_i(r_j)，及相应的C_i；通过多次重复计算，构建C_i的参考经验分布H。

更具体的，所述步骤S7具体为：

依据参考经验分布H，获得观测C_i对应的p_i值，以表示原假设为真时的强度,根据观测r_i和相应的p_i值，计算得到网格i的分类：

式中PP(正正)表示某一网格点和其周围网格点的r都显著为正；PN(正负)表示某一网格的r为正，但其周围网格的r值偏低或为负，ns(不显著)代表网格r_i及其周围网格r_j均不显著；NP(负正)代表负r_i被正的r_j围绕；最后一种情况NN(负负)表示的是负r_i被负r_j围绕；PP和NN表示r值具有较高的空间正相关性，提示存在区域集聚；而PN和NN则反映了r的空间分布存在异质性。

在具体实施过程中，本发明最后的网格分类结果具体为PP(显著正遥相关网格被相同符号遥相关的网格围绕)、PN(异常值)、ns(不显著网格)、NP(异常值)、PP(显著负遥相关网格被相同符号网格围绕)。

实施例2

更具体的，在实施例1的基础上，本实施例通过实验对方法的效果进行说明，以1982-2010全球季节网格降水数据及指标为例，计算厄尔尼诺-南方涛动(ENSO)与全球季节降水遥相关模式的空间自相关聚类，以显著性水平α为0.10为例，计算LISAAC，同时对比不考虑空间自相关的网格显著性分类结果和采用局部Moran指数计算的结果。

图2分别给出了采用三种方法计算的分类网格占比。其中图2a为利用本发明提供的LISAAC计算的网格分类。可以明显看到，虽然研究表明全球大部分区域季节降雨受到不同程度的ENSO遥相关作用，但大部分网格显示出不显著的遥相关作用，一方面可能是由于这种遥相关作用受到其它气象因素如融雪、土壤湿度和地形等的影响，一方面可能是由于ENSO信号在这些区域或季节本身偏弱。此外可以看到，夏季正遥相关作用所占网格比例低于其它三个季节，相反的，其负遥相关作用占比相对更高。相比而言，如果不考虑网格的空间信息，仅仅按照遥相关的显著性对网格进行分类，可以得到类似的结果(图2b)，即不显著遥相关作用占比最高。但是，一些特殊的空间异质性现象在图2b中无法得到体现。最后，图2c给出的局部Moran指数结果仅考虑了空间相关作用，而忽略遥相关本身的显著性，因此出现了和图2a-b截然不同的结果：不显著的网格占比在四个季节中都偏少。这是因为局部Moran指数考虑的是相比于全局的平均值的高低程度，且原假设不是遥相关是否显著，而是变量在空间上是否随机分布。因此，利用Moran指数进行分类虽然能够有效识别遥相关的空间特征，却难以区别显著的遥相关作用是否存在。

以上实验结果表明，本发明考虑了相关系数本身的标准化性质，同时结合变量的空间信息，实现对不同类型遥相关的空间聚类，结果便于对不同季节的遥相关程度进行横向对比。

附图中描述位置关系的用语仅用于示例性说明，不能理解为对本专利的限制；

显然，本发明的上述实施例仅仅是为清楚地说明本发明所作的举例，而并非是对本发明的实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动。这里无需也无法对所有的实施方式予以穷举。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明权利要求的保护范围之内。

Claims (5)

1.一种面向遥相关模式的空间自相关聚类方法，其特征在于，包括以下步骤：

S1：获取研究区域空间网格坐标信息，依据坐标信息计算空间权重矩阵；

S2：获取研究区域网格尺度降雨数据及相同时间范围的大尺度气象因子指标，得到降雨-气象指标时间序列；其中，记观测的气象因子指标序列为：

X＝[x_t]

记观测降雨序列为：

Y＝[y_t,i]

其中，x_t为第t年观测到的气象因子指标数值，X为x_t组成的集合；记网格总数为N，以i为网格索引；y_t,i为第t年在网格i观测到的降雨数值，Y为y_t,i组成的集合；

S3：依据所获得的降雨-气象指标时间序列，逐网格计算降雨-气象指标相关系数r；对于给定的网格i，由X和Y序列得到气象因子和降雨之间的相关系数：

式中，rx_t代表X序列中第t年气象因子指标数值在原始序列中从小到大的排序，ry_t代表Y序列中观测到的降雨数值在原始序列中从小到大的排序，而和/>则代表rx_t和ry_t的均值；相关系数r的值介于-1与1之间，当两个变量完全单调相关时，r的值为1或-1；对于网格i，所计算的相关系数可记为r_i，由此可得到相关系数组成的集合R＝[r_i]；

S4：根据相关系数r及空间权重矩阵，计算逐网格的空间自相关局部指标LISAAC；所述空间自相关局部指标LISAAC的计算公式具体为：

其中，w_i,j是用于将r_i与邻近网格r_j进行联系的权重系数；为网格i邻近网格的空间加权相关系数；C_i描述r_i和其邻近网格r_j之间的关系的强度；式中的权重系数是影响C_i的重要参数；

其中，为了使权重系数随距离增大而衰减，采用反距离权重法，以欧几里得法计算距离，并以该距离的平方的倒数作为两网格之间的距离权重系数：

式中d(i,j)为第i和j各网格之间的欧式距离；当d(i,j)增大，邻近网格相关系数r_j对C_i的影响越小，根据数据尺度预设距离阈值以减少计算量；

S5：重排步骤S2中的降雨及气象指标时间序列，获得新的降雨-气象指标时间序列，重复执行步骤S3-S4，直到达到预设的迭代次数n；

S6：根据步骤S4得到的n组随机LISAAC，构建参考经验分布H；其中，通过随机排列气象因子指标序列和第i或j个网格的观测降雨序列，重新计算r_i或r_j，及相应的C_i；通过多次重复计算，构建C_i的参考经验分布H；

S7：根据观测LISAAC、观测相关系数r以及参考经验分布H，得到在指定显著性水平下各网格的分类结果，用于对不同季节的遥相关程度进行横向对比。

2.根据权利要求1所述的一种面向遥相关模式的空间自相关聚类方法，其特征在于，在所述步骤S1的空间权重矩阵计算中，距离的表示方法为欧几里得距离。

3.根据权利要求1所述的一种面向遥相关模式的空间自相关聚类方法，其特征在于，在所述步骤S1的空间权重矩阵计算中，权重的表示方法为距离平方的倒数。

4.根据权利要求1所述的一种面向遥相关模式的空间自相关聚类方法，其特征在于，在所述步骤S5中，所述重排计算为直接随机重排每个网格的观测降雨-气象指标时间序列。

5.根据权利要求1所述的一种面向遥相关模式的空间自相关聚类方法，其特征在于，所述步骤S7具体为：

依据参考经验分布H，获得观测C_i对应的p_i值，以表示原假设为真时的强度,根据观测r_i和相应的p_i值，计算得到网格i的分类：

式中PP表示某一网格点和其周围网格点的r都显著为正；PN表示某一网格的r为正，但其周围网格的r值偏低或为负，ns代表网格r_i及其周围网格r_j均不显著；NP代表负r_i被正的r_j围绕；最后一种情况NN表示的是负r_i被负r_j围绕；PP和NN表示r值具有较高的空间正相关性，提示存在区域集聚；而PN和NN则反映了r的空间分布存在异质性。