時間序列趨勢相似性度量方法研究

譚章祿，王兆剛，胡 翰

中國礦業大學（北京）管理學院，北京 100083

1 引言

相似性度量是數據挖掘的基礎，針對時間序列數據的相似性度量已經成為時間序列數據挖掘的研究熱點之一[1-2]。其中，歐氏距離和動態時間彎曲（Dynamic Time Warping，DTW）是目前使用較多的時間序列相似性度量方法。但歐式距離只能處理等長度的時間序列，且無法識別變化趨勢，而動態時間彎曲距離雖然較好地克服了歐式距離的不足，但是計算復雜，時間復雜度較高，限制了其應用范圍[3]。

基于時間序列模式化的相似性度量方法，為時間序列的趨勢相似性度量提供了新的思路和可能，并取得了較多成果。

張海濤等[3]在以分段聚合近似實現數據降維的條件下，根據分段子序列的變化方向進行趨勢符號化，然后提出用于時間序列相似性度量的SMVT（Similarity Measurement of Variation-Trends）方法。劉慧婷等[4]通過經驗模態分解（Empirical Mode Decomposition，EMD）提取趨勢信息，分段后基于上升、保持、下降三種變化趨勢實現模式化轉化。肖瑞等[5]提出了基于時間序列變化趨勢的相似性度量方法和聚類方法，其中基于趨勢的相似性度量方法首先對時間序列進行區間劃分和區間內的趨勢形態判斷，生成短的趨勢符號序列，然后計算各趨勢符號的一階連接性指數，最后通過計算兩序列中各趨勢符號一階連接性指數的塔尼莫特系數完成相似性度量。王釗等[6]以上升、保持、下降的漲落模式保存原序列的趨勢變化信息，利用最長公共子序列算法計算漲落模式序列之間的形態相似性。基于局部變化方向實現序列的趨勢轉化，雖然可以有效描述短期趨勢，但其有序組合無法準確反映時間序列的長期變化方向，難以保留時間序列的整體趨勢信息，從而影響后續趨勢相似性的度量效果。

王達等[7]基于分段線性表示方法實現子序列劃分，在三元模式化的基礎上，提出模式距離概念，用以度量時間序列變化趨勢的相似性。董曉莉等[8]依據分段線性表示實現分段化處理，依據7種變化趨勢實現模式化轉化，提出基于形態的相似性度量方法。李正欣等[9-10]以擬合線段的傾斜角和時間跨度作為模式的描述方式，提出一種基于動態時間彎曲（DTW）的多元時間序列趨勢距離匹配方法，以斜率距離與時間跨度差距的線性綜合度量模式之間的差異。李海林等[11]基于SAX（Symbolic Aggregate Approximation）的均值符號化距離與分段導數序列DTW 距離的線性綜合，將動態時間彎曲與符號距離相結合來度量時間序列間的數值差異與形態差異距離。上述的趨勢相似性度量方法，均遵循“模式差異大，則數字距離大”的原則，但模式化的時間序列實際上是原始數據經過離散化處理后的一種符號化數據，尤其是上升、保持、下降等趨勢形態之間的差異應該是對等的，即上升與下降之間的距離和上升與保持之間的距離是一樣的，都是兩種不同的變化趨勢，而不是前者大于后者，統計距離難以準確度量衡量序列變化方向的差異，會在一定程度上影響趨勢相似性度量的準確性。

王燕等[12]在基于關鍵點實現時間序列分段的基礎上，對分段均值和斜率的符號化序列進行算術編碼，提出基于均值編碼距離和斜率編碼距離的分層歐氏距離的相似性度量方法，綜合考慮序列的統計距離和形態距離，達到序列整體趨勢匹配以及細節擬合的目標。一方面，分段均值無法反映局部形態，基于均值編碼距離實際上是數值差異，無法識別趨勢差異，其篩選結果會限制后續的相似性度量范圍；另一方面，斜率是對序列局部變化方向的描述，但斜率間的數值差異無法準確反映序列局部變化方向的異同，導致最終的相似性度量效果有限。

