基于LMD-CM-PCA的滾動軸承故障診斷方法*

付云驍， 賈利民， 秦 勇， 楊 杰

(1.北京交通大學(xué)軌道交通控制與安全國家重點(diǎn)實(shí)驗(yàn)室 北京，100044) (2.北京交通大學(xué)電氣工程學(xué)院 北京，100044) (3.北京交通大學(xué)北京市城市交通信息智能感知與服務(wù)工程技術(shù)研究中心 北京，100044)

?

基于LMD-CM-PCA的滾動軸承故障診斷方法*

付云驍1,2， 賈利民1,3， 秦 勇1,3， 楊 杰1

(1.北京交通大學(xué)軌道交通控制與安全國家重點(diǎn)實(shí)驗(yàn)室 北京，100044) (2.北京交通大學(xué)電氣工程學(xué)院 北京，100044) (3.北京交通大學(xué)北京市城市交通信息智能感知與服務(wù)工程技術(shù)研究中心 北京，100044)

為提高在非平穩(wěn)工況下對滾動軸承故障的直觀辨識能力，筆者提出基于LMD-CM-PCA的故障診斷方法。首先，對滾動軸承振動信號進(jìn)行局部均值分解(local mean decomposition，簡稱LMD)，提取乘積函數(shù)(product function，簡稱PF)矩陣；然后，計(jì)算PF矩陣與原振動信號的皮氏相關(guān)系數(shù)(pearson product-moment correlation coefficient，簡稱PPCC)，將PFs對應(yīng)的PPCC代入相關(guān)熵模型得到PF的相關(guān)熵矩陣(correntropy matrix，簡稱CM)，CM經(jīng)主成分分析(principal component analysis，簡稱PCA)進(jìn)行特征變換得到融合相關(guān)熵矩陣(integrated correntropy matrix,簡稱ICM)。分別在輕微和嚴(yán)重故障時(shí)，對滾動軸承不同工況下的振動樣本進(jìn)行交叉混合，并計(jì)算其ICM。結(jié)果證明，ICM在可視維度比傳統(tǒng)特征(如：能量矩和譜峭度)的融合特征更能隔離工況對故障可分性的干擾。LMD-CM-PCA方法為滾動軸承故障的直觀辨識提供了技術(shù)支持，在故障診斷方面具有良好的應(yīng)用前景。

局部均值分解； 融合相關(guān)熵矩陣； 主成分分析； 滾動軸承； 故障診斷； 可視化

引 言

滾動軸承是大型旋轉(zhuǎn)機(jī)械、關(guān)鍵電氣設(shè)備的力矩轉(zhuǎn)換部件，其健康狀態(tài)對設(shè)備安全具有重大影響。非平穩(wěn)運(yùn)行工況對滾動軸承故障診斷的置信度干擾明顯，僅僅依靠改進(jìn)模式識別方法不能從根本上解決此類問題，而優(yōu)秀的故障特征，才是提高診斷率的根本。若特征分布不能直觀可視，診斷結(jié)果也就不能簡單易懂，實(shí)時(shí)監(jiān)測的效率也隨之降低。但降低特征維度則影響特征的故障可分性。

振動信息采集法是目前旋轉(zhuǎn)機(jī)械故障診斷領(lǐng)域最常用的物理參量采集方法[1]?；谡駝有畔⒌慕?jīng)驗(yàn)?zāi)B(tài)分解(empirical mode decomposition,簡稱EMD)方法，已廣泛應(yīng)用在機(jī)電故障診斷領(lǐng)域[2-3]。然而EMD得到的非線性分量之間易產(chǎn)生明顯的模態(tài)混疊及端點(diǎn)效應(yīng)，使分量信號失真。局部均值分解(local mean decomposition,簡稱LMD)將信號分解為單分量的調(diào)頻-調(diào)幅信號，即乘積函數(shù)PF，能夠緩解經(jīng)驗(yàn)?zāi)B(tài)分解的欠包絡(luò)、過包絡(luò)和端點(diǎn)效應(yīng)等不足，且在滾動軸承振動分析中具有良好的應(yīng)用效果[4]。

