基于概率密度空間劃分的符號化時間序列分析及其在異常診斷中的應用

胡世杰， 錢宇寧， 嚴如強

(東南大學儀器科學與工程學院，江蘇 南京 210096)

引 言

在機械系統中異常檢測是指在給定的數據集中提取出特征，并利用提取出的特征來判斷其狀態是正常還是異常。隨著科學技術和現代工業的飛速發展，國民經濟的機械、能源、石化、運載和國防等行業的機械設備日趨大型化、高速化、集成化和自動化，這對中國經濟高速發展提供了有力保障。但是由于機械設備的故障失效引起的災難性事故屢有發生，若能準確及時識別機械系統運行當中異常狀態，對機械系統的安全運行，避免重大和災難性事故意義重大[1]。隨著對高質量、低能耗以及安全生產的需求，已經有多種信號處理方法引入到機械系統的異常檢測當中。其中，符號化時間序列分析由于具有計算效率高，高信噪比等特性被廣泛地應用在異常診斷當中。

符號化時間序列分析STSA起源于上世紀90年代中期，它是由符號動力學理論、混沌時間序列分析和信息論發展起來的一種新的信息分析方法。在符號化時間序列分析當中，最關鍵的一步是對原始時間序列的符號化。符號化即是對原始時間序列進行離散化，將有多種不同值的數據序列變為僅有幾個互不相同符號的序列。這一過程能夠保留原始時間序列當中的大尺度特征，從而降低動力學噪聲和測量噪聲的影響[2]；同時由于數據的符號化，其計算效率大大提高。目前，符號化方法可以分為兩類：1)基于值域的符號化方法；2)基于分布的符號化方法。基于值域的符號化方法主要是通過對時間序列值域的分析來進行符號化。例如，Asoky Ray提出了統一劃分符號化的方法[3]，該方法首先確定時間序列的最小值與最大值，然后將值域劃分為N個值域大小相同并且連續的區間。其中N為采用符號的個數。Rajagopalan在統一劃分符號化方法的基礎上提出了最大信息熵符號化方法[4]，該方法在符號化過程當中使每個符號出現的概率相同，從而保證在使用相同符號來符號化原始時間序列時系統擁有最大的信息熵。該方法可以通過以下步驟實現：首先將原始時間序列按照從小到大的順序排列，然后將整個序列劃分為N段長度相同的區間，最后將每個區間對應到一個符號即完成了符號化。這些基于值域的符號化方法(例如，等區間的統一劃分，等概率的最大信息熵劃分)在大多數情況下是有效的，但是機械系統的運行狀態與時間序列概率分布是密切相關的，符號化過程當中不考慮其分布將會導致符號化后的信息丟失。另一類是基于分布的符號化方法，這類方法基于時間序列符合高斯分布假設。例如，Lin假設時間序列符合高斯分布(N(0，1))，然后根據高斯分布確定“斷點”以得到若干等概率的區域[5]。然而實際應用中時間序列可能形成于一種未知的動力系統或者伴隨著不同水平的噪聲，這將導致時間序列實際分布與高斯分布相差甚遠。所以基于分布的符號化方法有一定的局限性，不是一種普遍適用的方法。受上述研究的啟發，本文提出一種稱為概率密度空間劃分的符號化方法。概率密度符號化方法結合了基于值域的符號化方法與基于分布的符號化方法的優點。該方法首先對原始時間序列進行統計分析并計算其概率密度圖，然后選擇時間序列的平均值作為概率密度圖的中心點，以中心點為對稱中心時間序列將被劃分為等概率的若干區域，之后將每個區域映射到一個符號。通過上述步驟，可以得到基于概率密度符號化方法的符號序列。為了對符號序列進行異常檢測，最后對符號化時間序列進行編碼。本文第二部分將介紹符號化時間序列分析的理論背景以及概率密度符號化方法的實現；第三部分通過實際軸承疲勞實驗數據驗證算法并與其他劃分方法的檢測效果進行比較，最后一部分得出結論。