綜上所述，在基于模式序列的趨勢相似性度量方面提出了較多的思路和方法，但度量效果并不十分理想，仍然存在較大的改進和提升空間。因此，本文在提出依據分段子序列的均值及其線性擬合函數的導數符號實現模式轉換的基礎上，以模式之間的異同性比較定義模式匹配距離，借鑒DTW方法的動態規劃原理，提出了一種動態模式匹配方法，并分析了該方法的特點。最后，運用實驗數據測試了該方法的趨勢相似性度量效果。

2 時間序列及其分段模式

2.1 時間序列

時間序列X是由n項與時間順序有關的數據記錄組成的元素的有序集合：

其中，(xi,ti)表示在ti時刻的值為xi，采集時間ti是嚴格增加的，間隔時間Δt=ti+1-ti通常相同，即ti+1-ti=ti+2-ti+1,因此一般將時間序列X簡記為。

2.2 時間序列分段線性表示

時間序列的模式特征是指時間序列的某種變化特征，通過提取時間序列的模式特征，將時間序列變換到模式空間，就得到了時間序列的模式表示[14]。時間序列可以用模式表示如下：

其中，f(w)是時間序列模式表示，e(t)是原序列與它的模式表示之間的誤差。將時間序列按時間分成多個子段，如k段，f(w)定義為連接子段兩端點的直線段，則時間序列的分段線性表示[14]為：

其中，fk(t,wk)表示連接時間序列分段點的線性函數，ek(t)是這段時間內時間序列與它的分段線性表示之間的誤差，tk,L和tk,R表示第k段直線的起始時刻與終止時刻，且t1,L=t1,tk-1,R=tk,L,tk,R=tn。

k的值取決于壓縮率[15]E，二者之間的關系如式（4）所示。

2.3 模式類型劃分

現有研究和方法在將時間序列轉化為模式序列時，均單純以模式作為分段子序列變化趨勢的離散化符號，這樣雖然可以有效表征局部形態，但可能導致以模式符號有序連接構成的模式序列，難以準確反映時間序列的整體趨勢。以最簡單的上升、保持、下降三元模式[7]為例，設第j段分段子序列的直線擬合函數的斜率為pj，若pj>0，則該分段的模式為上升，用+表示；若pj<0，則該分段的模式為下降，用-表示；若pj=0，則該分段的模式為保持，用0表示。如圖1所示。

圖1 整體趨勢不同的時間序列

圖1 中，時間序列數據X1和X2被劃分為7 個子序列，tj,R表示第j個分段子序列的結束時間，且j=1,2,…,7,tn=t7,R。

如圖1所示，圖中每個分段括號內為其對應的模式符號，其中+表示上升模式，-表示下降模式，0表示水平模式。從圖1中可以看出，X1和X2的模式序列完全一致，但實際上其整體趨勢是相反的。因此，時間序列模式轉化過程中，不僅要考慮分段子序列的變化方向，而且需包含其均值水平信息，才能使模式符號在反映局部形態的同時，其有序連接組合具備描述整體趨勢的能力。

基于上述分析，依據分段子序列的均值及其擬合直線的斜率符號，將時間序列的模式劃分為9種，如表1所示。

表1 模式符號對應表

表1中，xˉj表示第j段的均值，pj表示第j個分段擬合直線的斜率，xmin和xmax分別表示時間序列X的最小值和最大值，即xmin=min(x1,x2,…,xn) ，xmax=max(x1,x2,…,xn)，

依據上述定義，時間序列X在劃分為k段后，轉化為模式化數據，其中zj∈(A,B,C,D,E,F,G,H,I)。

3 動態模式匹配方法

3.1 模式匹配距離

時間序列的模式匹配距離是指，時間序列數據在分段模式化轉化后，兩種模式之間的距離，即：

從式（5）可以看出，模式匹配距離實際上是分段子序列的模式符號之間的異同性比較，符號相同則距離為0，符號不同則距離為1。

3.2 動態時間彎曲

