CN108537247A

CN108537247A - 一种时空多元水文时间序列相似性度量方法

Info

Publication number: CN108537247A
Application number: CN201810203059.9A
Authority: CN
Inventors: 冯钧; 郭涛; 杭婷婷; 李晓东; 朱跃龙
Original assignee: Hohai University HHU
Current assignee: Hohai University HHU
Priority date: 2018-03-13
Filing date: 2018-03-13
Publication date: 2018-09-14
Anticipated expiration: 2038-03-13
Also published as: CN108537247B

Abstract

本发明公开了一种时空多元水文时间序列相似性度量方法，首先对洪水的原始降雨数据进行栅格化，生成每个小时的雨量分布矩阵图。然后计算两个降雨分布矩阵序列之间的2D‑DTW距离，包括两个降雨分布矩阵的相似性计算方法以及对降雨分布矩阵序列的相似性度量方法。并利用得出的多个有关标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列之间的距离，判断出与标准模板洪水水文过程最相似的是哪一条或者哪几条测试洪水，可以输出与模板洪水水文过程最相似的几场测试洪水数据。

Description

一种时空多元水文时间序列相似性度量方法

技术领域

本发明涉及一种时空多元水文时间序列相似性度量方法，属于水文数据理技术领域。

背景技术

在过去的几十年中，时间序列数据的挖掘在研究领域引起了人们的极大关注。这些研究对许多领域有影响，从生物学，物理学，天文学，医学，纳米和股票市场分析等。序列模式的最早概念是由Agrawal和Srikant提出的，时间序列挖掘为在序列模式挖掘在关联挖掘中加入了时间属性，用以挖掘事务之间在时间上的顺序联系，其作用是能够从数据集中发现可以反映事务间联系和规律的一些模式，进而预测事务的将来发展趋势。时间序列相似模式挖掘主要涉及3个问题： 1)特征表示；2)相似模式度量；3)相似性搜索。本发明主要关注于时间序列的相似模式度量问题。

关于时间序列的相似性度量研究，国内外研究人员做了很多相关的研究工作，但绝大部分研究都是针对一元时间序列，针对多元时间序列的研究还不够成熟。现实世界中的状态刻画通常涉及多个变量，与一元时间序列相比，多元时间序列更具普遍性。水文学是一个数据密集型的研究领域，长期的观测和时间积累了大量的水文数据，水文时间序列的数据挖掘工作越来越受到重视。水文中通过多种途径来收集不同变量不同时间不同地点的数据，即水文时间序列具有多元的特点；考虑到水文时间序列中涉及空间分布概念，例如：暴雨中心等概念，即水文时间序列又具有空间属性，所以水文时间序列本质上是一种时空多元时间序列。

一元时间序列的相似性度量中，经典的时间序列相似性度量分为两大类，静态的锁步度量(lock-step measures)和动态的弹性度量(elastic measure)。锁步度量是时间序列“一对一”的比较；弹性度量则是允许时间序列“一对多”的比较。最常见的锁步度量即我们熟知的欧式距离(Euclidean Distance)。弹性度量方法主要包括：动态时间规整(Dynamic Time Wrapping)和基于编辑距离(edit distance)的度量。最新研究提出了一个新的时间序列相似度的度量：shapelets 分类器。shapelets分类器通过寻找能具有代表性的子形状来代表一类事物， shapelets是时间序列中能够最大限度地表示一个类别的子序列，能明显的规避噪点。

