長度分布懸殊的序列比對法

陸 成 剛

(浙江工業(yè)大學(xué) 理學(xué)院,杭州 310023)

1 引 言

動(dòng)態(tài)時(shí)間規(guī)整算法(Dynamic Time Warping,DTW)是一個(gè)經(jīng)典的時(shí)間序列比對算法,最早被應(yīng)用于短詞語音識別的研究[1].自DTW距離被引入到一般時(shí)間序列的模式分析領(lǐng)域后[2],它已經(jīng)成為經(jīng)典的模式相異(似)性度量技術(shù),并在近些年也屢獲關(guān)注和研究[3-7].對于DTW技術(shù)的研究主要分為以下幾類:

1)對于快速DTW計(jì)算方法的研究[7-12],例如考慮線性時(shí)間算法;

2)對于DTW有效改進(jìn)的研究,如導(dǎo)數(shù)DTW算法(簡稱DDTW)[13];

3)基于DTW距離的矢量均值中心的計(jì)算[14,15]、或多個(gè)序列之間的DTW比對[12],其研究目的是考慮將DTW用于均值聚類.

DTW算法是計(jì)算序列之間動(dòng)態(tài)對應(yīng)意義下的最小距離值,與常規(guī)的靜態(tài)對應(yīng)下的序列距離計(jì)算相比,該方法有一些突出的特點(diǎn):

1)它不僅能計(jì)算兩段長度相同的序列間距,而且能計(jì)算長度不同的序列間距;

2)雖然它不滿足傳統(tǒng)距離的三角不等式法則,但在處理時(shí)長伸縮變異的序列時(shí)有傳統(tǒng)距離不具備的特效.

例如,對于序列(a,b,b,a)和(a,a,b,a)如果依照傳統(tǒng)的距離(取范數(shù)為1的歐式距離)計(jì)算為|a-b|;采用DTW距離計(jì)算則為0.零DTW距離表明這兩個(gè)序列是相同序列,這種自動(dòng)“隱藏、忽略”足標(biāo)間隔時(shí)長不一的技術(shù)在工程上是非常有用的.但傳統(tǒng)的DTW距離也隱含著一種缺陷,例如考慮使用DTW距離計(jì)算以下三段直流信號的匹配情況:A=(0.1,0.1)、B=(0.25,0.25,0.25,0.25)、C=(0.28,0.28,0.28),信號A與B的DTW距離為0.6,其最優(yōu)匹配路徑的長度為4;而信號C與A的DTW距離0.54,最優(yōu)匹配路徑的長度為3.從我們直覺判斷信號B比信號C更接近A,但從DTW距離決定的相異性來看,反而是C比B更接近A.造成這個(gè)反直覺的效果是由于我們僅僅考慮了距離,而忽視了路徑長度,假如以DTW距離相對路徑長度的平均值而言,B比C更接近A.

雖然DTW算法處理長度不一致的序列間的距離計(jì)算是它的優(yōu)勢,但所謂“成也蕭何,敗也蕭何”,在參與比對的多對序列之間由于各對最優(yōu)路徑長度的不一致造成對“距離大小決定相異程度”準(zhǔn)則的干擾.上例兩條路徑長度只差1就顯示了這種干擾,一般來說,當(dāng)長度相差劇烈時(shí)這種干擾發(fā)生的概率更高.顯而易見,如果將DTW作內(nèi)置距離實(shí)現(xiàn)kNN最近鄰法的序列匹配時(shí),這種缺陷將造成致命的匹配錯(cuò)誤.由此本文設(shè)計(jì)研究一種祛除路徑長度因素的DTW距離,稱為均值DTW.后文第2節(jié)將建立均值DTW的組合優(yōu)化的目標(biāo)函數(shù),這是一個(gè)DTW的非平凡推廣,該節(jié)將重點(diǎn)分析實(shí)現(xiàn)優(yōu)化的算法;第3節(jié)將考慮受加窗限制時(shí)的均值DTW的處理方式,與常規(guī)DTW加窗形式有顯著的不同;第4節(jié)給出均值DTW和常規(guī)DTW的算法實(shí)驗(yàn)比較,最后第5節(jié)對整個(gè)工作進(jìn)行總結(jié).

2 均值DTW

假設(shè)兩段離散信號S=(s1,s2,…,sn)和H=(h1,h2,…,hm),設(shè)任意匹配路徑函數(shù)P=(p1,p2,…,pl),其中pi=(qi,ri),而1<=qi<=n,1<=ri<=m,DTW路徑函數(shù)滿足單調(diào)性、連續(xù)性、邊界值等性質(zhì),即:

1)對相鄰兩點(diǎn)pi+1=(qi+1,ri+1)和pi=(qi,ri)滿足qi+1≥qi且ri+1≥ri;

2)對相鄰兩點(diǎn)pi+1=(qi+1,ri+1)和pi=(qi,ri)滿足qi+1=qi或qi+1=qi+1且ri+1=ri或ri+1=ri+1,但qi+1=qi和ri+1=ri不同時(shí)滿足;

3)p0=(1,1)且pl=(n,m).

也可以對路徑函數(shù)提出窗口、斜率等限制要求.

常規(guī)DTW的缺陷是沒有考慮路徑長度平均的距離,從而可能導(dǎo)致兩個(gè)同類序列的距離大于兩個(gè)異類序列的距離.這就違背了“距離大小決定相異程度”的準(zhǔn)則.定理1從理論上揭示了這種可能.

定理1.設(shè)對于序列S1其采樣周期為T,先將S1使用sinc插值恢復(fù)成模擬信號,再使用T/n的采樣周期離散化形成序列S2,S2長度n倍于S1、但兩者仍屬同類序列,可證S1和S2的DTW距離隨著n的增加而無限制增加,勢必將大于S1與任一確定的其它類序列H的DTW距離.

圖1 信號重構(gòu)采樣后的DTW匹配Fig.1 Signal reconstruction resampled DTW matching

為克服常規(guī)DTW的缺陷,在滿足以上路徑條件1、2、和3的基礎(chǔ)上提出均值DTW的最優(yōu)路徑問題：

(1)

其中距離函數(shù)d可以是范數(shù)1的歐式距離,也可以是平方歐式距離.在DTW的原始文獻(xiàn)[1]中優(yōu)化目標(biāo)函數(shù)除以常數(shù)l,取l=n、l=m或l=m+n任定.因?yàn)槭浅?shù)所以之后的文獻(xiàn)中往往忽略這個(gè)值,然而既便除以常數(shù)仍難以保證祛除路徑長度的干擾,這從引言例子的l如取作信號A的長度便可知.目標(biāo)優(yōu)化問題(1)考慮各種路徑下的平均距離取優(yōu),由于路徑長度隨著路徑不同而不同,一般而言l不是常數(shù).通過簡單的組合計(jì)算可以得到路徑長度l的變動(dòng)范圍max(m,n)≤l≤m+n-1.

圖2 動(dòng)態(tài)路徑延展示意Fig.2 Dynamic path zigzag display

當(dāng)前點(diǎn)(i,j)延展到下一點(diǎn)時(shí)并不是直接取{(i+1,j+1),(i+1,j),(i,j+1)}中的距離最小值,而是將它們逐一考慮進(jìn)全局路徑的距離總量,而后取平均值中最小者對應(yīng)的點(diǎn).一個(gè)直接的思路是先定出一條完整路徑再計(jì)算它對應(yīng)的平均距離,并且枚舉所有的路徑就能找出最優(yōu)均值及其對應(yīng)的路徑.

2.1 DTW路徑數(shù)目

num(p0)=1

num((q,1))=num((q-1,1)),q=2,3,…,n;

num((1,r))=num((1,r-1)),r=2,3,…,m;

num((q,r))=num((q-1,r))+num((q,r-1))+num((q-1,r-1)),q=2,3,…,n且r=2,3,…,m;

最后輸出num((n,m))為路徑的數(shù)目.該算法成立的原理可以使用圖3作示意.

圖3 路徑計(jì)數(shù)原理Fig.3 Path counting principle

表1 DTW路徑數(shù)目的指數(shù)增長
Table 1 Exponential increase in the number of DTW paths

算例1234512(m+n)1020304045num((n,m))14625632.53E+102.57E+125.39E+131.8E+14

2.2 均值DTW解法

(2)

2.3 計(jì)算示例

如圖5示例同類序列間由于長度相差懸殊導(dǎo)致DTW距離反而比不同類序列大,同時(shí)均值DTW卻能抵制這種失效.

圖4 均值DTW算法Fig.4 Mean DTW algorithm

圖6例舉了兩個(gè)序列按照DTW、均值DTW匹配得到的路徑:細(xì)實(shí)線代表DTW方法、粗點(diǎn)虛線代表均值DTW方法.

圖5 DTW距離失效,而均值DTW有效的實(shí)例Fig.5 DTW distance failure,and the example of the mean DTW is valid

圖6 DTW和均值DTW計(jì)算示例Fig.6 DTW and mean DTW calculation example

