基于EEMD-Treelet和高斯過程的起重機齒輪故障診斷

李 勇，錢尼君，陳星斌，黃荔生，彭獻永

(1.江西省特種設備檢驗檢測研究院，南昌 330029；2.騰羿數據(上海)科技有限公司，上海 200240；3.上海交通大學 電子信息與電氣工程學院，上海 200240)

0 引言

齒輪傳動系統的主要機械振動通常歸因于齒輪故障，包括磨損/折斷的齒或不平衡的齒輪，并可能導致系統故障。眾所周知，對于不同的故障類型，其影響將以不同的頻率出現。準確而自動地評估其性能下降有兩個重要方面：特征提取和故障分類器。Randall[1]提供了有關齒輪診斷技術的綜述。目前，循環平穩性分析方法來診斷旋轉機械中的齒輪已經被廣泛采用。振動信號的頻譜分析也已經用以檢測齒輪中的剝落。

通常，在機器故障的情況下，測得的振動信號可能會出現不平穩的現象。因此，時頻分析方法，如小波變換方法和希爾伯特-黃氏變換(HHT)被用來觀察信號時變特性。Saravanan等[3]使用Morlet小波提取齒輪箱的故障特征作為SVM模型輸入。Wu and Chan[4]通過小波變換和人工神經網絡開發了齒輪故障識別模型。Peter等[5]利用小波分析和包絡線檢測進行軸承故障診斷。Peng等[6]比較了改進的希爾伯特-黃變換和小波變換在滾動軸承故障診斷中的性能差異。Lei和Zuo[7]用EEMD產生敏感的本征模態函數(intrinsic mode function,IMF)，用于旋轉機械的故障診斷。Fan和Zuo[8]利用本征模態函數提取機器故障特征。Yu等[9]用希爾伯特-黃變換的時頻熵方法進行齒輪故障診斷。Yang和Wu[10]用EEMD方法和PCA來診斷齒輪性能退化。

機器的故障信號通常顯示出很高的復雜性，許多研究人員利用機器學習來識別機械運行狀態故障。Samanta和Al-Balushi[11]提出了一種人工神經網絡，使用時域特征對滾動軸承進行故障診斷。Wang等[12]使用時間序列分析和神經網絡進行旋轉機械故障診斷。Shafei等[13]將神經網絡和模糊邏輯應用于齒輪的故障診斷。除了人工神經網絡，支持向量機(SVM，support vector machine)[2]和決策樹方法[14]也被用于不同的故障類型。Bordoloi和Tiwari[17]用SVM對齒輪頻域數據進行多故障分類。他們通過在4個不同的故障條件下檢查SVM模型的多類能力的性能。Rajeswari等[18]提出了一種基于遺傳算法和反向傳播神經網絡的齒輪不同狀態診斷混合框架。Sugumaran和Ramachandran[19]使用齒輪振動信號的時域統計和直方圖特征，作為SVM分類器輸入，找到了最佳故障特征。Zhang Ying等[2]提出了基于多種特征的故障診斷方法，他們用核主成分分析和粒子群優化支持向量機(來定位軸承故障。

高斯過程(GP)是一種基于核的非線性學習方法，可以處理帶有噪聲觀測值和系統不確定性的建模問題[20]。它已經被廣泛的應用于動作識別等。它具有嚴格的統計理論基礎和很好的適應復雜問題的能力。與ANN和SVM相比，它既保持良好性能又具有靈活的非參數推理能力的條件。在GP模型中，關鍵問題是估計一組適當的超參數，該參數會使邊際似然函數最小化。為了提高GP模型參數尋優的準確性[15]，用一種細菌覓食優化(BFO)算法[16]來優化GP模型的超參數。

為了提高起重機齒輪傳動系統故障識別性能，提出一種基于EEMD-Treelet的自動特征選擇方法，用于診斷齒輪故障。首先，采用EEMD方法來提取通常隱藏在振動信號中的與故障相關的特征。然后，通過Treelet變換對原始特征進行降維，降維后的特征可以表示齒輪的主要故障特征。所得的特征向量作為BFO-GP模型輸入，自動識別齒輪故障。

1 特征提取與降維

1.1 邊際譜

EEMD方法[3]的本征模式函數是單振蕩成分。因此，IMF瞬時頻率計算是可行的。這是建立時頻分布的必要和關鍵步驟。為了建立信號的時頻分布，第k個IMF可以用其包絡和載波部分表示為:

(1)

其中：包絡線Ak(t)表示隨時間變化的幅度項，可以通過復雜的三次樣條擬合計算得出，而被測信號x(t)的時頻分布H(ω,t)可以表示為：

(2)

此外，邊緣頻譜S(ω)可以通過積分沿時間軸的時頻分布表示如下：

(3)

其中：Hk(ω,t)表示第k個IMF的時頻分布，T表示信號長度，Ω表示所選IMF的索引集。

1.2 Treelets降維

Treelets降維算法是一個自適應的降維算法,它在數據每層結構上將兩列最相似的數據用一組逼近信號與一組細節信號來代替，采用PCA方法進行對這兩列數據進行分析[21]。

Treelets算法的結構可以表示為一個樹結構。圖1即為一個5層的樹結構信號：部分(主成分結構)代表的為逼近信號，(次成分信號)代表與主成分正交的細節信號。

圖1 五層Treelets結構

2 BFO-GP分類器

2.1 GP分類器

一個高斯過程可以通過其均值函數m(x)和它的協方差k(x,x')確定。它的超參數最佳估計可以在對數邊際似然函數最大情況下獲取。也就是說，可以通過對目標函數進行邊際化，獲得高斯過程超參數的概率，使用共軛梯度方法優化超參數。高斯過程的對數邊際似然函數定義為：

logp(y|x,Θ)=

(4)

超參數的選擇對于高斯過程分類器的整體分類性能至關重要。它的超參數{μ,σ,l1,l2,…,ld}可以在很大的值范圍內變化。通常的共軛梯度優化方法的性能或準確性取決于初始數值的選擇，并且相對于小擾動和局部最優而言可能不穩定。由于選擇超參數會極大地影響高斯過程模型的泛化性能，因此進化算法可能是一個可好的選擇。在這個工作中，細菌覓食優化算法被用作高斯過程超參數優化算法。細菌覓食優化算法的適應度函數采用式(7)中定義的對數邊際似然函數。

2.2 細菌覓食算法

BFO是一種自然啟發式優化算法[16]，它能夠找目標函數的最小值，并且不需要計算目標函數的梯度。BFO模仿了在實際細菌系統中觀察到的4個主要機制：趨化性，群聚，繁殖和遷徙，以解決非梯度優化問題。

1)趨向性操作：

大腸桿菌在整個覓食過程中有兩個基本運動：旋轉和游動。旋轉是為找一個新的方向而轉動，而游動是指保持方向不變的向前運動。BFO算法的趨向性操作就是對這兩種基本動作的模擬。通常，細菌會在食物豐盛或環境的酸堿性適中的區域中較多地游動，而在食物缺乏或環境的酸堿性偏高的區域則會較多地旋轉，即原地不動。

假設要找到目標函數J(θ)的最小值，其中θ∈Rp(即θ是實數p維向量)，并且沒有測量值或梯度▽J(θ)的數學解析描述。假設θi(j,k,l)代表第j個趨化性的第i個細菌,第k個繁殖和第l個消除-擴散事件。C(i)表示在翻滾指定的隨機方向上執行的步長。在計算趨化性中，細菌的運動可以表示為：

(5)

其中：Δ表示隨機方向上的向量，其元素位于[-1,1]。

2)群聚性操作：設目標函數的變化量如下：

