丁嘉鑫 王振亞 姚立綱 蔡永武

福州大學機械工程及自動化學院，福州，350116

0 引言

滾動軸承作為旋轉(zhuǎn)機械應用廣泛且容易損壞的零部件，對其進行故障診斷有著重要的理論和實際意義[1]。由于受負載、間隙、剛度、摩擦和沖擊等非線性因素的影響，滾動軸承振動信號通常表征為非平穩(wěn)和非線性等特性[2]。所以，眾多衡量機械動力學系統(tǒng)的非線性時間序列復雜性方法相繼被提出，并被應用于故障診斷領域，如近似熵、樣本熵、排列熵和模糊熵等[3-5]。其中，排列熵(permutation entropy, PE)無需考慮時間序列的數(shù)值大小，而是對相鄰樣本點進行對比分析，獲取相應特征信息，相較其他熵值方法更能捕獲序列的微弱變化，并且該算法具有理論簡單、抗噪能力強等優(yōu)勢，故在故障診斷領域應用較為廣泛[4-5]。但PE算法僅利用時間序列的序數(shù)結(jié)構(gòu)，忽視其幅值信息，因此，XIA等[6]在PE的基礎上提出了加權(quán)排列熵(weighted permutation entropy, WPE)；ZHOU等[7]將其應用于滾動軸承故障特征提取過程。與PE類似，WPE僅考慮單一尺度上時間序列的復雜性和動力學突變，忽視了其他尺度上的有用信息，而YIN等[8]將WPE與多尺度熵相結(jié)合，提出了多尺度加權(quán)排列熵(multiscale WPE, MWPE)。但將MWPE應用于滾動軸承特征提取過程，仍存在以下不足：①MWPE的熵值估計偏差會隨粗?；叨纫蜃拥脑龃蠖龃螅虎贛WPE粗?；^程忽略了其他粗?；蛄猩系挠杏眯畔?，從而影響熵值準確度。③MWPE進行粗粒化構(gòu)造時，利用均值處理方式在一定程度上會中和原始信號的動力學突變行為，影響特征提取結(jié)果。針對以上不足，需要研發(fā)一種廣義復合多尺度加權(quán)排列熵(generalized composite multiscale weighted permutation entropy, GCMWPE)新算法，通過采用廣義復合粗?；瘶?gòu)造方式，以此克服MWPE算法存在的不足，并將上述方法應用于滾動軸承故障特征提取過程。

滾動軸承故障診斷的本質(zhì)在于模式識別，支持向量機(support vector machine, SVM)在故障診斷領域應用廣泛，并且取得較好的識別結(jié)果[9]，但該算法性能易受懲罰因子和核函數(shù)參數(shù)的影響。有學者將粒子群優(yōu)化算法[10]、模擬退火算法[11]及人工魚群算法[12]應用于SVM參數(shù)尋優(yōu)過程，但上述尋優(yōu)方法易陷入局部最優(yōu)解并且尋優(yōu)耗時。本文采用一種新穎元啟發(fā)式優(yōu)化算法——天牛須搜索(BAS)優(yōu)化[13]對SVM參數(shù)進行尋優(yōu)，在此基礎上，提出了天牛須搜索優(yōu)化支持向量機(beetle antennae search based support vector machine, BAS-SVM)?；谏鲜隼碚?，建立一種GCMWPE、監(jiān)督等度規(guī)映射(supervised isometric mapping, S-Isomap)[14]和BAS-SVM相結(jié)合的滾動軸承故障診斷模型，將其應用于實驗數(shù)據(jù)分析過程。

1 廣義復合多尺度加權(quán)排列熵

1.1 加權(quán)排列熵算法

排列熵(PE)只考慮序列的順序結(jié)構(gòu)特征，忽略了幅度特性。XIA等[6]在PE的基礎上提出了加權(quán)排列熵(WPE)，具體算法流程如下。

(1)

式中，τ為時延；m為嵌入維數(shù)。

(2)計算出每個子序列的權(quán)重值wi：

(2)

(3)

