張延良， 張玉， 張偉濤

(1.河南理工大學物理與電子信息學院， 焦作 454150； 2.西安電子科技大學電子工程學院， 西安 710071)

盲源分離算法(blind source separation, BSS)是在源信號統(tǒng)計獨立，且源信號和傳輸信道等先驗知識均未知的情況下，僅根據(jù)觀測到的混合信號將各個源信號分離出來的技術[1-2]。盲源分離研究始于20世紀80年代，短短二十幾年間，其相關理論和實際應用快速發(fā)展[3]。目前如何將感興趣的信號從混合信號中有效地提取是信號處理領域的研究熱點和重點之一。

帶參考信號的盲源分離算法也被稱為約束獨立成分分析算法(constrained independent component analysis, CICA)[4]，是提取感興趣信號最有效的方法。許多情況下，源信號的特性并不是一無所知的。源信號的一些特征如載波頻率、調制方式等先驗信息是已知的，將這些先驗信息引入到獨立成分分析(independent component analysis, ICA)中可以用于感興趣信號的提取[5]。相對于ICA，CICA算法更適合實際應用。例如，在提取語音信號時可能只對某些信號感興趣，只需要根據(jù)感興趣信號的某些特征來進行提取，顯然CICA算法更有優(yōu)勢。

目前CICA算法已經被廣泛應用于機械故障診斷[6-7]、語音識別[8]、圖像處理[9]、通信[10]等領域。近年來，不少學者將CICA算法應用在生物醫(yī)學方面。王凱等[11]提出了一種實時約束獨立成分分析方法，在去除心電偽跡上得到了較好的效果；Shi等[12]提出了一種基于多目標優(yōu)化框架的約束時空ICA 方法，對感興趣信號的恢復能力在一定程度上得到了提高。但在這些應用中，提取信號與源信號之間還是存在較大誤差，因此研究如何減小提取信號與源信號之間的誤差具有重要意義。

為了CICA算法在其應用中有更好的提取效果，現(xiàn)提出直接將先驗信息耦合到基于負熵的目標函數(shù)中形成一個新的目標函數(shù)，再使用拉格朗日乘子法對新的目標函數(shù)優(yōu)化，使其更準確向感興趣信號方向收斂，以此來減小誤差。

1 盲源分離算法

盲源分離是用以分析多維數(shù)據(jù)的線性分析算法，試圖將一組隨機變量表示成統(tǒng)計獨立變量的線性組合[13]，它的數(shù)學模型為：假設L個源信號s(t)=[s1(t),s2(t),…,sL(t)]T，通過M個傳感器接收，得到觀測信號x(t)=[x1(t),x2(t),…,xM(t)]T，忽略噪聲后可得出源信號s(t)和觀測信號x(t)之間的關系表達式為

x(t)=As(t)

(1)

式(1)中：A為一個M×L的混合矩陣。原問題則轉化成根據(jù)觀測信號x(t)選定一個目標函數(shù)，使用優(yōu)化算法對該目標函數(shù)進行優(yōu)化得到最優(yōu)的分離矩陣w，可以得到分離信號y(t)與觀測信號x(t)之間的關系表達式為

y(t)=wTx(t)

(2)

把式(1)代入式(2)中，可得

y(t)=wTAs(t)

(3)

由式(3)可知，分離信號y(t)相當于是源信號是s(t)的一個有效估計，所以盲源分離算法的實質就是求一個最優(yōu)的分離矩陣w。

2 傳統(tǒng)的CICA

為提取感興趣信號，目標函數(shù)的選擇很關鍵，負熵可以反映信號的非高斯性大小，因此選擇快速ICA算法中基于負熵的目標函數(shù)，但是負熵的計算非常困難，常用非線性函數(shù)J(y)來近似：

J(y)=[E{G(y)}-E{G(v)}]2

(4)

式(4)中：y為分離的源信號；v為一個具有零均值和單位方差的高斯向量；G(·)為任何非線性函數(shù)；E(·)為非線性函數(shù)的期望值。

算法的目的是最大化目標函數(shù)，但是當目標函數(shù)達到最大值時每一個源信號都有被輸出的可能性，所以要提取感興趣信號可以將感興趣信號的先驗信息作為約束加入算法中，則CICA的框架定義為

maxJ(y)=[E{G(wTx)}-E{G(v)}]2

(5)

s.t.

(6)

式中：s.t.為約束條件，包括g(w)和h(w)兩種約束；ε(y,r)表示提取信號與參考信號的測量度，ε(y,r)=E{(y-r)2}；ξ為設置的閾值，閾值的設置應該注意，如果閾值的選擇過大可能會收斂到其他的源信號，相反如果閾值設置過小則收斂不到任何信號，因此閾值的設置需要謹慎。等式約束h(w)來確保目標函數(shù)J(y)和分離矩陣w有上界。

CICA算法中的參考信號應該具備以下兩個條件。

(1)具備感興趣信號的特性。這是參考信號存在的重要原因，例如，在提取胎兒心電信號時，參考信號只需要反映胎兒心電信號的某些特征即可，不用是全部特征。