(6)

其中:JCC(θ，P(j,k,l))表示目標函數J(θ)的目時變值，S是細菌總數，p是每個細菌中要優化的變量數量.θ=[θ1,θ2,…,θp]T是p維搜索域中的一個點。dattrac,wattrac,hrepel,wrepel是固定學習系數。

3)繁殖性操作：

生物進化過程一直服從達爾文進化準則，即“適者生存、優勝劣汰”。BFO算法執行一段時間的覓食過程后，部分尋找食物源能力弱的細菌會被自然淘汰，而為了維持種群規模不變，剩余的尋找食物能力強(適應度值低)的細菌會進行繁殖。在BFO算法中模擬這種現象稱為繁殖性操作。

(7)

4)遷徙性操作：

實際環境中的細菌所生活的局部區域可能會發生逐漸變化(如食物消耗殆盡)或者發生突如其來的變化(如溫度突然升高等)。這樣可能會導致生活在這個局部區域的細菌種群被遷徙到新的區域中去或者集體被外力殺死。在BFO算法中模擬這種現象稱為遷徙性操作。

2.3 BFO-GP 分類器

基于EEMD-Treelet和細菌覓食優化算法-高斯過程的齒輪故障診斷步驟如下：

步驟(1):數據準備，從原始振動信號中提取特征作為特征向量，標簽這些特征向量。使用Treelet消除原始特征向量的冗余；

步驟(2):GP分類初始化，使用常數均值函數和平方指數協方差核函數；

步驟(3):使用BFO訓練GP分類器：選擇不同的超參數選擇，使用式(4)計算相應的適應度函數值。設置和使用GP分類器：使用優化的超參數集來設置GP分類器。對于新的分類診斷問題，按照步驟(1)提取特征，并計算分類概率。

3 實驗結果與分析

