基于改進字典學(xué)習(xí)的單通道振動信號盲源分離算法

余 路， 曲建嶺， 高 峰， 田沿平， 郭超然， 李 儼

(1.海軍航空大學(xué)青島校區(qū) 航空儀電控制工程與指揮系，青島 266041；2.西北工業(yè)大學(xué) 自動化學(xué)院，西安 710072)

盲源分離(Blind Source Separation,BSS)是機械設(shè)備進行復(fù)合故障診斷的一種有效方法。目前在機械設(shè)備的盲分離方法中，通常假設(shè)觀測信號的個數(shù)不小于源信號的個數(shù)，然而這在實際的故障診斷實踐中并不總是成立。受限于機械設(shè)備關(guān)鍵部件的安裝位置和傳感器的布置成本，觀測信號個數(shù)小于源信號個數(shù)的欠定盲源分離問題更加普遍，而單通道盲源分離(Single Channel Blind Source Separation，SCBSS)問題可以看作是欠定盲源分離問題中的特例。由于觀測信號個數(shù)越少意味著獲得的可用信息越少，對于源信號的恢復(fù)越困難，因此，單通道盲源分離問題一直以來是盲源分離問題研究的熱點和難點問題[1-2]。

多年來，已有眾多學(xué)者對此問題展開了研究。一些學(xué)者采用經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition，EMD)及其變種理論的方法展開研究。例如，文獻[3]利用EMD算法增加觀測信號個數(shù)，將SCBSS問題轉(zhuǎn)化為正定或超定BSS問題，再使用成熟的FastICA算法求解。為了抑制EMD算法的端點效應(yīng)問題，文獻[4]采用EEMD算法(Ensemble EMD,EEMD)和基于四階累積量的特征矩陣聯(lián)合近似對角化算法(JADE，Joint Approximative Diagonalization of Eigenmatrix)并成功應(yīng)用于地聲信號的分離。這些算法雖然在一定條件下都可以較準確的恢復(fù)源信號，但EMD算法本質(zhì)的端點效應(yīng)及模態(tài)混疊問題并沒有得到根本的解決，同時，這些方法改變了信號的結(jié)構(gòu)特征，受人工經(jīng)驗的因素影響較大，因此給這類算法的廣泛應(yīng)用帶來了限制。還有部分學(xué)者考慮觀測信號的稀疏性問題，文獻[5]和文獻[6]采用時頻分析和聚類的方法解決欠定盲源分離問題，但對于SCBSS問題，線性稀疏性的條件常常不能滿足，而更復(fù)雜的非線性混合問題則需要更大的運算量和新的求解算法，目前在這方面的研究還尚在起步階段。

基于以上分析，本文采用了另外一種思路求解SCBSS問題，即將移不變字典學(xué)習(xí)算法加入到SCBSS求解問題的框架中。

具有移不變特性的稀疏編碼算法最早由Smith在2006年的《Nature》中提出[7]，其本質(zhì)是將信號以多個基函數(shù)及其相應(yīng)系數(shù)卷積的形式表示，自適應(yīng)的學(xué)習(xí)蘊含在信號中的移不變基函數(shù)，從而實現(xiàn)對原信號的稀疏表達。對于旋轉(zhuǎn)機械來說，反應(yīng)信號特征的分量往往都是重復(fù)出現(xiàn)的，比如故障軸承的沖擊，齒輪的磨損及旋翼轉(zhuǎn)速異常等[8]。在恒定工況下，不同激振源的信號特征相似，不同激振源的信號有所不同。利用移不變字典學(xué)習(xí)算法(Shift Invariant Dictionary Learning，SIDL)求解觀測信號中的移不變基函數(shù)，引入改進的模糊C均值聚類算法自適應(yīng)地將估計源信號個數(shù)的同時，實現(xiàn)對各個移不變基函數(shù)時域重構(gòu)移不變分量進行聚類，實現(xiàn)同類源信號的聚類和不同類源信號的區(qū)分，從而實現(xiàn)單通道振動信號的盲源分離。

1 SIDL學(xué)習(xí)算法

1.1 SIDL學(xué)習(xí)算法簡介

SIDL與標準稀疏編碼不同，其將輸入信號y表示成各個基函數(shù)與其系數(shù)的卷積

*xj+ε

(1)

式中，“*”表示卷積算子，dj∈Rm為基函數(shù)，M為基函數(shù)個數(shù)，xj∈Rl為稀疏系數(shù)，m為基函數(shù)長度，l為信號的分段長度，n為原始信號長度，且三者滿足m

1.2 系數(shù)求解

求解稀疏系數(shù)的過程相當(dāng)于把所有基函數(shù)固定，求解每個基函數(shù)對應(yīng)的系數(shù)