目前，多元时间序列的相似模式度量方法主要包括：Euclidean距离、动态时间弯曲距离(DTW)、奇异值分解方法(SVD)、基于点分布特征(PD)的方法等。 1)Euclidean距离。Euclidean距离是最常用的一种距离度量，其计算复杂度不高，与序列长度成线性关系。然而，它对时间序列的各种形变比较敏感，且要求时间序列的长度必须相等。2)DTW距离。该算法基于动态规划(DP)的思想，它能够度量不等长时间序列的相似性，并且支持时间轴的弯曲和伸缩等形变，因此在时间序列相似模式度量中被广泛采用。3)SVD方法。它将时间序列中的变量视为随机变量，将各个时刻的记录值视为随机变量的样本点，以相关系数矩阵作为特征表示的基础，用扩展Frobenius范数来度量两个多元时间序列的相似性。SVD 方法能够体现变量间的相互关系，支持不同长度多元时间序列的相似模式度量。 4)PD方法。PD方法在三维空间中描述多元时间序列，并抽取局部重要点作为特征。然后根据局部重要点的统计分布，构建一个特征模式向量，并采用 Euclidean距离度量特征模式向量之间的差异。PD方法能够有效地刻画多元时间序列的形状特征，且支持不等长序列的相似模式度量，对小规模的多元时间序列能够得到相对较好的匹配效果。

传统的水文时间序列的相似性挖掘主要通过计算每个水文特征维度的相似性以及通过降维思想来度量整体时间序列的相似性。李士进针对水文时间序列的特点和当前多元时间序列相似性的分析研究的不足，提出了一种新的基于 BORDA计数法多元水文时间序列相似性度量计算方法，该方法首先根据多元水文时间序列数据复杂性进行主成分变换降维，然后利用线性分段斜率距离分别计算单一序列的相似性,最后采用基于BORDA计数法的投票方法获得多元水文时间序列整体相似性。陈志飞对于中小河流的水文时间序列采用基于符号化的方法进行数据挖掘，涉及水文中影响洪水过程的重要因素：降雨量、面雨量、土壤含水、流量、暴雨中心等多个水文特征维度。水文中通过多种途径来收集不同变量不同时间不同地点的数据，即水文时间序列本质上是一种时空多元时间序列。传统的多元时间序列相似性计算方法并不能很好对带有时空特征的水文时间序列进行度量。RishaadAbdoola和Barend van Wyk在大气湍流图像的处理中提出了2D-DTW算法对多帧图像取平均进行图像模糊的修正，涉及到多帧图像相似性问题，多帧图像可以清晰的表达时空概念，这与水文时空多元时间序列有很大相似性，但是缺点也很明显，领域的适用性是一个亟待解决的问题。

DTW思想：

假设现在有一个标准的参考模板R，是一个M维的向量，即R＝{R(1)，R (2)，……，R(m)，……，R(M)}，每个分量可以是一个数或者是一个更小的向量。现在有一个测试的模板T，是一个N维向量，即T＝{T(1)，T(2)，……， T(n)，……，T(N)}同样每个分量可以是一个数或者是一个更小的向量，注意M不一定等于N，但是每个分量的维数应该相同。

由于M不一定等于N，现在要计算R和T的相似度，就不能用以前的欧式距离等类似的度量方法了。DTW就是为了解决这个问题而产生的。

首先R中的一个分量R(m)和T中的一个分量T(n)的维数是相同的，它们之间可以计算相似度(即距离)。在运用DTW前，首先计算R的每一个分量和T中的每一个分量之间的距离，形成一个M*N的矩阵。(为了方便，行数用将标准模板的维数M，列数为待测模板的维数N)。然后举例说明计算步骤。

这个例子中假设标准模板R为字母ABCDEF(6个)，测试模板T为1234(4 个)。R和T中各元素之间的距离已经给出，如图1所示。

既然是模板匹配，所以各分量的先后匹配顺序已经确定了，虽然不是一一对应的。现在目的是要计算出测试模板T和标准模板R之间的距离。因为2个模板的长度不同，所以其对应匹配的关系有很多种，需要找出其中距离最短的那条匹配路径。现假设例子中满足如下的约束：当从一个方格((i-1,j-1)或者(i-1,j)或者 (i,j-1))中到下一个方格(i,j)，如果是横着或者竖着的话其距离为d(i,j)，如果是斜着对角线过来的则是2d(i,j)；其约束条件如下公式所示：

