SE-KSD 優化的FRFT-VMD 軸承故障診斷方法

韋明輝 , ，江麗霞 , ，涂鳳秒 , ，姜蓬勃 ,

（1.西南石油大學 機電工程學院，成都 610500；2.石油天然氣裝備技術四川省科技資源共享服務平臺，成都 610500）

滾動軸承對機車等旋轉器械是不可或缺的重要部件之一，但由于其工作環境致使其具有易損壞的特點，而滾動軸承的運行狀況對于設備有著巨大影響，任何一個小故障都可能影響設備的正常運行從而給國家以及人民帶來巨大損失。若能在滾動軸承失效之前便對其進行故障診斷，不僅能最大限度減少損失，還能降低設備的維修成本[1]。因此，對滾動軸承進行故障診斷具有必要性。然而，采集到的滾動軸承故障信號?；煊写罅扛蓴_信號，使得所采集到的信號具有非平穩、非線性等特點。由于干擾信號的存在，導致在進行故障診斷時，無法從中快速辨別故障特征信號，使故障診斷的難度上升，時間增加。如何從這些信號中提取所需的故障特征信息成為許多學者研究的問題[2]。

Huang 等[3]提出了經驗模態分解方法（Empirical mode decomposition，EMD），該方法根據信號在時間尺度上的局部特征結構，自適應地提取反應信號本質特征的固有模態分量，可該算法缺少理論支撐，會出現端點效應和模態混疊現象影響信號處理結果。吳振華等[4]將噪聲輔助信號分析方法引入到EMD中，提出了集成經驗模態分解（EEMD），雖然該方法在抑制模態混疊方面起到了一定作用，可效果不徹底，且會產生多個虛假分量，而且噪聲的加入會導致重構誤差增加。Dragomiretskiy 等[5]提出了變分模態分解（Variational mode decomposition，VMD），一種全新的非遞歸自適應信號處理方法，近年來，該方法成為了一個研究熱點，但都只是在整數域上的改進。由于算法本身的局限性，過分解和欠分解的問題無可避免。

本文結合上述文獻，提出了一種基于樣本熵（Sample entropy，SE）和峭度均方差（Kurtosis standard deviation，KSD）優化的分數階變分模態分解（Fractional variational mode decomposition，FRFT-VMD）和隨機森林（Random forest，RF）分類器相結合的軸承故障診斷方法，首先，通過搜尋最小SE 值確定分數階Fourier 變換的最佳階數，并對信號進行最佳階數分數傅里葉變換，其次，通過峭度均方差準則來確定變分模態分解的最優參數，將變換后的信號進行變分模態分解，之后，計算分解后信號的峭度和脈沖因子作為特征向量，最后，將特征向量輸入至隨機森林分類器進行故障識別分類。通過與不同參數的FRFT-VMD分解所得各模態分量對比，本文方法所得的模態分量具有更多故障特征信息，與VMD-RF 的故障診斷方法比，本文作者所提出的基于SE-KSD 優化的FRFT-VMD-RF 方法所得故障特征信息聚類性更強，具有更高的識別準確率。

1 分數階Fourier 變換

1.1 分數階Fourier 變換定義

分數階傅里葉變換作為傅里葉變換在分數域空間上更為廣義的存在方式，使得信號可以在以時間軸旋轉到任意角度坐標軸上進行表示，在一定程度上對時域信息和頻域信息進行融合，進而將其某部分特征進行突出顯示[6]。

傳統的傅里葉變換被定義為存在于信號空間中的連續線性算子ξ，其特征方程為

式（1）中，傳統傅里葉變換所對應的特征值為λn=e-jnπ/2，其 特 征 函 數 為Hermite-Gauss 函 數Hn(t)e-t2/2，其中Hn(t)為n階Hermite 多項式，表達式為

從式（1）看出，隨機信號 ψn的傅里葉變換等價于其自身與復數 ψn之積。接下來將引入一個定義，該定義可由分數階傅里葉變換基本定義進行直接推導得到。

定義1 設 ψn表示為Hermite-Gauss 函數，作為普通Fourier 變換中特征值為 λn的特征函數，并可用作構成有限信號空間的特征函數。然后，分數階傅里葉變換可以定義為線性和滿足

根據上述表達式，有限能量函數xn可以展開為Fourier 變換特征函數的線性疊加，其表達式為

式（5）被稱為分數階Fourier 變換核的頻譜展開，Hermite-Gauss 函數 ψn須得滿足：

式中： α=pπ/2 ， 當p=1時，分數階Fourier 變換變為傳統Fourier 變換。同時用核函數的形式表示分數階Fourier 變換，即：

式中：Kp（u,t）為分數階傅里葉變換的核函數，且同時可以看出以p為變量的參數是以4 為周期的，因此F4n和F4n±2分別相當于恒等算子I和奇偶算子P。 當p=0時，經分數階傅里葉變換后得到的是原信號。其逆變換形式為