1 符號化時間序列分析與概率密度符號化方法

時間序列的符號化分析分為兩步[6,7]：先將時間序列轉化為符號序列，再對符號序列進行統計分析。

為了詳細介紹算法流程，先定義以下變量：xn={x(1),x(2),…,x(n)}表示直接從傳感器獲得的時間序列。而原始時間序列xn的符號化表示為sn={s(1),s(2),…,s(n)}。其中sn通過將xn劃分為q(q≥2)個不相交的區域并且將每個區域映射到一個符號s(i)∈{0,1，2，… ,q-1}=S，其中符號集S是有限個符號的集合，集合的大小為q。符號化過程所采用的空間劃分方法對接下來的符號序列統計分析影響十分重大。在下一小節，將會詳細介紹文中提出的概率密度符號化方法。

要將時間序列符號化，首先要選擇符號集S的大小q，這是目前還有待研究的一個問題。若q太小將導致“粗粒化”過程當中原始時間序列有用信息丟失；若q太大將會導致計算效率大大降低。文獻[3～5]選擇q=2，3，4。為了闡述符號化的原理，在這里選取q=4簡化說明

(1)

式中Ci為對xn進行空間劃分得到的不相交集合。通過上述閾值函數，時間序列xn被轉化符號序列sn。

在將獲取的時間序列轉化為符號序列之后，為了提取符號序列當中的特征信息需要對其編碼[8]。首先選擇一個標準長度L(L≥2)，L個連續的符號組成一個字，每個字被編碼成qL進制，這樣就形成了新的編碼序列。圖1為L=4,q=4時符號序列的編碼示意圖。

圖1 編碼示意圖(L=4, q=4)

符號序列當中所含有的特征信息可以通過對編碼序列進行統計學分析提取出來。對編碼序列的分析方法有很多種，最常見的方法有信息熵法、標準差法。香農熵是通過計算編碼序列的復雜度來確定系統的運行狀態。香農熵越大則表示編碼序列復雜度大，進而可以確定系統的不穩定性較大[9]。標準差法通過計算編碼序列的標準差來衡量系統的運行狀態，當機械系統運行正常時其編碼序列的標準差較小，同時它會隨著機械系統故障程度的加深而增大[10]。鑒于標準差法計算效率較高，本算法當中選取標準差法來對編碼序列進行統計分析。

綜上所述，符號化時間序列分析異常檢測方法流程如圖2所示。

圖2 符號化時間序列分析異常檢測算法流程圖

具體步驟如下所示：

(1)符號化：確定符號化當中符號集大小q，并選取劃分方法將原始時間序列xn進行符號化。為了對比概率密度劃分方法與其他符號化方法，本文中統一選取q=4。

(2)編碼：確定編碼字長L，并對符號序列sn進行編碼。

(3)統計學分析：對編碼序列進行統計學分析，本文中選取標準差為特征量。

(4)異常檢測：將統計學分析中獲得的表征系統運行狀況的特征量與正常時系統運行特征量進行對比，從而判斷當前系統運行狀況。

2 概率密度劃分方法

為了對原始時間序列進行概率密度劃分，可以假設原始時間序列為連續型隨機變量。然后對其進行數理統計分析，可以得到其概率密度函數f(x) 。密度函數f(x) 有如下性質：

(a)f(x)≥0；

