基于EEMD 和PCA 滾動軸承性能退化指標的提取方法

肖順根， 馬善紅， 宋萌萌， 孔慶光

(1.寧德師范學院 物理與電氣工程系，福建 寧德352100;2. 浙江理工大學 機械與自動控制學院，浙江杭州310018)

滾動軸承是旋轉(zhuǎn)機械中應用最為廣泛的零部件之一，具有承受載荷的作用和相對運動的功能，其運行狀態(tài)直接關系到設備整機性能及壽命［1］。為了避免災難性事故的發(fā)生，有效地實施狀態(tài)基維修，越來越迫切需要了解滾動軸承運轉(zhuǎn)過程中各階段的健康狀態(tài)即運行過程中的衰退狀況［2］。為此，如何提取從開始運行直至失效全過程的滾動軸承性能退化指標尤為重要。目前，表征性能退化指標的特征域有時域、頻域和時頻域3 種。其中，時域和頻域因計算簡單且不涉及參數(shù)選取問題，從而應用較廣泛。徐小力等［3］利用振動烈度等時域指標對大型旋轉(zhuǎn)機械設備的運行狀態(tài)進行跟蹤。Shichang等［4］通過提取時域和頻域共17 個特征域信息作為滾動軸承壽命預測模型的輸入向量，并獲得了較好的剩余壽命預測效果。

近年來，隨著小波分析、經(jīng)驗模態(tài)分解(Empirical Mode Decomposition，EMD)等非線性信號處理方法的出現(xiàn)，時頻域特征指標也被廣泛采用。文獻［5-7］利用小波包變換方法或EMD 技術提取各頻帶能量及其總能量來組成時頻域特征向量，較好地實現(xiàn)滾動軸承的性能退化評估。上述研究成果在某種程度上確實為滾動軸承的性能退化趨勢預測取得了一定效果，但是，這些成果的性能退化指標基本上是基于某一種特征域建立。實際上，滾動軸承的運行狀況是復雜的，僅依靠時域、頻域或時頻域中某一種特征域的退化指標分析其健康狀況往往是不夠全面的。

基于此，文中在文獻［8］的基礎上，通過對軸承全壽命振動信號進行集成經(jīng)驗模態(tài)分解(Ensemble Empirical Mode Decomposition，EEMD)運算，得到一組本征模態(tài)函數(shù)(Intrinsic Mode Function，IMF)分量，利用相關系數(shù)準則對本征模態(tài)函數(shù)分量進行篩選并重構(gòu)。對重構(gòu)信號分別提取時域、頻域和時頻域3 種特征向量，利用主成分分析(Principal Component Analysis，PCA)法加權融合這3 個特征域，從而構(gòu)建能夠全面表征滾動軸承全壽命過程的性能退化指標。

1 EEMD 原理

1.1 EEMD 算法

WU［9］等提出的EEMD 是EMD 的一種改進算法，能夠克服EMD 算法產(chǎn)生非白噪聲干擾所引起的模態(tài)混疊現(xiàn)象。其思想是利用白噪聲的頻譜均勻性，將噪聲數(shù)據(jù)隨機分布到合適的時間點上，有效地抑制了信號奇異點的影響，又利用白噪聲的零均值性能，經(jīng)過多次平均運算后有效地抵消了噪聲的影響，其本質(zhì)是一種疊加高斯白噪聲的多次經(jīng)驗模式分解［10］。

EEMD 算法分解步驟［11］如下:

1)給初始信號x(t)分別加入N 次均值為零，幅值標準差為常數(shù)的高斯白噪聲ωi(t)，構(gòu)建一個總體信號xi(t)，即

式中，ωi(t)的強度取決于高斯白噪聲標準差與初始信號標準差之比Rstd，i = 1 ～N。

2)對各個信號xi(t)分別進行EMD 分解，獲得l 個本征模態(tài)函數(shù)IMF 分量及1 個殘余項ri(t):

式中，yij(t)是第i 次加入高斯白噪聲ωi(t)經(jīng)EMD分解后得到的第j 個IMF 分量。