由圖1 可知，在時頻面上，將 （t,ω）按逆時針方向作 α角度的旋轉，變換到 （u,v）域，便可以得到新的信號表示。因此，利用分數階傅里葉變換，通過改變角度 α 或 階數p，可以觀察到信號從時域到頻域的所有特征信息。

圖1 時頻平面上分數階傅里葉變換旋轉示意圖Fig.1 Schematic diagram of fractional Fourier transform rotation on the time-frequency plane

1.2 分數階Fourier 變換階次選擇

在分數階Fourier 域可更好解決時頻域難以解決的數據可分性差等問題，但階次的選擇對其影響較大。

本文采用樣本熵[7]為指標，度量數據在分數階域的集中程度[8]。樣本熵是用來評價混亂度的指標[9]，樣本熵值越大，表明包含的有效信息越少，反之表明表明包含的有效信息越多[10]。

樣本熵算法步驟如下：

1）令原始數據 {Xi}={x1,x2,···xn} ，長度為N，根據原始信號重構一個m維的向量，即

2）定義dij為xi與xj的距離，且為兩者對應元素差值的絕對值的最大值，即

3）計算dij小于相似容限r的數目及此數目與距離dij的總數N-m-1 的比值，用Bmi(r)表示，即

4）計算Bmi(r)的平均值

5）對維數m+ 1 重復上述過程1）到5），得到(r) ，便可得到得到Bm+1（r）

6）樣本熵的定義為

當N為有限數時，式（14）可表示為

2 變分模態分解

2.1 變分模態分解定義

VMD 算法核心思想是認為原始信號每個模態的絕大多數特征內容都是緊緊圍繞在某一中心頻率周圍的，通過維納濾波、頻率混合和希爾伯特變換的方式將各子模態的帶寬問題分解為多個約束問題，再利用一系列迭代求解的方法得到中心頻率的最優解[11]。相較于傳統EMD 方法相比，VMD 的數學推導基礎更加堅實，具有更好的收斂表現和更出色的魯棒性。

綜上所述，VMD 最后是將給定信號f(t)分解成k個固有模態分量（IMF），并且保證分解后的IMF分量的總帶寬和最小。該算法實現步驟如下；首先通過希爾伯特變換得到各個IMF 分量的解析信號，并且計算出相應的單邊頻譜[12]，即

之后預先估計各個IMF 分量的中心頻率，通過頻率混合使得各個IMF 的頻譜調制到對應的基頻帶上，即

最后由高斯平滑指標計算得到調制信號梯度平方L2的L2范數，據此估計每個IMF 分量的帶寬，通過約束方程來使得IMF 的帶寬的和為最小，該方程表達式為

式中： dt為對函數求時間t的偏導數； {uk}為信號f(t) 分 解得到的k個IMF 分量， {uk}={u1,u2,u3,···,uk}；{ωk} 為 各IMF 的中心頻率， {ωk}={ω1,ω2,ω3,···,ωk}。引入含有二次懲罰因子 α和拉格朗日乘數算子 λ的增廣拉格朗日方程，使得上述變分約束方程轉變為了非約束性變分方程，從而得到該模型的最優解：

因為懲罰因子和拉格朗日乘數算子的增加，從而使得模型在高斯白噪聲存在的情況下具有更好的重構精度和更加嚴格的約束條件。通過交替方向乘子法對式（19）中的、及 λn+1進行交替更新，從而得到非約束變分方程的鞍點。

利用埃爾米特實信號的對稱性，將式（21）轉換為非負頻率的半空間積分方程，同時進行二次優化可得的更新方程為

同理根據帶寬之和最小的約束條件，再利用傅里葉等距變換將最小問題轉換到頻域中，得到與 λn+1交替更新方程：

VMD 采用遞歸循環的方式對原始信號進行層層剝離，有效避免了EMD 模態混疊的缺點，具有自適應的能力[13-14]。

2.2 變分模態分解算法參數選擇

變分模態分解算法雖然克服了傳統經驗模式分解及其改進方法的缺點，但分解前需要設定分解層數K和懲罰因子α，參數的選擇對分解結果的影響很大[15]。本文采用峭度均方差準則來對參數進行尋優[16]。峭度是反映隨機變量分布的數值統計值，對滾動軸承振動信號類的沖擊信號非常敏感[17]。均方差可以反映數據組中個體之間的離散程度，衡量分布程度的結果。通過峭度均方差準則則可反應VMD 分解后不同模態分量之間的差異性，峭度均方差值越大，表明不同模態分量間差異越大，反之則表明不同模態分量間差異越小。

3 軸承故障特征提取

