基于VMD-HMM的滾動軸承磨損狀態識別

李奕江， 張金萍， 李允公

(1.沈陽化工大學 機械工程學院，沈陽 110142；2.東北大學 機械工程與自動化學院，沈陽 110004)

滾動軸承是旋轉機械的重要組成部分，在數控機床和機器人等運動精度要求較高的領域，其運行狀況直接關系到設備的精度和可靠性。磨損是滾動軸承失效形式中較為常見的一種，隨著軸承長時間的工作，摩擦導致軸承內部間隙增大，粗糙度增加，軸承的運轉精度因而下降。在精密機械中所廣泛應用的軸承經常是因為磨損失效而導致設備的運轉精度降低。相關研究表明百分之三十的旋轉機械故障是由軸承的失效引起的[1]。因此，對滾動軸承進行故障診斷和磨損狀態監測的研究非常重要。

在滾動軸承磨損失效方面，大量學者進行了相關研究，文獻[2]綜述了聲發射技術在滾動軸承磨損研究上的進展；文獻[3]在預測軸承磨損壽命上使用了灰色系統理論，得到了較為精確的預測結果；文獻[4]仿真驗證了圓柱滾子軸承的不均勻磨損特性；然而上述文獻對滾動軸承的磨損狀態分類研究少有涉及，也沒有找到可以分類的特征值。

隱馬爾可夫模型具有強大的動態時間處理能力，可以對軸承的磨損狀態進行識別與診斷[5]；在以往的軸承故障診斷研究中，其信號分解常常采用的方法有：小波變換，希爾伯特變換、經驗模式分解(Empirical Mode Decomposition,EMD)、總體經驗模式分解(Ensemble Empirical Mode Decomposition,EEMD)等。目前，文獻[6]提出了基于二進制小波包變換的軸承故障診斷方法，能夠準確過濾出故障瞬態沖擊成分；文獻[7]提出了基于EMD的無量綱索引方法，可以實現旋轉機械的故障診斷；文獻[8]提出一種基于EEMD的多變量多尺度統計過程監測方法，開發的方法能夠有效地提取軸承故障特征，得到可靠的軸承故障檢測和診斷策略。上述方法在軸承故障診斷上取得了良好的效果，但也存在著一些局限性，小波分析中的小波基難以選擇，需要較平穩的信號且抑制脈沖干擾能力較低。EMD的算法效率低下，存在模式混疊現象，EEMD雖然抑制了模式混疊現象，但算法效率低的弊端依然存在。由于現場采集的軸承振動信號是非線性非平穩的，往往存在窄帶脈沖干擾和隨機信號，大量的干擾會掩蓋軸承的故障特征信號，嚴重影響了軸承磨損狀態的識別。

變分模態分解(Variational Mode Decomposition, VMD)是一種誕生于2014年的新的非遞歸、自適應的數據分解方法[9]，它是一種新的數據分解算法，分解層數少算法效率高，沒有模態混疊和端點效應[10]，且分解后所得的各層IMF(Intrinsic Mode Function)中包含著豐富的故障信息。因此本文從滾動軸承不同磨損狀態的振動信號入手，充分利用VMD分解對非平穩信號的處理能力和HMM(Hidden Markov Model)模型極好的模式識別功能，將VMD分解和能量熵與HMM模型結合起來，應用到滾動軸承的磨損狀態識別中，最后通過實驗驗證了該方法的準確性與穩定性。

1 變分模態分解

變分模態分解在“變分”過程中使用二次懲罰因子α和拉格朗日懲罰算子λ(t)，使有約束變分轉化為無約束變分，然后通過交替方向乘子法迭代搜尋變分模型最優解，來不斷更新每個模式分量的帶寬和頻率中心，最終自適應的實現信號的有效分解。

1.1 VMD分解原理

VMD算法分解所得的本征模式分量(Local Mean Decomposition,IMF)與EMD和LMD分解的IMF不同，它被重新定義為一個調幅-調頻信號，數學表達式為

uk(t)=Ak(t)cos[φk(t)]

(1)

式中：Ak(t)為uk(t)的瞬時幅值，且Ak(t)≥0;φk(t)為相位，它為非負函數，且φk(t)≥0；ωk(t)為uk(t)的瞬時頻率。

(2)

式中:uk(t)可以認為是一個諧波信號，它的幅值為Ak(t)、頻率為ωk(t)。

VMD算法分解信號x(t)后，得到K個離散的的模態uk(t)(k∈1,2,…,K)，它在頻域中具有稀疏特性，通過以下方法估計各模態信號頻率帶寬：

(1)通過對模態函數uk(t)進行Hilbert變換，得到邊際譜。

*uk(t)

(3)