3)對所有的IMF 分量yij(t)進行均值運算，消除多次加入ωi(t)對真正的IMF 影響，最終得到了IMF 分量yj(t)和殘余項r(t)，如式(3)與式(4)所示。

式中，yj(t)為EEMD 對信號x(t)分解后得到的第j個IMF 分量。因此，被EEMD 分解的信號x(t)可由l個yj(t)和1 個r(t)構(gòu)成，如式(5)所示。

1.2 有效IMF 分量提取

初始信號經(jīng)EEMD 分解后，由于分解誤差和插值誤差等多重因素，使分解結(jié)果中IMF 分量易產(chǎn)生虛假分量。如果虛假分量存在特征域中，顯然會使提取的退化性能指標產(chǎn)生較大誤差。因此，將虛假分量剔除掉是十分必要的?？紤]虛假分量與初始信號的關聯(lián)性較小，因此可通過引入相關系數(shù)來辨別虛假分量。相關系數(shù)的數(shù)學表達式如下:

式中，λj(j = 1，…l)為所有IMF 分量與原始信號的相關系數(shù)，x 和ˉx 分別表示為初始信號及其平均值，y和ˉy 分別表示為IMF 分量信號及其平均值，n 為振動信號采集數(shù)據(jù)樣本點數(shù)。

隨后，設置辨別虛假分量閾值μ 值，若IMF 分量與初始信號的相關系數(shù)小于μ 值則將其剔除，最終將大于μ 值的IMF 分量保留下來并重構(gòu)，從而獲得有效的初始信號，達到降噪目的。

2 提取多特征域信息

采集的全壽命滾動軸承振動信號信息非常龐大，不利于數(shù)據(jù)分析。為此，非常有必要對采集的數(shù)據(jù)進行特征提取處理，以簡化后續(xù)數(shù)據(jù)分析工作量。另外，采集的振動信息僅是一維的時域信息，而軸承的運行狀況包含了軸承的多種形式，多種特征集合，多種狀況征兆信息，若僅用一維的時域信息很難全面表征滾動軸承運行特征，亟需從時域、頻域和時頻域3 種特征域來獲取軸承的全面信息。

2.1 時域特征提取

時域特征指標包括有量綱和無量綱。在反映滾動軸承性能退化趨勢上，有16 種時域特征指標均能較好地表征，其中有量綱包含均值、方差、最大值、最小值、方根幅值、絕對平均幅值、歪度、峭度、均方根和峰峰值10 種，無量綱包含波形指標、峰值指標、脈沖指標、裕度指標、峭度指標和偏斜度指標6 種。關于16 種時域特征指標的數(shù)學表達式詳見相關文獻［8］。

2.2 頻域特征提取

滾動軸承性能退化頻域特征指標主要選擇與頻譜相關的信息，如故障位置的頻率及相應的幅值等。常提取的頻域特征指標反映頻域振動能量大小、頻譜的分散或者集中程度和主頻帶位置的變化3 種特征，各指標的數(shù)學表達式如式(7)至式(19)所示。其中p1(均值頻率)反映頻域振動能量大小，p2(標準差)、p3(特征1)、p4(特征2)、p6(特征3)、p10(特征6)、p11(特征7)、p12(特征8)、p13(特征9)反映頻譜的分散或者集中程度，p5(頻率中心)、p7(均方根頻率)、p8(特征4)、p9(特征5)反映主頻帶位置的變化。

2.3 時頻域特征提取

獲得時域和頻域特征集的確能較好地反映全壽命性能退化狀況，但時域和頻域特征主要表征線性信號特征，對于表征非線性信號呈現(xiàn)能力不足，而軸承運行過程較為復雜，其振動信號具有非線性特征。因此，僅采用時域和頻域兩種特征還不足以全面反映軸承性能退化趨勢。時頻域是表征非線性信號的一種特征指標，若與時域和頻域特征共同構(gòu)成集合則能全面反映軸承全壽命過程的運行狀態(tài)。目前，小波包分析是一種時間窗和頻率窗都可以改變的時頻分析方法，可同時對信號進行低頻和高頻分解，實現(xiàn)非常精細的刻畫，從而應用非常廣泛?；诖?，文中采用小波包分解并提取各頻段能量作為軸承性能退化的時頻域特征集，具體提取方法:

1)對預處理信號進行EEMD 分解、篩選與重構(gòu)，然后選用db3 小波基函數(shù)進行3 層小波包分解得到8 個頻段S3i(3 表示分解層數(shù);i = 0，1，2，…，7，代表8 個頻段位置)。

2)分別計算各頻段的能量值，如式(20)所示:

式中xi，j為第i 個頻段信號的離散點幅值，n 為采樣點數(shù)。

3 基于PCA 加權融合多特征域的性能退化指標

上述分析的時域、頻域和時頻域特征均有諸多指標構(gòu)成，若用單純的某一指標來全面表征滾動軸承的衰退狀況顯然不合適。另外，因軸承全壽命的運行狀態(tài)較為復雜，不同域特征集分布于空間不同的維度中，為此由不同維度的特征來表征軸承的運行特征也不合適，可見尋找眾多域特征信息進行加權融合方法是十分必要的。目前，用于信息特征融合的方法主要有PCA 和流行學習算法，其中流行學習算法計算較為復雜，相比PCA 算法而言其穩(wěn)定性不足。另外，流行學習算法在高維信息融合過程中參數(shù)選擇對于信息融合的效果有很大影響。因此，文中選擇穩(wěn)定性較高的PCA 加權融合滾動軸承時域、頻域和時頻域多特征的性能退化指標集。PCA方法本質(zhì)是通過坐標變換，實現(xiàn)用少數(shù)的主成分描述高維的數(shù)據(jù)，并將數(shù)據(jù)集轉(zhuǎn)化為由維數(shù)較少的“有效”特征成分來表示，并且不減少原始數(shù)據(jù)所包含的信息內(nèi)容，使其在統(tǒng)計意義下達到方差最優(yōu)的目的。PCA 加權融合步驟如下:

1)提取滾動軸承的時域、頻域和時頻域的退化性能指標，將其構(gòu)成數(shù)據(jù)矩陣X:

2)對矩陣X 標準化處理:

式中

(j = 1，2，…，p)

3)計算退化性能指標的相關系數(shù)矩陣R:

式中，假定原始數(shù)據(jù)標準化后仍用X 表示，則經(jīng)標準化處理后數(shù)據(jù)的相關系數(shù)為

4)用雅克比方法求相關系數(shù)矩陣R 的特征值(λ1，λ2，…，λp)和相應的特征向量:

ai= (ai1，ai2，…，aip)，i = 1，2，…，p

5)選擇貢獻率在85% 以上的主成分:

根據(jù)PCA 建立的衰退性能指標不難看出，PCA可以在最小化的數(shù)據(jù)情況下，將多個特征信息的冗余數(shù)據(jù)予以消除，并完整地反映滾動軸承的退化特征。

4 實驗分析

4.1 實驗裝置

研究滾動軸承的性能退化過程，必須要有軸承全壽命周期的振動數(shù)據(jù)，這與傳統(tǒng)的故障診斷研究存在本質(zhì)區(qū)別。文中的實驗數(shù)據(jù)［12］來自Cincinnati大學，其實驗裝置如圖1 所示。裝置中通過電動機帶動軸上的4 個滾動軸承轉(zhuǎn)動，其中每個滾動軸承分別在X 和Y 方向上安裝ICP 加速度振動傳感器。為了實現(xiàn)軸承加速疲勞壽命試驗，在每個滾動軸承上均施加徑向載荷26 670 N，軸承從開始運行至失效共持續(xù)7 d。實驗中，電機轉(zhuǎn)速恒為2 000 r/min，采樣頻率為20 kHz，共獲得984 個文本(每個文本視為1 個數(shù)據(jù)點)，每個文本包含20 480 個采樣數(shù)據(jù)，每個文本之間有10 min 間隔。通過持續(xù)7 d 采集，共采集20 152 320 個數(shù)據(jù)，如圖2(每隔10 個數(shù)據(jù)進行抽 取得到此圖形)所示。