動態時間彎曲（DTW）是一種通過彎曲時間軸來更好地對時間序列形態進行匹配映射的相似性度量方法。它最早被應用于處理語音數據，后來Berndt等人將它用于度量時間序列的相似性。從此，DTW 在時間序列數據挖掘領域中得到廣泛的應用[16]。

DTW 不僅可以度量長度相等的時間序列，也可以對不等長的時間序列進行相似性度量，且對時間序列的突變點或異常點不敏感，比較適用于此類數據的度量，可以實現異步相似性比較[16]。

假設有兩個時間序列Q和U，且Q={q1,q2,…,qn}和U={u1,u2,…,um}，那么兩個時間序列數據點之間形成的距離矩陣Dn×m={d(i,j)}n×m，其中 1 ≤i≤n，1 ≤j≤m。d(i,j)的值由qi和uj之間的歐氏距離的平方來確定，即d(i,j)=(qi-uj)2。也就是說，矩陣D存儲了兩個時間序列不同時間點上數據之間的距離[17]。

如圖2 所示[17]，圖中的每個方格相當于D中元素值，那么DTW 就是從該矩陣中找到一條連續的路徑H={h1,h2,…,hs}，使得路徑上的元素值相加之和最小，同時這條路徑必須滿足以下三個條件，即邊界限制、連續性和單調性[17]。

圖2 動態時間彎曲路徑

在矩陣D中，滿足以上三個條件的路徑有很多，但只需要一條路徑作為動態時間彎曲距離，即：

最優路徑的查找方法是通過動態規劃來實現的，構造一個累計矩陣R={r(i,j)}n×m來記錄從起始位置到結束位置的最短路徑，且

其中，r(0,0)=0，r(i,0)=r(0,j)=∞。

最終兩個時間序列的動態時間彎曲距離可由累計距離表示，即LDTW(Q,U)=r(n,m)。由上述算法可以知道，實現長度分別為n和m的兩個時間序列之間的動態時間彎曲距離的時間復雜度為O(nm)[17]。

DTW 不要求數據點一一匹配，支持時間軸的伸縮和彎曲，可以有效度量不等長序列的相似性[18]，尤其對時間序列在時間軸上的形狀扭曲有非常優秀的辨識能力[5]。但DTW 需要計算基于歐氏距離的代價矩陣及其最優彎曲路徑，因此對數值變化敏感，計算復雜度較高，限制了其在大規模高維度時間序列數據中的應用。

時間序列的趨勢相似性度量建立在分段模式化基礎上，衡量時間序列趨勢變化的相似性，尤其是整體趨勢的相似性，因此要求相似性度量方法對短期局部的噪聲具備較好的抗干擾能力。DTW方法通過彎曲時間軸實現序列的異步相似性度量，使其不僅能夠根據時間序列的形態度量相似性，而且異步相似性度量可以有效避免短期局部的噪聲對序列趨勢相似性的干擾。因此，DTW通過彎曲時間軸的異步匹配原理有利于衡量時間序列的趨勢相似性。此外，以模式序列為對象，構建基于模式匹配距離的代價矩陣，在降低時間復雜度的同時，可以有效彌補DTW對數值變化敏感的缺陷，有利于度量大規模高維度時間序列的趨勢相似性。

3.3 動態模式匹配

動態模式匹配方法的主要思想，就是遵循動態時間彎曲距離的計算過程，并以兩條模式序列的模式之間的匹配距離構建距離矩陣D，最終的累計距離即為兩條模式序列的動態模式匹配（Dynamic Pattern Matching，DPM）距離。

假設兩個時間序列數據，經過分段模式化后形成兩條趨勢序列Y4和Y5,Y4={z41,z42,…,z4j,…,z4k},Y5={z51,z52,…,z5i,…,z5v}，Y4和Y5并不一定等長。兩條序列Y4和Y5的基元形態之間形成距離矩陣Dk×v={d(j,i)}k×v，其中 1 ≤j≤k,1 ≤i≤v,d(j,i)的值為z4j和z5i之間的模式匹配距離，即d(j,i)=d(z4j,z5i)，以模式匹配距離代替動態時間彎曲距離中計算時間序列數據點之間的歐氏距離，形成模式匹配距離矩陣Dk×v={d(j,i)}k×v={d(z4j,z5i)}k×v。