(2)調整指數修正，將每個估計的uk(t)中心頻帶調制到相對應的基頻帶。

(4)

(3)計算解調信號梯度的平方L2范數。

(5)

在上述求解的基礎上構造VMD變分約束模型，如式(5)所示

(6)

為了求解式(5)的變分約束模型的最優解，通過引入二次懲罰因子α和拉格朗日乘法算子λ(t)，即可將約束模型轉化為式(6)所示的非約束模型問題，該增廣拉格朗日函數如式(7)所示

L({uk},{ωk},λ)=

(7)

不斷地使用ADMM法迭代搜索找到拉格朗日函數L的極小值點，原始信號x(t)被最終分解為K個模態函數uk(t)。

從以上VMD原理可以看出，VMD通過頻域上的更新將傅里葉逆變化轉換到時域上，相比于經驗模態分解和局部均值分解，VMD具有更好的數學模型和理論且避免了模態混疊和端點效應，在機械振動信號和心電圖信號的處理上有良好的效果。文獻[11]將VMD分解的自適應形態學的特征提取方法并與粒子群算法相結合。通過仿真和試驗分析，將其與EMD的自適應形態學、包絡解調方法比較，結果表明VMD具有較高的精度，分解層數少效率高，克服了EMD的缺點，最終有效的提取了滾動軸承故障信號的特征頻率。

2 HMM模型

隱馬爾可夫模型(HMM)是一個概率模型，經常用來描述一個系統的隱性狀態之間的轉移和表現概率,已經在語音識別方面取得了顯著的成果[12]。HMM是一個雙重隨機的過程，在該過程中隱狀態必須滿足馬爾可夫性且隱狀態能夠大概被估計。隱馬爾可夫模型中的每個狀態不能被直接觀測到但是可以通過觀測值和隨機過程進行間接的感知。每個觀測值可以通過概率密度分布表現為它所對應的狀態。由于其在時間序列上有著強大的分類處理能力，適合對滾動軸承的磨損狀態演變進行分類。

一個離散隱馬爾科夫模型一般由以下幾個參數來表達:

(1)馬爾可夫鏈隱狀態數N。N個狀態為S={S1,S2,…,SN}，t時刻的狀態為qt,qt∈S。

(2)在各個隱狀態下的可觀測值數M。M個觀測值為L={L2,L2,…,LM}，t時刻的觀測值為Ot,Ot∈L。

(3)隱狀態間狀態轉移概率矩陣

A={αij},αij=(qt+1=Sj|qt=Si),1≤i,j≤N

(4)觀測值概率矩陣

B={βj(k)},βj(k)=P(Ot=Lk|qt=Sj),

1≤j≤N,1≤k≤M

(5)初始狀態概率分布矢量

π={πi},πi=P(q1=Si),1≤i≤N

該離散HMM模型由N,M,A,B和π 5個參數表示，一般簡記為λ=(π,A,B)。

3 VMD-HMM磨損狀態識別過程

隨著滾動軸承磨損狀態的不斷改變，軸承振動信號經VMD分解后的IMF分量也會有相應的改變。應用HMM在滾動軸承故障診斷中已經取得了一定的成果，而VMD分解在自適應性、正交性上有很大的優勢且避免了EMD和LMD等的模態混疊和端點效應等缺點。因此本文提出了基于VMD-HMM的滾動軸承磨損狀態識別算法，圖1所示為利用VMD-HMM進行狀態識別的主要步驟。

(1)隨機選取滾動軸承磨合時期、正常使用、失效前期和失效期4種振動信號各20組共80組來作為訓練信號，對每組信號進行VMD分解，根據能量熵定義計算每層IMF分量的能量熵參數，提取滾動軸承振動信號各層IMF分量能量熵構成能量熵序列共80組，然后將能量熵序列矢量量化后利用此序列輸入到HMM模型中使用Baum-Welch算法進行訓練，得到對應的各類磨損狀態分類，形成磨損狀態模型庫。

(2)模式識別，對待診斷信號(剩余的80組)進行VMD分解并計算各層IMF能量熵，得到的能量熵序列矢量量化后，將此能量熵序列輸入到已經訓練好的4種磨損狀態模型庫中，HMM模型將會輸出對數似然概率lnP(O∣λ)，所輸出的最大的對數似然概率所對應的磨損狀態即為該待診斷信號的磨損狀態，從而實現軸承磨損狀態的識別。

圖1 VMD-HMM模型磨損狀態識別過程Fig.1 Patter recognition procedure with VMD-HMM model

4 仿真實驗

為了驗證本文所使用的VMD算法的有效性，利用式(8)所示的仿真信號進行實驗，并與EMD算法進行比較。

