基于時(shí)變自回歸參數(shù)模型的滾動(dòng)軸承智能故障診斷

李健寶 彭 濤

1.湖南工業(yè)大學(xué)，株洲，412008 2.中南大學(xué)，長沙，410083

0 引言

滾動(dòng)軸承是旋轉(zhuǎn)機(jī)械中應(yīng)用最廣泛的零部件之一，軸承工作狀況的好壞決定著機(jī)器能否正常工作。統(tǒng)計(jì)顯示，旋轉(zhuǎn)機(jī)械設(shè)備的功能失效有30%是由軸承故障引起的，因此，對滾動(dòng)軸承進(jìn)行故障診斷具有十分重要的意義。

軸承運(yùn)行時(shí)的振動(dòng)信號是典型的非線性非平穩(wěn)時(shí)間序列，很難用一個(gè)完全確定的數(shù)學(xué)函數(shù)來表達(dá)。因而對軸承進(jìn)行故障診斷常通過提取振動(dòng)信號的特征參數(shù)并建立其與運(yùn)行狀態(tài)之間的關(guān)系來實(shí)現(xiàn)。時(shí)間序列分析（time series analysis）通過將觀測數(shù)據(jù)擬合為一個(gè)參數(shù)模型，實(shí)現(xiàn)對系統(tǒng)動(dòng)態(tài)特征與內(nèi)在結(jié)構(gòu)關(guān)系的近似描述［1］。常用的時(shí)間序列參數(shù)模型，如自回歸（autoregressive，AR）模型、滑動(dòng)平均（moving average，MA）模型和自回歸滑動(dòng)平均（autoregressive moving average，ARMA）模型等，都是在假定數(shù)據(jù)序列為平穩(wěn)的條件下建立的，在實(shí)際應(yīng)用中存在一定的局限性。時(shí)變參數(shù)模型因其具有較高的非平穩(wěn)信號時(shí)頻分布分辨率而受到普遍關(guān)注［2－6］。時(shí)變自回歸 （time－varying autoregressive，TVAR）模型是目前應(yīng)用最多的一種時(shí)變參數(shù)模型，如文獻(xiàn)［5］通過對滾動(dòng)軸承振動(dòng)信號的TVAR模型求取的時(shí)頻譜進(jìn)行奇異值分解，提取奇異值作為特征參數(shù)，文獻(xiàn)［6］對轉(zhuǎn)子系統(tǒng)振動(dòng)信號的TVAR模型，提取其基函數(shù)的組合權(quán)值作為特征參數(shù)，實(shí)現(xiàn)了對旋轉(zhuǎn)機(jī)械的狀態(tài)監(jiān)測和故障診斷。目前，直接提取振動(dòng)信號TVAR模型的時(shí)變參數(shù)作為反映機(jī)械設(shè)備運(yùn)行狀態(tài)特征參數(shù)的研究還不多見。

通過對軸承振動(dòng)信號的TVAR模型的時(shí)變參數(shù)進(jìn)行大量實(shí)驗(yàn)研究分析，發(fā)現(xiàn)時(shí)變參數(shù)能有效利用信號的時(shí)頻分布信息，較好地表征非平穩(wěn)信號的動(dòng)態(tài)特征。為此，本文提出一種基于時(shí)變自回歸參數(shù)模型的滾動(dòng)軸承智能故障診斷方法。首先對滾動(dòng)軸承的振動(dòng)信號建立時(shí)變自回歸參數(shù)模型，提取時(shí)變參數(shù)的均值作為反映軸承運(yùn)行狀態(tài)的特征參數(shù)，然后采用支持向量機(jī)（support vector machines，SVM）分類器實(shí)現(xiàn)對滾動(dòng)軸承的智能故障診斷。

1 TVAR建模原理

對離散時(shí)間序列｛x（t），t＝1，2，…，N｝可建立如下TVAR模型：

式中，ai（t）為時(shí)變參數(shù)；p 為模型的階次；ε（t）為均值為零、方差為σ2的白噪聲。

時(shí)變參數(shù)ai（t）可用一組基函數(shù)的線性組合來表示［7］，即