圖1 滾動軸承全壽命數(shù)據(jù)采集裝置Fig.1 Lifetime data acquisition device of the rolling bearing

圖2 軸承運行的初始信號Fig.2 Initial signal of the bearing running

4.2 EEMD 預處理初始信號

采集的信號一般都有隱含的背景噪聲，若不處理必定會對性能退化指標提取產(chǎn)生很大影響。文中采用EEMD 方法對初始信號進行分解、篩選和重構(gòu)，從而達到降噪的目的。具體的實施過程為:對初始信號進行EEMD 計算，其中EEMD 算法中設置Rstd =0.2，添加高斯白噪聲次數(shù)為1 000，分解得到如圖3 所示的15 個IMF 分量和1 個殘余分量，然后利用式(6)所示的相關系數(shù)準則分別計算15 個IMF分量與初始信號的相關程度，各IMF 分量的相關系數(shù)如表1 所示。根據(jù)經(jīng)驗公式［13］，設定μ = 0.1，若IMF 分量相關系數(shù)小于閾值，則表明該IMF 分量與初始信號的相關系數(shù)過小，因此可以將該IMF 分量視為偽分量給予直接剔除。由表1 不難看出，IMF7～IMF15均小于0.1，從而將它們直接剔除，剩余的IMF1～IMF6再進行重構(gòu)后得到新的初始信號。

圖3 EEMD 分解軸承信號情況Fig.3 EEMD decomposed bearing signal

表1 各IMF 分量的相關系數(shù)Tab.1 Correlation coefficient of each IMF component

4.3 建立加權融合多域特征的性能指標

根據(jù)時域指標數(shù)學表達式，分別建立16 個時域特征指標，如圖4 所示。同理，根據(jù)式(7)至式(19)分別建立13 個頻域特征指標，如圖5 所示;根據(jù)式(20)分別建立小波包8 個節(jié)點的能量時頻域特征指標，如圖6 所示。由圖4 ～6 不同域的性能退化指標曲線可以看出，有些指標對于性能退化敏感度很弱，如時域指標的均值、偏斜度和峭度，頻域特征指標8，時頻域指標的S34、S35 和S37 頻帶均是等到軸承快失效時才體現(xiàn);有些指標過于信息嘈雜，如時域指標的峰值指標、脈沖指標、裕度指標、頻域指標特征1 和特征2。對于敏感度太弱和信息過于嘈雜的性能指標都不適宜用于表征滾動軸承的性能退化狀態(tài)，因此應將這些指標給予剔除，而其它保留下來的性能指標則通過PCA 方法進行加權融合，從而得到既代表軸承的退化趨勢，又使維數(shù)得到約簡的性能指標。圖7 為PCA 加權融合各域性能指標的第1 主成分，不難看出，第1 主成分與各域單個的性能指標相比，性能指標得到改善:其一，信息不冗余，非常純凈;其二，提升了軸承退化的敏感度，PCA 加權融合之前，各域單一指標大部分在600 點左右才體現(xiàn)軸承退化趨勢，加權后的性能指標在490 點左右就顯現(xiàn)出軸承衰退趨勢。可見，利用PCA 加權融合各域特征指標使冗余信息更純凈，能夠全面表征軸承性能退化趨勢，為后續(xù)剩余壽命預測研究奠定基礎。

圖4 16 個時域特征指標Fig.4 Feature indexes of sixteen time domain

圖5 13 個頻域特征指標Fig.5 Domain feature indexes of thirteen frequency

圖6 8 個時頻域特征指標Fig.6 Domain feature indexs of eight time -frequency

圖7 PCA 加權融合多域特征的第一主成分Fig.7 First principal component of the PCA weighted fusion multi-domain feature

5 結(jié) 語

1)采用EEMD 分解，結(jié)合相關系數(shù)準則，提取相關程度高的IMF 分量，并重構(gòu)有效IMF 各分量信息，從而實現(xiàn)了信號降噪。

