基于貝葉斯序列分割的高維時間序列在線分類算法

尹 君 張建業 李德高 景 康 周 平

1(國家電網新疆電力有限公司烏魯木齊供電公司 新疆 烏魯木齊 830000)2(國家電網新疆電力有限公司 新疆 烏魯木齊 830002)3(新疆信息產業有限責任公司 新疆 烏魯木齊 830026)

0 引 言

時間序列數據現已成為許多行業和工程領域中一種重要的數據形式，對時間序列進行在線挖掘分析具有極大的價值[1]。時間序列之間往往為非對齊的形式，所以基于歐氏距離的傳統分類算法無法實現理想的效果。研究人員提出動態時間規整(Dynamic Time Warping,DTW)算法[2]解決不對準的時間序列相似性度量問題，但基于DTW的相似性度量無法度量時間序列串聯結構的階段間差異。文獻[3]針對該問題提出了重要的shapelet方法，并得到了廣泛的關注和應用，也實現了很高的分類精度，但shapelet類的方法存在時間復雜度高的問題。雖然許多研究人員設計了shapelet的加速算法[4-5]，但是時間復雜度依然較高。

基于概率密度的方法[6]是另一種有效的時間序列分類算法，其時間復雜度較低，能夠實現在線的時間序列分類。此類方法[7]使用密度估計算法評估時間序列之間的相似性，實現快速的在線分類處理。密度估計的準確性是此類時間序列分類算法的關鍵部分，核密度估計(Kernel Density Estimation,KDE)[8]是最為常用的一種方法，但該方法無法應用于高維數據，而其他的非參數化密度估計方法[9]對高維數據的時間效率較低，難以滿足在線密度估計的要求。

動態時間規整解決了時間序列的不對準問題，對低維度數據流的效果較好，但是高維時間序列包含豐富的時空信息，動態時間規整則忽略了這些時空信息。本文將高維時間序列投影至重建相位空間，保留高維時間序列的時空信息，然后在重建相位空間完成密度估計和相似性度量，以期提高高維時間序列的分類準確率。此外，近期一些研究人員將貝葉斯序列分割技術應用于高維數據的密度估計問題，證明該技術對于高維數據的計算效率較高。受此啟發，本文將貝葉斯序列分割技術應用于時間序列的在線密度估計模型，以期對高維時間序列進行快速、準確的密度估計。

1 基于貝葉斯序列分割的密度估計

1.1 重建相位空間

如果重建空間和原空間的動態拓撲相同，那么該空間稱為重建相位空間(Reconstructed Phase Space,RPS)[10]。本文采用時間延遲嵌入方法將時間序列觀察投影到RPS，給定一個時間序列xm(m=1,2,…,N)，將xm投影到RPS的結果為：

xn=[xnxn+τxn+2τ…xn+(d-2)τxn+(d-1)τ]

(1)

式中:n=1,2,…,(N-(d-1)τ)，d為嵌入維度，τ為時間延遲。xn的完整時間序列可表示為：

(2)

式中:矩陣X的一行(向量)表示相位空間的一個點xn。

時間延遲嵌入方法采用滑動窗口訪問時間序列的數據，嵌入維度d對應窗口的大小，時間延遲參數τ決定了下一次采樣的步長。采用假近鄰法調節參數d，采用最小互信息法調節參數τ。假近鄰法把在d+1維空間距離遠，但在d維空間的近鄰點定義為假近鄰，選擇假近鄰的時間長度低于閾值的維度作為參數d的值。

使用最小互信息函數調節時間延遲參數τ，互信息函數定義為：

(3)

式中：p(·)為概率分布函數。

式(3)評估了兩個窗口Xt和Xt+τ之間的依賴性，即量化了Xt和Xt+τ之間的共享信息量。最小互信息函數的思想是選擇第一次出現兩個窗口互信息最小化的τ值，此時滑動窗口之間的依賴性最小。

1.2 貝葉斯序列分割

貝葉斯序列分割方法建立一個多維的直方圖，再不斷地對樣本空間進行二分類處理。給定一個由N個樣本構成的D維數據集X，將樣本空間逐漸分割為若干子區域。在經過若干次的分割處理之后，每個子區域的密度可粗略計算為該子區域的數據點數量和總數據點的比例。每次分割序列，嘗試M種不同的分割方式，由此可提高密度估計的準確率。

考慮一個二維的樣本空間。第一次分割(j=1)產生兩個子區域，后續的分割方案(j>1)記為gj={cut2,cut3,…,cutj-1}，樣本空間共有j-1個子區域，設子區域p的空間體積為vp，數據點數量為np。第j次分割共有(j-1)×D種可能的分割方式，基于一個概率函數隨機選擇一種分割方式。分割的結束條件為獲得了最優的分割結果，或者達到預設的最多分割次數。假設經過t次分割獲得了最優的分割結果，子區域1≤p≤t的概率密度計算為np/(Nvp)。