式中，m為基函數(shù)維數(shù)；aij為基函數(shù)的組合權(quán)值；gj（t）為一組基函數(shù)。

常用的基函數(shù)有時(shí)間多項(xiàng)式基函數(shù)、傅里葉基函數(shù)、離散余弦基函數(shù)、勒讓德多項(xiàng)式基函數(shù)和離散長球序列基函數(shù)等。由于軸承的振動(dòng)信號可近似看作是循環(huán)周期信號［8-9］，而傅里葉基函數(shù)比較適于周期變化情況［10］，因此，本文采用傅里葉基函數(shù)，其表達(dá)形式為

將式（2）代入式（1）可得

進(jìn)而可以求得模型殘差ε（t）的方差σ2的最小二乘估計(jì)值：

利用式（9）求解模型參數(shù)時(shí)存在矩陣求逆的問題，當(dāng)矩陣較大時(shí)，求解式（9）需要很大的存儲空間和很長的時(shí)間?？梢岳眠f推算法來求解式（9）。令則有如下遞推算法［11］：

估計(jì)出參數(shù)a10，…，a1m，…，ap0，…，apm后，就可根據(jù)式（2）求出模型參數(shù)ai（t）在各個(gè)時(shí)刻的值。

在實(shí)際應(yīng)用中，模型的階次通常是未知的，需要根據(jù)觀測數(shù)據(jù)來適當(dāng)?shù)嘏卸?。AIC準(zhǔn)則（akaike information criterion）是目前應(yīng)用比較廣泛的定階方法，通過使平均對數(shù)似然函數(shù)為最大或Kullback信息量最小來確定AR模型的最佳階次，其定義如下［12］：

2 基于TVAR模型的智能故障診斷方法

本文提出的滾動(dòng)軸承智能故障診斷方法流程見圖1，首先對滾動(dòng)軸承運(yùn)行時(shí)的振動(dòng)信號建立時(shí)變自回歸參數(shù)模型，用模型參數(shù)的均值構(gòu)建特征向量來表示第i個(gè)樣本信號（i＝1，2，…，N），然后采用支持向量機(jī)實(shí)現(xiàn)狀態(tài)的智能識別。

圖1 故障診斷流程圖

具體實(shí)現(xiàn)過程如下：

（1）分別在軸承正常、內(nèi)環(huán)故障、滾動(dòng)體故障和外環(huán)故障狀態(tài)下，按一定的采樣頻率進(jìn)行采樣，獲得一定數(shù)量的振動(dòng)信號作為樣本（1024個(gè)采樣點(diǎn)為一個(gè)樣本）。

（2）對來自傳感器的含有大量噪聲的原始信號進(jìn)行低通濾波等預(yù)處理，得到待分析信號x（t），利用式（12）確定p、m，然后由式（11）估計(jì)出參數(shù)A，建立時(shí)變自回歸參數(shù)模型。

（4）利用訓(xùn)練樣本對SVM分類器進(jìn)行訓(xùn)練，SVM的建立詳見文獻(xiàn)［13］。

（5）用訓(xùn)練好的SVM多類分類器對測試樣本進(jìn)行故障模式識別。

3 實(shí)驗(yàn)分析

3.1 實(shí)驗(yàn) Ⅰ

本部分采用的數(shù)據(jù)全部來自美國Case Western Reserve大學(xué)滾動(dòng)軸承數(shù)據(jù)中心［14］。測試平臺如圖2所示，由1492W（2馬力）三相交流電機(jī)（左）、轉(zhuǎn)矩傳感器（中）、測力計(jì)（右）和電子控制裝置（未顯示）組成。電機(jī)軸由測試軸承支撐，通過放電加工（electro－discharge machining，EDM）技術(shù)在測試軸承中植入單一局部故障缺陷，故障點(diǎn)直徑分別為177.8μm（7mil）、355.6μm（14mil）和533.4μm（21mil），深度為279.4μm，測試軸承為6025－2RS JEM SKF型深溝軸承。

圖2 測試臺

