一種盲源提取方法及其在滾動軸承故障特征提取中的應(yīng)用

王宏超, 陳 進(jìn), 董廣明

(上海交通大學(xué)機(jī)械系統(tǒng)與振動國家重點實驗室, 上海 200240)

引 言

由于各種干擾源和噪聲的影響，來自現(xiàn)場傳感器的軸承故障信號非常復(fù)雜，直接對這些信號進(jìn)行故障特征提取是非常困難的。盲源提取(Blind source extraction,BSE)技術(shù)的出現(xiàn)為這個問題的解決提供了新的技術(shù)和手段。BSE技術(shù)是近10年來新興并引起廣泛重視的一類信號處理方法，有時也叫盲信號抽取，它源于盲源分離(Blind source separation, BSS)技術(shù)。BSS技術(shù)在語音、通訊和醫(yī)學(xué)工程等領(lǐng)域都有著成功的應(yīng)用[1～3]，但其在旋轉(zhuǎn)機(jī)械故障診斷中的應(yīng)用有著較大的局限性，其原因可歸結(jié)為如下[4]：瞬時混合模型的幅值不確定性和順序不確定性；適用于振動信號的多通道卷積混合盲分離算法迄今仍沒有令人滿意的解決方案；源物理相關(guān)性的不確定性；源數(shù)的估計存在偶然性；盲源分離要求系統(tǒng)具有可逆性，但實際工程中很多機(jī)械系統(tǒng)是不可逆的。BSE相對于BSS有以下優(yōu)點[5]：充分利用被提取對象的先驗信息；僅僅提取感興趣的信號；無需估計源的數(shù)目；計算量??；特別適用于源數(shù)較多而感興趣信號較少的情況(如機(jī)械系統(tǒng))。綜上所述，BSE在某種程度上彌補(bǔ)了BSS的缺點。在既有的BSE技術(shù)中，約束獨立成分分析最具代表性，其最先由Wei Lu等人于2000年提出[5]。隨后Wei Lu分別將CICA用于圖像盲分離與醫(yī)學(xué)信號盲分離[6，7]，相對于傳統(tǒng)獨立成分分析(Independent component analysis, ICA)方法取得了不錯的效果。Zhang Z L將CICA方法加以改進(jìn)，提出了基于形態(tài)的CICA方法并將其用于微弱腦電信號的盲源提取[8]，通過仿真合成信號和實際腦電信號驗證了所提出方法的正確性及實際應(yīng)用價值。在旋轉(zhuǎn)機(jī)械故障診斷方面，CICA也取得了一定的應(yīng)用[9]，在文獻(xiàn)[9]中Wang Z Y等將CICA用于滾動軸承的故障診斷中，并通過仿真信號和滾動軸承全壽命加速疲勞實驗信號驗證了CICA在旋轉(zhuǎn)機(jī)械盲信號提取中的有效性。然而，上述CICA在盲源提取中的成功應(yīng)用都需要以下兩個前提條件：用于構(gòu)建參考信號的目標(biāo)源信號基本周期的精確估計；正確參考信號的構(gòu)建。如文獻(xiàn)[9]所驗證，不正確參考信號的選取會造成誤判或漏判的結(jié)果。此外，用于構(gòu)建參考信號的目標(biāo)源信號基本周期的估計值也會對提取結(jié)果造成很大的影響。然而在滾動軸承的實際工程應(yīng)用中由于安裝精度、制造誤差及滾子相對于滾道的隨機(jī)滑動等影響，都會造成目標(biāo)源故障信號的實際基本周期與理論計算基本周期的偏差。這樣就限制了CICA方法在旋轉(zhuǎn)機(jī)械故障診斷中的應(yīng)用。在文獻(xiàn)[10]中，Barros等提出了一種基于最佳估計延時的能快速提取周期信號的盲提取方法。但所述方法本質(zhì)上還是以目標(biāo)信號的基本周期或基本周期的整數(shù)倍作為最佳延時。本文首先用基于二階統(tǒng)計量的自相關(guān)函數(shù)算法估計目標(biāo)源信號基本周期[10]，再與基于高階統(tǒng)計量的固定點算法相結(jié)合的盲提取方法成功用于復(fù)雜運(yùn)行環(huán)境下滾動軸承故障信號的盲提取[11,12]，此種方法在某種意義上來講是對文獻(xiàn)[10]所述方法的改進(jìn)。此外，此種方法相對于CICA方法具有較強(qiáng)的魯棒性：只需要大致估計目標(biāo)信號的基本周期并將其作為所述方法相關(guān)步驟的基本周期即可。同文獻(xiàn)[9]，將所述方法用于滾動軸承全壽命加速疲勞實驗故障信號的盲提取，并與文獻(xiàn)[9]CICA在滾動軸承全壽命加速疲勞故障信號盲提取中的應(yīng)用作以比較，得出前者相對于后者具有更強(qiáng)的目標(biāo)源信號基本周期估計誤差的容錯性。