表2 DTW、均值DTW對比
Table 2 Comparison of DTW,Mean DTW

距離平均距離路徑長度DTW43.710.70562均值DTW46.857980.65080572

由表2可知從平均距離的角度看均值DTW小于常規(guī)DTW,而從(總值)距離的角度看DTW小于均值DTW.這是自然的、且符合這兩個(gè)概念的設(shè)計(jì)初衷的.最后討論均值DTW路徑長度與常規(guī)DTW路徑長度的關(guān)系.

定理2.對序列S=(s1,s2,…,sn)和H=(h1,h2,…,hm)它們的DTW距離對應(yīng)的路徑長度必不長于均值DTW對應(yīng)的路徑長度.

表2給出的例子顯示均值DTW的路徑長度大于DTW距離的路徑長度,定理2說明這是必然的.

3 均值DTW的加窗機(jī)制

圖7 均值DTW算法加窗限制Fig.7 Mean DTW algorithm windowing limit

圖8 “緊貼”窗口區(qū)域的外部點(diǎn)l值更新原理Fig.8 Principle of external point l value update in the “close to”window area

圖9給出圖6對應(yīng)的比較實(shí)例在應(yīng)用了加窗機(jī)制后的情形,窗口大小為10.表3給出該情形的距離值的對比.

圖9 DTW和均值DTW的加窗計(jì)算結(jié)果Fig.9 Windowing calculation results for DTW and mean DTW

表3與表2對比可見無論是常規(guī)DTW還是均值DTW,加窗處理會導(dǎo)致對應(yīng)的平均距離值有所升高,這是因?yàn)榧哟昂蟮钠ヅ浞秶s小了的緣故.此外,在加窗情形下就平均距離而言均值DTW依然是低于DTW,且路徑長度長于DTW.

表3 加窗計(jì)算的對比
Table 3 Comparison of windowing calculations

距離平均距離路徑長度DTW54.861010.97966156均值DTW60.499550.9166666

4 實(shí)驗(yàn)與結(jié)果

實(shí)驗(yàn)分析數(shù)據(jù)使用電離輻射時(shí)間序列,它們是取自于核電站周邊的環(huán)境輻射監(jiān)測儀.核電站環(huán)境輻射檢測數(shù)值通常和降水、風(fēng)向等天氣因素有關(guān),這是由于地球大氣層微粒吸收源自太空的高能粒子的電離輻射能量,通過雨、雪、風(fēng)等各種氣象狀況聚集于地表.所以研究核電站的環(huán)境背景電離輻射離不開對于降水導(dǎo)致的輻射劑量異常的分析和識別,是精確分析、檢測由核電泄漏事故引發(fā)的劑量異常的基礎(chǔ).核輻射監(jiān)測時(shí)間序列分析需要實(shí)時(shí)檢測輻射劑量波峰、以及識別波峰的具體類型(不同的氣象氣候所致、或者核泄漏事故等).此外,核輻射檢測儀器故障自動(dòng)重啟、分布式監(jiān)測網(wǎng)絡(luò)數(shù)據(jù)傳輸中斷等設(shè)備異常會造成核輻射時(shí)間序列夾雜一定的隨機(jī)外點(diǎn)或噪聲.

本文采用的算法有層次聚類、kNN分類,它們都可以使用DTW或均值DTW作內(nèi)置距離,數(shù)據(jù)是已經(jīng)人工標(biāo)注的序列,因而便于驗(yàn)證算法效果.聚類分類是模式識別的核心問題,常用的方法有k均值聚類、模糊c均值聚類、層次聚類、kNN最近鄰法分類、PCA、譜聚類、DBSCAN以及深度學(xué)習(xí)神經(jīng)網(wǎng)絡(luò),等等.如何判斷模式之間的相似(異)性又是聚類分類的核心問題,DTW距離是聚類分類、模板匹配的經(jīng)典的技術(shù)之一.對均值DTW算法與常規(guī)DTW作比較,分別植入層次聚類、和kNN最近鄰法進(jìn)行標(biāo)注數(shù)據(jù)集的聚類分類測試.