(4)計算時間序列X的加權(quán)排列熵WPE值：

(4)

1.2 多尺度加權(quán)排列熵算法

多尺度加權(quán)排列熵(MWPE)克服了WPE單一尺度分析的不足，能夠從多個尺度全面表征時間序列復雜性，具體過程如下：

(1)對時間序列X進行粗?；幚?，得到粗粒序列y(s)=(y(s)(j))：

(5)

式中，s為尺度因子。

(2)計算不同尺度因子下粗?；蛄衴(s)的WPE值：

MWPE(X,m,τ,s)=WPE(y(s),m,τ)

(6)

圖1 MWPE粗粒化構(gòu)造方式(s=2,3)

1.3 廣義復合多尺度加權(quán)排列熵算法

圖2 GCMWPE算法流程

針對MWPE存在的不足，本文作出以下改進：①采用復合粗?；瘶?gòu)造方式，考慮同一尺度下多個粗?；瘯r間序列的加權(quán)排列熵值，以此抑制由粗?；瘯r間序列變短而導致的熵值突變，得到更為精準的熵值特征。②將粗?；^程的均值計算變更為方差計算，避免中和原始信號的動力學突變現(xiàn)象的發(fā)生。上述改進即為廣義復合多尺度加權(quán)排列熵(GCMWPE)，流程如圖2所示，具體步驟如下。

(7)

j,k,s∈Z+

(3)均化同一尺度下多個WPE值，即可得到s尺度下的GCMWPE值，對應表達式如下：

(8)

GCMWPE算法需人為設定以下4個參數(shù)：樣本長度N，嵌入維數(shù)m，尺度因子s和時延τ。其中，s的取值目前尚無選定標準，通常設置為s>10，本文設定s=20。時間序列長度應滿足N>200 s，故本文選取N=4096。m對GCMWPE具有一定影響，若m過小，相空間重構(gòu)的向量包含較少信息，算法將無法有效監(jiān)測序列的動力學突變；反之，若m過大，相空間重構(gòu)的向量會忽略序列的細微變化，并且增加運算時間。通常設定m的取值范圍為[4, 7]。時延τ對熵值計算影響較小，一般設定為τ=1[7]。

1.4 實驗驗證

為探究嵌入維數(shù)m對所提GCMWPE算法的影響，將其應用于凱斯西儲大學軸承數(shù)據(jù)實驗分析[15]。實驗軸承為6205-2RS深溝球軸承，利用電火花技術(shù)在軸承表面進行單點故障加工，其中，故障直徑為0.3356 mm，故障深度為0.2794 mm。本次實驗設置電機轉(zhuǎn)速為1797 r/min，負載為0，在采樣頻率為12 kHz的條件下，分別采集軸承正常和具有外圈故障振動信號各20組，每組信號包含4096個采樣點。

將不同嵌入維數(shù)m(分別取4, 5, 6, 7)下的GCMWPE應用于兩種狀態(tài)軸承數(shù)據(jù)分析過程，兩種狀態(tài)的熵值均值曲線以及對應的線性擬合分別如圖3a和圖3b所示。由圖3可知：①當嵌入維數(shù)m較小(為4或5)時，GCMWPE熵值擬合線較為平緩，無法體現(xiàn)多尺度的優(yōu)勢；而當嵌入維數(shù)m較大(為7)時，相空間重構(gòu)的向量會忽略序列的細微變化，導致兩種狀態(tài)熵值曲線較為接近，無法有效區(qū)分故障類型(兩種狀態(tài)下軸承的GCMWPE熵值均值曲線對比如圖4所示)。因此，本文設定m=6。②同一m下，第一個尺度上，正常狀態(tài)的熵值高于外圈故障狀態(tài)。原因在于，當軸承處于正常狀態(tài)時，振動信號波動較為隨機，信號的無規(guī)則性較高，自相似性較低，故熵值較大；而當軸承出現(xiàn)局部故障時，振動信號波動出現(xiàn)一定規(guī)律性，信號的規(guī)則性和自相似性較高，故熵值較小,因此，GCMWPE能夠監(jiān)測軸承是否發(fā)生故障。