振動(dòng)信號通過加速度傳感器采集，采樣頻率為48kHz。在不同電機(jī)負(fù)載／轉(zhuǎn)速工況工作條件（0、746W、1492W 和2238W）、4種不同故障類型（正常狀態(tài)、滾動(dòng)體故障、內(nèi)環(huán)故障、外環(huán)故障）和3種不同故障程度情況下記錄振動(dòng)加速度信號數(shù)據(jù)，共獲得40組數(shù)據(jù)（正常狀態(tài)無故障大小，只需采集不同工況下的4組數(shù)據(jù)），并將每組數(shù)據(jù)以1024個(gè)采樣點(diǎn)為一個(gè)樣本構(gòu)成一組樣本集。按相同負(fù)載／轉(zhuǎn)速、相同故障程度分類，可將40組樣本集分成12個(gè)數(shù)據(jù)集，每個(gè)數(shù)據(jù)集由在相同負(fù)載／轉(zhuǎn)速和相同故障點(diǎn)直徑條件下的4種不同故障類型樣本集組成，如表1所示，表中，N表示正常狀態(tài)，I表示內(nèi)環(huán)故障，B表示滾動(dòng)體故障，O表示外環(huán)故障（下同）。以相同負(fù)載／轉(zhuǎn)速、相同故障類型分類（不包含正常情況），亦可將36組樣本集分成12個(gè)數(shù)據(jù)集，每個(gè)數(shù)據(jù)集由在相同負(fù)載／轉(zhuǎn)速和相同故障類型條件下的3種不同故障點(diǎn)直徑的樣本集組成，如表2所示，表中，F(xiàn)D7表示故障直徑為177.8μm，F(xiàn)D14 表示故障直徑為355.6μm，F(xiàn)D21表示故障直徑為533.4μm（下同）。按第2節(jié)所述方法對各數(shù)據(jù)集進(jìn)行仿真實(shí)驗(yàn)研究，在此以數(shù)據(jù)集D007＿0為例，對具體過程加以詳細(xì)說明。

表1 樣本集分類描述（一）

表2 樣本集分類描述（二）

選取4種故障類型的樣本各一個(gè)，進(jìn)行模型階次（p×m）的確定。得到正常狀態(tài)的最佳模型階次為19×4，內(nèi)環(huán)故障的最佳模型階次為18×3，滾動(dòng)體故障的最佳模型階次為18×3，外環(huán)故障的最佳模型階次為20×3。為便于進(jìn)行SVM分類，需對數(shù)據(jù)集中的各樣本建立相同階次的TVAR模型，為此，選取最佳模型階次為18×3。然后求出數(shù)據(jù)集中各樣本的時(shí)變參數(shù)，圖3所示為數(shù)據(jù)集D007＿0中每類故障各一個(gè)樣本的部分時(shí)變參數(shù)情況?？梢?，不同故障類型的時(shí)變參數(shù)之間具有較大的區(qū)分度，如正常和內(nèi)環(huán)故障的時(shí)變參數(shù)a1（t）與滾動(dòng)體和外環(huán)故障的時(shí)變參數(shù)a1（t）之間有較大的差異等。因此，時(shí)變自回歸模型的參數(shù)可以用作反映軸承運(yùn)行狀態(tài)的特征參數(shù)。

