李 勇，錢尼君，陳星斌，黃荔生，彭獻(xiàn)永

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

0 引言

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

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

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

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

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

1 特征提取與降維

1.1 邊際譜

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

(1)

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

(2)

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

(3)

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

1.2 Treelets降維

Treelets降維算法是一個(gè)自適應(yīng)的降維算法,它在數(shù)據(jù)每層結(jié)構(gòu)上將兩列最相似的數(shù)據(jù)用一組逼近信號與一組細(xì)節(jié)信號來代替，采用PCA方法進(jìn)行對這兩列數(shù)據(jù)進(jìn)行分析[21]。

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

圖1 五層Treelets結(jié)構(gòu)

2 BFO-GP分類器

2.1 GP分類器

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

logp(y|x,Θ)=

(4)

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

2.2 細(xì)菌覓食算法

BFO是一種自然啟發(fā)式優(yōu)化算法[16]，它能夠找目標(biāo)函數(shù)的最小值，并且不需要計(jì)算目標(biāo)函數(shù)的梯度。BFO模仿了在實(shí)際細(xì)菌系統(tǒng)中觀察到的4個(gè)主要機(jī)制：趨化性，群聚，繁殖和遷徙，以解決非梯度優(yōu)化問題。

1)趨向性操作：

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

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

(5)

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

2)群聚性操作：設(shè)目標(biāo)函數(shù)的變化量如下：

(6)

其中:JCC(θ，P(j,k,l))表示目標(biāo)函數(shù)J(θ)的目時(shí)變值，S是細(xì)菌總數(shù)，p是每個(gè)細(xì)菌中要優(yōu)化的變量數(shù)量.θ=[θ1,θ2,…,θp]T是p維搜索域中的一個(gè)點(diǎn)。dattrac,wattrac,hrepel,wrepel是固定學(xué)習(xí)系數(shù)。

3)繁殖性操作：

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

(7)

4)遷徙性操作：

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

2.3 BFO-GP 分類器

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

步驟(1):數(shù)據(jù)準(zhǔn)備，從原始振動(dòng)信號中提取特征作為特征向量，標(biāo)簽這些特征向量。使用Treelet消除原始特征向量的冗余；

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

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

3 實(shí)驗(yàn)結(jié)果與分析

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

表1 實(shí)驗(yàn)齒輪故障類型

3.1 特征提取

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

圖2 案例2中振動(dòng)信號的EEMD分解產(chǎn)生的IMF

表2 每個(gè)IMF的瞬時(shí)頻率分布的平均值

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

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

步驟3：沿時(shí)間軸對步驟2中的頻譜進(jìn)行積分，得到S的希爾伯特邊際包絡(luò)譜。

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

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

根據(jù)以上分析，從振動(dòng)信號中提取的特征如表3所示。

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

圖3 特征空間中的齒輪數(shù)據(jù)分布

表3 每個(gè)IMFS的瞬時(shí)頻率分布的平均值

BFO-GP輸出特征空間中每個(gè)數(shù)據(jù)點(diǎn)的分類概率。以0.5概率輪廓線作為決策線，將特征空間分為兩類。圖4顯示了齒輪診斷結(jié)果。顯然，在使用Treelet進(jìn)行特征降維之后，正常齒輪和退化齒輪數(shù)據(jù)集在特征空間中得到很好的分類，并且彼此之間有著很好的區(qū)分。分類結(jié)果是令人滿意的?；贐FO-GP的方法是適合齒輪故障診斷。

圖4 齒輪診斷結(jié)果

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

表6 不同案例的診斷結(jié)果

4 結(jié)束語

基于EEMD-Treelet-GP的齒輪故障診斷結(jié)果發(fā)現(xiàn)EEMD-Treelet算法能夠?qū)W習(xí)到更多振動(dòng)信號精細(xì)特征。由EEMD-HT獲得齒輪振動(dòng)信號的希爾伯特頻譜，更多齒輪故障特征被提取。通過Treelet變換對這些特征進(jìn)行篩選，獲得線性不相關(guān)特征。這些特征輸入進(jìn)BFO-GP模型，實(shí)現(xiàn)齒輪故障診斷與退化檢測。