(a) 正常信號分析結(jié)果

圖4 GCMWPE對兩種信號分析結(jié)果

2 BAS算法和BAS-SVM

2.1 BAS算法

BAS算法通過模擬自然界天牛根據(jù)觸角接受的氣味強度進行搜索食物，從而獲取全局最優(yōu)解，具體步驟如下。

(1)初始化參數(shù)。初始兩須距離p(0)，初始步長δ(0)，迭代次數(shù)T，初始位置為u(0)。

(2)計算兩須坐標。左右須坐標分別如下：

(9)

p(t)=0.95p(t-1)+0.01

(10)

(11)

其中，sgn(·)為符號函數(shù)，δ(t+1)為t+1次迭代時的移動步長，其表達式為

δ(t+1)=0.95δ(t)

(12)

(4)判斷是否達到最大迭代次數(shù)，若滿足，則計算終止；否則，繼續(xù)循環(huán)。

2.2 BAS-SVM

針對支持向量機(SVM)算法性能易受懲罰因子c和核函數(shù)參數(shù)g影響的問題，本文研發(fā)一種天牛須搜索優(yōu)化支持向量機(BAS-SVM)新算法，流程如圖5所示，具體過程如下。

圖5 BAS-SVM流程

(1)輸入訓練集與測試集樣本，并對兩個樣本集分別進行歸一化處理。初始化BAS和SVM參數(shù)，包括初始兩須距離p(0)，初始步長δ(0)，迭代次數(shù)T，初始位置為u(0)；參數(shù)c和g的取值范圍為[0.001, 100]，選用徑向基函數(shù)。

(2)利用式(9)計算出天牛左右兩須的坐標。

(3)計算出天牛左右兩須的氣味強度。其中，以訓練樣本三折交叉驗證后的平均識別率為其適應度函數(shù)值。

(4)利用變步長法(即式(11))確定天牛下一步的位置。

(5)判斷是否滿足最大迭代次數(shù)，若滿足，則計算終止；否則，繼續(xù)循環(huán)。

(6)輸出天牛搜索食物的位置(即為最優(yōu)c和g)，并將最優(yōu)參數(shù)建立SVM預測模型，并利用該模型對測試集樣本進行測試。

3 滾動軸承故障診斷實驗分析

3.1 故障診斷模型

本文提出一種廣義復合多尺度加權(quán)排列熵(GCMWPE)、監(jiān)督等度規(guī)映射(S-Isomap)和BAS-SVM相結(jié)合的滾動軸承故障診斷方法，診斷流程如圖6所示，具體步驟如下。

圖6 滾動軸承故障診斷流程

(1)信號采集。在一定采樣頻率fs下，利用加速度傳感器分別采集滾動軸承不同狀態(tài)下的振動加速度信號，并將其分為測試和訓練樣本集。

(2)計算高維故障特征。利用GCMWPE算法對訓練樣本與測試樣本信號進行熵值特征提取，合并成高維故障特征集。

(3)降維處理。由于GCMWPE高維特征集存在冗余特征，影響最終識別效果，故利用S-Isomap流行學習算法對該特征集進行維數(shù)約簡，獲取低維、易于區(qū)分故障類型的敏感特征集。

(4)故障識別。利用BAS-SVM分類器對GCMWPE+S-Isomap低維故障特征集進行訓練與測試，診斷出測試樣本各故障類型。

3.2 故障診斷實例

為模擬滾動軸承實際工作條件，利用Spectra Quest公司研發(fā)的動力傳動系統(tǒng)故障診斷實驗臺進行實驗數(shù)據(jù)的采集，實驗平臺如圖7所示。本次實驗中，輸入軸轉(zhuǎn)速為20 Hz，負載為0，在采樣頻率3000 Hz下，利用加速度傳感器分別采集滾動軸承4種狀態(tài)振動加速度信號各100組，即正常(normal, NOR)、外圈故障(outer race fault, ORF)、內(nèi)圈故障(inner race fault, IRF)和滾動體故障(ball fault, BF)。每組信號包含4096個采樣點，4種狀態(tài)共計400組樣本信號，對應時域波形如圖8所示。其中，每種狀態(tài)隨機選取20組樣本作為訓練樣本，剩余80組樣本作為測試樣本，4種狀態(tài)共計80組訓練樣本、320組測試樣本。