如果直接將時(shí)變參數(shù)作為表征軸承運(yùn)行狀態(tài)的特征向量輸入到SVM，則一個(gè)樣本的輸入即為18×1024＝18 432維的高維向量，不利于SVM的處理，且影響分類效率。為此，在求出軸承振動(dòng)信號的TVAR模型參數(shù)的基礎(chǔ)上，需進(jìn)一步研究時(shí)變參數(shù)的特性來獲取低維特征參數(shù)，用以表征軸承的運(yùn)行狀態(tài)。本文在建立軸承振動(dòng)信號的TVAR模型后，分別對每一個(gè)時(shí)變參數(shù)求平均，這樣18個(gè)時(shí)變參數(shù)可獲得18個(gè)平均值，則特征參數(shù)的維數(shù)降至18維。在大量實(shí)驗(yàn)研究分析的基礎(chǔ)上發(fā)現(xiàn)，在軸承的不同運(yùn)行狀態(tài)下，對應(yīng)的TVAR模型參數(shù)的均值亦有著較大的可分性。以數(shù)據(jù)集D007＿0為例，各樣本的部分TVAR模型參數(shù)的均值如圖4所示，顯然，不同故障類型所對應(yīng)TVAR模型參數(shù)的均值之間仍具有較大的區(qū)分度。選取表2中的數(shù)據(jù)集DINN＿2進(jìn)行分析，結(jié)果如圖5所示，可見，不同故障程度所對應(yīng)TVAR模型參數(shù)的均值之間同樣具有較大的區(qū)分度。大量實(shí)驗(yàn)研究分析表明，時(shí)變自回歸模型參數(shù)的均值可以較好地刻畫運(yùn)行在不同負(fù)載條件下不同故障類型、不同故障程度軸承的動(dòng)態(tài)特性，因而可以用作描述軸承振動(dòng)信號的特征參數(shù)。

圖3 單個(gè)樣本的部分時(shí)變參數(shù)

圖4 D007＿0各樣本的部分時(shí)變參數(shù)均值

圖5 DINN＿2各樣本的部分時(shí)變參數(shù)的均值

因此，對軸承的振動(dòng)信號建立18×3階TVAR模型，采用模型參數(shù)的均值作為反應(yīng)軸承運(yùn)行狀態(tài)的特征參數(shù)，即構(gòu)建特征向量用以反映軸承的運(yùn)行狀態(tài)。由特征向量Xi構(gòu)造出訓(xùn)練樣本集和測試樣本集利用訓(xùn)練樣本對SVM分類器進(jìn)行訓(xùn)練，并用訓(xùn)練好的SVM多類分類器對測試樣本進(jìn)行故障模式識別。

對表1、表2中的數(shù)據(jù)集分別進(jìn)行仿真實(shí)驗(yàn)，結(jié)果如表3和表4所示，各狀態(tài)均有較高的分類識別率。仿真實(shí)驗(yàn)結(jié)果說明，采用本文所提的特征提取方法可以有效地提取出反映軸承運(yùn)行狀態(tài)的特征參數(shù)，從而實(shí)現(xiàn)滾動(dòng)軸承故障的有效診斷。將本文方法與文獻(xiàn)［6］所提方法進(jìn)行對比實(shí)驗(yàn)（均采用SVM作為分類器），限于篇幅，此處只列出各數(shù)據(jù)集整體分類識別率，分別如表3和表4中文獻(xiàn)［6］方法一欄所示，可知，采用本文方法所獲得的整體分類識別率均高于文獻(xiàn)［6］所述方法，結(jié)果進(jìn)一步說明了本文所提方法的有效性。

表3 測試結(jié)果（一）

表4 測試結(jié)果（二）

3.2 實(shí)驗(yàn)Ⅱ

本實(shí)驗(yàn)是在筆者所在實(shí)驗(yàn)室的QPZZ－Ⅱ型旋轉(zhuǎn)機(jī)械故障平臺上進(jìn)行的，裝置如圖6所示，軸承型號為N205型滾珠軸承。分別在不同負(fù)載、轉(zhuǎn)速下通過安裝在軸承座外殼上的加速度傳感器采集軸承運(yùn)行時(shí)的振動(dòng)信號，采樣頻率為20kHz，模擬的故障分別為外環(huán)、內(nèi)環(huán)、滾動(dòng)體出現(xiàn)細(xì)小裂紋的情況。然后采用第2節(jié)所述的方法對其進(jìn)行故障診斷，同樣以1024個(gè)采樣點(diǎn)為一個(gè)樣本，得到各狀態(tài)的訓(xùn)練樣本和待識別（測試）樣本數(shù)分別為20和480，診斷結(jié)果如表5所示。由表5可知，在不同的故障發(fā)生形式下，采用本文所提的方法對各種故障類型進(jìn)行分析亦具有較高的分類識別率，說明本文所提方法具有一定的實(shí)用性。

圖6 QPZZ－Ⅱ旋轉(zhuǎn)機(jī)械故障模擬平臺