(2)參考信號是非高斯性信號或者可以轉換為非高斯性信號[14]。ICA算法以數(shù)據(jù)的非高斯性為基礎，將目標函數(shù)的非高斯性進行最大化來尋找最優(yōu)分離矩陣[15]。CICA算法是在ICA算法上發(fā)展起來的，因此CICA算法可以看成是對一個以負熵為目標函數(shù)進行約束優(yōu)化的問題，然后將對比函數(shù)的非高斯性進行最大化來尋找最優(yōu)分離矩陣。在ICA算法中加入?yún)⒖夹盘柺且龑惴ǖ氖諗糠较?。因此參考信號是非高斯性信號或者可以轉換為非高斯性信號才會改變算法的收斂方向。

參考信號的設計有許多種，在進行實驗時脈沖法是使用最多的一種方法。脈沖法是將感興趣信號的先驗信息的脈沖信號作為參考信號，這也就說明了脈沖法滿足第一個條件。此外，脈沖信號是方波序列本身具有非高斯性，因此滿足第二個條件。感興趣信號與參考信號之間的測量度越小，那么分離結果就越精確。源信號與它自身的接近度是最大的，但往往源信號是未知的，因此只能利用其先驗信息近似的設計參考信號。

3 改進的CICA

(7)

s.t.

(8)

通過引入松弛因子z可以將不等式約束g(w)≤0轉化為等式g(w)+z2=0。利用拉格朗日乘子法來解決上述的約束優(yōu)化問題，該框架下的增廣拉格朗日函數(shù)L(w,μ,λ)為

(9)

利用類牛頓算法可以得到分離矩陣w的學習算法為

w←w-η(L2)-1Γ1

(10)

式(10)中：η為學習率；L2為Hessian矩陣；Γ1為L(w,μ,λ)對w的一階導數(shù)。為了對算法的進一步簡化，可以將L2近似為Γ2RXX-1，那么近似后的分離矩陣w學習算法可以定義為

(11)

式(10)中：RXX為觀測信號的協(xié)方差矩陣；Γ1、Γ2分別為L(w,μ,λ)對w的一階導數(shù)和二階導數(shù)：

Γ1=2aρE{x(t)g′(y)}-(2b+μ)E{x(t)(y-

r)}-λE{x(t)y}

(12)

Γ2=2aρE{g″(y)}-2b-μ-λ

(13)

式中：a=1/[ε(y,r)+1]，b=[E{G(y)}-E{G(v)}]2/[ε(y,r)+1]2，G(y)為任意的非線性函數(shù)，取值G(y)=lg[cosh(y)]；g′(y)為G(y)的一階導數(shù)，g′(y)=tanh(y)；g″(y)表示G(y)的二階導數(shù)，g″(y)=1-tanh2(y)；ρ=E{G(y)}-E{G(v)}。

基于梯度下降算法可以得到拉格朗日乘子的迭代公式為

μ=max{0,μ+γg(w)}

(14)

λ=λ+γh(w)

(15)

對改進的帶參考信號盲源分離算法步驟如下。

(1)對觀測信號進行去均值、去白化得到處理后的數(shù)據(jù)。

(2)根據(jù)感興趣信號的先驗知識來設計參考信號r。

(3)用均勻分布隨機數(shù)初始化分離矩陣w并對w進行標準化。

(4)初始化拉格朗日乘子λ和μ。

4 實驗分析

4.1 仿真實驗

為展現(xiàn)算法的有效性，對改進算法進行仿真實驗，選取5個相互獨立的信號作為實驗的源信號，源信號的函數(shù)為

(16)

式(16)中：f1=50 Hz；f2=610 Hz；f3=540 Hz；s1、s2和s3為3個周期信號。s4為均值是0標準差是1的標準正態(tài)分布隨機數(shù)；s5為[-1，1]的均勻分布隨機數(shù)；s4和s5統(tǒng)稱為噪聲信號；N表示噪聲函數(shù)的取樣長度；k的取值為1～N。本實驗的采樣點為500(201～700)。在提取期望源信號時所運用到的參數(shù)設置為：拉格朗日乘子λ、μ和懲罰因子γ的初始值為1，閾值ξ為1.75，學習率η=1，最大的迭代步長為200。

通過兩組仿真實驗(提取源信號s1和提取源信號s2)對本文提出的改進算法與文獻[16]提出的算法進行性能比較。利用信噪比和均方誤差來衡量其性能，衡量標準如下。

(1)信噪比。若源信號與其對應的分離信號分別為si和yi，則兩者之間的信噪比(signal-noise ratio, SNR)表示為

(17)

根據(jù)信噪比的定義可以判斷SNR越大表明信號功率越大，說明算法的分離性能越好。

(2)均方誤差。均方誤差(mean-square error, MSE)是衡量兩組數(shù)據(jù)的偏差程度，yi為分離信號，si為源信號，N為采樣點數(shù)，衡量標準為