(2)

1.3 字典學(xué)習(xí)

在該階段中，固定系數(shù)x求解字典矩陣D，則原優(yōu)化問題轉(zhuǎn)化為l2約束優(yōu)化問題

(3)

式中，N為輸入信號個數(shù)。直接求解上式計算量巨大，本文在字典學(xué)習(xí)階段采用時頻域轉(zhuǎn)化的思想，將時域計算量較大的卷積運算轉(zhuǎn)換到頻域，由頻域的乘積運算代替時域的卷積運算，使運算復(fù)雜度由O(N2)降低為O(NlgN)，提高了運算效率。由Parseval定理可知，式(3)可轉(zhuǎn)化為如下優(yōu)化問題[13]

(4)

(5)

其中，λ∈RM為對偶變量，單位向量I∈RM，并且：

(6)

2 模糊C均值聚類算法

模糊C均值(Fuzzy C-Means，F(xiàn)CM)聚類算法已被廣泛應(yīng)用于數(shù)據(jù)挖掘和信號分析[14]，與傳統(tǒng)的基于硬閾值聚類的算法不同，F(xiàn)CM聚類算法利用模糊集理論實現(xiàn)對數(shù)據(jù)的軟劃分。文獻[15]指出，對于振動信號這類非平穩(wěn)信號，采用軟閾值比硬閾值劃分往往能取得更佳的效果，因此本文采用模糊C均值聚類算法實現(xiàn)移不變分量的聚類。

給定數(shù)據(jù)集X=[x1,x2,…,xn]為s維空間Rs的有限觀測樣本，采用FCM算法針對每個數(shù)據(jù)點產(chǎn)生模糊隸屬關(guān)系uij=ui(xj)∈[0,1]。聚類中心由采樣點的數(shù)值及其隸屬化目標函數(shù)的最小值確定

(7)

(8)

(9)

3 基于改進字典學(xué)習(xí)的單通道振動信號盲源分離算法

利用SIDL自學(xué)習(xí)觀測信號中的移不變基函數(shù)，將各個基函數(shù)進行時域重構(gòu)得到移不變分量(Shift Invariant Component，SIC)，將各個SIC進行分段并提取關(guān)鍵點，將FCM算法對各SIC進行自適應(yīng)聚類，以聚類后源信號的相關(guān)度最小作為對源信號數(shù)目的最佳估計。輸出聚類后的SIC作為最終分離的源信號。

算法流程圖如圖1所示。本文算法的具體實現(xiàn)步驟如下。

圖1 本文算法簡要流程圖

步驟1 將觀測信號進行分段和平滑預(yù)處理。

步驟2 利用SIDL算法自學(xué)習(xí)信號中的移不變基函數(shù)。

步驟3 將各移不變基函數(shù)進行時域重構(gòu)得到各個移不變分量(SIC)。

步驟4 求解各SIC的包絡(luò)譜。

步驟5 將各SIC的包絡(luò)譜進行分段并執(zhí)行局部最大值檢測，得到反應(yīng)包絡(luò)譜概略波形特征的峰值點。

步驟6 將反映波形特征的峰值進行模糊C均值聚類，聚類個數(shù)依次從2到SIC的總個數(shù)，依次計算不同聚類數(shù)下的恢復(fù)的源信號相關(guān)系數(shù)平均值。

步驟7 以相關(guān)系數(shù)平均值最低作為對源信號個數(shù)的估計，并對此時的SIC進行合并，得到分離的源信號。

4 實驗驗證

4.1 實驗設(shè)置

在大型設(shè)備中，齒輪和軸承是重要的機械部件。當(dāng)齒輪和軸承發(fā)生故障時，往往給機械設(shè)備的健康性能帶來損失。因此，本文根據(jù)齒輪與軸承的故障模型，分別構(gòu)造了故障齒輪信號和故障軸承信號的仿真模型[17-18]，同時加以白噪聲。

當(dāng)齒輪發(fā)生故障時，往往表現(xiàn)出幅值調(diào)制與相位調(diào)制的現(xiàn)象，其振動信號模型為

(10)

式中，Am(t)和θm分別為第m次諧波的幅值和相位，f1為齒輪嚙合頻率。在本例中取M=5，Am(t)=2.5[sin(60πt)+2]/m，θm=mπ/4，f1=340 Hz。

軸承發(fā)生故障時往往伴隨著振動信號的沖擊響應(yīng)，因此故障軸承的模型可以簡化為

y2=a(t)e-ptsin(2πf2t+φ)

(11)