基于時(shí)域、頻域以及時(shí)頻域的故障特征已經(jīng)被用于滾動軸承恒工況的故障診斷[5]，時(shí)頻域特征[6-7]雖然包含更全面的時(shí)頻故障信息，但也增加了冗余量，弱化了故障可辨性。近幾年特征提取技術(shù)開始采用信息理論作為指導(dǎo)[8]，此后發(fā)展的相關(guān)熵[9]屬于廣義相關(guān)性函數(shù)，可作為隨機(jī)變量間局部相似性的度量。雖然相關(guān)熵已用于實(shí)現(xiàn)時(shí)間序列建模、非線性分析、相關(guān)濾波等信號處理算法中[10]，然而尚未在故障診斷領(lǐng)域得到推廣[11]。此外，多維特征中的冗余信息，需用降維變換解決上述問題。目前數(shù)據(jù)降維已應(yīng)用于不同研究領(lǐng)域[12]，然而針對原始維度不高且計(jì)算求簡的情形，最有效的方法就是主成分分析PCA。PCA將原特征矩陣進(jìn)行線性組合并最大化樣本方差，得到相互獨(dú)立的新特征矩陣。依照方差大小順序可以解析出對應(yīng)的主要特征。

1 特征提取方法

提取融合相關(guān)熵矩陣大致分為3個(gè)步驟。首先對信號進(jìn)行LMD，得到多個(gè)PF；然后計(jì)算PF與原始信號的相關(guān)系數(shù)，將相關(guān)系數(shù)帶入相關(guān)熵公式得到CM；最后經(jīng)過PCA變換后，以方差貢獻(xiàn)率為判據(jù)，篩分主要向量構(gòu)成融合特征矩陣，即ICM。

1.1 局部均值分解

LMD是自適應(yīng)時(shí)頻分析方法，依據(jù)信號自身特點(diǎn)確定信號在不同頻域尺度的分辨率，優(yōu)化信號分析。對非線性和非平穩(wěn)的旋轉(zhuǎn)機(jī)械故障信號具有很好的泛化分析能力。LMD的本質(zhì)是自適應(yīng)地將非線性、非平穩(wěn)信號按頻率遞減順序分離出調(diào)頻信號和包絡(luò)信號，兩者相乘得到一個(gè)瞬時(shí)頻率具有物理意義的PF向量，循環(huán)迭代得到全部PF向量和一個(gè)單調(diào)的調(diào)頻分量。若令ci(t)表示第i個(gè)PF分量，則PF分量矩陣(PF matrix,簡稱PFM)可由C=(c1(t),c2(t),…,cn(t))′表示，若每個(gè)c均是m維行向量，則第i個(gè)PF向量有ci(t)=(ci1,ci2,…,cim)。令單調(diào)分量為un(t)，其結(jié)構(gòu)模型見式(1)。LMD的計(jì)算過程見文獻(xiàn)[4]。

(1)

在k個(gè)狀態(tài)下分別采集l個(gè)s(t)，則得到三維k×l×m的樣本矩陣(signal matrix,簡稱SM)S。對S進(jìn)行LMD后，得到四維k×l×n×m的PFM。

1.2 相關(guān)熵矩陣

相關(guān)熵(correntropy)是通過核函數(shù)將數(shù)據(jù)從非線性空間投影到高維的線性空間，得到高維空間數(shù)據(jù)的線性相關(guān)信息。以核寬度為σ的Gauss核為底，設(shè)m維向量X與Y的聯(lián)合概率為PXY，可得二維離散變量的相關(guān)熵估計(jì)為

(2)