在實際應用中，很難預知數據的實際密度，為了解決該問題，原貝葉斯序列分割算法定義了分區的評分指標。設一個分區為p，p含有j個子區域，分區的評分方法為：

(4)

式中:nk為子區域k的數據量;Vk為區域的體積;參數α和β為常量;D(·)為狄利克雷分布，其參數為(α,α,…,α)，D(·)作為分區的后驗分布。參數β是分區先驗分布(exp(-j))的相關參數。

貝葉斯序列分割技術對于高維稀疏數據的密度估計準確率也存在不足之處，本文使用copula變換提高對高維數據的處理效果。將每個維度的邊際密度估計與copula變換的聯合密度估計的乘積作為最終的密度估計結果。為copula變換空間的每個維度設立邊界[0,1]。

1.3 概率密度估計方法

假設一個D維數據集共有N個數據實例，原貝葉斯序列分割算法將全部數據集作為樣本空間Ω，然后將Ω分為若干的子區域，基于每個子區域的數據量和體積估計子區域的密度。

為了減少計算復雜度，使用貝葉斯序列分割技術估計每個子區域的密度。本文的貝葉斯序列分割算法在訓練階段首先將樣本空間均勻分為B個子分區，每個分區b視為原樣本空間的一個近似，記為Ω(b)(b=1,2,…,B)，每個分區包含L=N/B個數據實例。使用貝葉斯序列分割技術獨立處理每個分區，最終獲得B個子分區的集合及其相應的密度，該集合包含了B個概率密度的估計。

(5)

根據Sklar定理[11]，任意的多元分布均可以轉換為帶變量邊際分布的形式，將一個有限維的聯合分布分解為它的邊緣分布和一個表示結構關系的copula函數，copula函數描述了變量間的相關性和一致性。

(1) 估計邊際密度。首先估計數據的邊際密度，使用邊際密度獲得累積分布函數(Cumulative Distribution Function,CDF)，使用邊際CDF在copula空間內構建一個多維的密度分區，獲得一個均勻邊際分布的D維樣本空間，記為[0,1]D。

基于邊際密度和copula變換空間的密度，可獲得測試數據z的總密度：

(6)

式中:fd為邊際密度;Fd為對應的邊際CDF;對copula相關性進行求導可計算出c。

(2) 維度對齊。上文對B個數據分區進行了不同維度的擴展，所以獲得的樣本空間大小不等。本文將CDF的范圍限定為[0,1]，copula變換域不存在空間不等的問題，因此，計算B個分區的平均值僅需要對齊B個樣本空間。因為所有的邊際密度均為一維空間，所以設計了高效的維度對齊方法。對齊方法的步驟為：

Step1在B個分區中搜索最小數據擴展和最大數據擴展。

Step2為B個分區的所有擴展設立相同的邊界。

Step3擴展每個分區密度的開始部分和結尾部分，與設立的邊界對齊。

Step4重新計算修改后分區的邊際密度，分區的數據點數量保持不變。

Step5使用更新的邊際密度計算新的CDF。

圖1 維度對齊方法的示意圖

最終使用邊際密度和copula變換密度計算每個分區的密度，將每個分區的密度代入式(5)產生最終的概率密度估計函數。

1.4 概率密度函數的度量方法

信息領域存在多個常用的多樣性距離度量方法，如Kullback-Leibler divergence，但大多數方法為非對稱方法，且計算效率低，多樣性度量方法(Integrated Squared Error,ISE)是其中計算效率較高的一種方法，本文采用ISE計算兩個概率密度函數之間的距離，ISE的計算方法為：

(7)

式中：p和q表示兩個概率密度函數，p和q越接近，ISE(p,q)則越接近0。

2 時間序列在線分類算法設計

本文方法將時間序列觀察表示為概率密度函數，利用K近鄰(K Nearest Neighbor,KNN)模型將概率密度函數在線分類。首先，通過時間延遲嵌入將時間序列數據投影到重建相位空間，圖2(b)所示是重建相位空間的實例圖。然后，采用本文基于貝葉斯序列分割方法基于重建相位空間的觀察估計其概率密度函數，如圖2(c)所示。之后，使用積分平方誤差計算概率密度函數間的相似性，建立所有概率密度函數的相似性矩陣。最終，使用KNN算法將時間序列分類。

圖2 時間序列的處理實例圖