式中，a(t)為調(diào)幅幅值，p為衰減系數(shù)，f2為軸承固有頻率，φ為相位。在本例中取a(t)=20[sin(60πt)+2]，φ=0，f2=3 400 Hz，T=0.01 s，采樣頻率fs=12 kHz。觀測信號可以表示為y=y1+y2+yn，其中yn為高斯白噪聲，使觀測信號與原始信號的信噪比為10 dB，源信號與觀測信號如圖2所示(為了更好的展示細節(jié)，各個信號取1 000個點，下同)。應(yīng)用本文提出的字典學(xué)習(xí)算法，設(shè)置參數(shù)如表1所示。

圖2 源信號及加噪混合信號

參數(shù)數(shù)值基函數(shù)長度50基函數(shù)個數(shù)6子信號長度256子信號重疊率50%稀疏度160迭代次數(shù)60

4.2 仿真實驗

經(jīng)過60次迭代后算法收斂，學(xué)習(xí)到的基函數(shù)及其對應(yīng)的移不變分量(SIC)如圖3所示，從基函數(shù)中可以隱約看出，故障軸承的沖擊響應(yīng)與故障齒輪的調(diào)制特性，但從對應(yīng)的SIC中可以比較明顯的看出兩種故障已經(jīng)實現(xiàn)了區(qū)分。為了以數(shù)值計算的方式而不是主觀臆斷作為兩種故障區(qū)分的依據(jù)，考慮將各個SIC轉(zhuǎn)化至頻域，采用包絡(luò)譜對其進行進一步分析。

由于采樣點數(shù)以及采樣頻率的原因，原始包絡(luò)譜呈現(xiàn)出凌亂復(fù)雜的波形，不利于后期的處理分析。因此本文提出了局部最大值檢測的處理策略，即將包絡(luò)譜進行不重疊分段處理，以每一區(qū)間段內(nèi)的最大值點作為關(guān)鍵點，將所有關(guān)鍵點平滑連接形成了易于理解與處理的包絡(luò)譜波形特征，同時波形的峰值關(guān)鍵點可以作為接下來模糊C均值聚類的特征向量加以分析處理，如圖4所示。考慮到信號源個數(shù)的未知性，本文采取自適應(yīng)的源信號數(shù)目估計，即將聚類個數(shù)從2開始逐一增加到預(yù)先設(shè)定的基函數(shù)個數(shù)，在每個聚類數(shù)目下執(zhí)行模糊C均值聚類得到源信號的估計，若達到理想的分離效果，則分離的源信號之間的皮爾遜相關(guān)系數(shù)的均值達到最小。在本例中，聚類數(shù)目與分離的源信號相關(guān)系數(shù)均值(corrmean)的關(guān)系如圖5所示，可以看出，當(dāng)取值為2時達到最小值，這也符合本例的實際情況。

圖3 學(xué)習(xí)到的基函數(shù)及其對應(yīng)的時域移不變分量

圖4 各個移不變分量的包絡(luò)譜(實線)及局部最大值檢測(虛線)的結(jié)果

Fig.4 Envelope spectrums of SICs(line) and results of local maximum detection(dotted line)

圖5 聚類數(shù)與corrmean關(guān)系

聚類數(shù)目為2時分離的源信號的時域波形如圖6所示。從圖中可以看出，分離效果良好，保留了故障軸承和故障齒輪的主要特征信息，只是在波形的幅值上與原始信號存在一定的出入。

圖6 本文算法分離結(jié)果

作為對比，分別以文獻[3]和文獻[4]的方法進行實驗，實驗結(jié)果分別如圖7和圖8所示，可以看出并沒有達到良好的分離效果。為了定量衡量本文算法與其他算法的分離效果，采用峰值信噪比(SNR)和均方差(MSE)作為指標參數(shù)，不同算法的實驗結(jié)果如表2所示。可以看出，本文算法在三者中保持最高信噪比的同時具有最低的均方差，表明本文算法較好的保持了源信號的特征。

圖7 EMD+ICA法分離結(jié)果

圖8 EEMD+JADE法分離結(jié)果

4.3 實測航空發(fā)動機振動信號實驗

為驗證本文算法在實際裝備中的可行性，將其應(yīng)用于某型直升機振動數(shù)據(jù)采集與處理。眾所周知，直升機的飛行安全至關(guān)重要，一旦機體發(fā)生故障將會給機組人員的生命和裝備造成損失，因此常利用振動信號傳感器采集機體各位置處的振動信號作為故障診斷的依據(jù)，但受限于直升機各部件的緊湊性，單個傳感器采集到的信號可能是多個機體振源的疊加，因此有必要對其進行分離作進一步的判斷。直升機部分振動傳感器的安裝如圖9所示。

表2不同算法分離效果對比

Tab.2Comparisonofseparatingresultsofdifferentmethods