圖7 滾動軸承故障診斷實驗平臺

圖8 時域波形

3.2.1特征提取

首先，利用GCMWPE對振動信號進行熵值特征提取，構(gòu)建原始高維特征集。為驗證所提方法有效性，將其與GMWPE(僅采用廣義粗?；瘶?gòu)造方式，未采用復合粗?；瘶?gòu)造方式，即式(7)中k=1，其余與GCMWPE相同)、MWPE和復合多尺度加權(quán)排列熵(composite multiscale weighted permutation entropy, CMWPE)[16]三種算法進行比較。4種算法在滾動軸承不同狀態(tài)下的熵值均值曲線見圖9。其中，設置對比方法的尺度因子s=20，時延τ=1，嵌入維數(shù)m=6。

(a) MWPE和CMWPE均值曲線

由圖9可知：①與MWPE和CMWPE相比，本文所提的GMWPE和GCMWPE方法提取得熵值均值曲線較為平滑，并且能夠有效將4類樣本區(qū)分開，體現(xiàn)出廣義粗粒化構(gòu)造方式的優(yōu)越性。②就起始尺度而言，MWPE和CMWPE得到的滾動軸承4種狀態(tài)中滾動體狀態(tài)熵值最大，而GMWPE和GCMWPE得到的滾動軸承4種狀態(tài)中正常狀態(tài)熵值最大。就實際工況而言，當滾動軸承處于正常狀態(tài)時，振動信號波動較為隨機，信號的無規(guī)則性較高，自相似性較低，故熵值較大；而當軸承出現(xiàn)局部故障時，振動信號波動出現(xiàn)一定規(guī)律性，信號的規(guī)則性和自相似性較高，故熵值較小，故相比MWPE和CMWPE，GMWPE和GCMWPE算法更適用于滾動軸承的特征提取過程。③GMWPE和GCMWPE算法對滾動軸承每種狀態(tài)分析得到的熵值均值曲線較為接近，但隨機選取1組信號進行分析時(圖10)，GMWPE提取的熵值曲線較GCMWPE波動較大，表明采用復合粗?；?即考慮同一尺度多個粗粒化序列)的GCMWPE提取的熵值較采用原始粗?；?即僅考慮同一尺度下一個粗?；蛄?的GMWPE算法更為穩(wěn)定；同理，CMWPE和MWPE對隨機選取的熵值分析結(jié)果(圖10b)也能得到相同的結(jié)論。上述分析驗證了將所提GCMWPE應用于滾動軸承故障特征提取的可行性和優(yōu)越性。

(a) 隨機1組信號的GMWPE和GCMWPE熵值曲線

3.2.2維數(shù)約簡

GCMWPE提取的故障特征集具有高維、非線性和冗余等特性，直接利用BAS-SVM進行故障識別，會增加識別時間甚至影響識別效果，因此，本文采用S-Isomap算法對其進行降維處理，提取出易于區(qū)分故障類型的低維特征集。同時，為了進行對比，利用該流形學習算法對MWPE、CMWPE和GMWPE進行維數(shù)約簡。S-Isomap算法對4種特征集的降維結(jié)果如圖11所示。其中，通過交叉驗證的方式確定各參數(shù)(表1)。

表1 S-Isomap參數(shù)設置

(a) S-Isomap對MWPE降維結(jié)果