因為重建相位空間能夠表示非線性動態時間序列數據的時間模式，將混沌不規則的時間序列映射到重建相位空間能夠增強時間序列的信息量，有利于后期的密度估計和距離度量處理。

2.1 建立相似性矩陣階段

算法1為建立相似性矩陣的算法。輸入參數包括：時間序列的觀察訓練集Tser[·],時間延遲嵌入方法的參數d和τ，以及密度估計的參數M，算法1采用核密度估計函數，M為核帶寬參數，本文設M≥1。首先，運用時間延遲嵌入方法將時間序列觀察s轉化為重建相位空間sRPS，然后，基于貝葉斯序列分割的密度估計方法計算sRPS數據點的概率密度函數Tpdf[]。最終，輸出所有時間序列的概率密度函數集Tpdf[]。

算法1建立相似性矩陣的算法

輸入：Tser[],d,τ,h。

輸出：Tpdf[]。

1.i=0；

2.forsinTser[] do

3.sRPS= delay_embed(s,d,τ);

//延遲嵌入

4.Tpdf[i] = density_est(sRPS,M);

//估計密度

5.i++;

6.end for

2.2 分類階段

算法2為時間序列在線分類算法。輸入參數包括：時間序列觀察s,時間序列觀測的概率密度函數集Tpdf,時間延遲嵌入的模型參數d和τ,密度估計的參數M,KNN的近鄰數量k。

首先，將時間序列觀察s轉化為概率密度函數，運用時間延遲嵌入方法處理s。然后，采用密度估計算法計算sRPS的概率密度函數。使用積分平方誤差計算目標時間序列概率密度函數和訓練集概率密度函數之間的距離。最終，運用KNN預測目標時間序列觀察s的分類。

算法2時間序列在線分類算法

輸入：s,Tpdf[],d,τ,h,k。

輸出：cl。

1.sRPS= delay_embed(s,d,τ);

//延遲嵌入

2.pdf= density_est(sRPS,h);

//估計密度,h=0.1

3.i=0;

4.forpinTpdf[] do

3.DISE[i]=ISE(pdf,p);

//計算積分平方誤差

4.i++;

5.end for

6.cl= KNN(DISE[·],k)；

2.3 算法的計算復雜度分析

將長度為N的時間序列轉化為d維重建相位空間的程序中，需要運行(N-(d-1)τ)×d次的時間延遲嵌入。因為訓練集的M個時間序列均需要該處理，所以共需要M×(N-(d-1)τ)×d次的時間延遲嵌入。

在時間序列的分類程序中，需要計算全部訓練時間序列的積分平方誤差，該過程包含兩個循環體。平方誤差的計算次數等于時間序列重建相位空間矩陣的行數，即(N-(d-1)τ)。計算每個測試時間序列和訓練時間序列間積分平方誤差的復雜度為(N-(d-1)τ)2×d。

最終，時間序列分類的總復雜度為O(M×(N-(d-1)τ)2×d)，其中:M為訓練集的時間序列數量;N為copula變換一維空間的維度。因此本文分類算法對于時間序列的維度具有魯棒性。

3 仿真實驗和結果分析

實驗環境為Intel i7-3820 CPU，主頻為3.6 GHz，內存為32 GB。基于MATLAB編程實現實驗中的所有算法。

3.1 實驗數據集

從UCR時間序列分類數據集[12]中選擇7個維度高于200的數據集，評測本文方法對于高維時間序列的分類性能。首先使用z-score將7個數據集歸一化處理，然后將數據集分為訓練集和測試集。表1為實驗數據集的基本屬性，7個數據集來自于不同的領域。這7個數據集被許多時間序列分類文獻所采用，因此便于完成對比實驗和分析。

表1 實驗數據集的基本屬性

3.2 密度估計實驗

首先從兩個角度評估本文基于貝葉斯序列分割的密度估計算法性能，所考慮的性能指標為密度估計誤差和計算時間。

(1) 密度估計的準確性。為了觀察密度估計的細節信息，該組實驗采用了人工合成的數據集，隨機生成服從多元高斯分布的合成數據集，數據集共有400個數據，維度為16，前2個維度的數據服從三峰值的正態分布，第3個維度的數據服從單峰的正態分布，4～64維的數據服從雙峰值的正態分布。對于不同分區大小分別測試密度估計的性能。采用KL散度(Kullback Leibler Divergence,KLD)評估密度估計算法的準確率，圖3為密度估計的KLD結果，圖中分別將每個分區的數據數量設為10、40和80。結果顯示，分區的數據量越大，KLD的性能越好，當數據量大于100時，密度估計的準確性較好。