按照動態時間彎曲距離的求解過程，最優路徑的查找方法通過動態規劃來實現，以累計矩陣RDPM={rDPM(j,i)}k×v記錄從起始位置到結束位置的最短路徑。

其中，rDPM(0,0)=0,rDPM(j,0)=rDPM(0,i)=∞。

最終，兩條趨勢序列之間的動態模式匹配距離LDPM(Y4,Y5) 可由累計距離表示，即LDPM(Y4,Y5)=rDPM(k,v)。

3.4 算法分析

DPM方法以模式序列間的動態模式匹配距離來衡量原始時序數據間的趨勢相似性，因此，模式序列完整保留原始數據的局部形態與整體趨勢信息，是該方法能夠準確度量時間序列數據趨勢相似性的前提和基礎。

作為對模式形態間的異同性比較，模式匹配距離不遵循傳統的“模式差異大，則數字距離大”的原則，而是將模式符號看作定性數據，對不同形態之間的差異等同化對待，而沒有等級之分與大小之別，且計算過程較為簡單，計算量小。

時間序列的分段模式化過程，依據分段子序列的均值及其線性擬合函數的導數符號將時間序列轉化為模式序列，不僅降低了序列維度，而且可以在濾除噪聲平滑序列的同時，保留序列趨勢特征，從而為趨勢相似性度量奠定基礎。此外，DPM 方法基于模式匹配距離的代價矩陣，尋找最優路徑的過程，通過彎曲序列軸允許異步異同性比較，可以避免短期局部的突變過程和異常情況對趨勢相似性衡量的干擾。

該方法借鑒了DTW 方法的計算原理，因此繼承了其優良特性，不僅適用于序列等長的情況，而且也適用于序列不等長的情況。此外，該方法主要用以度量時間序列的趨勢相似性，尤其是整體趨勢的相似性，通過異步相似性度量模式序列的趨勢相似性，因此不要求序列之間對齊。該方法需要根據序列的變化幅度劃分區間，因此主要適用于離線類時間序列數據的相似性度量。而且，劃分的區間數量可以根據數據的變化幅度與復雜程度來選擇，如當數據值的變化幅度較大，或者各分段線性擬合函數的斜率差異較大時，可以設定較多的模式類型。

根據計算過程，容易分析得到該方法的計算時間效率。長度為n和m的時間序列進行平均分段的時間復雜度分別為O(n)和O(m)，對長度為k和v的模式序列進行DPM距離度量需要的時間復雜度為O(kv)，因此整個算法過程的時間復雜度近似為O(m+n+kv)。因為分段數目k和v通常小于原始序列的長度，所以DPM的時間復雜度O(m+n+kv)要小于DTW的時間復雜度O(mn)，且壓縮率越大，二者的差距越大，即DPM 的計算時間效率越高。

4 實驗分析

為了驗證動態模式匹配方法的可行性及有效性，選取數據集進行實驗對比分析。

實驗運行環境為 Intel?Xeon?CPUE5-2650 v4@2.20 GHz（2 處理器）、64 GB 內存、512 GB & 7200 rpm SSD 和Microsoft Windows 7 操作系統，開發工具為Matlab2014a。

4.1 實驗數據簡介

使用UCI 知識發現數據庫檔案（http：//archive.ics.uci.edu/ml/datasets.html）中的控制圖數據，選擇其中Normal類的前5個樣本數據，Increasing trend類的前5個樣本數據，Decreasing trend類的前5個樣本數據，共15個時間序列樣本數據。每個樣本數據等長，均為60 維。各類數據如圖3所示。

圖3 各類時間序列數據

以選擇的15 個時間序列樣本數據作為實驗數據集，其中 Normal 類數據的樣本序號依次為 1、2、3、4、5，Increasing 類數據的樣本序號依次為 6、7、8、9、10，Decreasing類數據的樣本序號依次為11、12、13、14、15。

4.2 實驗方案