1 問題的提出

假設(shè)有一多維觀測信號矢量x(t)可表示為

x(t)=As(t)+n(t)

(1)

式中s(t)=[s1(t),s2(t),…,sn(t)]T是n×1滿足統(tǒng)計獨立性的源信號矢量，x(t)=[x1(t),x2(t),…,xm(t)]T是m×1的混合信號矢量，矩陣A={aij}是由混合參數(shù)aij構(gòu)成的m×n階的混合矩陣。n(t)=[n1(t),n2(t),…nm(t)]T是與源信號矢量統(tǒng)計獨立的加性噪聲矢量。

(2)

式中 矩陣C=(ci,j)n×n被稱為混合-分離復(fù)合矩陣，該矩陣的每一行與每一列有且僅有一個非零元素。瞬時混合模型的混合和分離過程如圖1所示。

圖1 瞬時混合和分離模型示意圖

根據(jù)傳感器數(shù)m和源數(shù)n的關(guān)系，瞬時ICA模型可劃分為3種：

當(dāng)m=n時，此時混合矩陣A是方陣，稱為平方ICA模型；

當(dāng)m

當(dāng)m>n時，稱為超定ICA模型。

本文是基于平方ICA模型展開討論的。

此外，為提高計算速度及簡化本文算法，首先對觀測信號x(t)進(jìn)行白化處理，即

(3)

式中V稱為白化矩陣。以下所述所有觀測信號x(t)均為白化信號。

對于要提取的目標(biāo)源信號si，假設(shè)對一特定整數(shù)τ*有以下關(guān)系式成立

(4)

式中sj為其他非目標(biāo)源信號，k代表時間，τ*即為最佳時間延時[10]。

J(W)=E{y(k)y(k-τ*)}=WTE{x(k)x(k-τ*)T}W

(5)

在約束條件‖W‖=1下，最大化公式(5)所示的目標(biāo)函數(shù)就可以提取出目標(biāo)函數(shù)si。因為對于目標(biāo)函數(shù)si，其關(guān)于延時為τ*的延時自相關(guān)函數(shù)是一個大于0的正值，而其他非目標(biāo)源函數(shù)sj關(guān)于τ*的延時自相關(guān)函數(shù)都將為0。

忽略Rx(τ*)與Rx(τ*)T的微小差距，利用標(biāo)準(zhǔn)梯度算法，可得出文獻(xiàn)[10]的算法

(6)

式中Rx(τ*)=E{x(k)x(k-τ*)T}。盡管文獻(xiàn)[11]所述方法具有算法簡潔、計算速度快等優(yōu)點，但是其具有以下缺點：

(1)即使目標(biāo)源信號si關(guān)于τ*的時延自相關(guān)函數(shù)E{si(k)si(k-τ*)}>0，但不能保證其他所有非目標(biāo)源信號sj關(guān)于τ*的時延自相關(guān)E{sj(k)sj(k-τ*)}=0(j≠i)；

(2)即使所有源信號嚴(yán)格獨立不相關(guān)，但是在實際計算中由于用源信號有限點數(shù)的算術(shù)平均來代替源信號的數(shù)學(xué)期望，就會造成源信號的自相關(guān)計算值不為0[15,16]，即：即使有下式成立