皮氏相關(guān)系數(shù)PPCC作為信號結(jié)構(gòu)線性相關(guān)測度，可度量兩個(gè)變量的相關(guān)程度，其值不會隨變量尺度及順序的變化而改變，故用PPCC作為相關(guān)熵置信系數(shù)，可度量更加完備的相關(guān)信息。若ci(t)與s(t)的PPCC用γ表示，將γ帶入式(2)，則m維c(t)和s(t)的相關(guān)熵為

(3)

于是，s(t)的相關(guān)熵計(jì)算如下

(4)

由于CM數(shù)值較大，為方便計(jì)算，利用式(5)對CM進(jìn)行歸一化處理

(5)

將式(4)和(5)推廣到PFM，令Δ=k×l，Ω=Δ×n，則相關(guān)熵矩陣(CM)可由式(6)求出

(6)

其中：CΔ，SΔ表示二維k×l的樣本矩陣。

為了便于進(jìn)一步計(jì)算的表述，需要對φΩ部分轉(zhuǎn)置，令Γ=(k×l)′×n，則φΓ為所求CM。

1.3 仿真測試

以式(7)所示仿真信號為例證明相關(guān)熵的優(yōu)勢。故障軸承的振動信號具有短時(shí)高頻沖擊的特性，這里以信號0.5cos(60πt1)為正常仿真信號，添加信噪比為1.4 dB的白噪聲wn(t)得到含噪原始信號x0。將短時(shí)高頻信號添加至x0得到合成故障仿真信號x1，時(shí)域分布見圖1。

(7)

圖1 仿真信號時(shí)域分布Fig.1 Simulated signal time domain distribution

圖2為故障仿真信號x1的頻譜，f0為固有頻率，f1和f2為微弱沖擊倍頻。再對x0和x1分別進(jìn)行LMD，并計(jì)算PF分量的相關(guān)熵和能量矩，結(jié)果如圖3所示。從圖中可知，相關(guān)熵分布相較于能量矩分布明顯不同，因此可知相關(guān)熵比能量矩在應(yīng)對局部沖擊時(shí)更敏感。

圖2 仿真信號頻域分布Fig.2 Simulated signal frequency domain distribution

圖3 仿真信號的相關(guān)熵與能量矩對比Fig.3 Correntropy and energy moment comparison of simulated signals

1.4 融合特征

為了消除相關(guān)熵矩陣的冗余分量，同時(shí)實(shí)現(xiàn)故障診斷的可視化，需要對φΓ進(jìn)行降維融合處理。采用PCA通過減小n實(shí)現(xiàn)降維。PCA能保持樣本矩陣的總方差不變，使變換后的特征按統(tǒng)計(jì)方差大小降序排列，且每個(gè)樣本之間均保持相互獨(dú)立。對相關(guān)熵矩陣而言，令CM降維后的融合矩陣為

(8)

其中：φl×k×p為滾動軸承故障辨識的特征全集，p∈{1,2,3}表示可視化的特征維度。

用距離測度構(gòu)造特征可分性判據(jù)，首先計(jì)算類內(nèi)及類間散度

(9)

則可由SW和SB構(gòu)造如下可分性判據(jù)

(10)

其中：Tr為矩陣的跡；J1，J2∈(0,1)，J1越大，J2越小則可分性越好，反之可分性越差；J3∈(-∞,+∞)越大則可分性越好，越小可分性越差。

J1，J2和J3在任何非奇異變換下保持不變。

2 案例分析實(shí)驗(yàn)

2.1 實(shí)驗(yàn)流程

工況交叉的滾動軸承故障辨識目的是檢驗(yàn)ICM故障辨識能力的同時(shí)，檢驗(yàn)其對工況的抗擾性和對故障程度的魯棒性。