(18)

源信號s1到源信號s5如圖1所示。

圖1 源信號Fig.1 Source signal

將5個獨立的源信號進行混合，混合矩陣A為[-1，1]的均勻分布隨機數(shù)，混合后的信號如圖2所示。

圖2 混合信號Fig.2 Mixed signal

本組仿真試驗運用改進的算法來提取源信號s2，設計參考信號時采用第2節(jié)所提到的脈沖法。一般情況下，感興趣信號的周期、初始相位和脈沖寬度是脈沖法所需要的先驗信息。根據(jù)源信號s2的先驗信息設計參考信號r1周期為1/f2，初始相位設置為1，脈沖的寬度為8，參考信號r1和提取信號y1如圖3(a)所示。

為了驗證參考信號對算法分離性能的影響，可以改變參考信號的某些值以觀其效果，如改變參考信號的初始相位或者脈沖的寬度，設計參考信號r2周期為1/f2，初始相位設置為3，脈沖的寬度為8，參考信號r3周期為1/f2，初始相位設置為1，脈沖的寬度為10，參考信號r2、r3和提取信號y2、y3如圖3(b)和圖3(c)所示。

由圖3可以看到，利用參考信號r2和參考信號r3分離出的信號和真實源信號都有所不同?？梢钥闯龈淖儏⒖夹盘柕某跏枷辔缓兔}沖寬度會影響算法的分離效果，只有根據(jù)源信號的先驗信息設計的參考信號才能正確提取感興趣信號。

本組仿真實驗進行1 000次最后求均值。文獻[16]改進的算法和本文改進算法性能比較如表1所示。

由表1可知，本文改進算法在均方誤差和信噪比上都有提高，可以進一步說明提出的算法有效。

為了試驗的完備性，用改進的算法分離源信號s1，參考信號的設計同樣采用第二節(jié)所提到的脈沖法，則參考信號周期為1/f1，初始相位設置為2，脈沖的寬度為100。r1為參考信號，y1為提取信號，s1為源信號，提取效果如圖4所示。

圖3 提取信號s2效果圖Fig.3 The effect picture of extract signal s2

表1 性能評估

同樣，為了衡量改進算法分離源信號s1的誤差，利用信噪比和均方誤差來計算，本組仿真實驗進行1 000次最后求均值，文獻[16]改進的算法和本文改進算法性能比較如表2所示。

由表2可知，本文改進的算法在均方誤差和信噪比上都稍有提高，進一步的說明該算法的有效性。

圖4 提取信號s1效果圖Fig.4 The effect picture of extract signal s1

表2 性能評估

4.2 滾動軸承故障診斷實驗

以滾動軸承故障診斷為例，采用改進的CICA方法對其驗證，利用如圖5所示的試驗臺對滾動軸承外圈單一故障進行試驗，試驗軸承安裝在主軸上，實驗機主軸與變頻電機用聯(lián)軸器連接，通過調節(jié)變頻電機轉速來模擬主軸旋轉。

圖5 試驗臺Fig.5 Test bench

試驗軸承型號為D276126NQ1U的雙半內圈三點接觸球軸承，其軸承外圈內徑為142.9 mm，外徑為190.0 mm，寬度為33 mm，球個數(shù)為10，球徑為24.6 mm，節(jié)徑為166.45 mm，采樣頻率為20 000 Hz，采樣時間為10 s。根據(jù)軸承外圈故障頻率計算公式在1 000 r/min轉速下的外圈故障頻率為127 Hz。

傳感器采集的外圈混合信號如圖6所示，從圖6中可以看到采集到的信號雜亂無章，看不到任何有效的故障特征。根據(jù)滾動軸承的理論外圈故障頻率(f=127 Hz)，現(xiàn)產生一方波序列作為參考信號。將參考信號和觀測信號同時作用于改進的CICA算法中，可得到提取信號。參考信號r與提取信號y如圖7所示。

為了驗證提取信號的正確性，對提取的故障信號進行Hilbert包絡譜分析，結果如圖8所示，從圖8可以看到，在125～130 Hz處幅值較高，計算的軸承外圈故障特征頻率127 Hz在其區(qū)間，因此可以說明CICA算法準確地提取了故障特征。

圖6 外圈故障采集信號Fig.6 Outer ring fault acquisition signal

圖7 參考信號r與提取信號yFig.7 Reference signalrand extracted signaly

圖8 提取信號包絡圖Fig.8 Extract signal envelope

5 結論

為減小CICA算法中源信號和提取信號之間的誤差，提出了一種新的CICA算法。該算法在標準對比函數(shù)中，耦合了含有先驗信息的測量度函數(shù)，以此得到了新的目標函數(shù)；然后在新目標函數(shù)上采用拉格朗日乘子法，最終有效地獲得了最優(yōu)分離矩陣。在仿真實驗中，與其他算法相比，該算法具有較小的誤差。在滾動軸承故障診斷實驗中，該算法正確提取了滾動軸承外圈故障特征，進一步說明了該算法的有效性。