一种基于气象数据的监测指标时空延迟相关性计算方法

摘要

本发明公开了一种基于气象数据的监测指标时空延迟相关性计算方法。本方法为：1)采集各地的气象数据和监测指标数据；2)对步骤1)的采集数据进行数据对齐，保留同一时刻同时存在的监测指标数据和气象数据；3)从步骤2)处理后的采集数据中分别提取监测地点A、B的设定时间段T内的监测指标数据及气象数据；4)保持地点A的监测指标数据X不变，地点B的监测指标数据Y向前多取t

说明书

技术领域

本发明属于气象学的大数据应用领域，具体涉及一种针对不同地理位置的两个地点，计算此两地监测指标在时空维度上延迟相关性的方法。

背景技术

相关性分析方法是对两个或两个以上的变量进行一系列的统计分析，以考察变量之间所具有的线性关系。通俗地讲，相关性分析的目的是探索当一个变量变化时，另一个变量如何变化，这种相互关系需要用相关性分析这样的统计方法来进行定量的描述。在相关性的计算中，“相关系数r”则是两个变量之间相关程度的定量指标。在分析过程中，不仅要考察相关程度的大小，还需要考虑数据之间相关的方向，即r存在正负值，正值代表正相关，负值则代表负相关，所以r的取值为-1到1，其绝对值越大代表相关程度越高。

在实际分析中，相关性的大小不能直接以“相关系数”下定论，即不能从相关系数r直接定论两个变量的相关性，从而需要引进了“显著性水平α”和假设检验中常用到的p值。“显著性水平α”是用于判断被检验的统计量与假设之间是否差异明显，确切地讲，在相关性分析中，显著性检验是指样本的相关系数与总体相关系数的差异。而假设检验中的p值代表概率，反映着某一事件发生的可能性大小，在这里则代表“样本间差异是由抽样误差导致的”这一事件的的概率。在实际操作中，一般认为显著性水平α＝0.05是统计学意义的边界线，根据假设检验方法所得到的p值，通常认为如果p≤0.05，则说明相关系数r可以接受。值的指出的是，相关关系不能直接混淆为因果关系，如果两个变量之间出现了很强的相关性且通过显著性检验，也无法直接说明这两个变量存在着因果的关系，而只能说明两个变量之间存在较多的相同特征。

时间序列的延迟相关性分析是指，两个序列在同一时间段内的数据不具有最大的相关性，而在某一序列往后延迟多个时刻后，两个序列重叠时间段内的数据具有最大的相关性。随着延迟数的增加，两序列重叠时间段逐渐变短，若重叠时间段太短，则所计算出的相关性也没有实际应用价值，所以，一般设置序列最大可延迟数为序列长度的一半。在实际生产中，延迟相关性被运用到了多种领域，例如气象学分析，空气质量分析以及金融市场分析等多方面，延迟相关性分析对实际生产生活起到了很优秀的辅助决策作用。

空间延迟相关性分析是指针对不同的地理位置地点，它们的某些指标例如空气质量，受到其他空间性因素例如风速的影响，从而使得这些地点的空气质量在时间上呈现出某些延迟相关的特性。

对于两个时间序列的延迟性分析，业界提出了多种科学的探测方法；然而对于不同时空下的多点空间延迟相关性分析，最主要面临的问题是如何确定相关性计算的窗口，特别是在处理真实的时空数据时，时间和空间的特性都需要在分析过程中被考虑。

发明内容

本发明针对上述问题，提出了一种针对气象和监测指标数据进行时空延迟相关性分析的方法。本方法不仅考虑了监测指标时间序列数据在时间上的相关性特性，同时加入考虑了两地在空间上的特性，从而确定了一种新颖的延迟相关性计算窗口的定义方法，然后再进一步进行延迟相关性分析。

在本发明中，监测指标指各种在空间分布中的监测站点所监测的各种时序数据，例如PM2.5,PM10,氮氧化物等。

具体来说，本发明采用的技术方案为：

一种基于气象数据的监测指标时空延迟相关性的计算方法，主要步骤包括：

1)采集各地的气象数据和监测指标数据；