s(t)=x1(t)+x2(t)+x3(t)+n(t)

(8)

x1(t)=0.6sin(2π×10t)

(9)

x2(t)=0.4sin(2π×50t)

(10)

x3(t)=0.2sin(2π×100t)

(11)

式中：x1(t),x2(t),x3(t)為幅值和頻率均不一樣的正弦信號；n(t)為均值為0、方差為1的高斯白噪聲。

設置仿真信號的采樣頻率fs=1 000，采樣點數N=500，仿真信號的時域波形圖如圖2所示，然后分別對仿真信號進行EMD,VMD分解，分解后得到的結果如圖3～圖6所示。

圖2 仿真信號時域波形Fig.2 Time-domain waveform of simulation signals

圖3 仿真信號EMD分解結果Fig.3 The result of simulation signal after EMD

圖4 EMD分解結果幅頻圖Fig.4 Amplitude spectrum of EMD decomposition results

圖5 仿真信號VMD分解結果Fig.5 The result of the simulation signal after VMD

圖6 VMD分解結果幅頻圖Fig.6 Amplitude spectrum of VMD decomposition results

由圖3和圖4可以看出，EMD分解后的信號波形不規則，局部出現失真；IMF2出現了50 Hz和100 Hz兩種頻率成分，出現了模式混疊；EMD分解需要較多的分解層數才能達到好的分解效果。與圖5和圖6進行對比，VMD分解后所有的正弦頻率被完整分解出來，IMF1對應x1(t)、IMF2對應x2(t)、IMF3對應x3(t)、IMF4為噪聲信號，且波形和幅值均與原信號一致；各IMF之間沒有模態混疊現象；VMD所需的分解層數較少，且效果更好；VMD分解有更好的噪聲魯棒性。通過對比可以發現本文采用VMD進行信號分解可以有效提取出信號的特征頻率且能抑制白噪聲。

5 實驗與分析

為了驗證以上方法的有效性和準確性，進行了滾動軸承磨損狀態全壽命周期實驗。實驗裝置如圖7所示，實測軸承來自某機械加工廠的一臺線切割機床導輪軸承，信號采集與檢測系統包括北京東方振動與噪聲研究所得INV3018C智能信號處理分析儀、INV9822型ICP加速度傳感器及DASP信號采集分析系統。

圖7 實驗測試情況Fig.7 Experimental test conditions

由于加工環境造成該線切割導輪軸承經常因磨損失效導致導輪精度降低，所以采集該線切割導輪軸承的全壽命周期振動信號，該軸承經過磨合，正常工作，性能退化直至磨損實效，最后對其進行更換處理。在導輪總成上安裝傳感器，從軸承開始更換時開始，每兩天采集一次數據，持續60天采集了30組數據，采樣頻率19 682 Hz。

5.1 VMD分解與能量熵特征提取

根據軸承的性能，按時間順序將其服役分為4個階段：磨合階段、正常階段、早期磨損階段、磨損失效階段。磨合階段的軸承剛開始使用，正常階段的軸承可以正常使用，早期磨損階段下的軸承有輕微的磨損，但是不影響使用，磨損失效階段下軸承不能再繼續使用，應進行更換。圖8為測試軸承在正常狀態和磨損失效狀態下的時域振動波形圖。將采集到的滾動軸承全壽命周期數據按軸承磨損狀態分為4組，每種狀態各分為40組數據共160組，對每組數據進行VMD分解，以軸承正常狀態為例，選取最優分解層數為6層，懲罰因子α=500，其中一組正常狀態的軸承數據分解結果如圖9所示。圖10為該組數據的EMD分解結果。

(a)正常狀態

(b)磨損狀態圖8 實測信號Fig.8 Measured signal

從圖10的分解結果可以看出IMF5,IMF6,IMF7三者之間存在模態混疊，同時EMD的分解層數較多，對接下來的數據處理會造成不便。反觀圖9的VMD分解結果，結果IMF1～IMF6由低頻到高頻排列，可以看出各IMF之間不存在模態混疊，且分解層數比EMD少，更便于后期處理數據。實驗結果同時驗證了VMD分解的優勢。

然后將160組樣本信號進行VMD分解，得到960組本征模態函數，根據能量熵的定義求取各組本征模態函數(IMF)的能量熵，IMF1～IMF6的特征參數為F1～F6，本文給出每種磨損狀態各3組能量熵如表1所示。

圖9 正常狀態軸承信號VMD分解結果Fig.9 Normal state bearing signal VMD decomposition result

圖10 正常狀態軸承信號EMD分解結果Fig.10 Normal state bearing signal EMD decomposition result

5.2 HMM模型訓練與磨損狀態識別