2)建立時域、頻域和時頻域共37 個滾動軸承性能退化指標，剔除不適合表征退化趨勢的指標后，采用PCA 方法加權融合其它指標，得到了更純凈且敏感度更高的退化性能指標。

3)獲得的退化性能指標為滾動軸承的剩余壽命評估奠定了良好的基礎。

［1］艾樹峰.基于雙樹復小波變換的軸承故障診斷研究［J］.中國機械工程，2011，22(20):2446-2451.AI Shunfeng. Research on bearing fault diagnosis based on dual-tree complex wavelet transform［J］. China Mechanical Engineering，2011，22(20):2446-2451.(in Chinese)

［2］徐東.球軸承疲勞剩余壽命分析與預測方法研究［D］.長沙:國防科技大學，2011.

［3］徐小力，王紅軍.大型旋轉(zhuǎn)機械運行狀態(tài)趨勢預測［M］.北京:科學技術出版社，2011.

［4］Shichang D，Jun L，Lifeng X.Degradation process prediction for rotational machinery based on hybrid intelligent model［J］.Robot Comput Integr Manuf，2012，28(2):190-207.

［5］肖文斌，陳進，周宇，等.小波包變換和隱馬爾可夫模型在軸承性能退化評估中的應用［J］. 振動與沖擊，2011，30(8):32-35.

XIAO Wenbin，CHEN Jin，ZHOU Yu，et al. Wavelet packet transform and hidden Markov model based bearing performance degradation assessment［J］.Journal of Vibration and Shock，2011，30(8):32-35.(in Chinese)

［6］李巍華，李靜，張紹輝.連續(xù)隱半馬爾科夫模型在軸承性能退化評估中的應用［J］.振動工程學報，2014，27(4):613-620.

LI Weihua，LI Jing，ZHANG Shaohui. Application of continuous hidden semi-Markov model in bearing performance degradation assessment［J］.Journal of Vibration Engineering，2014，27(4):613-620.(in Chinese)

［7］王立，于重重，施彥，等.基于振動信號分析的電機性能退化特征提?。跩］.計算機仿真，2014，31(4):416-420.

WANG Li，YU Chongchong ，SHI Yan，et al.Degradation characteristics extraction of electric motor properties based on vibration signal analysis［J］.Computer Simulation，2014，31(4):416-420.(in Chinese)

［8］董紹江.基于優(yōu)化支持向量機的空間滾動軸承壽命預測方法研究［D］.重慶:重慶大學，2012.

［9］WU Z H，HUANG N E.Ensemble empirical mode decomposition:a noise-assisted data analysis method［J］.Advances in Adaptive Data Analysis，2009，1(1):1 -41.

［10］余發(fā)軍，周鳳星.基于EEMD 和自相關函數(shù)特性的自適應降噪方法［J］.計算機應用研究，2015，32(1):206-209.

YU Fajun，ZHOU Fengxing.Adaptive de-noising method based on EEMD and autocorrelation function property［J］. Application Research of Computers，2015，32(1):206-209.(in Chinese)

［11］陳仁祥，湯寶平，呂中亮.基于相關系數(shù)的EEMD 轉(zhuǎn)子振動信號降噪方法［J］.振動、測試與診斷，2012，32(4):542-546.

CHEN Renxiang，TANG Baoping，LU Zhongliang.Ensemble empirical mode decomposition de-noising method based on correlation coefficients for vibration signal of rotor system［J］.Journal of Vibration，Measurement and Diagnosis，2012，32(4):542-546.(in Chinese)

［12］Lee J，QIU H，YU G.Rexnord technical services“bearing data set”［R］.Moffett Field，CA:IMS，University of Cincinnati，1990.

［13］邊杰，王平，梅慶. EEMD 結(jié)合能量特征和小波降噪的軸承故障診斷［J］. 廣西大學學報:自然科學版，2014，39(6):1206-1211.

BIAN Jie，WANG Ping，MEI Qing. Fault diagnosis of bearings by using EEMD combined with energy feature and wavelet denoising［J］.Journal of Guangxi University:Nat Sci Ed，2014，39(6):1206-1211.(in Chinese)