圖3 密度估計的KLD結果

(2) 密度估計的時間性能。統計了估計B個塊的密度所需的總時間，每個分區的大小為L=N/B。處理N個數據的總時間計算為：

(8)

式中:tov為計算平均密度的時間，可忽略不計;tL(N)為處理N個實例的總時間;L為分區的大小;B為分區的數量;tb(L)為處理第b個分區的時間。

圖4為密度估計的時間結果，分區的數據量越大，處理時間越長。但本文算法對不同數據量的處理時間幾乎為常量，因此本文方法同時適用于穩態數據流和非穩態數據流。

圖4 密度估計的時間結果

3.3 時間序列分類實驗

目前主流的時間序列在線分類算法主要包括基于距離(基于密度)的分類方法和基于深度學習的分類方法兩種類型，本文算法分別和這兩種類型的分類方法作比較，深入評估本文方法的有效性。

1) 基于距離的時間序列分類方法。本文方法是一種基于距離的時間序列在線分類算法，首先選擇4個經典的方法作為對比方法。

(1) 高斯混合模型和重建相位空間結合的分類方法(GMMRPS)[13]。該方法首先將時間序列投影到重建相位空間，然后利用高斯混合模型建模數據，再采用最大期望算法對時間序列分類。

(2) K近鄰和歐氏距離結合的分類方法(KNNED)[14]。該方法將時間序列表示為t維空間的一個向量，采用歐氏距離度量測試時間序列和訓練時間序列之間的距離，從而對測試樣本進行實時分類。該算法易于實現，且計算效率較高，但對于序列不對準較為敏感。

(3) K近鄰和動態時間規整結合的分類方法(1NNDTW)[15]。該方法與KNNED較為相似，不同之處主要在于采用動態時間規整表示時間序列。

(4) 基于動態時間規整的快速分類算法(SDTW)[16]。該方法設計了時間序列的質量評價方法，并對低質量的部分時間序列提前剪枝，從而實現加速分類的目標。

圖5為基于距離分類方法的分類準確率結果，可看出，GMMRPS和KNNED對于高維數據集的準確率均較低，GMMRPS的高斯混合模型對于高維數據的度量效果較差，KNNED的歐氏距離對高維數據的度量效果也較差。KNNDTW和SDTW兩種基于動態時間規則的分類方法實現了較高的準確率，但隨著維度的提高，這兩種方法的分類準確率呈現明顯的下降趨勢。本文方法對于7個數據集均實現了較為理想的分類準確性，并且對數據維度顯示出明顯的魯棒性。

圖5 基于距離時間序列分類方法的準確率結果

2) 基于深度學習的時間序列分類方法。深度學習技術是近期性能極好的一種學習方法，選擇3個經典的深度學習方法作為對比方法。

(1) 基于多層感知機的分類方法(MLP)[17]。該方法的神經網絡包含三個隱層，每層包含500個神經元，采用ReLU激活函數，采用softmax作為輸出層。

(2) 基于全卷積神經網絡的分類方法(FCN)[18]。該方法的神經網絡包含3個隱層，濾波器數量為128 256 128，采用全局池化機制，采用softmax作為輸出層。

(3) 基于殘差神經網絡的分類方法(resnet)[19]。該方法的神經網絡包含3個殘差塊，每個殘差塊包含3個隱層神經元 ，采用全局池化機制，采用softmax作為輸出層。

圖6為基于深度學習分類方法的分類準確率結果，總體而言，基于深度學習的方法優于基于距離的方法。多層感知機對于高維時間序列的準確率較低，FCN和resnet兩種基于深度學習的分類方法實現了較高的準確率，但隨著維度升高，這兩種方法的分類準確率呈現明顯的下降趨勢。本文方法對于7個數據集均實現了較為理想的分類準確性，并且對數據維度顯示出明顯的魯棒性。

圖6 基于深度學習時間序列分類方法的準確率結果

4 結 語

重建相位空間能夠表示非線性動態時間序列數據的時間模式，將混沌不規則的時間序列映射到重建相位空間能夠增強時間序列的信息量，有利于后期的密度估計和距離度量處理。貝葉斯序列分割技術對于數據的維度具有魯棒性，本文將貝葉斯序列分割技術應用于時間序列的在線密度估計模型，對高維時間序列進行快速、準確的密度估計。在基于多組高維數據集上進行仿真實驗，本文方法的時間性能和分類準確率均對時間序列的維度具有魯棒性，并且實現了較好的分類準確率。目前本文僅考慮了常規的高維時間序列問題，未來將研究本文方法在演化高維數據流和混沌高維時間序列等問題上的應用，擴大本文方法的應用價值。