上述式子中，P(a

在確定概率密度函數f(x) 之后，選擇時間序列的平均值xave作為中心點，可以按照下式分別計算右概率函數Fr(x) 和左概率函數Fl(x)：

(2)

(3)

式中Fr(x) 表示時間序列當中取值在中心點xave到xave+x之間的概率。同理Fl(x)表示時間序列取值在(xave-x,xave]之間的概率。

接下來需要將原始時間序列進行空間劃分, 即將其劃分為q個區域。由于每個區域是連續的，所以其可以由兩端“斷點”來確定。“斷點”是一系列點的集合，例如，概率密度劃分P=[P1，…，Pi，…，Pq-2Pq-1]，其中“斷點”Pi到“斷點”Pi+1為一個劃分區域并且該區域的概率為1/q。根據符號集q的奇偶性，“斷點”由兩種不同的方式確定。如果q為偶數，按照表1來確定“斷點”。如果q為大于2的奇數，其“斷點”的確定方式如表2所示。

表1 偶數符號化斷點

表2 奇數符號化斷點

得到“斷點”之后，空間劃分的q個區域也就確定了。通過類似式(1)當中的閾值函數，原始時間序列將被轉化為符號化時間序列。如圖3所示為對實際機械振動信號進行概率密度劃分的實例。其中圖3(a)為原始信號，圖3(b)為概率密度分布圖。按照上述劃分方法，信號被劃分為4(q=4)個區間,每個區間分別表示為符號“0”, “1”,“2” ,“3”。

圖3 概率密度劃分示例

3 實驗驗證

為了在實際實驗中檢驗本劃分方法的有效性，將基于概率密度劃分的符號化方法應用在實際軸承故障數據上進行實驗，并對結果進行分析。

軸承疲勞實驗數據來自美國智能維護系統中心(IMS)[11],實驗裝置如圖4所示。一個軸上安裝了4套Rexnord ZA-2115 雙列滾子軸承，每列滾子數量為16，滾子組節圓直徑為75.501 mm，滾子直徑為8.470 74 mm，接觸角為15.17°。軸的轉速保持2 000 r/min恒定不變，通過彈簧裝置在軸上加載6 000 lb(2 721.554 kg)的徑向載荷。所有軸承潤滑固定，并且每個軸承座都安裝2個PCB加速度傳感器用來采集軸承的振動數據。振動信號由NI公司DAQCard-6062E數據采集卡每隔10 min采集一次，采樣長度為20 480個點，采樣頻率為20 kHz。試驗臺中的4套軸承從2月12日11∶16∶18運行至2月19日06∶22∶39，一共采集到984個文件數據。在疲勞實驗結束時，軸承1檢測到外圈故障。

圖4 軸承實驗裝置

對軸承1的振動數據進行基于概率密度劃分的時間序列異常診斷。首先對軸承1的振動數據進行統計學分析，從而得到時間序列的概率密度分布。然后確定符號集q的大小。采用較大的q值劃分之后的符號序列中含有更多的細節信息，但是這也會導致由噪聲引起的錯誤符號增多，同時計算量也大大提高。采用較小的q值劃分，則可能導致信息丟失[12]。合理的選取q值是很重要的，但是q值的選取不是本文主要研究內容。在本實驗當中，參照以往文獻經驗性地選取符號集q=4。在完成對符號序列編碼之后，每個文件編碼序列的標準差如圖5所示。編碼序列的標準差越大，則其偏離正常狀態越遠。由于疲勞試驗前期(前540個點)軸承處于健康狀態，所以其編碼序列標準差比較平穩并且處于較低水平。在541點時編碼序列標準差檢測到劇烈的增長，可以判斷軸承故障在此時開始。編碼序列的標準差隨著故障程度的加深也逐漸增大。

為了檢驗基于概率密度劃分符號序列分析算法的優劣，進行了對比實驗2。在對比試驗2中不采用符號序列分析方法，直接計算原始時間序列的標準差，其實驗結果如圖6所示。從圖中可以看到在741點時標準差有劇烈增長，由異常檢測原理可以判斷該點為故障的起始點。與實驗1的結果對比，檢測到異常起始點滯后了200個點。本次對比實驗可以證明基于概率密度劃分的符號化方法對于系統中的異常更加敏感。

圖5 基于概率密度劃分符號序列分析結果

圖6 原始時間序列標準差分析結果

對比實驗3使用傳統統一劃分的STSA算法與本文方法進行比較。實驗結果如圖7所示。通過曲線可以觀察到，當軸承處于健康狀態時，曲線較為平穩。但是當異常發生時曲線發生劇烈的抖動，從而導致人們無法判斷異常起始點。從本實驗可以看到劃分方法對于基于符號化時間序列分析異常診斷的影響，同時此實驗結果也證明概率密度符號化方法相比統一劃分方法更好。

圖7 基于統一劃分符號序列分析結果

4 結 論

在符號化時間序列分析當中，生成符號化序列這一步至關重要。本文提出了一種新穎的符號化方法，該方法利用原始時間序列的概率密度分布來對其劃分，具有適用性強、符號化效果好的特點。同時將基于概率密度劃分的STSA應用于實際軸承系統的異常診斷中。對比實驗結果表明，相比直接標準差分析，本文提出的基于概率密度劃分的符號化方法對于軸承的異常信息更加敏感，能夠更早地發現故障。這對于故障的預防、軸承維護、安全生產具有重大意義。通過對比實驗進一步發現，概率密度符號化方法相比于傳統的統一劃分在軸承的異常診斷當中更加有效、穩定。

本文當中也存在一些不足之處，例如未對符號集q的大小以及編碼序列長度L進行細致的研究，只是根據以往文獻經驗性地取值。然而這些取值對與符號化時間序列分析具有很大的影響，這些問題需要更加深入的研究和驗證。

參考文獻：

[1] 何正嘉, 陳進, 王太勇. 機械故障診斷理論及應用[M].北京：高等教育出版社，2010.He Z J, Chen Jin , Wang T Y. Theories and Applications of Machinery Fault Diagnostics[M]. Beijing:Higher Education Press, 2010.

[2] 郭勁松, 衛武迪. 基于稀疏符號時間序列分析的轉子碰摩故障早期檢測方法[J]. 振動與沖擊, 2008, 27(12): 148—150.Guo J S, Wei W D.Methods for identification and early detection of rub-impact fault in rotors based sparse symbolic time series analysis[J]. Journal of Vibration and Shock, 2008, 27(12): 148—150.

[3] Ray A, Symbolic dynamics analysis of complex for anomaly detection[J]. Signal Processing, 2004, 84(7):1 115—1 130.

[4] Rajagopalan V, Ray A. Symbolic time series analysis via wavelet-based partitioning[J]. Signal Processing, 2006, 86:3 309—3 320.

[5] Lin J, Keogh E, Lonardi S, et al. A symbolic representation of time series with implications for streaming algorithms[A].Proceedings of the 8th ACM SIGMOD workshop on Research Issues in Data Mining and Knowledge Discovery[C]. ACM, 2003: 2—11.

[6] 王妍, 徐偉. Lorenz 系統中時間序列的相空間重構方法與特性[J]. 振動工程學報, 2006, 19(2): 277—282.Wang Y, Xu W. The methods and performance of phase spacereconstruction for the time series in Lorenz system[J].Journal of Vibration Engineering, 2006, 19(2): 277—282.

[7] Daw C S, Finney C E A, Tracy E R. A review of symbolic analysis of experimental data[J]. Review of Scientific Instruments, 2003, 74(2): 915—930.

[8] Das G, Lin K I, Mannila H, et al. Rule discovery from time series[J]. Knowledge Discovery and Data Mining, 1998,98:16—22.

[9] Tang X Z, Tracy E R, Boozer A D, et al. Symbol sequence statistics in noisy chaotic signal reconstruction[J]. Physical Review E, 1995, 51(5): 3 871—3 889.

[10] Bishop C M. Pattern recognition and machine learning[M]. New York: Springer, 2006.

[11] Qiu H,Lee J,Lin J,et al. Robust performance degradation assessment methods for enhanced rolling element bearing prognostics[J]．Advanced Engineering Informatics,2003,17(3):127—140.

[12] Finney C E A, Nguyen K, Daw C S, et al. Symbol-sequence statistics for monitoring fluidization[J]. ASME HEAT TRANSFER DIV PUBL HTD, 1998, 361: 405—412.