HMM建模時采用VMD分解得到的6層IMF能量熵組成特征向量F={F1,F2,F3,F4,F5,F6}作為觀測序列，初始概率π，初始狀態轉移矩陣A，初始觀測值概率矩陣B均由隨機函數取得，得到160×6特征矩陣。

在HMM的建模與訓練時，把磨合狀態、正常狀態、早期磨損狀態、磨損失效狀態代表4種隱狀態，記作λ1,λ2,λ3,λ4,HMM模型的觀測序列為上述提取的6維能量熵特征序列；選擇80組特征向量(每種狀態20組)作為訓練的樣本，Baum-Welch算法訓練HMM模型需要輸入的特征為離散值，所以對每組特征向量進行矢量量化。HMM訓練過程中，迭代次數增加最大似然估計也隨之增加，最后模型達到收斂。訓練結束后得到了4種狀態所對應的HMM模型，圖11給出了4種狀態HMM的訓練曲線。所有模型的迭代次數在30次時達到收斂，收斂較快。

表1 滾動軸承各磨損狀態各層IMF能量熵Tab.1 IMF energy entropy of different wear state of roller bearings

經過上述訓練后得到4種磨損狀態HMM模型，將剩余的80組能量熵特征向量(每種狀態20組)作為測試樣本。測試前同樣對特征向量進行矢量量化，將矢量量化后的能量熵序列分別輸入到個狀態的HMM模型中，每個模型都會輸出一個對數似然概率lnp(O|λ)，這個對數似然概率表示觀測值序列與所對應的HMM模型的相似程度，當該估計值越大時觀測值與該HMM模型越相似，識別算法采用Viterbi算法，圖12給出了每種測試樣本在各HMM模型上的測試結果，通過對數似然概率可以準確識別出軸承當前所處的磨損狀態。圖12(a)中第10組數據分類錯誤，可能是因為該組磨合信號與正常信號有些相似。由圖12(b)可知，正常狀態的能量熵特征序列輸入到該狀態的HMM模型上時，其對數似然概率最大。從圖12還可以看出，隨著磨損狀態的加劇4種狀態的觀測序列在同一模型上輸出的對數似然概率也變大，滾動軸承4種磨損狀態的識別成功率較高。實驗結果表明VMD-HMM滾動軸承磨損狀態識別模型可以成功的識別出滾動軸承的磨損狀態，算法效率較高，識別準確率高，穩定性好。

圖11 HMM訓練曲線Fig.11 Curves of HMM training

(a)磨合狀態數據在各HMM模型上的輸出值

為了比較驗證該VMD-HMM模型的識別效果，采用EMD-HMM模型、諧波小波樣本熵與HMM模型對相同的訓練樣本進行處理，最后進行識別結果的對比。EMD-HMM模型采用文獻[13]介紹的方法，諧波小波樣本熵與HMM模型結合的方法采用文獻[14]中介紹的方法，最終得到的識別結果如表2所示。

由表2結果分析，在較少的樣本情況下，EMD-HMM、諧波小波樣本熵結合HMM和VMD-HMM對磨損狀態的識別率分別為90%,95%和98.75%，EMD-HMM模型的各狀態識別成功率在3種模型中最低，諧波小波樣本熵結合HMM模型的識別成功率高于EMD-HMM模型，VMD-HMM模型的各狀態識別成功率和綜合識別成功率在3種模型中最高。這表明VMD-HMM模型可以有更好的識別率且其性能更為優越，可以有效的解決軸承磨損故障的分類問題。

表2 3種模型識別正確率對比Tab.2 Recognition accuracy of three models

6 結 論

鑒于VMD分解和HMM的優點，本文提出了一種基于VMD分解和HMM模型的滾動軸承磨損狀態識別方法。最后，為了驗證所提出方法的可行性與有效性，實際記錄了一位于易受磨損部位的的滾動軸承60天全壽命周期振動信號。將該組信號分為四類磨損狀態，并訓練了對應的4種HMM模型，最終用待測數據實現了滾動軸承磨損狀態的識別，并與EMD-HMM模型和諧波小波樣本熵結合HMM模型進行了對比。實驗結果表明：

(1)VMD分解適合滾動軸承磨損的非平穩信號處理，具有很強的識別能力，分解出的各層本征模態函數沒有模式混疊和端點效應，可以反映出各個頻段的大致信息。

(2)基于VMD分解與能量熵結合，可以捕捉到不同磨損狀態間所隱含的信息，從而可以有效的提取出軸承磨損狀態的特征向量。

(3)HMM模型對軸承4種磨損狀態進行訓練，成功地對軸承磨損狀態進行了識別，識別率較高，非常適合軸承磨損這種隨時間變化，信號非平穩，特征重復性不佳的故障。