基于經驗模態分解、多尺度熵算法和支持向量機的滾動軸承故障診斷方法

□ 張文哲 □ 張為民， □ 林文波

1.同濟大學機械與能源工程學院 上海 201804

2.同濟大學中德學院 上海 201804

1 研究背景

作為機械設備中的關鍵零部件，滾動軸承是主要故障源［1］。與其它機械零部件相比，滾動軸承壽命的離散性很大，對工況的依賴性也很大，有的軸承已經遠超設計壽命卻依然良好工作，有的軸承遠未達到設計壽命卻出現各種故障，可見，滾動軸承的故障診斷方法是機械故障診斷中的重點研究方向，在監測軸承狀態、及時發現故障等方面起著重要作用。好的故障診斷方法可以提高機械設備的維修效能，從而顯著地提高經濟效益。據統計，在旋轉機械故障中，由滾動軸承引起的占30%左右；在感應電機故障中，由滾動軸承引起的占40%左右；在齒輪箱各類故障中，軸承故障率占20%左右，且僅次于齒輪［2］。

滾動軸承的振動信號通常是非平穩、非線性的，一般將其轉化為平穩信號處理，但有一定局限性。另一方面，用時頻分析提取滾動軸承振動信號故障特征的方法越來越多地應用于軸承故障診斷中。常用的時頻分析方法有經驗模態分解（EMD）［3-4］、多尺度熵（SME）算法［5］和小波包分解［6-7］等，這些方法的優缺點比較見表1。

表1 時頻分析方法對比

完整的滾動軸承故障診斷過程包括信號采集、特征提取、狀態識別、故障分析、決策干預等五個基本環節，其中的關鍵環節是特征提取和狀態識別。常見的狀態識別方法有人工神經網絡［8］、支持向量機（SVM）［9］和隱馬爾可夫模型［10］等，各方法的優缺點比較見表2。

表2 狀態識別方法對比

在滾動軸承出現早期故障時，故障信號淹沒在噪聲信號中。為了提高信噪比，筆者首先利用小波包變換對滾動軸承振動信號進行降噪，然后對重構信號進行EMD，根據本征模態函數（IMF）與降噪后信號的相關因子選擇值較大的IMF，計算篩選后的IMF分量多尺度樣本熵，確認最佳尺度，并將該尺度下的IMF樣本熵輸入SVM，對滾動軸承故障進行分類，這一算法的流程如圖1所示。

▲圖1 算法流程

2 基于小波包變換的信號降噪

小波分析能滿足低頻部分頻域分辨率高、高頻部分時域分辨率高的要求，但小波分析僅對分解之后的低頻部分進行再分解，未考慮高頻部分。相比之下，小波包分析在每層都能同時對信號的低頻部分和高頻部分進行分解，克服了小波分析不能細分高頻部分的缺點，應用更為廣泛。

2.1 小波包算法

定義函數庫｛un（t）｝為由尺度函數 u0（t）=φ（t）所確定的小波包。

式中：u2n（t）為尺度函數φ （t）在平方可積函數空間L2（R）的推廣；u2n+1（t）為小波函數 ψ（t）在 L2（R）的推廣；hk可看作低通濾波器組因子；gk可看作高通濾波器組因子。

式（1）中小波包的尺度固定不變，則將｛2-j/2un（2-jtk）|n=0，1，2，...k，j∈Z｝作為由多尺度函數 φ（t）導出的小波庫，j表示尺度，k表示位移，n表示振蕩因子。

2.2 小波包降噪步驟

小波包分解對信號特征具有局部細化的能力，可以降低噪聲對信號的干擾，提高信噪比。筆者基于小波包變換進行信號降噪，通常有以下幾個步驟［11］：① 選擇一個合適的小波基，確定小波分解層數n，然后對信號進行n層小波包分解；② 對于給定熵值，計算最佳小波包分解樹，確定最優小波包分解基；③選擇適當的閾值對最佳小波包分解樹各個分解尺度下的小波因子進行閾值量化處理；④ 對小波包分解的第n層的尺度因子和閾值量化處理后的小波因子進行重構。

以上步驟中，閾值大小及閾值量化處理方法比較關鍵，直接關系到小波包降噪的質量。

3 信號EMD

EMD將復雜信號分解為若干IMF分量之和，每個IMF分量在任意時刻都只有唯一頻率成分，能夠反映信號內部的振蕩模式。IMF分量必須滿足以下兩個條件：①在時域范圍內，信號序列局部極值點和過零點個數相等或者差值為1；② 在任意時刻，由信號局部極大值和極小值確定的上下包絡線的平均值為0。

篩分信號需要基于如下假設：任意復雜信號均由若干獨立的IMF分量組成；所有IMF分量均具有數量相等的極值點和零點，或者數量相差為1；由極值點形成的上下包絡線關于時間軸局部對稱。

