最近阅读了一篇paper"Searching and Mining Trillions of Time Series Subsequences under Dynamic Time Warping" 在此做一个Summary小节
该论文主要讨论了一种如何在大数量级(Trillions)的时间序列数据上应用基于DTW的搜索算法,来找出相似度高的子串。
一、基础概念的梳理:
1.什么是时间序列?
时间序列是一种带有顺序关系的数据串,经常应用于,语音识别,脑电波分析,气温分析,股票数据分等领域,文本也能算是一种时间序列。
时间列序列T定义如下:
其中
2. 什么是DTW(Dynamic Time Wrapping)?
说到DTW必须得先讲一讲“欧几里得距离”,那么什么是欧几里得距离呢?两个点之间的欧几里得距离定义为两个点的所有维度坐标差的平方之和,再次基础上再求一个根号。他是一种用来衡量在平面空间中坐标点之间的“平面距离”,当然这里的平面不限定于2维,可以扩展到高维。
如果要计算两个序列的相似度,最直观的计算方法是就算两个序列之间所有对应点的欧几里得距离。假设有两个序列Q和C,并且他们的长度都为n,那么他们之间的欧几里得距离计算如下:
![a04f5c807fc406f0ac5f602b3b6131d0.png](https://i-blog.csdnimg.cn/blog_migrate/6073ef32aafe3f20a490bf5e4ddb92ee.jpeg)
可是如果两个序列的长度不同怎么办,又或者如果两个时间序列之间存在一个误差(offset)的情况?这种基于欧几里得距离的计算就会存在很大的误差。因此有学者提出了一种基于warp(扭曲)对应关系的距离计算Dynamic Time Warping。
DTW距离定义为,两个序列经过扭曲(warp)以后最小的欧几里得距离,是一种基于动态规划的计算。DTW定义如下:
其中<>表示一个空的序列, First(C) 表示C序列的第一个元素, Rest(C)表示C序列中出去第一个元素的所其他元素组成的子序列,Last(C)表示,C序列中的最后一个元素。他是一种基于动态规划的算法。其中Dbase可以使用任意的Lp函数来计算他们之间的距离。Lp定义如下:
当P值取2的时候,就是使用欧几里得距离,
直观上DTW就是在一个|C| x |Q|的矩阵上寻找最优短距离的问题,这方面很类似计算字符串的最短编辑距离的问题:
![610a40b59df185260cf349a0f60651d0.png](https://i-blog.csdnimg.cn/blog_migrate/1874861af8cff4e973f5fa5c0d830a85.jpeg)
红色的轨迹就是经过DTW计算出来两个时间序列对应的最短距离。
那么DTW和欧几里得距离有什么区别呢?可以看下面这幅图:
![4349a70ae91672d7b4769bc6a4d3e7c1.png](https://i-blog.csdnimg.cn/blog_migrate/809ba34569bb2e2a7abb1b76eccd6683.jpeg)
可以看到欧几里得距离要求两个序列长度必须相等,且是垂直方向的一一对应;而DTW距离允许一对多,并且对应关系不是完全垂直方向上的, 可以有一个offset倾斜(这个在计算DTW的时候可以设置一个最大允许的offset偏差的大小)
3. 什么是Lower Bounding计算?
可以看出每一次DTW的计算是
实际上我们不需要每次都完整计算两个序列之间的Lower Bounding,因为有序列的匹配明显是不符合最优解的情况的,我们可以对他进行排除(论文里叫Early Abandoning),这个什么意思呢?
可以这样来理解,比如我们开始对序列Q进行匹配,先计算好一次DTW距离,存在一个叫“best-so-far”的变量里,表示目前为止最好(最小)的DTW距离,在后面的计算中,我们使用一种叫Lower Bounding的计算,并且这种计算一定满足Lower Bounding距离 < DTW距离,而且这种Lower Bounding的复杂度是远低于DTW距离的(通常复杂度是
![c461e3af977224a1a1d3739a6dbc7c6a.png](https://i-blog.csdnimg.cn/blog_migrate/821d467bb2a6b62d7fc8bc8e0392363b.jpeg)
关于Lower Bounding的计算方法常用的有如下:
1)Lower Bounding Kim 简称 LB_kim(由学者Sang-Wook Kim在2001年提出)。
定义如下:
就是分别计算序列C和Q中的第一个元素值,最后一个元素值,最大值和最小值。然后对他们分别取距离并取最大的那个值那么:
2) 还有一种叫LB_Keogh的方法(由学者Keogh在2004年提出)
定义如下:
![a4bb6916cfc71e1a406f8cd12f92a2cd.png](https://i-blog.csdnimg.cn/blog_migrate/086dba1da99d76d651bef6106befd781.jpeg)
其中U和L的定义如下:
r是自己定义的一个窗口大小,这里的U和L直观的理解是在原始序列的周围包装另一个“封套”类似下图的效果:
![11be7cfa34066cf36d2179892522f3e7.png](https://i-blog.csdnimg.cn/blog_migrate/dc8dc07b82d7a0f09cc31820edfdaf62.jpeg)
计算Q和C之间的LB_keogh实际上等效于下图的阴影部分的面积:
![8baf00ddcc5818aa04bb9ecbc8bbe3c0.png](https://i-blog.csdnimg.cn/blog_migrate/9e10257fb325807c8d3e06190f74d39b.jpeg)
以上就是这篇paper中涉及到的一些基础知识。
二、这篇论文主要解决了什么问题?
基于时间序列的算法,比如计算序列相似度,或者进行挖掘的时候,通常需要计算序列或者子序列的相似度。而相似度计算中比较经典的算法是Dynamic Time Wrapping(他是一种基于动态规划和欧几里得距离来计算两个序列的相似度,有点类似于计算两个字符串的编辑距离)。但是,当序列的数量级很大的时候,欧几里得距离的计算就成为了一种性能瓶颈。例如,当序列的长度达到10^12级别以后,对序列数据的索引也成为了一种困难,因为序列将无法被完全加载到内存中。而且通常为了保证计算的精确度还需要对序列数据进行标准化操作(叫z-score Normalize),但标准化操作需要加载完整个序列并计算均值和标准差,对于大的时间序列这也是一种困难。
虽然目前已经有一些方案能处理百万级别长度的时间序列,但是他们是基于近似计算的,而该论文讨论的主要范畴是:基于DTW的“精确”计算时间序列的相似度,并且能处理的数据长度是目前主流论文长度的5~6倍。
作者从以下几个方面一步一步提出了对传统的DTW算法进行了优化:
1)欧几里得距离计算的优化
时间序列的相似度加算中通常需要计算欧几里得距离,部分一中介绍的传统欧几里得距离需要进行平方求和以后再计算平方根,其实平方根这一步骤是不影响最终的结果的,这里做了一个省略。
2)分层Lower Bounding的计算
下图是常用的两种Lower Bounding的计算:
![74abc8675f9a3c1dabb4128b6f09ffad.png](https://i-blog.csdnimg.cn/blog_migrate/d6077cc2bec583a8ad7a08f88764d11d.png)
上边左是一个简化的LB_kim的计算, LB_Kim从序列中抽取了First, Last, Greatest, Least4个特征,需要
上边右图是LB_keogh计算,需要复杂度O(N).
作者在计算的时候采用了分层Lower Bounding的策略,也就是先计算LB_kim如果发现其大于best-so-far值则直接排除,接着计算LB_keogh,如果发现它大于best-so-far也直接排除当前匹配,最后才进行DTW的计算。根据作者提供的源码进行了测试(实际上最终匹配的所有串中,50%以上的计算可以通过LB_kim过滤; 而LB_keogh占了30%;LB_keoghEC (才能考第4部分)占了10%,最后DTW计算只占了5%左右)
其中:
3)Early Abandoning of ED 和 LB_keogh
如果在计算欧几里得距离或者LB_keogh ,如果发现对当前序列点之间的求和超过best-so-far值,则停止计算,因为我们已经计算过更好的(low bound)
![4c4dd97d29a3a11d7d6bb0a9b190e402.png](https://i-blog.csdnimg.cn/blog_migrate/4f0cab34acb7f41d5a4f1baa581fe89f.jpeg)
4)增量计算DTW和LB_keogh
假如计算了LB_keogh发现他小于best-so-far值,需要重新计算当前两个序列的DTW,可以采用以下优化。
![1d662c9eb85d73d159325c07a56eb983.png](https://i-blog.csdnimg.cn/blog_migrate/5bedfbe385747f02ad6a6c165d4efd16.jpeg)
这里采用了一个K分割线,其中:
以上公式可以作为整个DTW(Q,C)的Lower Bound,叫做LB_keogh2,如果这个lower Bound超过当前的best-so-far值,则可以停止计算,这种方法避免大量重复计算。
5)多核心优化,这个算法可以将被被查询的序列分割成多“段”然后开启多线程的方式下并行搜索,但作者提供的源码中似乎并没有采用这种“优化”而且在没开启多线程的模式下,该算法的优化程度已经远好多其他同类的算法。
6)在计算Z-Normalization的同时计算欧几里得距离(或者LB_keogh)。这也意味着,如果我们提前中止的话(LB值大于best-so-far值),我们也提前中止了normalization步骤。
先进行一次扫描计算均值和方差:
计算子序列的时候可以采用以下公式从总序列的均值中得出:
7)总的算法伪代码如下:
![348a3038e1fc5cd18c9eb5e8d5f7d745.png](https://i-blog.csdnimg.cn/blog_migrate/a87be221f6d62f142975d4696e462cb9.jpeg)
8)对于早期终止计算顺序的优化:
![c179f16e897ef3f949d250c900b8d862.png](https://i-blog.csdnimg.cn/blog_migrate/c74a770bfd4eadfa013f4bab54939a2d.jpeg)
左边的图,如果从序列的开头开始计算,需要计算到第九个元素才能终止计算。
右图采用了一种特殊的计算顺序只需要计算5个元素就能终止计算。
该策略的方案是:对于查询的序列Q进行先Z-Normalization化,然后对所有元素的绝对值进行排序,返回降序排列的下标,后续的LB计算会基于这个下标的排列顺序。一个直观的理解:因为Qi要和许多Ci对比,对于一个进过Z-normalize化的序列,Ci会变成一个均值为0的高斯分布。所以距离均值最远的点将会给距离产生最大的贡献值,因此选择从这些点开始计算。
9)反向优化LB_keogh中定义的C和Q的职责
![59721d8b3a73b36462e2050f19ee6842.png](https://i-blog.csdnimg.cn/blog_migrate/6b6642ab6279a7b6a89e5091b43ddbec.png)
左图是LBkeogh算法中定义的{U,L}根查询序列Q来计算,而我们可以对该定义进行转换,将{U,L}计算与被查询的序列C上,叫做LB_keoghEC,这样计算出来的LB值又会更接近于DTW,能作为另一层的Lower Bounding优化。
10)对于各种Lower Bound算法的对比:
![6988060a2a9b46a1afaadfdc476facaf.png](https://i-blog.csdnimg.cn/blog_migrate/9b2eed78625663ac6c3416a6b44a593f.jpeg)
最后采用了分层的LBKimFC, LBkeoghEQ, LBKeoghEC, Earlyabondoning_DTW, DTW来进行逐层筛选。
三、试验结果和启发:
最后经过优化的DTW算法,这里称为UCR_DTW和UCR_ED和一些其他目前最佳算法的对比:
![4b9c14dcec834b3ce682720434c5172f.png](https://i-blog.csdnimg.cn/blog_migrate/a45ac44f0104551ec5f4bd44a586dd97.jpeg)
可以看出作者提出的各种优化思路大大提升了该算法的性能,明显优于当前最优实践。
这是一种在搜索匹配算法中,通过一个叫Lower Bounding的思想来进行逐层筛选,排除了很多不相关的计算(也叫剪枝)而达到优化的目的。作者的“分层计算”的思想,通过更简单的计算作为一个“过滤”只有通过“过滤”的数据样本才能进入下一层进行更“复杂”的计算,这其实避免了很多直接进行DTW复杂计算带来的性能开销。
![de18c7f21b5e7b9d6c534fed7989e694.png](https://i-blog.csdnimg.cn/blog_migrate/1ce7f2a39c74e1a82fce1e0fb95f4b20.jpeg)
在实际计算过程中复杂度最高的DTW在不同参数下计算的次数很低。各种Lower Bounding计算进行了大量的优化。
作者对欧几里得距离计算的理解很深,直接去除了对结果不影响但计算量不菲的平方根。
计算数据的顺序进行了排序的思想,认为数据应该服从“高斯分布”,这样先计算权重值更大的点,能高优先应用“早期排除(early abandoning)”的思想。
这篇算法的重点是讲解各种分层的Low Bounding和Early Abandoning,其实稍微改造一下就可以应用在多个匹配的应用场景。比如我要找K个最接近的匹配?(这在搜索引擎中可能使用最广)可以将算法中的best-so-far进行一下改造,加入一个最大堆maxHeap size大小为k,那么maxHeap的堆顶元素就是当前k个候选中距离最大的点,那么只要把堆顶元素作为当前的best-so-far值进行对比,继续应用early abondong逻辑,如果发现比top值小的匹配,则计算DTW,如果DTW也比top值还小,就可以将当前元素入堆,然后pop出当前的top值,作为后续的“best-so-far”来计算。这样最后heap中的元素就是候选的最接近的k个匹配。
稍微改造了一下作者的代码,用maxHeap来维护topK,发现在增加K值以后DTW的计算量确实会稍微增加,但是Low Bounding的优化作用仍然有效。参考值top1, top10, top50 分别是0.935秒, 1.056秒,1.183秒,low Bounding的优化仍然减少了大量的DTW计算。修改后的代码参考[6]
![56abd3cb4676bb52a550a332a37002da.png](https://i-blog.csdnimg.cn/blog_migrate/11182f15919eccc89c31fe366f5bae7c.jpeg)
四、这篇文章有哪些贡献?
这篇文章的技术可以应用在大数据时代的信息搜索上,例如现在搜索引擎中还没有比较好的能基于视频的搜索,这种基于时间序列的的高效算法,加上计算机性能的提升,说不定在未来会出现可以基于“视频”,“声音”的搜索引擎的应用。
在文本处理领域,这种技术也可以用于文本搜索上,在大的语料库中进行文本语料的匹配,比如在医疗领域对临床病历的信息检索而后续进行医学研究等等。
Reference:
[1] Thanawin Rakthanmanon, Bilson Campana. 2012. Searching and Mining Trillions of Time Series Subsequences under Dynamic Time Warping. KDD'12, August 12-16
[2] H. Ding, G. Trajcevski, P. Scheuermann, X. Wang, and E. J. Keogh. 2008. Querying and mining of time series data: experimental comparison of representations and distance measures. PVLDB 1, 2, 1542-52.
[3] Eamonn Keogh, Chotirat Ann Ratanamahatana. 2004. Exact indexing of dynamic time warping. Springer-Verlag London Ltd.
[4] Sang-Wook Kim, 2001. An Index-Based Approach for Similarity Search Supporting Time Warping in Large Sequence Databases. ICDE 2001 Paper-ID:207
[5] http://www.cs.ucr.edu/~eamonn/UCRsuite.html
[6] UCR_DTW top-k implementation
https://github.com/sesiria/Algs/blob/master/Lib/TimeSeries/UCR_Suite/UCR_DTW_k.cppgithub.com