本文提出一種基于樣本熵和峭度均方差優化的分數階變分模態分解的故障特征提取新方法。該方法的核心思想是通過樣本熵搜尋分數階傅里葉變換的最優階次，把原始數據中可分性差的數據映射到合適的分數階空間，再通過峭度均方差準則得到變分模態最優參數后在分數域對其進行變分模態分解。

實驗數據由CWRU 電氣實驗室軸承數據中心提供，本文選取的軸承不同狀態數據如表1 所示。

表1 軸承不同狀態數據Tab.1 Different status data of bearings

在故障特征提取部分，選用外圈故障這一組數據，該數據故障特征頻率為107.365 Hz，基于分數階Fourier變換的變分模態分解算法流程圖如圖2 所示。

圖2 基于分數階Fourier 變換的變分模態分解算法流程Fig.2 Procedures of variational modal decomposition algorithm based on fractional fourier transform

步驟1 計算分數階階次P從0.01 依次增加0.01 至2 的樣本熵值，樣本熵最小值處的階次即為分數階傅里葉最優階次，對數據進行最佳階數分數階Fourier 變換。

步驟2 計算α從1 000 依次增加100 到10 000，K從2 依次增加1 到11 的峭度均方差值，進行變分模態分解的分解層數K和懲罰因子α尋優，對分數階Fourier 變換后的數據進行變分模態分解。

步驟3 對變分模態分解后的每個變分模態分量進行分數階Fourier 逆變換。

步驟4 對每個變分模態分量進行包絡譜分析。

計算數據在0 ～ 2 階次分數階Fourier 變換時的樣本熵值，由圖3 可知在X＝93 時，樣本熵值最小，所以分數階Fourier 最優階次為0.93。

圖3 分數階Fourier 最優階次Fig.3 Fractional Fourier optimal order

計算最優階次分數階傅里葉變換后懲罰因子α 由1 000 依次遞增100 到10 000，分解層數K由2 依次遞增1 到11 的峭度均方差值，由圖4 可知，懲罰因子α增加了15 次，即為2 500，分解層數K增加了3 次，即為4，由此可知變分模態分解最優參數α和K分別為2 500 和4。圖5 為基于SE-KSD 優化的FRFT-VMD 分解后的各模態分量。

圖4 懲罰因子和分解層Fig.4 Penalty factor and decomposition layer

圖5 分數階VMD 分解各分量Fig.5 Fractional order VMD decomposition of each component

從圖6 可看出，基于SE-KSD 優化的FRFT-VMD分解所得各模態分量中包含更多更明顯的故障特征信息。對比不同參數的FRFT-VMD 所提取的特征頻率的幅值，本文所提出的故障特征提取方法所得的特征頻率及其相應的倍頻的幅值明顯大于其他參數FRFT-VMD 的幅值，表明本文所提方法可以為軸承故障診斷提供更加有效的故障特征信息。

圖6 不同參數的FRFT-VMD 各分量包絡譜圖Fig.6 Envelope spectra of FRFT-VMD components with different parameters

4 實驗驗證

4.1 實驗流程

為了進一步驗證該方法在軸承故障診斷中的作用，將采用西儲大學數據庫數據以及自己搭建的實驗平臺所測數據進行驗證，實驗流程如圖7 所示。

圖7 實驗流程圖Fig.7 Experimental flow chart

4.2 數據庫故障診斷

采用的實驗數據來自西儲大學軸承數據庫，選用的軸承不同狀態的數據如表1 所示。軸承各狀態下原始時域圖如圖8 所示。

圖8 軸承各狀態下原始時域圖Fig.8 Original time domains diagram of bearings in different states

為了說明基于SE-KSD 優化的FRFT-VMD 方法在軸承故障診斷中的作用，本文將峭度和脈沖因子作為故障特征向量，對基于SE-KSD 優化的FRFT-VMD 所提取的特征向量和基于VMD 所提取的特征向量進行K值聚類分析，從圖9 聚類結果圖可知，基于SE-KSD 優化的FRFT-VMD 所提取的特征向量具有較強的聚類性，不同故障之間的區分度更高，而基于VMD 所提取的特征向量聚類效果較差，故障類型之間區分不清且較分散。

為了進一步驗證所提方法的有效性，本文采用峭度和脈沖因子作為故障特征向量，使用隨機森林分類器[18-19]進行故障診斷。首先將不同類型的故障分為10 000 個數據為一組的10 組數據，總共分為40 組數據，對這40 組數據分別進行基于SE-KSD優化的FRFT-VMD 分解，再對所得變分模態分量提取峭度和脈沖因子作為故障特征向量，每個故障類型將會得到40 組峭度和脈沖因子特征向量，將各種故障類型中32 組作為隨機森林分類器的訓練集，8 組作為測試集數據對整體模型的分類準確率進行驗證，并將各軸承狀態進行編號，如表2 所示，設置隨機森林分類器中的決策樹數目為800，隨后進行故障診斷。將相同數據使用變分模態分解后，提取相應故障特征向量，同樣使用隨機森林進行故障診斷。