由圖11可知：①MWPE和CMWPE降維結(jié)果中，內(nèi)圈故障和外圈故障狀態(tài)的樣本出現(xiàn)了混疊現(xiàn)象，并且4類樣本較為分散；而GMWPE和GCMWPE的降維結(jié)果中，可以將4類樣本基本分離開，并且每類樣本擁有較好的聚類效果，以此進一步驗證廣義粗?；瘶?gòu)造方式的優(yōu)勢。②與GMWPE相比，GCMWPE的降維結(jié)果中可以將4類樣本完全分離開，未出現(xiàn)GMWPE中個別樣本遠離聚類中心現(xiàn)象，表明采用復合粗?；瘶?gòu)造方式的GCMWPE算法能夠提取出更為穩(wěn)定的熵值特征信息。③S-Isomap對GCMWPE特征集的降維結(jié)果中，能夠有效將4類樣本完全區(qū)分開，4類樣本聚集性較好，降維效果最佳，這說明本文所提GCMWPE與S-Isomap相結(jié)合的特征提取方式能夠有效提取出易于區(qū)分滾動軸承故障特征信息的低維、敏感特征集。

3.2.3模式識別

為量化上述4種降維后特征提取效果，分別將其輸入BAS-SVM分類器中進行診斷識別，識別結(jié)果和混淆矩陣如圖12所示。其中，設置BAS-SVM算法中初始兩須距離p(0)= 2，初始步長δ(0)=4，迭代次數(shù)T=100。

由圖12可知：①BAS-SVM分類器對GCMWPE(CMWPE)+S-Isomap的平均識別率比GMWPE(MWPE)+S-Isomap高出1.25%(0.625%)，驗證了復合粗?；膬?yōu)越性。②BAS-SVM分類器對GCMWPE(GMWPE)+S-Isomap的平均識別率比CMWPE(MWPE)+S-Isomap高出2.5%(1.875%)，體現(xiàn)了廣義粗?；膬?yōu)勢。③本文所提GCMWPE+S-Isomap特征集的故障識別率達到100%，能有效、精準地識別出滾動軸承故障類型。

為驗證所提BAS-SVM分類器的優(yōu)越性和高效性，將其與粒子群優(yōu)化支持向量機(particle swarm optimization SVM, PSO-SVM)、模擬退火優(yōu)化支持向量機(simulated annealing SVM, SA-SVM)和人工魚群優(yōu)化算法支持向量機(artificial fish swarm algorithm SVM, AFSA-SVM)三種分類器進行對比，4種分類器(參數(shù)設置見表2)對上述4種降維后特征集的識別結(jié)果如圖13所示，對應識別時間見表3。其中，本次實驗各算法編譯于MATLAB 2017軟件，運行環(huán)境如下：Intel (R) Core (TM) i7-8550 CPU, 16 GB RAM, Windows 10系統(tǒng)。

表2 3540Fe基粉末化學成分

圖13 滾動軸承故障診斷模型

表3 平均識別時間

由圖13和表3可知：①就每種特征集而言，BAS-SVM均具有較好的故障識別效果，以此驗證BAS-SVM分類器在模式識別方面的優(yōu)越性。②PSO-SVM、SA-SVM和AFSA-SVM分類器對4種特征集的平均識別時間分別是BAS-SVM分類器用時的8.039倍、16.670倍和19.655倍，體現(xiàn)了BAS-SVM分類器的高效性。③所提GCMWPE+S-Isomap方法提取的故障特征在每種分類器下的故障識別率均為100%，進一步驗證所提故障診斷方法的優(yōu)越性。

4 結(jié)論

(1)提出一種衡量時間序列復雜性新算法——GCMWPE，該算法克服了MWPE粗?；^程的不足。故障診斷實驗結(jié)果表明，所提方法優(yōu)于MWPE、CMWPE和GMWPE。

(2)研發(fā)一種BAS-SVM模式識別新算法，滾動軸承實驗分析結(jié)果表明，BAS-SVM分類器的識別效果和識別速度優(yōu)于現(xiàn)有的PSO-SVM、SA-SVM和AFSA-SVM。

(3)建立一種GCMWPE、S-Isomap和BAS-SVM相結(jié)合的滾動軸承故障診斷新方法。實驗分析結(jié)果表明，該方法能夠有效提取出軸承的各種故障特征信息，識別率達到100%。