其中g(i,j)表示2个模板都从起始分量逐次匹配，已经到了M中的i分量和T中的j分量，并且匹配到此步是2个模板之间的距离。并且都是在前一次匹配的结果上加d(i,j)或者2d(i,j)，然后取最小值。所以将所有的匹配步骤标注后如图2所示。

例如g(1,1)＝4,当然前提都假设是g(0,0)＝0，就是说g(1,1)＝g(0,0)+2d(1,1)＝0+2*2＝4.g(2,2)＝9是一样的道理。首先如果从g(1,2)来算的话g(2,2)＝g(1,2)+d(2,2) ＝5+4＝9，因为是竖着上去的。

如果从g(2,1)来算的话是g(2,2)＝g(2,1)+d(2,2)＝7+4＝11，因为是横着往右走的。如果从g(1,1)来算的话，g(2,2)＝g(1,1)+2*d(2,2)＝4+2*4＝12，因为是斜着过去的。综上所述，取最小值为9，所有g(2,2)＝9，当然在这之前要计算出g(1,1),g(2,1) ,g(1,2)，因此计算g(i,j)也是有一定顺序的。其基本顺序可以体现在如图3所示，计算了第一排，其中每一个箭头表示最小值来源的那个方向。当计算了第二排后的结果如图4所示，最后计算的结果如图5所示。得到答案，即2个模板直接的距离为26，还可以通过回溯找到最短距离的路径，通过箭头方向反推回去。如图6所示。

2D-DTW算法思想：图片数据实际上是一个一个像素点组成，可以认为一张图片也就是一张像素矩阵。如果有两帧图像，假设每帧有n*m个像素点，第二帧图像仅仅是第一帧图像发生了少量的偏移，其实他们两是同一张图像，如何确定这两张图像是同一张图像呢？也即两帧图片之间的相似度如何计算？定义如下相似性矩阵：D_row

其中D₁₁代表第一帧图像的第一行像素点向量与第二帧图像的第一行像素点向量之间的距离，同理D₁₂代表第一帧图像的第一行像素点向量与第二帧图像的第二行像素点向量之间的距离，这里的距离使用DTW进行计算，最后使用动态规划计算D_row矩阵的最短累积距离得到两帧图像的相似程度。

发明内容

发明目的：为了解决现有的水文多元时间序列的相似性度量缺少空间维度的问题，本发明提供一种时空多元水文时间序列相似性度量方法。

技术方案：一种时空多元水文时间序列相似性度量方法，依次包含以下步骤：

(1)对洪水的原始降雨数据进行栅格化；

(2)通过步骤1生成的栅格化矩阵数据利用2D-DTW算法计算两个降雨矩阵序列中相互两帧矩阵中相互两行行向量的距离；

(3)使用步骤2计算出的两帧矩阵的相互两行行向量的距离生成两帧矩阵的距离矩阵

(4)使用DTW算法计算步骤3的两帧矩阵的距离矩阵的最短累积距离

(5)使用步骤4中所有计算过的多帧最短累积距离，生成标准模板降雨分布矩阵序列T_R与测试模板降雨分布矩阵T_T的距离矩阵

(6)对步骤5生成的距离矩阵使用DTW算法进行计算，得到最小累计距离；

(7)通过步骤6计算的最小累计距离的大小，进行寻找相似降雨序列模式，并根据降雨序列与洪水水文序列的对应关系输出相似洪水数据。

有益效果：本发明提供的一种时空多元水文时间序列相似性度量方法，由于将原始洪水的降雨数据格式化生成矩阵数据，体现多元洪水数据的时空分布特征；由于自定义了解决时空序列相似度度量的新算法，相似模式输出模块可以对降雨分布矩阵序列进行相似度的计算工作；对降雨分布矩阵序列相似性结果进行整理，输出相似洪水模式的结果(降雨分布矩阵序列与洪水水文序列的具有对应关系，不同的降雨分布导致不同的洪水水文过程)，完成时空多元水文时间序列相似性的度量。

附图说明

图1 DTW计算标准模板R与测试模板T最短累计距离过程示意图一；