實(shí)驗(yàn)采用4907條人工標(biāo)注的共分為4個(gè)類別的核輻射波峰片段,它們被分為三個(gè)集合:集合A含有2324條4類別波峰段,序列段長的平均值為274.9、長度分布標(biāo)準(zhǔn)差為174.1;集合B含有1974條4類別波峰段,段長平均為194.4、標(biāo)準(zhǔn)差為5.7;集合C含有609條4類別波峰段,段長平均為186.1、標(biāo)準(zhǔn)差為178.3.此外集合C與A、B不同,幾乎每一條都含有隨機(jī)外點(diǎn)的成份,而集合A、C有一個(gè)共同點(diǎn)就是信號長度分布的懸殊程度較B更為劇烈,這容易從它們的長度標(biāo)準(zhǔn)差數(shù)值得知.實(shí)際上在核輻射檢測中一方面由于劑量上升異常所持續(xù)的時(shí)間長度確實(shí)不一,有時(shí)候差別甚為巨大(例如不同時(shí)間長度的持續(xù)降雨,等等);另外一方面,由于波峰提取階段所用技術(shù)的限制導(dǎo)致有些波峰段提取顯得較為“緊湊”、有些則含有較長的兩端低幅度延伸部分.圖10示例了A、B集合的少部分4種類別的波峰段,圖11示例了C集合的帶外點(diǎn)的波峰段,外點(diǎn)發(fā)生比較隨機(jī),幅值可能突然上升,也有可能突然下降.

圖10 集合A、B少部分的波峰段示例,分別以“▽+xo” 標(biāo)記4種類型Fig.10 Example of a peak segment with a small number of sets A and B.Mark 4 types with "▽+xo" respectively

圖11 集合C的4種含外點(diǎn)的波峰段類型示例,其中3個(gè)外 點(diǎn)是突降型,另外一個(gè)突升型Fig.11 Example of four types of peak segments with outlier in Set C,three of which are dip-type and the other is a sub-type

圖12給出了對集合A、B作層次聚類所得的計(jì)數(shù)矩陣(nij)4×4的三維柱狀圖,對集合B在波峰長度變化不懸殊時(shí)(長度分布標(biāo)準(zhǔn)差較小),常規(guī)DTW和均值DTW都能取到比較高的精度;當(dāng)對集合A進(jìn)行聚類時(shí)由于波峰長度相差較為懸殊(長度分布標(biāo)準(zhǔn)差較大),這時(shí)路徑長度的不一致對DTW距離值的干擾較大,此時(shí)常規(guī)DTW近乎完全失效,而均值DTW還能取到較好的效果.表4列出圖12計(jì)數(shù)矩陣計(jì)算的F值.

圖12 集合A、B進(jìn)行層次聚類的計(jì)數(shù)矩陣三維顯示Fig.12 Three-dimensional display of the counting matrix for hierarchical clustering of sets A and B

由表4可見在數(shù)據(jù)集A上使用均值DTW至少比使用DTW(對類別1而言)提升了精度達(dá)50%以上.

表4 兩類數(shù)據(jù)集的算法效果F值
Table 4 Algorithmic effectFvalues of two types of data sets

表5 對集合C進(jìn)行kNN近鄰法分類的定量分析
Table 5 Quantitative analysis ofknearest neighbor classification for set C

從表5看出在集合C所含具有較懸殊的信號長度差別的動(dòng)態(tài)特性下,將常規(guī)DTW和均值DTW用于kNN模式分類時(shí),后者的識別精度顯著地優(yōu)于前者.

5 總 結(jié)

對于模式聚類分類而言,在研究強(qiáng)有力的分類邏輯和機(jī)制中,離不開對數(shù)據(jù)統(tǒng)計(jì)規(guī)律的建模,同時(shí)對數(shù)據(jù)尋找有效的相似(異)性度量技術(shù)也是一個(gè)繞不開的課題.DTW算法因?qū)r(shí)間序列時(shí)長節(jié)奏的伸縮性有極佳的適應(yīng)性而被廣泛使用,常規(guī)DTW算法常常忽略路徑長度不一致所導(dǎo)致的對距離值的干擾,這在一定程度上抑制了DTW對長度相差懸殊的序列數(shù)據(jù)集的適用性.該文的研究工作從單位路徑長度的距離值的組合目標(biāo)優(yōu)化出發(fā),給出了捕捉最小平均距離的DTW對齊方式,并將此用作時(shí)間序列相似(異)性的度量方法.該方法較常規(guī)DTW在聚類分類上有顯著的精度改善.它和常規(guī)DTW有等價(jià)的計(jì)算復(fù)雜度,只在內(nèi)存消耗上比常規(guī)DTW增加一倍.將來可以考慮在最優(yōu)路徑搜索中融合自動(dòng)變系數(shù)的二階IIR濾波機(jī)制,以避免存儲路徑長度矩陣的內(nèi)存占用,并探討這樣處理給高含噪量的時(shí)間序列動(dòng)態(tài)匹配帶來的積極作用.