在多組恒定轉(zhuǎn)速、負(fù)載工況下進(jìn)行滾動軸承實(shí)驗(yàn)驗(yàn)證融合相關(guān)熵矩陣的可行性。本研究采用美國西儲大學(xué)提供的滾動軸承實(shí)驗(yàn)數(shù)據(jù)，軸承為6205-2RS JEM SKF型深溝球軸承，這里按故障程度分兩組實(shí)驗(yàn)進(jìn)行驗(yàn)證，采用放電加工技術(shù)得到故障深度為0.177 8和0.533 4 mm的故障情形(以兩種情形分別代表故障的輕微情形和嚴(yán)重情形)，驗(yàn)證算法的有效性和魯棒性。每組工況下狀態(tài)辨識實(shí)驗(yàn)流程如圖4所示。圖中l(wèi)，m和n的含義與前文一致。

圖4 故障特征提取流程圖Fig.4 Flow diagram of fault feature extraction

1) 整合樣本矩陣。電機(jī)運(yùn)行參數(shù)為恒速穩(wěn)定載荷，分別采集轉(zhuǎn)速為1 730，1 750，1 772和1 797 r/min，且各自載荷為3，2，1，0 hp的軸承振動加速度信號，并以A組、B組、C組和D組命名4組工況參數(shù)。采樣頻率為12 kHz，每個(gè)樣本長度為1 024。每組工況各有472個(gè)信號樣本，將A，B，C，D的樣本按照組序依次排列，所有樣本輕微故障和嚴(yán)重故障下各自的全樣本矩陣維數(shù)分別為k×l×m=4×1 888×1 024。

2) 局部均值分解。對SM1和SM2進(jìn)行LMD，每個(gè)樣本提取前4個(gè)PF分量，則得到k×l×n×m=4×1 888×4×1 024的PFM1和PFM2。圖5所示為A組內(nèi)圈故障信號及內(nèi)圈輕微故障和嚴(yán)重故障的PF分量圖。

3) 特征提取。按照圖4所示的流程計(jì)算k×

l×n=4×1 888×4的CM1和CM2，將4類故障數(shù)據(jù)整合為一個(gè)樣本矩陣，輸入PCA。通過PCA的空間變換得到k×l×n=4×1 888×2的ICM1和ICM2。將常用的時(shí)頻域特征能量矩[2]及頻域特征譜峭度[13]進(jìn)行融合處理，得到融合能量矩(integrated energy moment,簡稱IEM)和融合譜峭度(integrated spectral kurtosis,簡稱ISK)。依據(jù)圖4所示流程，將相關(guān)熵算子改為能量矩算子及譜峭度算子提取輕微故障及嚴(yán)重故障時(shí)的IEM及ISK。

圖5 A組工況下原信號及內(nèi)圈故障樣本PF分量Fig.5 Raw signal and PF components of inner race fault sample in A group

2.2 結(jié)果分析

圖6為在不同工況下的輕微故障的特征樣本直方圖。x軸為不同工況下的融合特征維度；y軸為故障類型；z為特征幅值，其中，{I，II，III，IV}代表{正常，內(nèi)圈故障，外圈故障，滾子故障}。從x軸方向比較，在同故障類型下的ISK幅值比IEM和ICM具有明顯變化，說明IEM和ICM的故障聚類性對工況變化不如敏感ISK。從y軸方向比較，IEM具有明顯一致性，而ISK和ICM在各工況下均具有明顯差異，說明ISK和ICM的故障可分性更直觀。而進(jìn)一步觀察圖7可知，ICM比IEM和ISK具有更大的類間間距，和明顯的故障聚類性和故障可分性。因?yàn)镾K受軸承轉(zhuǎn)速影響，不同轉(zhuǎn)頻得到的ISK值差異很大，因此ISK故障可分性很差。ICM則包含了IEM忽略的在不同故障時(shí)的振幅差異信息，因此IEM可分性不如ICM。綜上所述，ICM在輕微故障時(shí)對工況干擾不敏感，而IEM和ISK的故障可分性則受工況影響明顯。