(7)

但也極有可能下式成立

(8)

綜上所述，文獻(xiàn)[10]所述基于二階統(tǒng)計量的盲提取方法受上述兩點限制，其盲提取性能的穩(wěn)定性及精確性就受到了很大的影響，甚至?xí)崛〕龇悄繕?biāo)源信號。由于本文的研究是基于瞬時線性ICA基礎(chǔ)上的，也就是假設(shè)源信號在物理意義上是相互獨立的，所以可以用高于二階統(tǒng)計量的更高階的統(tǒng)計量去改善文獻(xiàn)[10]所述的方法。以下是所述方法的介紹。

2 方法的提出

所用的方法大概分為3個步驟，如圖2所示。

圖2 所述方法流程圖

各個步驟的具體計算過程如下：

(1)目標(biāo)源信號基本周期的估計：當(dāng)滾動軸承發(fā)生故障時，其故障基本周期可以由理論計算頻率的倒數(shù)得出；

(2)初始目標(biāo)源信號及分離矩陣計算：

式(5)所示的目標(biāo)函數(shù)

(9)

(10)

(11)

由式(10)可得出式(11)的相應(yīng)算法如下

(12)

(3)初始分離矩陣W的進(jìn)一步優(yōu)化：

(13)

為提高算法對野點和沖擊噪聲的魯棒性，可以采用文獻(xiàn)[11]的改進(jìn)固定點算法

(14)

式(13)或(14)收斂時，即可得到最佳分離矩陣W。以下的仿真及實驗結(jié)果均取式(13)，(14)二者之間效果最好的。

圖3 原始信號

圖4 觀測信號

圖5 不同盲提取方法得到的結(jié)果

3 仿 真

圖3為4個信號，其中圖3中的s1為準(zhǔn)周期信號，s2為正弦信號，s3為余弦信號，s4為隨機(jī)高斯白噪聲。在Matlab中任意產(chǎn)生一隨機(jī)4×4階的矩陣A對4個信號進(jìn)行混合。圖4為混合后的觀測信號。提取目標(biāo)信號是s1。圖5為用不同方法得到盲抽取結(jié)果。其中圖5中的y1為文獻(xiàn)[10]所述方法的盲抽取結(jié)果，得到了錯誤的抽取結(jié)果；圖5中y2為文獻(xiàn)[15]所述方法的提取結(jié)果；圖5中的y3為本文所述方法的提取結(jié)果。雖然直觀上看y2與y3要都為正確的目標(biāo)源信號s1?？梢杂孟率絹韺μ崛〗Y(jié)果y2和y3的盲抽取精度進(jìn)行比較

(15)

將s1和y2，s1和y3分別代入式(15)得到的PI值分別為22.4和47.6 dB，由此，所述方法相對于文獻(xiàn)[15]所提出的方法具有更高的提取精度。

圖6 基于不同參考信號的CICA抽取結(jié)果

圖6是CICA方法用不同參考信號對觀測信號的盲抽取結(jié)果。其中參考信號r1為正確的參考信號，將其與觀測信號一起輸入到CICA算法中，y11為抽取的結(jié)果。從y11可以看出CICA在正確輸入?yún)⒖夹盘柕那闆r下可以很好地抽取出目標(biāo)信號s1；不改變參考信號r1的方波周期，而改變方波的寬度即得到參考信號r2，將r2與觀測信號一起輸入到CICA算法中，y22為抽取的結(jié)果，可見抽取不出目標(biāo)信號。同樣不改變r1方波的寬度，而改變方波的周期得到參考信號r3，將r3與觀測信號一起輸入到CICA算法中，y33為抽取的結(jié)果，可見仍抽取不出目標(biāo)信號。由此，CICA算法的能否成功應(yīng)用在某種程度上取決于參考信號的能否正確設(shè)計，其使用受到了很大的限制。而分別將r1，r3的周期作為所述方法的基本估計周期(或基本周期的整數(shù)倍)，均能提取出目標(biāo)源信號s1。限于篇幅限制在此不再給出抽取結(jié)果。