2)根据1)中采集的数据，对数据进行预处理工作，剔除空数据和错误数据，以及剔除所有没有成对的监测指标数据和气象数据的时刻下的数据；一般情况下，监测指标数据与气象数据在采集数据时，具有不同的时间间隔，如指标数据每隔一小时采集一次数据，而气象数据可能每隔两个小时采集一次数据，这样情况下，需要将两类数据对齐，即只取同一采集时刻具有监测指标数据和气象数据的采集数据。

3)确定一个将要分析的时间段T，其长度为t₀；

4)确定两个地点A和B，取出两个地点在时间段T内所对应的监测指标数据X和Y以及气象数据(如风速数据)；

5)第一个地点A的监测指标数据X不变，第二个地点B的监测指标数据Y向前多取t₁个时刻的监测指标数据，t₁的计算方法为A地气象指标影响B的时间(如风速数据，则t₁为两地距离除以B地当时风速)，同时再向后多取t₂个时刻的监测指标数据，t₂的计算方法为A地被B地气象指标影响的时间(如风速数据，则t₂为两地距离除以A地当时风速)；

6)此时，计算所取的B点监测指标数据Y对应的时间长度为(t₀+t₁+t₂)，B点监测指标数据Y从数据起始，以一个时刻为单位向后依此进行窗口移动，总计移动d个时刻，d取值范围为[0,t₁+t₂]。为保证A、B两组数据在计算时具有相同数目的时序数据值，A、B两点在一次计算时时间范围(窗口长度)一直为检测指标数据X的长度t₀；

7)在6)每次窗口移动时，计算Y移动后窗口内的数据与未移动的X数据之间的延迟相关性，相关性计算采用皮尔森相关系数计算方法，取循环移动后相关系数最大值和对应的移动数值d_max作为两地的监测指标的延迟相关性。在该步骤中，窗口是指计算时所取的时间范围t₀，窗口滑动，是指起始时刻从数据起始，到(t₁+t₂)，窗口每次向后进行一个时刻的移动。

8)根据t₁-d_max的值确定A、B两地之间的影响关系。正值的t₁-d_max代表A地在延迟|t₁-d_max|个时刻后与B地相关，负值的t₁-d_max代表B地在延迟|t₁-d_max|个时刻后与A地相关。

与现有方法相比，本发明的积极效果为：

现有的方法在计算延迟相关性时，没有考虑监测指标数据的空间特性，本方法扩展了相关性分析方法到真实时空数据的处理，加入考虑了数据的空间特性，即将两地之间互相影响的时间，做为计算窗口滑动的范围。同时，现有的方法在计算延迟相关性时，只计算重叠时间段内数据的相关性，重叠时间段内的数据数量随着延迟移动逐渐减小，当数据量较小时计算出的相关性没有实际意义，然而本方法每次计算的时间窗口的大小不变，一直为原始时间段长度，从而计算出的相关性更准确，多个相关性数值之间更具有可比较性(原因：并非只计算重叠时间段，而是每次窗口滑动时，均计算相同时间长度的相关性)。此外，普通延迟相关性计算方法，只可单向计算两地的延迟相关关系，然而本方向可以双向计算出两地的相关关系，即，可以通过最优相关系数对应的窗口移动的时刻数d，来确定是A地延迟于B地还是B地延迟于A地。

附图说明

图1是本发明方法流程图。

具体实施方式

为使本发明的上述目的、特征和优点能够更加明显和清晰，下面给出实施例，对本发明做进一步说明。在本实施中，我们选择风速作为气象数据，选择空气质量作为监测指标数据。

首先采集数据各地的气象小时数据和空气质量的小时数据，数据格式(示例，真实数据可能与该格式不同，但都具有时间点及对应的检测数值)分别如下所示。

同时采集发布数据的气象观测站和空气质量观测站的地理坐标(经纬度)，(示例，真实数据可能与该格式不同，但都具有站点标识及对应的地理位置)数据格式分别如下所示；

所有数据以bson数据形式存储在MongoDB数据库中，本发明的处理流程如图1所示。