實驗中對時間序列數據進行等長分段直線擬合，即tj,R-tj,L=tj+1,R-tj+1,L=c,c的值與壓縮率有關，c=為了避免無法等長分段的情況，最后一段包含剩余時間序列數據點，即tk,R-tk,L=n-c×(k-1)-1，其余分段子序列的長度皆為c。

圖4 方案1的聚類系統樹圖

圖5 方案2的聚類系統樹圖

圖6 方案3壓縮率為50%時的聚類系統樹圖

對選擇的數據使用層次聚類方法進行歸類分析，對比歐氏距離、動態時間彎曲距離、動態模式匹配距離之間的相似性度量效果，層次聚類方式均選擇average法。

方案1 使用原始數據，以歐氏距離作為原始數據樣本之間距離的計算方法。

方案2 使用DTW作為原始數據樣本之間距離的計算方法，對原始數據進行聚類。

方案3 在壓縮率分別為50%、80%的情況下，針對分段模式轉化后的趨勢序列，使用DPM 計算趨勢序列之間的距離，分析兩種壓縮率條件下的聚類效果。

方案4 分析DPM對不等長模式序列的趨勢相似性度量效果，其中樣本4、5、9、10、14、15的壓縮率為75%，其余樣本的壓縮率為80%，且不進行對齊處理。

4.3 實驗結果

實驗結果如圖4～圖8 所示，其中圖4 為方案1 的聚類譜系圖，圖5為方案2的聚類譜系圖。

由圖6、圖7可知，動態模式匹配距離的應用效果較好，在壓縮數據的條件下依然能準確衡量數據間的趨勢相似性。

將上述結果及聚類過程的時間消耗予以統計匯總，結果如表2所示。

圖7 方案3壓縮率為80%時的聚類系統樹圖

圖8 方案4不等長模式序列的聚類系統樹圖

表2 實驗結果

由表2 可知，動態模式匹配距離在壓縮率分別為50%、80%的條件下，聚類準確率均為100%，且其時間消耗隨著壓縮率的增大而極大縮短，遠遠小于動態時間彎曲的時間消耗，即在保證相似性度量準確率的前提下有效縮短了時間消耗。

方案4 為DPM 方法用于不等長模式序列的趨勢相似性度量。由圖8 可知，在模式序列不等長的情況下，只有樣本5 被錯誤地劃分為Increasing 類，其余樣本均實現正確的類別劃分，聚類準確率為93%，而且由表2可知，其時間消耗較低，遠小于DTW 方法的時間消耗。因此，DPM 方法對不等長序列具有較好的度量效果和性能。

5 結束語

本文在實現時間序列分段線性表示的基礎上，依據分段子序列的均值及線性擬合函數的導數符號，實現模式化轉化，并以模式間的異同性比較定義模式匹配距離，借鑒DTW方法的原理，構建了一種動態模式匹配方法，用以度量時間序列的趨勢相似性。實驗數據分析表明，該方法不僅能夠在不同的壓縮率下準確度量序列間的趨勢相似性，且相對于DTW方法，其時間消耗大幅降低。

因此，在分段模式化過程中，不僅要考慮局部趨勢方向，而且應該將分段均值水平作為模式劃分的標準之一，有利于模式在準確描述局部形態的同時，其有序連接組合能夠完整保留時間序列的整體趨勢信息，從而為趨勢相似性度量奠定良好的基礎。另一方面，作為對分段子序列的符號化表示，應該將模式作為定性數據，不同的模式代表不同的局部形態，模式之間并沒有等級高低之分，不同模式之間的差異沒有大小之別。以時間序列的模式化為基礎，挖掘序列頻繁模式與關聯規則，具有重要價值。

由于動態模式匹配方法可以度量不等長序列之間的趨勢相似性，而序列不等長實際上是存在數據缺失情況的一種表現形式，數據存在不同程度的缺失也是不同領域的時序數據挖掘研究所面臨的客觀現實情況，因此，對于不等長序列的趨勢相似性度量效果，尤其是該方法的相似性度量效果與數據缺失比例之間的關系，值得進一步的深入研究。