4 實 驗

同樣采用文獻(xiàn)[9]中的滾動軸承全壽命加速疲勞實驗數(shù)據(jù)。實驗具體實施步驟及數(shù)據(jù)采集參數(shù)均可見參考文獻(xiàn)[9]。選用和文獻(xiàn)[9]同一個故障軸承進(jìn)行分析，文獻(xiàn)[9]給出了故障軸承的全壽命周期的均方根值(Root mean square, RMS)值曲線圖。本文給出故障軸承全壽命周期的峭度指標(biāo)曲線圖，如圖7所示，因為峭度指標(biāo)對滾動軸承發(fā)生故障時的沖擊特征更為敏感。

圖7 實驗軸承全壽命周期的峭度值

對第2 304組數(shù)據(jù)進(jìn)行分析。將3個加速度傳感器采集到的信號作為觀測信號，其時域圖如圖8所示(注：文獻(xiàn)[9]的時域圖縱坐標(biāo)單位為加速度，而本文的縱坐標(biāo)單位為電壓)。

圖9為圖8所示信號相對應(yīng)的包絡(luò)解調(diào)譜，從圖9中均無法很好地提取出滾動軸承內(nèi)圈故障通過頻率fi=246 Hz及其諧頻：圖9(a)雖然能大致提取出內(nèi)圈故障通過頻率，但其諧頻及調(diào)制頻率即轉(zhuǎn)頻提取效果并不好；圖9(b)，(c)的提取結(jié)果會造成誤判或漏判。

將本文所述方法用于上述觀測信號，提取出的信號的時域圖及其包絡(luò)解調(diào)譜圖如圖10所示，由包絡(luò)解調(diào)譜可以看出其很好的提取出內(nèi)圈故障通過頻率fi=246 Hz及其諧頻。在用所述方法對觀測信號進(jìn)行處理時，估計目標(biāo)信號的基本周期設(shè)置為T=(1/fi)×fs=104(用點數(shù)表示，可參考文獻(xiàn)[9])。

圖11是基于不同參考信號的CICA盲提取結(jié)果，其中圖11(a)的r1是正確參考信號：周期為T=104(采用點數(shù)表示，可參考文獻(xiàn)[9])；方波寬度為20(點數(shù))，圖11(b)為盲抽取的信號，圖11(c)為11(b)的包絡(luò)解調(diào)譜，由此看見，CICA很好地提取了內(nèi)圈故障通過頻率；改變圖11(a)的周期，即T′=106，方波寬度不變，得到參考信號r2如圖11(d)；圖11(e)為盲提取結(jié)果的時域圖，圖11(f)為圖11(e)的包絡(luò)解調(diào)譜，提取不出內(nèi)圈故障頻率；同樣，不改變圖11(a)的周期，而改變方波寬度(方波寬度由20設(shè)置為30)即得到參考信號r3如圖11(g)；圖11(h)，(i)分別為盲提取結(jié)果的時域圖及包絡(luò)解調(diào)譜，同樣提取不出內(nèi)圈故障特征頻率?？梢?，CICA能否正確地提取故障特征頻率，對精確參考信號的設(shè)計有著苛刻的要求(參考信號的周期及方波形狀的微小改變都會對CICA提取結(jié)果造成很大的影響)。在工程實際應(yīng)用中，由于設(shè)計精度、安裝誤差及滾動體的隨機(jī)滑動都會造成滾動軸承理論計算故障頻率與實際故障頻率的誤差。而本文所述的方法只需要大致估計目標(biāo)源信號的基本周期即可。圖12為將T=106作為目標(biāo)源信號的估計周期，用本文所述方法得到盲提取信號的時域圖及其相應(yīng)的包絡(luò)解調(diào)譜，仍可以提取出內(nèi)圈故障特征頻率。說明了本文所述方法相對于CICA方法有較強(qiáng)的目標(biāo)源信號周期誤差的容錯性。

圖8 實驗觀測信號