對于滾動軸承降噪信號 x（t），EMD流程如下。

（1）找出x（t）的所有局部極值點，利用三次樣條插值，由局部極大值和局部極小值分別構成上下包絡線。

（2）將上下包絡函數的均值記為m1（t），將原始信號與上下包絡函數均值序列之差記為 h1（t），即 h1（t）=x（t）-m1（t），若 h1滿足 IMF 分量的條件，則 h1（t）是原始信號的第一個IMF分量，否則用h1（t）代替上式中的x（t），重復 k 次得到 h1k（t）=h1（k-1）（t）-m1k（t）。 直到 h1k（t）滿足條件，令 c1（t）=h1k（t），則 c1（t）為 x（t）的第一個IMF分量。

（3） 在 x（t）中剔除 c1（t），得到余項 r1（t），將其作為原始信號，重復前兩個步驟，得到第二個IMF分量c2（t）。 如此循環 n 次，直到rn（t）呈單調趨勢或者|rn（t）|很小時停止循環，得到n個IMF分量，循環過程如下：

因此，原始信號可以表示為：

式中：rn（t）為殘余項，表示信號的平均趨勢。

4 MSE算法

作為非線性動力學參數，MSE被廣泛用于滾動軸承故障診斷中，取得了不錯的效果。

假設一組原始信號 X=｛x1，x2，x3，...，xN｝，信號長度為N，其MSE算法過程如下。

（1）對原始信號進行粗粒化處理。給定嵌入維數m 與相似容限 r，構建粗粒序列

式（7）中τ∈N，稱為尺度因子，每個粗粒序列的長度為原始信號的長度與尺度因子的商，因此當τ=1時，就是原始信號。

式中：L為粗粒序列的長度。

（4）對于任意一個m維向量，都有L-m個距離。給定相似容限 r，對于每一個 i，統計 dij＜r的個數 n，計算其與距離總數的比值，記為，則有：

求其對所有 i的均值，記為 Bτ，m（r），則有：

（5） 將維數增加到 m+1，重復步驟（2）～（4），可以得到均值 Bτ，m+1（r）。

（6）計算原始信號的樣本熵：

當原始信號長度有限時，樣本熵可以表示為：

樣本熵隨著尺度因子τ的變化而變化，原始信號的MSE為不同尺度下樣本熵的集合，其值為：

由式（14）可以看出，MSE的值與嵌入維數m、相似容限r、尺度因子τ及數據長度N相關。m越大，所需 N 越多。 r可以取（0.1～0.25）SD，SD為原始數據標準偏差。 筆者取 τ=15，m=2 048，r=0.2SD。

5 基于SVM的故障分類

SVM是基于兩分類問題的一種學習算法，其基本思想是最大限度地分開兩類訓練樣本，即對于一個兩類問題的訓練樣本，構造一個分類超平面使分類間隔達到最大。 給定一組訓練樣本 U=｛（xi，yi）|（xi，yi）∈Rd×R，xi∈x，yi∈y，i=1，2，3，...，l｝， 其中 x 為訓練樣本，y∈｛1-，+1｝，即樣本分為兩類，l表示樣本總數，d表示訓練樣本中每一個樣本向量的維數。若在l×d維空間中，有超平面wTx+b=0，使樣本均能被正確分開，且樣本兩邊數據間隔最大，那么該平面為最優超平面，w為加權向量，b為偏移向量。

對應上述條件的優化問題如下：

引入拉格朗日函數，將約束條件融入目標函數中，上述問題轉化為：

αi為拉格朗日算子， 當 xi為支持向量時，αi＞0，否則 αi=0。

式（16）的對偶問題為：

式（17）是線性可分條件下尋求最優超平面的對偶算法，其優點為：①計算復雜度不再取決于樣本空間維數，而僅取決于支持向量數；②可以引入核函數，進而推廣到非線性分類問題。

當在低維空間不能利用超平面對樣本進行分類時，可以將非線性的樣本點映射到高維空間，然后在高維特征空間構造出最優超平面。為避免維數爆炸，基于相關泛函理論，引入核函數 K（xi，yi），對任意函數 φ（x）≠0，若滿足 Mercer條件，則 K（xi，yi）表示對應高維空間樣本點的內積。Mercer條件為：

定義非線性映射 φ：Rn→H，若有 K（xi，xj）=＜φ（xi），φ（xi）＞，則非線性SVM的對偶問題轉化為：

對于任意測試樣本x，分類函數為：

標準的SVM是二分類器，需要結合多個SVM實現多分類。筆者采用一對一算法，即n類樣本點需要構造［n（n-1）］/2 個 SVM，可以大大減少重復訓練，適用于多樣本分類。

6 算例分析