實驗設備由起重機齒輪箱組成，模擬不同齒輪故障的工作條件。安裝在驅動軸上的1號齒輪與安裝在從動軸上的2號齒輪嚙合。驅動軸連接到由伺服電動機驅動的撓性聯軸器。兩個軸由四套軸承支撐。機械振動由加速度計測量。表1顯示了此實驗中6種齒輪故障情況。驅動馬達的旋速為1 500 RPM，驅動軸和從動軸的轉動頻率分別為25 Hz和11.3 Hz。根據軸的轉速和相應的齒輪嚙合頻率，將數據采集裝置的采樣率設置為6 400 Hz。此實驗中，每個實驗案例有50組振動信號，總共記錄了300組信號。

表1 實驗齒輪故障類型

3.1 特征提取

這里以案例3的信號為例進行說明。將所有振動信號用EEMD分解為本征模態函數 (C1,C2,…,Cn)，如圖2所示。然后，計算每個IMF的瞬時頻率，表2列出了每個瞬時值的平均值。結果顯示本征模態函數C7和C8與主動軸和從動軸的轉動頻率有關，分別為25 Hz和11.3 Hz。研究結果也顯示嚙合頻率575 Hz位于本征模態函數C2和C3的帶寬之內。由于在嚙合頻率附近的振動信號分量通常由與機器故障有關，因此提取C2和C3相關的特征，更具實際意義。為了提取與故障有關的特征，對C2和C3進行包絡分析如下。

圖2 案例2中振動信號的EEMD分解產生的IMF

表2 每個IMF的瞬時頻率分布的平均值

步驟1：將與嚙合頻率相關的C2和C3表示為S，提取S信號的包絡信號Envs。

步驟2：在Envs上應用EEMD和希爾伯特變換，得到S的希爾伯特時頻分布頻譜。

步驟3：沿時間軸對步驟2中的頻譜進行積分，得到S的希爾伯特邊際包絡譜。

步驟4：利用三次樣條擬合希爾伯特邊際包絡譜的包絡信號。

包絡分析結果顯示案例3的包絡信號在0～50 Hz范圍內比案例2的能量集中度更高，這是因為案例3的齒輪故障水平比案例2更為嚴重。齒輪故障在25 Hz(1號齒輪的旋轉頻率)處觀察到峰值。還要注意的是，驅動齒輪轉速峰值的大小取決于齒輪損壞的程度，這意味著嚴重磨損的牙齒(如案例3)的包絡信號具有比輕微磨損的情況更集中的能量。另外，由于斷齒的劣化程度比輕度和嚴重磨損的牙齒更嚴重，因此在25 Hz的情況下斷齒(案例4)的包絡譜峰值高于案例2和案例3。2號齒輪的主要特征頻率是11.3 Hz。在齒輪的邊緣希爾伯特包絡譜中觀察到嚴重不平衡情況的能量集中度高于輕微不平衡情況和正常情況(11.3 Hz)的能量集中。這是因為嚴重不平衡情況下的故障水平比輕微失衡和正常情況下要嚴重。

根據以上分析，從振動信號中提取的特征如表3所示。

從案例1-5中隨機選擇30個振動信號，形成樣本集，計算表3的特征，獲得21維振動信號狀態特征。Treelet變換算法層次設置為2，將這21維特征降為為2維。圖3顯示了經Treelet變換后齒輪特征數據分布。顯然，Treelet方法可以將正常數據與故障數據分開。每種案例的數據都在特征空間中較好地聚集。僅有一對不平衡數據特征存在重疊。但是，BFO-GP仍然可以獲得95%的診斷精度，這是令人滿意的。對于嚴重磨損的齒輪和斷齒的齒輪，也發現了相同的情況。

圖3 特征空間中的齒輪數據分布

表3 每個IMFS的瞬時頻率分布的平均值

BFO-GP輸出特征空間中每個數據點的分類概率。以0.5概率輪廓線作為決策線，將特征空間分為兩類。圖4顯示了齒輪診斷結果。顯然，在使用Treelet進行特征降維之后，正常齒輪和退化齒輪數據集在特征空間中得到很好的分類，并且彼此之間有著很好的區分。分類結果是令人滿意的。基于BFO-GP的方法是適合齒輪故障診斷。

圖4 齒輪診斷結果

進一步，C1-C55案例的故障診斷其結果如表6所示。顯然，BFO-GP的診斷精度高于其它模型如GA-GP，PSO-GP和BFO-SVM，但是BFO-GP消耗更多的運行時間成本。測試結果顯示對于GP的超參數優化能力，BFO比GA、PSO等其它算法表現得更好。

表6 不同案例的診斷結果

4 結束語

基于EEMD-Treelet-GP的齒輪故障診斷結果發現EEMD-Treelet算法能夠學習到更多振動信號精細特征。由EEMD-HT獲得齒輪振動信號的希爾伯特頻譜，更多齒輪故障特征被提取。通過Treelet變換對這些特征進行篩選，獲得線性不相關特征。這些特征輸入進BFO-GP模型，實現齒輪故障診斷與退化檢測。