圖9 實驗觀測信號的包絡(luò)解調(diào)譜

圖10 所述方法的盲提取結(jié)果

圖13為文獻(xiàn)[15]所述方法的提取效果，從其包絡(luò)解調(diào)譜中無法提取中滾動軸承內(nèi)圈故障通過頻率。

圖11 用不同參考信號CICA的提取結(jié)果

圖12 用所述方法以T=106作為估計基本周期的提取結(jié)果

圖13 用文獻(xiàn)[15]所述方法的盲提取結(jié)果

5 結(jié) 論

本文將一種盲提取方法用于復(fù)雜運(yùn)行環(huán)境下滾動軸承故障信號的盲提取中，取得了不錯的結(jié)果。并通過仿真驗證了所述方法相對于其他盲提取方法有較高的精度；通過仿真和實驗驗證了所述方法相對于CICA方法具有目標(biāo)源信號基本周期誤差容錯性的優(yōu)點；此外，所述方法只需要估計目標(biāo)源信號的基本周期即可，而CICA方法不但需要目標(biāo)源信號基本周期的精確估計，而且還需要構(gòu)建準(zhǔn)確的參考信號，其使用相對于所述方法具有較大的局限性。

參考文獻(xiàn)：

[1] Hyvarinen A, Oja E. Independent component analysis:algorithms and applications[J]. Neural Networks， 2000, 12:411—430.

[2] Choi S J, Cichocki A, Park H M, et al. Blind source separation and independent component analysis:A review[J]. Neural Information Processing-letters and Reviews, 2005, 6(1):1—57.

[3] James C J, Gibson O J. Temporally constrained ICA:An application to artifact rejection in electromagnetic brain signal analysis[J]. IEEE Transaction on Biomedical Engineering, 2003, 50(1):1 108—1 116.

[4] Antoni J. Blind separation of vibration components:rinciples and demonstrations[J]. Mechanical Systems and Signal Processing 2005, 19:1 166—1 180.

[5] Lu W, Rajapakse J C. Constrained independent component analysis, in:Advance in Neural Information Processing Systems, vol.13(NIPS 2000)[M]. MIT Press, Cambridge, MA, 2000:570—576.

[6] Lu W, Rajapakse J C. Approach and applications of constrained ICA[J]. IEEE Transaction on Neural Networks, 2005, 16(1):203—212.

[7] Lu W, Rajapakse J C. ICA with reference[J]. Neurocomputing, 2006, 69:2 244—2 257.

[8] Zhang Z L. Morphologically constrained ICA for extracting weak temporally correlated signals[J]. Neurocomputing, 2008, 71:1 669—1 679.

[9] Wang Z Y, Chen J, Dong G M, et al. Constrained independent component analysis and its application to machine fault diagnosis[J]. Mechanical Systems and Signal Processing, 2011, 25:2 501—2 512.

[10] Barros A K, Cichocki A. Extraction of specific signals with temporal structure[J]. Neural Computation, 2001, 13(9):1 995—2 003.

[11] Cichocki A, Amari S I. Adaptive Blind Signal and Image Processing. Learning Algorithms and Applications[M]. New York:John Wiley & Sons,2002.

[12] Hyvarinen A, Oja E. A fast fixed-point algorithm for independent component analysis[J]. Neural Computation, 1997, 9(7):1 483—1 492.

[13] Sabri K, Badaoui M E, Guillet F, et al. A frequency domain-based approach for blind MIMO system identification using second-order cyclic statistics[J]. Signal Processing, 2009, 89:77—86.

[14] Belouchrani A, Meraim K A, Cardoso J F, et al. A blind source separation technique using second-order statistics[J]. IEEE Transactions on Signal Processing, 1997, 45(2):434—444.

[15] Zhang Z L, Zhang Y. Robust extraction of specific signals with temporal structure[J]. Neurocomputing, 2006, 69:888—893.

[16] Bermejo S. Finite sample effects in higher order statistics contrast functions for sequential blind source separation[J]. IEEE Signal Processing Letters,2005, 12(6):481—484.