圖8為嚴(yán)重故障的特征樣本直方圖。從x軸方向比較，可發(fā)現(xiàn)ISK在同故障下的幅值變化很大，例如滾子故障的軸承轉(zhuǎn)速越快，負(fù)載越輕，ISK2幅值越小；而IEM和ICM沒有明顯變化，說明IEM和ICM的聚類特性對工況變化不敏感。y軸方向比較，C工況下的IEM具有明顯一致性，而ISK和ICM在各工況下均具有明顯的差異性，說明ISK和ICM具有更靈敏的故障可分性。進(jìn)一步觀察圖9可知，ICM比IEM和ISK具有更大的類間間距，和明顯的故障聚類性和故障可分性。綜上所述，ICM在嚴(yán)重故障時(shí)對工況干擾不敏感，而IEM和ISK的故障可分性則受工況影響明顯。

IEM，ISK和ICM的可分性判據(jù)見表1，可從J1判斷三者的聚類性相當(dāng)；在輕微和嚴(yán)重故障時(shí)，從J2可知ICM的總體類間間距均最大；從J3可知ICM的可分性均最強(qiáng)。因此，ICM更適于作滾動軸承故障特征。

圖6 輕微故障的融合特征直方圖Fig.6 Specimen chromatogram of integrated features in slight fault situation

圖7 輕微故障的融合特征分布Fig.7 Integrated feature distribution in slight fault situation

圖8 嚴(yán)重故障的融合特征直方圖Fig.8 Specimen chromatogram of integrated feature in serious fault situation

圖9 嚴(yán)重故障的融合特征分布Fig.9 Integrated feature distribution in serious fault situation

表1 故障特征有效性判據(jù)

3 結(jié)束語

筆者提出了基于LMD-ICM-PCA的滾動軸承故障魯棒辨識方法，證實(shí)了該方法有效地隔離工況變化對軸承故障辨識的影響。通過在輕微故障與嚴(yán)重故障兩種情形下，對交叉工況的滾動軸承故障數(shù)據(jù)進(jìn)行實(shí)驗(yàn)驗(yàn)證。對比融合能量矩(IEM)及融合譜峭度(ISK)的空間分布可知，ICM比IEM和ISK的故障可分性以及故障聚類性都要好，且ICM對工況干擾具有明顯的魯棒性，而IEM和ISK不具備這一特性。

盡早地做出準(zhǔn)確預(yù)警，對提高預(yù)警滾動軸承安全隱患的效率具有重要意義。而驗(yàn)證不同損傷程度時(shí)的識別能力，亦有助于進(jìn)一步研究軸承性能的退化規(guī)律，對軸承運(yùn)行的可靠性分析有重要的理論參考價(jià)值。

[1] Chen Fafa, Tang Baoping, Song Tao, et al. Multi-fault diagnosis study on roller bearing based on multi-kernel support vector machine with chaotic particle swarm optimization [J]. Measurement, 2014, 47(1): 576-590.

[2] Bin Guangfu, Gao Jinji, Li Xuejun, et al. Early fault diagnosis of rotating machinery based on wavelet packets—empirical mode decomposition feature extraction and neural network [J]. Mechanical Systems and Signal Processing, 2012, 27: 696-711.

[3] 鄭近德, 潘海洋, 張俊, 等. APEEMD及其在碰摩故障診斷中的應(yīng)用 [J]. 振動、測試與診斷, 2016, 36(2): 257-263.

Zheng Jinde, Pan Haiyang, Zhang Jun, et al. Adaptive partly ensemble empirical mode decomposition and its application for rotor rubbing fault diagnosis [J]. Journal of Vibration, Measurement & Diagnosis, 2016, 36(2): 257-263. (in Chinese)

[4] Han Minghong, Pan Jiali. A fault diagnosis method combined with LMD, sample entropy and energy ratio for roller bearings [J]. Measurement, 2015, 76: 7-19.