对数据进行预处理和清洗工作。由于在之后计算时需要同一时刻同时存在气象数据和空气质量数据，所以需要清除掉同一时刻下气象数据不存在的空气质量数据，以及同一时刻下空气质量数据不存在的气象数据。并且，还需要清除掉空数据以及错误由于观测站设备原因发布的错误数据。

数据清洗完毕后，选定两个地点A和B，以及确定一个需要分析的时间段T，分别取出A和B两地在此时间段内的空气质量PM2.5小时数据，形成两个时间序列X和Y。

本方法特别加入考虑了空气质量数据在空间上的传播特性，所以引入了气象数据(风速)对空气质量的时空特性进行研究。首先设置A地的空气质量数据X不变，而对在B地的Y序列，在其前部额外多取T开始时刻前的t₁个时刻数据，t₁的计算方式为A、B两点的直线距离除以B地在T开始时刻全天的平均风速，数学表示如下：

$t_1=\frac{dis\tan ce(A,B)}{windSpeed\left(B\right)}$

同时，在Y序列的后部额外多取T结束时刻后的t₂个数据，t₂的计算方式为A、B两点直线距离除以A地在T结束时刻全天的平均风速，数学表示如下：

$t_2=\frac{dis\tan ce(A,B)}{windSpeed\left(A\right)}$

此时，Y序列的长度比X序列多出的时刻数据个数为：

t＝t₁+t₂

一般地，对于计算两个普通序列的相关系数，最常采用的是“皮尔森相关系数”计算方法，也被称作“皮尔森积差相关系数”，其定义如下：

$r_{XY}=\frac{{\Sigma }_{i=1}^N(x_i-\bar{x})(y_i-\bar{y})}{\sqrt{{\Sigma }_{i=1}^N(x_i-\bar{x}){\Sigma }_{i=1}^N(y_i-\bar{y})}}$

其中r_XY的取值范围是-1到1，正负号代表相关的方向是正相关或负相关，其绝对值越大代表相关程度越高，和分别代表序列X和Y的平均值，x_i和y_i分别为X和Y中第i个值。两个时序数据的延迟相关性计算的一般思路为，一个序列依次后移，计算两个序列重叠部分的数据的相关性，定义如下：

$r_{XY}\left(d\right)=\frac{{\Sigma }_{i=1}^{N_d}(x_{i+d}-\bar{x})(y_i-\bar{y})}{\sqrt{{\Sigma }_{i=1}^{N_d}(x_{i+d}-\bar{x}){\Sigma }_{i=1}^{N_d}(y_i-\bar{y})}}$

$\bar{x}=\frac{1}{N_d}\underset{i=1}{\overset{N_d}{\Sigma }}{x}_{i+d},\bar{y}=\frac{1}{N_d}\underset{i=1}{\overset{N_d}{\Sigma }}{y}_i,N_d=N-d$

其中d为Y序列往后移动的时间间隔数，其取值范围业界一般规定其最大值为序列长度的一半，而N_d为移动后X和Y两序列的公共部分长度；此种针对普通时序数据的延迟相关性计算方法。然而普通的延迟相关性计算方法并不适用于真实时空数据的延迟相关性计算，而本方法特别考虑了空气质量数据在空间上的特性，空间特性由气象数据中的风速和地理距离给出。本方法主要思想为：首先通过风速和两地距离计算出Y序列额外需要向前多取的时刻数t₁和向后多取的时刻数t₂，分别取出这些时刻对应的监测指标数据并添加在Y中，之后再对X和Y两个时间序列数据进行“延迟相关性”计算，计算方法为：设置Y上的窗口大小等于X的长度，将Y的窗口向后移动d个时间间隔数，即Y往后延迟一定的时间间隔数d后，从Y延迟后的起始点开始，依次取出与X等长的数据，即取出了Y的窗口内的所有数据，然后再计算两个序列的相关系数。本方法的数学定义如下：

$r_{XY}\left(d\right)=\frac{{\Sigma }_{i=1}^{N_X}(x_i-\bar{x})(y_{i+d}-\bar{y})}{\sqrt{{\Sigma }_{i=1}^{N_X}(x_i-\bar{x}){\Sigma }_{i=1}^{N_X}(y_{i+d}-\bar{y})}}$