表2 各狀態軸承對應編號Tab.2 Corresponding numbers of bearings in each state

兩種方法訓練得到的隨機森林結構圖如圖10所示。通過測試集的測試，結果顯示兩種方法的分類結果圖如圖11 所示，從圖中可以看出，基于SE-KSD優化的FRFT-VMD 特征提取方法所提取的不同故障類型間的特征向量可以更好的進行識別并分類，分類結果圖顯示預測結果與期望值基本吻合，而基于VMD 特征提取方法的分類結果中，預測結果與期望值之間有較大出入。

圖10 隨機森林結構圖Fig.10 Random forest structure

圖11 分類結果圖Fig.11 Classification results

兩種方法分類結果如表3 和表4 所示。由表3和表4 可知，基于SE-KSD 優化的FRFT-VMD-RF故障診斷方法所得的分類總的準確率結果為96.875%，且對于編號1、2、3 的故障類型診斷準確率為100%。使用基于VMD-RF 故障診斷方法所得的分類總的準確率結果為78.1%。實驗結果顯示，基于SE-KSD 優化的FRFT-VMD 分解所得到的特征向量具有更好的區分度，且訓練得到的分類器模型具有更好的泛化能力。

表3 基于SE-KSD 優化的FRFT-VMD-RF 分類結果匯總Tab.3 Summary of FRFT-VMD-RF classification results based on SE-KSD optimization

表4 基于VMD-RF 分類結果匯總Tab.4 Summary of VMD-RF classification results

4.3 實測數據故障診斷

本節由自行搭建的滾動軸承實驗平臺所測數據進行故障診斷實驗，滾動軸承實驗平臺由旋轉機械和數據采集系統兩部分組成，如圖12 所示。旋轉機械部分由正常/缺陷滾動軸承和三相鼠籠式電動機組成，軸承安裝在三相鼠籠式電動機的風扇端，轉速通過光電傳感器測得。數據采集系統由安裝于電機風扇端的ICP 加速度傳感器和INV 3060A 型24 位網絡分布式采集儀（含采集內嵌服務軟件Coinv DASP V10）組成。軸承損傷方式采用EDM。在內圈、外圈和滾動體上分別加工一條深度和寬度為0.2 mm 的直線，模擬內圈故障、外圈故障和滾動體故障。采樣頻率為12 kHz。實驗平臺如圖13 所示。

圖12 實驗平臺框圖Fig.12 Experimental platform's block diagram

圖13 實驗平臺圖Fig.13 Experimental platform

通過實驗平臺獲取了4 種狀態下振動數據，長度均為119 000 點左右。將每種狀態的振動數據進行截取分組，每組10 000 個采樣點，每種狀態20 組，共得到80 組數據，對這80 組數據分別進行基于SE-KSD 優化的FRFT-VMD 分解，再對所得變分模態分量提取峭度和脈沖因子作為故障特征向量，每個故障類型將會得到80 組峭度和脈沖因子特征向量，將各種故障類型中56 組作為隨機森林分類器的訓練集，24 組作為測試集數據對整體模型的分類準確率進行驗證，設置隨機森林分類器中的決策樹數目為800，訓練得到的隨機森林結構圖如圖14 所示，通過測試集的測試，結果顯示分類準確率為90.625%，分類結果如圖15 所示。

圖14 實測數據隨機森林結構圖Fig.14 Random forest structure of measured data

圖15 實測數據基于SE-KSD 優化的FRFT-VMD-RF 分類結果Fig.15 FRFT-VMD-RF classification results of measured data based on SE-KSD optimization

通過診斷結果得知，本文提出的基于SE-KSD優化的FRFT-VMD-RF 軸承故障診斷方法不僅對于數據庫的數據能夠進行故障診斷的準確率較高，將該方法應用于實測振動數據的故障診斷，也能對故障類型進行準確率較高的分類識別。

5 結論

1）對基于SE-KSD 優化的FRFT-VMD 處理后的信號與不同參數的FRFT-VMD 處理后的信號進行包絡譜分析對比結果表明，與不同參數的FRFTVMD 相比，該方法處理后的信號能夠獲得更為明顯以及更多的故障特征信號。

2）與單獨VMD 提取的特征向量進行K值聚類對比，基于SE-KSD 優化的FRFT-VMD 提取的特征向量聚類性更強，不同故障之間區分度更大。為故障分類提供了優良的特征向量。

3）進一步的從基于SE-KSD 優化的FRFT-VMD分解后的信號中提取峭度和脈沖因子作為特征向量，送入隨機森林分類器進行故障診斷，與單獨的變分模態分解方法相比，本文提出的方法故障診斷的準確率高出許多，達到了96.875%，且用于實測數據的故障診斷準確率也較高。