方法y1y2SNRMSESNRMSEEMD+ICA(文獻[3])2.88.353.55.38EEMD+JADE(文獻[4])5.23.954.34.46本文算法13.51.0514.21.21

將位于右發(fā)動機前部的振動傳感器采集到的信號作為研究對象，其時域波形和頻域包絡(luò)譜如圖10所示，從圖中無法看出明顯的特征。設(shè)置基函數(shù)個數(shù)為4，其余參數(shù)與表1相同，應(yīng)用本文的字典學(xué)習(xí)算法學(xué)習(xí)到的基函數(shù)如圖11所示，對其時域重構(gòu)的移不變分量(SIC)進行包絡(luò)譜分析及局域最大值檢測，實驗結(jié)果如圖12所示。從圖中可以看出SIC1與SIC3的包絡(luò)譜略有相似，SIC2與SIC4的包絡(luò)譜也有相似之處，為了驗證猜想，應(yīng)用模糊C均值聚類算法，聚類個數(shù)從2到4，發(fā)現(xiàn)在聚類個數(shù)為2時達到corrmean最小值，因此判斷采集到的振動信號包含2個信號源。

(a) 觀測信號時域波形

(b) 觀測信號包絡(luò)圖

圖11 學(xué)習(xí)到的基函數(shù)

圖12 各移不變分量的包絡(luò)譜及局部最大值檢測

Fig.12 Envelope spectrums of SICs and local maximum detection

分離的時域信號如圖13所示，對各個分離信號求取包絡(luò)譜如圖14所示，可以看出分離信號1為慢車時右發(fā)動機轉(zhuǎn)子的轉(zhuǎn)速fm(293 Hz)，分離信號2為減速器齒輪振動信號(fc，2fc，3fc，4fc)。因此該振動傳感器采集到的信號已成功得到分離，一旦某部位發(fā)生故障異常振動時可以準確定位到具體的位置。

(a) 分離信號1時域波形

(b) 分離信號2時域波形

(a) 分離信號1包絡(luò)圖

(b) 分離信號2包絡(luò)圖

作為對比，分別將文獻[3]和文獻[4]的方法應(yīng)用于上述實驗數(shù)據(jù)，分離后信號頻譜分別如圖15和16所示。

(a) 分離信號1包絡(luò)圖

(b) 分離信號2包絡(luò)圖

Fig.15 Envelope spectrum of separated signals by paper[3]

(a) 分離信號1包絡(luò)圖

(b) 分離信號2包絡(luò)圖

Fig.16 Envelope spectrum of separated signals by paper[4]

通過分離信號的包絡(luò)譜對比可以看出，文獻[3]和文獻[4]方法均一定程度上存在兩種信號源混疊的現(xiàn)象。例如，文獻[3]中的分離信號1中含有減速器齒輪振動信號的各階譜，但也部分含有右發(fā)動機轉(zhuǎn)子轉(zhuǎn)速的頻率，給分離信號的識別增加了困難。同理，文獻[4]方法的分離結(jié)果也存在類似的問題，表明分離效果較本文算法存在差距。

取三種算法各運行20次的平均運算時間結(jié)果進行對比，實驗軟件環(huán)境配置采用Linux系統(tǒng)下Matlab2014a，硬件環(huán)境采用i5HQ6300處理器，8 GB內(nèi)存。實驗結(jié)果如表3所示。

表3 不同算法運算時間對比

不可否認，文獻[3]和文獻[4]算法相比于本文算法具有更快的運算速度，本文算法的時間花費主要在于移不變基函數(shù)的迭代求解，但通過迭代求解的移不變基函數(shù)可以更好地學(xué)習(xí)到觀測信號中存在的本質(zhì)特征，因此可以獲得更佳的信號分離效果。

5 結(jié) 論

本文提出了一種采用改進字典學(xué)習(xí)對單通道振動信號進行特征自學(xué)習(xí)和盲源分離的解決方案。移不變字典學(xué)習(xí)算法能夠自學(xué)習(xí)反映信號本質(zhì)特征的移不變基函數(shù)，解決了特征提取依賴專家知識與人工提取的困難；以局部最大值檢測算法提取各個移不變分量的包絡(luò)譜關(guān)鍵點，降低了后期聚類運算的效率和準確性，而且更有利于反映頻域信號的波形特征；利用自適應(yīng)模糊C均值聚類算法實現(xiàn)對源信號的恢復(fù)，解決了源數(shù)估計和信號恢復(fù)的問題。由仿真和實測數(shù)據(jù)的實驗結(jié)果可知，本文提出的算法具有較高的源信號恢復(fù)精度，有利于對單通道振動信號的進一步分析，具有較強的工程應(yīng)用推廣價值。