表5 驗(yàn)證實(shí)驗(yàn)結(jié)果

4 結(jié)束語

通過對滾動(dòng)軸承的振動(dòng)信號建立時(shí)變自回歸參數(shù)模型，提出了一種用模型的時(shí)變參數(shù)均值構(gòu)建表征軸承運(yùn)行狀態(tài)特征參數(shù)的方法，經(jīng)支持向量機(jī)分類器對所提取特征進(jìn)行故障診斷與分類，實(shí)現(xiàn)了滾動(dòng)軸承故障的智能診斷。利用不同運(yùn)行狀態(tài)下的數(shù)據(jù)進(jìn)行仿真實(shí)驗(yàn)，結(jié)果證明了本文所提方法的有效性。

［1］ 吳今培，孫德山.現(xiàn)代數(shù)據(jù)分析［M］.北京：機(jī)械工業(yè)出版社，2006.

［2］ Chen Z S，Yang Y M，Hu Z，et al.Detecting and Predicting Early Faults of Complex Rotating Machinery Based on Cyclostationary Time Series Model［J］.Journal of Vibration and Acoustics，2006，128（5）：666－671.

［3］ Zhan Y M，Mechefske C K.Robust Detection of Gearbox Deteriorateion Using Compromised Autoregressive Modeling and Kolmogorov－Smirnov Test Statistic－PartⅠ：Compromised Autoregressive Modeling with the Aid of Hypothesis Tests and Simulation Analysis［J］.Mechanical System and Signal Processing，2007，21：1953－1982.

［4］ 王勝春，韓捷，李志農(nóng)，等.基于TVAR的自適應(yīng)時(shí)頻分析及在故障診斷中的應(yīng)用［J］.軸承，2007（6）：28－31.

［5］ Wang G F，Luo Z G，Qin X D，et al.Fault Identification and Classification of Rolling Element Bearing Based on Time－varying Autoregressive Spectrum［J］.Mechanical System and Signal Processing，2008，22：934－947.

［6］ 張龍，熊國良，柳和生，等.基于時(shí)變自回歸模型與支持向量機(jī)的旋轉(zhuǎn)機(jī)械故障診斷方法［J］.中國電機(jī)工程學(xué)報(bào)，2007，27（9）：99－103.

［7］ Pally R K，Louis Beex A A.Modeling of Time－varying Instantaneous Frequency in a Finitely Correlated Environment［C］／／In 16th International Conference on Digital Signal Processing.Santorini，2009：1165－1170.

［8］ McCormick A C，Nandi A K.Cyclostationarity in Rotating Machine Vibrations［J］.Mechanical System and Signal Processing.1998，12：225－242.

［9］ Capdessus C，Sidahmed M，Lacoume J L.Cyclostationarity Processes：Application in Gear Faults Early Diagnosis［J］.Mechanical System and Signal Processing，2000，14：371－385.

［10］ 張海勇，李勘.非平穩(wěn)隨機(jī)信號的參數(shù)模型分析方法［J］.系統(tǒng)工程與電子技術(shù)，2003，25（3）：386－390.

［11］ Girault J M，Ossant F.Time－varying Autoregressive Spectral Estimation for Ultrasound Attenuation in Tissue Characterization［J］.IEEE Transactions on Ultrasonics，F(xiàn)erroelectrics，and Frequency Control，1998，45（3）：650－659.

［12］ Zhan Y M，Jardine A K S.Adaptive Autoregressive Modeling of Non－stationary Vibration Signals Under Distinct Gear States.Part 1：Modeling［J］.Journal of Sound and Vibration，2005，286：429－450.

［13］ 楊正友，彭濤.基于振動(dòng)信號分析和支持向量機(jī)的滾動(dòng)軸承故障診斷［J］.湖南工業(yè)大學(xué)學(xué)報(bào)，2009，23（1）：96－99.

［14］ Case Western Reserve University.Bearing Data Center［EB／OL］.［2010－02－10］.http：／／www.eecs.cwru.edu／laboratory／bearing／welcome＿overview.htm.