图2 DTW计算标准模板R与测试模板T最短累计距离过程示意图二；

图3 DTW计算标准模板R与测试模板T最短累计距离过程示意图三；

图4 DTW计算标准模板R与测试模板T最短累计距离过程示意图四；

图5 DTW计算标准模板R与测试模板T最短累计距离过程示意图五；

图6 DTW计算标准模板R与测试模板T最短累计距离过程示意图六；

图7为本发明方法的数据流图；

图8为本发明方法中的流域栅格化降雨分布矩阵；

图9为本发明方法中的流域栅格化降雨分布矩阵序列；

图10为本发明方法中的2D-DTW距离计算模块计算示意图。

具体实施方式

下面结合具体实施例，进一步阐明本发明，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

如图7所示，时空多元水文时间序列相似性度量方法，依次包含以下步骤：

(1)对洪水的原始降雨数据进行栅格化；

(6)对步骤5生成的距离矩阵使用DTW算法进行计算，得到最小累计距离

(7)通过步骤6计算的最小累计距离进行寻找相似降雨序列模式，并根据降雨序列与洪水水文序列的对应关系输出相似洪水数据。

原始降雨数据处理模块用于对洪水的原始降雨数据进行栅格化，生成每个小时的雨量分布矩阵图。根据该流域的雨量站点获取的时段降雨量数据，将该流域降雨数据分别填入该流域的栅格化矩阵中，不在流域范围内的地区作补零处理，这样即生成了时段降雨的栅格化矩阵。由于一场洪水历时会有多个小时，对应的时段降雨的栅格化矩阵会有多个，这样就形成了一场洪水的栅格化降雨分布矩阵序列。

降雨数据是洪水水文过程的演变的重要决定因素，降雨量的累积可以间接反映土壤含水量程度，降雨分布可以间接反映出面雨量的大小及暴雨中心的位置，所以在二维平面的降雨时间序列本质上是一个涵盖多个变量维度的时空多元时间序列。如图8，某流域栅格化之后的形状为灰色部分单元格所示，一个单元格边长代表1KM实际距离，栅格化之后的流域图为10行乘以8列的矩阵。根据该流域的雨量站点获取的时段降雨量数据，将该流域降雨数据分别填入该流域的栅格化矩阵中，不在流域范围内的地区作补零处理，这样即生成了时段降雨的栅格化矩阵。由于一场洪水历时会有多个小时，对应的时段降雨的栅格化矩阵会有多个，这样就形成了一场洪水的栅格化降雨分布矩阵序列，如图9所示。

计算两个降雨分布矩阵序列之间的2D-DTW距离。首先从两个降雨分布矩阵的相似性度量方法开始描述，之后在对降雨分布矩阵序列的相似性度量方法进行描述。

1.降雨分布矩阵的相似性度量方法

两个降雨分布矩阵M_R和M_T，M_R为标准模板降雨分布矩阵，M_T为测试模板降雨分布矩阵，都为n行m列的降雨分布矩阵，如何定义他们之间的相似度，可以类比于2D-DTW的定义。

定义1：两个降雨分布矩阵的距离矩阵

构成的矩阵为n行n列的相似性矩阵。其中D_R1T1代表标准降雨分布矩阵M_R的第一行行向量M_R1与测试模板降雨分布矩阵M_T的第一行行向量M_T1之间的距离(这里的行向量为1*m的矩阵)，同理D_R1T2代表标准降雨分布矩阵M_R第一行行向量M_R1与模板降雨分布矩阵M_T的第二行行向量M_T2之间的距离，这里的距离使用DTW距离进行计算。这里降雨分布矩阵的行向量在水文中代表为将流域分为n块长条区域，每块长条面积为1*m(km²)，对两个不同时间的长条区域降雨情况进行相似度的计算。图8为10行8列的降雨分布矩阵(m＝10,n＝8)，两个降雨分布矩阵构成的矩阵即为10行10列的相似性矩阵，其中的行向量为1*8的矩阵。