[5] Henr′iquez P, Alonso J B, Ferrer M A, et al. Review of automatic fault diagnosis systems using audio and vibration signals [J]. IEEE Transactions on Systems Man and Cybernetics: Systems, 2014, 44(5): 642-652.

[6] 張?jiān)茝?qiáng), 張培林, 吳定海, 等. 基于CSLBP的軸承信號時(shí)頻特征提取方法 [J]. 振動、測試與診斷, 2016, 36(1): 22-27.

Zhang Yunqiang, Zhang Peilin, Wu Dinghai, et al. Time-frequency feature extraction method based on cslbp for bearing signals [J]. Journal of Vibration, Measurement & Diagnosis, 2016, 36(1): 22-27. (in Chinese)

[7] 付云驍, 賈利民, 季常煦，等. 基于多維振動特征的滾動軸承故障診斷方法[J]. 噪聲與振動控制, 2014, 34(3): 165-169.

Fu Yunxiao, Jia Limin, Ji Changxu, et al. Fault diagnosis method of rolling bearings based on multi-dimensional vibration features[J]. Noise and Vibration Control, 2014, 34(3): 165-169. (in Chinese)

[8] 王奉濤, 陳守海, 閆達(dá)文, 等. 基于流形-奇異值熵的滾動軸承故障特征提取[J]. 振動、測試與診斷, 2016, 36(2): 288-294.

Wang Fengtao, Chen Shouhai, Yan Dawen. et al. Fault feature extraction method for rolling bearing based on manifold and singular values entropy [J]. Journal of Vibration, Measurement & Diagnosis, 2016, 36(2): 288-294. (in Chinese)

[9] Santamaria I, Pokharel P P, Principe J C. Generalized correlation function: definition, properties, and application to blind equalization[J]. IEEE Transactions on Signal Processing, 2006, 54(6): 2187-2197.

[10]Hassan M, Terrien J, Marque C, et al. Comparison between approximate entropy, correntropy and time reversibility: application to uterine electromyogram signals[J]. Medical Engineering & Physics, 2011, 33(8):980-986.

[11]付云驍, 賈利民, 秦勇, 等. 基于乘積函數(shù)相關(guān)熵的滾動軸承故障辨識方法 [J]. 應(yīng)用基礎(chǔ)與工程科學(xué)學(xué)報(bào), 2016, 24(2): 333-343.

Fu Yunxiao, Jia Limin, Qin Yong, et al. Roller element bearing fault identification method based on product function correntropy [J]. Journal of Basic Science and Engineering, 2016, 24(2): 333-343. (in Chinese)

[12]李巍華, 林龍, 單外平. 基于廣義S變換與雙向2DPCA的軸承故障診斷 [J]. 振動、測試與診斷, 2015, 35(3):499-506.

Li Weihua, Lin Long, Shan Waiping. Bearing fault diagnosis based on generalized S-Transform and directional 2DPCA[J]. Journal of Vibration, Measurement & Diagnosis, 2015, 35(3): 499-506. (in Chinese)

[13]Wang Yanxue, Xiang Jiawei, Markert R, et al. Spectral kurtosis for fault detection, diagnosis and prognostics of rotating machines: A review with applications [J]. Mechanical System and Signal Processing, 2015, 66-67: 679-698.

*科技部科技支撐計(jì)劃資助項(xiàng)目(116B300011)；軌道交通控制與安全國家重點(diǎn)實(shí)驗(yàn)室自主研究課題項(xiàng)目(116K00100)

2016-03-24；

2016-06-03

10.16450/j.cnki.issn.1004-6801.2017.02.006

TH133.33； U279.3+23

付云驍，男，1990年3月生，博士生。主要研究方向?yàn)檐壍懒熊囍悄馨踩珯z測。曾發(fā)表《基于多維振動特征的滾動軸承故障診斷方法》(《噪聲與振動控制》2014年第34卷第3期)等論文。 E-mail：yunxiaof2012@163.com