$\bar{x}=\frac{1}{N_X}\underset{i=1}{\overset{N_X}{\Sigma }}{x}_i,\bar{y}=\frac{1}{N_X}\underset{i=1}{\overset{N_X}{\Sigma }}{y}_{i+d}$

其中N_X为X序列的长度，d为Y序列往后移动的时间间隔数，取值范围为[0，t]，x_i为X中第i个值，y_i+d为Y中的第i+d个值(即在窗口移动d个时间间隔数后，Y在窗口内的第i个值)。在计算完毕所有d值对应的延迟相关性后，将r_XY(d)最大值记录，其所对应的d记作d_max，分别作为A、B两地间的最优相关系数和最优延迟数。在本方法中，根据t₁-d_max的值可确定A、B两地之间的影响关系。正值的t₁-d_max代表A地在延迟|t₁-d_max|个时刻后与B地相关，负值的t₁-d_max代表B地在延迟|t₁-d_max|个时刻后与A地相关。

下面给出一个具体的算例：

1)选取两个地点A(位于石家庄市)和B(位于北京市)，两地点对应的空气质量观测站点的信息如下：

2)选择要分析的时间段T为：2015年6月1日00:00点至2015年6月3日00:00点，数据以小时为时间间隔，T的长度为t₀＝48；

3)从数据库中取出A和B两地在时间段T内的空气质量数据集X和Y(以PM2.5数值为例)，分别为：

X：{130 146 137 122 122 81 90 90 90 99 117 116 134 134 161 136 149 90 76 53 55 61 68 79 63 62 59 57 52 42 33 31 23 27 30 29 31 32 35 39 40 44 56 81 81 69 74 56}，

Y：{88 68 90 93 93 123 119 119 86 84 45 25 59 59 102 141 153 73 53 84 74 76 58 3 3 9 7 3 3 6 4 5 3 6 3 9 5 14 158 36 39 34 32 30 34 28 30}；

4)查询位于B地的气象观测站发布的2015年6月1日00:00点的气象数据(以风速为例)，此时B地风速为8.43km/h；查询A地的气象观测站2015年6月3日00:00点的气象数据(以风速为例)，此时A地风速为9.89km/h；

5)根据A、B两地经纬度计算得到两地空间距离为253.90km，可以计算出Y需要向前多取t₁＝253.90/8.43≈30个时刻的空气质量数据，以及Y需要向后多取t₂＝253.90/9.89≈26个时刻的空气质量数据。取出这些数据后，Y一共含有104个数据，此时Y变成了：

Y：{60 67 72 72 44 46 47 59 62 57 57 45 42 41 44 44 53 54 54 59 52 53 72 78 75 89 92 78 86 87 88 68 90 93 93 123 119 119 86 84 45 25 59 59 102 141 153 73 53 84 74 76 58 3 3 9 7 3 3 6 4 5 3 6 3 9 5 14 158 36 39 34 32 30 34 28 30 33 32 39 38 23 23 22 18 17 10 19 21 24 27 28 33 38 43 50 47 45 43 40 48 46 44}；

6)设置Y上的窗口长度等于X的长度48，将窗口依次向后移动d个时间间隔，并根据本发明提出的延迟相关性计算公式计算此时Y窗口中的数据和X的数据的时序相关性r_XY(d)，d的取值范围为[0,30+26]即[0,56]。

7)根据6)可以分别计算得到57个相关性数据，记录57个数据中的最优延迟相关性结果为r_XY(33)＝0.785，则最优延迟数为d_max＝33。根据t₁-d_max＝30-33＝-3可以得出结论：B地在延迟3个小时后与A地的空气质量最相关。

在本方法中，克服了普通方法只能利用重叠时间段计算两序列的延迟相关性，随时延迟的依次增大，重叠时间段内的数据量越来越小，当参与计算相关性的数据量过小的时候，容易导致所计算出的相关性不再具有实际意义。同时，本方法可以通过t₁-d_max双向判定两地之间的空气质量的延迟相关关系，而不再局限于被移动序列的选定。

以上实施例仅用以说明本发明的技术方案而非对其进行限制，本领域的普通技术人员可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明的精神和范围，本发明的保护范围应以权利要求所述为准。