定义2：距离矩阵的距离

为矩阵中的每条规整路径对应的累积距离，取最小的累积距离作为对应距离矩阵的距离，也即降雨分布矩阵M_R与降雨分布矩阵M_T之间的距离。水文中代表两个时间点该流域降雨情况的相似程度。

2.降雨分布矩阵序列的相似性度量

两个降雨分布矩阵序列的相似性度量方法如下：

定义3:标准模板降雨分布矩阵序列：

T_R＝{M_R1-t1,M_R2-t2,M_R3-t3,......M_Rn-tn}。

定义4:测试模板降雨分布矩阵序列：

T_T＝{M_T1-t1,M_T2-t2,M_T3-t3,......M_Tn-tn}。

T_R为标准模板降雨分布矩阵序列，由n个标准模板降雨分布矩阵组成。 M_R1-t1表示t1时刻标准模板降雨分布矩阵为M_R1，同理M_R2-t2表示t2时刻标准模板降雨分布矩阵为M_R2。

T_T为标准模板降雨分布矩阵序列，由n个标准模板降雨分布矩阵组成。M_T1-t1表示t1时刻标准模板降雨分布矩阵为M_T1，同理M_T2-t2表示t2时刻标准模板降雨分布矩阵为M_T2。

定义5:降雨分布矩阵序列标准模板与测试模板分布矩阵序列的距离：

为标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的距离。

DTW(T_R,R_T)为使用DTW算法计算标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的距离。两个模板序列在水文中代表两场不同的降雨时间序列，DTW(T_R,R_T)表示两场降雨情况的相似程度。

定义6:标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的距离矩阵

构成的矩阵为n行n列的矩阵。其中代表标准降雨分布矩阵序列T_R的t1时间点的矩阵与测试模板降雨分布矩阵序列T_T的t1 时间点的矩阵之间的距离，同理代表标准降雨分布矩阵序列 T_R的t1时间点的矩阵与测试模板降雨分布矩阵序列T_T的t2时间点的矩阵之间的距离这里的距离使用定义2中的公式进行计算。

3.2 D-DTW距离计算模块算法过程

下面通过示意图10来详细说明算法流程：

2D-DTW距离计算模块对降雨分布矩阵序列相似度的计算分为两层，第一层为利用DTW算法进行降雨分布矩阵序列相似度的计算，第二层利用2D-DTW算法进行降雨分布矩阵之间的相似度计算。

第一层：图中假定标准模板降雨分布矩阵序列和测试模板降雨分布矩阵序列的时间序列为t1-t4，也即总数为4的降雨分布矩阵序列。因为所以需要计算降雨分布矩阵序列T_R与T_T的DTW距离，类似一维DTW算法，需要构造一个n*m的矩阵(这里的n代表标准模板降雨分布序列的长度，m代表测试模板降雨分布序列的长度)，这里示意图的降雨分布矩阵序列长度都为4，所以需要构造一个4*4的矩阵，如图10中(a)所示。d(1,1) 表示标准模板降雨分布矩阵序列t1时刻的降雨分布矩阵M_R1-t1和测试模板降雨分布矩阵序列t1时刻的降雨分布矩阵M_T1-t1之间的距离。类似的d(2,4)表示标准模板降雨分布矩阵序列t2时刻的降雨分布矩阵M_R2-t2和测试模板降雨分布矩阵序列t4时刻的降雨分布矩阵M_T4-t4之间的距离。两降雨分布矩阵之间的距离使用定义2中的公式进行计算。得到计算好的4*4矩阵之后求得最小累积距离的规整路径，该累计距离即为标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的相似程度。

第二层：图10中(b)所示矩阵为第一层d(2,4)距离计算的展开描述。两个降雨分布矩阵的距离d'(1,1)代表降雨分布矩阵R2 的第一行行向量R2-row1与降雨分布矩阵T4的第一行行向量T4-row1的DTW 距离。两个降雨分布矩阵的距离矩阵由定义1计算生成。