數據取自凱斯西儲大學軸承數據庫［12］，實驗臺架由電動機、扭矩傳感器、譯碼器、測功計及電子控制器組成，待檢測的滾動軸承支承電動機轉軸安裝在電動機驅動端，實驗臺如圖2所示。

測試滾動軸承選用斯凱孚6205-2RS深溝球軸承，結構參數見表3，實驗臺工況見表4。根據滾動軸承理論、故障特征頻率計算公式可以得到不同部件的故障特征頻率，取軸承正常狀態、內圈故障、外圈故障及滾動體故障四種類型各45組樣本，其中30組為訓練集，其余為測試集，每組樣本長度為2 048，故障均為單點電蝕，理論故障特征頻率及不同位置的故障尺寸見表5。

表3 6205-2RS軸承結構參數

▲圖2 實驗臺示意圖

表4 實驗臺工況

表5 軸承故障數據

采集滾動軸承振動信號的過程中，會混有大量噪聲，不利于提取故障特征，因此需要對試驗數據進行小波包降噪處理。小波基選擇正交小波db10，考慮到故障頻率范圍及不同層對應的尺度函數與小波函數中心頻率，小波包分解層數選擇7。由于采用默認閾值會過濾掉部分有用信號，因此最終選取閾值比默認閾值小4。以軸承內圈故障為例，降噪前后時域信號如圖3所示。相比原始信號，降噪效果明顯，但是無法直接判斷出故障信號。

對180組降噪后的時域信號進行EMD，根據IMF分量與降噪信號的相關因子選擇包含故障信息的前3個IMF分量，基于不同故障類型進行希爾伯特包絡分析，如圖4、圖5、圖6所示。從內圈故障、外圈故障及滾動體故障信號的包絡譜中都能得到實際故障特征頻率，依次為 164.1 Hz、105.5 Hz、140.6 Hz，與理論頻率相近。在圖4中出現了內圈故障特征頻率的二倍頻到五倍頻，以及滾動軸承轉頻的二倍頻58.59 Hz，并且在內圈故障特征頻率的邊頻帶可以找到調制頻率接近滾動軸承轉頻。在圖5中出現了外圈故障特征頻率的二倍頻到七倍頻，以及滾動軸承轉頻29.3 Hz與其二倍頻。在圖6中可找到滾動體的故障特征頻率140.6 Hz及外圈故障特征頻率105.5 Hz，但故障特征的幅值不明顯，因此通過EMD與包絡譜分析并不能較好判斷滾動體故障，僅能較好地分析內圈或外圈故障。

▲圖3 內圈故障信號波形

▲圖4 內圈故障信號包絡譜

▲圖5 外圈故障信號包絡譜

▲圖6 滾動體故障信號包絡譜

取訓練集中3個樣本為例，計算前8階IMF分量與重構信號的相關因子，對不同IMF分量按故障類型取樣本平均值，見表6。剔除相關因子小于0.01的IMF分量，選取前5個IMF分量，計算不同IMF分量的多尺度樣本熵，找到能區分故障狀態及故障類型的最佳尺度因子τ，綜合比較發現τ=1時各狀態都能明顯區分，IMF1分量MSE如圖7所示。

將τ=1時180組降噪信號經EMD的前5個IMF分量的樣本熵作為故障診斷特征向量，對其進行歸一化預處理，見表7。交叉選取120組數據作為訓練集，將特征向量輸入SVM，SVM類型選擇C-SVC，核函數類型為高斯徑向基核函數，對影響分類質量的重要參數懲罰因子及核參數進行參數尋優，通過網格參數尋優法進行選取，得懲罰因子取0.000 98，核參數取2。

為驗證模型的泛化性能，將剩余的60組數據作為測試集，得到故障診斷結果，見表8。可以看出，訓練后的SVM模型對四種狀態的診斷率為100%，沒有出現過擬合現象。

表8 故障診斷結果

7 結論

為了有效識別滾動軸承故障類型，提高診斷準確率，筆者提出基于EMD、MSE算法和SVM的故障診斷方法。

表6 IMF分量相關因子

表7 故障診斷特征向量

▲圖7 IMF1分量MSE

（1）采用db10小波基對原始信號進行7層分解，利用閾值量化處理小波因子，與原始信號比較，確認重構信號故障特征更為明顯，更有利于特征提取。

（2）針對滾動軸承振動信號的非線性和非平穩性，采用EMD和MSE算法將降噪后的信號分解為若干IMF分量，對前3個IMF分量進行希爾伯特包絡分析，可以得到內圈、外圈的故障頻率，但滾動體故障不能很好分辨，因此要對包含故障信息的IMF分量MSE作進一步分析。

（3）選擇最佳尺度因子，將篩選后的IMF分量樣本熵作為故障特征向量，輸入SVM進行訓練，得到故障診斷模型。經過試驗，這一模型對滾動軸承故障診斷準確率為100%，具有很好的滾動軸承故障識別效果。