Claims

1.一种时空多元水文时间序列相似性度量方法，其特征在于，依次包含以下步骤：

(1)对洪水的原始降雨数据进行栅格化；

2.如权利要求1所述的时空多元水文时间序列相似性度量方法，其特征在于，步骤(1)对洪水的原始降雨数据进行栅格化，生成每个小时的雨量分布矩阵图；根据流域的雨量站点获取的时段降雨量数据，将该流域降雨数据分别填入该流域的栅格化矩阵中，不在流域范围内的地区作补零处理，这样即生成了时段降雨的栅格化矩阵；由于一场洪水历时会有多个小时，对应的时段降雨的栅格化矩阵会有多个，这样就形成了一场洪水的栅格化降雨分布矩阵序列。

3.如权利要求1所述的时空多元水文时间序列相似性度量方法，其特征在于，设两个降雨分布矩阵M_R和M_T，M_R为标准模板降雨分布矩阵，M_T为测试模板降雨分布矩阵，都为n行m列的降雨分布矩阵，两个降雨分布矩阵的距离矩阵

构成的矩阵为n行n列的相似性矩阵；其中D_R1T1代表标准降雨分布矩阵M_R的第一行行向量M_R1与测试模板降雨分布矩阵M_T的第一行行向量M_T1之间的距离，同理D_R1T2代表标准降雨分布矩阵M_R第一行行向量M_R1与模板降雨分布矩阵M_T的第二行行向量M_T2之间的距离，这里的距离使用DTW距离进行计算；这里降雨分布矩阵的行向量在水文中代表为将流域分为n块长条区域，每块长条面积为1*m(km²)，对两个不同时间的长条区域降雨情况进行相似度的计算；

两帧矩阵的距离矩阵的距离为：

为矩阵中的每条规整路径对应的累积距离，取最小的累积距离作为对应距离矩阵的距离，也即降雨分布矩阵M_R与降雨分布矩阵M_T之间的距离；水文中代表两个时间点该流域降雨情况的相似程度。

4.如权利要求1所述的时空多元水文时间序列相似性度量方法，其特征在于，定义标准模板降雨分布矩阵序列T_R及测试模板降雨分布矩阵序列T_T：

T_R＝{M_R1-t1,M_R2-t2,M_R3-t3,......M_Rn-tn}；

T_T＝{M_T1-t1,M_T2-t2,M_T3-t3,......M_Tn-tn}；

T_R为标准模板降雨分布矩阵序列，由n个标准模板降雨分布矩阵组成；M_R1-t1表示t1时刻标准模板降雨分布矩阵为M_R1，同理M_R2-t2表示t2时刻标准模板降雨分布矩阵为M_R2；

T_T为标准模板降雨分布矩阵序列，由n个标准模板降雨分布矩阵组成；M_T1-t1表示t1时刻标准模板降雨分布矩阵为M_T1，同理M_T2-t2表示t2时刻标准模板降雨分布矩阵为M_T2。

5.如权利要求1所述的时空多元水文时间序列相似性度量方法，其特征在于，定义降雨分布矩阵序列标准模板与测试模板分布矩阵序列的距离：

为标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的距离；

DTW(T_R,R_T)为使用DTW算法计算标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的距离；两个模板序列在水文中代表两场不同的降雨时间序列，DTW(T_R,R_T)表示两场降雨情况的相似程度。

6.如权利要求1所述的时空多元水文时间序列相似性度量方法，其特征在于，定义标准模板降雨分布矩阵序列与测试模板降雨分布矩阵序列的距离矩阵：

构成的矩阵为n行n列的矩阵；其中代表标准降雨分布矩阵序列T_R的t1时间点的矩阵与测试模板降雨分布矩阵序列T_T的t1时间点的矩阵之间的距离，同理代表标准降雨分布矩阵序列T_R的t1时间点的矩阵与测试模板降雨分布矩阵序列T_T的t2时间点的矩阵之间的距离这里的距离使用公式进行计算。