付永慶，郭慧，，蘇東林，劉焱

(1.哈爾濱工程大學 信息與通信工程學院，黑龍江哈爾濱150001;2.北京航空航天大學電磁兼容技術研究所，北京100191)

盲源分離(blind signal separation，BSS)起源于雞尾酒問題［1］，可在源信號和傳輸信道參數(shù)均未知的情況下僅根據(jù)源信號的統(tǒng)計獨立性和觀測信號來恢復源信號，被廣泛應用于生物醫(yī)學［2］、地質(zhì)信號處理和分析［3］、圖像處理［4］和移動通信［5］等領域。根據(jù)觀測信號和源信號的數(shù)目關系，將盲源分離問題分為正定、超定和欠定3種情況。上述情況都以信源數(shù)目的正確估計為前提實現(xiàn)信號的分離。目前，在信源數(shù)目估計領域中，以前兩種情況下的盲源分離為主，基于高階累積量擴展了很多成熟的算法，而欠定盲源分離(undetermined blind signal separation，UBSS)方面的研究尚少。目前欠定盲源分離中信源數(shù)目的估計大多基于稀疏混合信號的線性聚類特點，如:2001年Bofill利用勢函數(shù)的局部最大值的個數(shù)對信源數(shù)目進行估計［6］;2008年譚北海等直接利用概率統(tǒng)計知識，通過劃分區(qū)域并統(tǒng)計觀測點在不同區(qū)域的概率分布實現(xiàn)信源數(shù)目估計［7］;2009年張燁等提出采用拉普拉斯勢函數(shù)判斷其局部最大值達到源數(shù)估計的目的［8］。但這些方法都存在一些問題，如算法的抗噪聲性和信號的稀疏性易受異常值的影響從而造成誤判。針對上述問題，提出一種Hough加窗法，利用信號的稀疏特性，先確定經(jīng)平滑后的聚類區(qū)域，減小異常值的影響，再尋找區(qū)域最大值，避免陷入局部最大，在此基礎上對混合矩陣進行估計。針對信號分離問題，提出一種無約束分離算法，采用“內(nèi)點法”選擇合適的混合矩陣列向量作為初始迭代值，修正了現(xiàn)有算法對混合矩陣列向量的選取標準。

1 稀疏欠定盲源分離模型

通常欠定盲源分離的線性瞬時混合模型寫作:

式中:X(t)為m個觀測信號矢量，A為m×n維的混合矩陣，S(t)為n個源信號矢量，t=1，2，…，T為觀測點時刻，當m＜n時即為欠定盲源分離。

為了實現(xiàn)欠定情況下的盲源分離，必須對源信號提出進一步的約束條件，假設源信號是稀疏信號，即源信號在絕大多數(shù)采樣點的取值為零或接近于零，只有少數(shù)采樣點的取值遠離零［9］。假設在采樣時刻t，只有s1(t)起主導作用，則有

觀測信號可以看作:在以x1(t)…xm(t)為坐標軸的m維空間中，所有以源信號s1(t)為主導的采樣時刻確定的一條直線，其斜率取決于混合矩陣A的第一個列向量 [a11，…，am1]T。如果混合矩陣A的任意m×m的子矩陣都滿足可逆性，則每個源信號可確定一條直線。因此，源信號的數(shù)目估計就轉(zhuǎn)化為對觀測空間中直線的條數(shù)估計。

實際應用中，很多信號具有稀疏特性，或者可以通過適當變換(如短時傅里葉變換，STFT)［10］、Gabor變換［11］、小波變換［12］等)使得信號在變換域中滿足稀疏特性。

2 Hough加窗法

2.1 Hough 變換

Hough變換原用于圖像的匹配，將原圖中給定形狀的所有點都集中到變換空間中的某些位置，形成峰值點。基于該思想，把觀測信號空間中直線的檢測問題轉(zhuǎn)化為尋找變換空間中峰值點的問題。

觀測數(shù)據(jù)可看作m維空間中經(jīng)過零點的一條直線［13］，據(jù)此構(gòu)造轉(zhuǎn)換關系，將觀測信號空間中的直線轉(zhuǎn)換成變換空間中的角度變換量［14］:

式中:k=1，…，π/h，表示高維直方圖的分區(qū)數(shù)目，h是Hough變換的量化步長，k越大，即h越小，變換量的分類精度越高;在變換量的每一個分類中，i=1，…，m－1，包含了變換量的所有維數(shù)，[·]表示比·小的最大整數(shù)。則落入高維直方圖某一區(qū)域的變換量數(shù)目為:

2.2 加窗法

得到變換量的直方圖后，若直接檢測峰值點容易陷入局部最大值。為此對直方圖做加窗處理。由于m維的觀測信號經(jīng)Hough變換后變?yōu)閙－1維變換量，故對其加m－1維的超方體窗。加窗的目的是平滑直方圖中的異常值，得到聚類區(qū)域。加窗后重新計算變換量的分類統(tǒng)計數(shù)為

式中:j=0，…，[π/hd]－[L/d]+1表示重新計算后的變換量分類數(shù)，L表示超方體窗口的邊長，d表示超方體窗口沿觀測空間每個坐標軸方向的移動步長表示權重，由落入直方圖某一區(qū)域中的變換量數(shù)占變換量總數(shù)的比率獲得。

對新的直方圖搜索峰值即可得到源數(shù)的估計。峰值對應的角度變換估計量與混合矩陣元素有式(7)的關系，由此得到混合矩陣的估計值:

式中:ami表示混合矩陣第m行i列的元素。

3 無約束分離算法

3.1 無約束問題

文獻［2］中分離混合信號有如下的約束條件:

由于式(1)可寫為

式中:Am稱作混合矩陣Am的滿秩因子，Amn－m稱作混合矩陣Am的余式因子。又有

式中:sign(·)表示取·的符號。因此，可以將式(8)帶約束的優(yōu)化問題轉(zhuǎn)化為上述無約束優(yōu)化問題，要滿足，令J(Sn－m)對Sn－m求導，即

3.2 初始迭代值選取

初始迭代值的選取直接影響算法的迭代次數(shù)和收斂速度，因此，需要合理選擇初始值?，F(xiàn)有算法通過選擇觀測點與混合矩陣列向量夾角最小的前m列作為分離矩陣，該方法并不適用于任何情況。

圖1 二維矢量分解示意圖Fig.1 Schematic diagram of 2D vector decomposition

如圖1所示，以兩個混合信號為例，在單位化對稱化后的散點圖中，混合信號x與列向量a1、a2的夾角最小，但是其分解矢量的長度和遠大于混合信號x在列向量a2、a3上的分解矢量。這是因為當兩個列向量分布于混合信號的同一側(cè)時，有一個分解矢量必在其中一列向量的反方向上，此時，分解矢量與觀測信號的夾角并不是最小。

根據(jù)散點圖分布采用“內(nèi)點法”選取初始迭代值。首先，令混合矩陣單位化和對稱化，且使混合矩陣第一行和混合信號第一行數(shù)據(jù)符號為正。然后，去處混合矩陣和混合信號第一行數(shù)據(jù)后繪制散點圖，以混合信號為原點，混合矩陣列向量分布在其四周。從n個點中選取m個點，這m個點圍成的多邊形應能包含混合信號點。從滿足要求的組合中選擇與混合信號夾角和最小的作為滿秩因子Am，否則按照最小夾角法選取Am，其余列向量作為余式因子An－m，Sm則由式(10)得到，令Sn－m=0。

源信號越稀疏，J'(Sn－m)值越小，在迭代過程中，分離信號逐漸逼近源信號S，稀疏性增強，J'(Sn－m)逐漸減小，當J'(S(k+1)n－m≥J'()時迭代終止，第k次迭代計算得到的即為分離信號。

4 仿真實驗

4.1 算法性能評價標準

4.1.1 泛化交擾誤差

采用泛化交擾誤差(generalized crosstalking error，GCE)評價混合矩陣的估計精度。定義如下:

式中:Π表示所有n×n維可逆矩陣組成的集合，表示與一個尺度矩陣和置換矩陣的乘積，以消除的幅度不確定性和排序不確定性。當且僅當和A完全等價時，Err(A，)=0。

4.1.2 信干比

信號分離后第i個源信號si和分離信號的信干比(signal to interference radio，SIR)定義為:

盲源分離的幅度具有不確定性，在計算信干比前需先確定源信號si和分離信號的幅度因子λi:

通過式(15)使得si和的幅度盡量一致。在該度量下，信干比SIR越大，說明分離信號越接近源信號。

4.1.3 互相關系數(shù)

互相關系數(shù)(cross correlation coefficient，CCC)描述兩個信號的相似性，定義為

式中:sj為源信號，為還原信號。在該度量下，ξij=1和0分別表示與sj完全相似和不相似。

在上述3個標準中，泛化交擾誤差和互相關系數(shù)都是與特定值做比較(GCE值與0比較，CCC值與1比較)，可以體現(xiàn)單獨算法的有效性，但因其計算值很小，在多算法性能比較時不利于突出算法的差異，因此采用信干比來衡量不同算法的性能。

4.2 實驗及結(jié)果分析

4.2.1 實驗一:驗證Hough加窗法的性能

為了驗證信號稀疏度對算法的影響，需構(gòu)建稀疏度可度量的稀疏信號。滿足廣義Gaussian分布的信號，其概率密度函數(shù)為:

據(jù)此構(gòu)建信噪比為40 dB，α＜1.5的稀疏度不同的廣義Gaussian分布的隨機信號和α=1，信噪比不同的信號，采樣值T=5 000，隨機選取2×6的混合矩陣得到混合信號，分別進行100次蒙特卡洛運算，分析算法的正確率，如圖2所示。

圖2 源數(shù)估計的正確性Fig.2 Accuracy of the source number estimation

Hough加窗法在稀疏性度量值為1.2及以上時，算法的正確率較低，在稀疏性度量值小于1.1和信噪比大于2 dB時算法正確率分別可達81%和85%以上，在稀疏性度量值小于0.8和信噪比大于15 dB時算法正確率達到100%?？梢娝惴▽π盘柕南∈杳舾行暂^低，且具有較好的抗噪聲性能，即使在低信噪比下也能以較高的正確率估計出信源數(shù)目。

源數(shù)估計的高正確率保證了算法的穩(wěn)定性。選取源數(shù)估計正確的情況估計混合矩陣，計算泛化交擾誤差取平均值，并與勢函數(shù)法、K均值法進行比較，得到不同算法下噪聲和稀疏度對混合矩陣估計性能的影響圖，如圖3和圖4所示。

從兩幅圖中看出，無論是在低信噪比下還是低稀疏度下，Hough加窗法得到的混合矩陣的泛化交擾誤差都很小，趨近于0，尤其是在稀疏度為1.2以下和信噪比為8 dB以上時，混合矩陣估計精度更高，性能明顯優(yōu)于其他算法。勢函數(shù)法直接對累計量進行峰值搜索，無法平滑奇異值對峰值檢測的影響，因此性能不如Hough加窗法，K均值法由于受初始點選取的影響較大，故估計性能最差。

圖3 噪聲對混合矩陣估計精度的影響Fig.3 The influence of noiseon the accuracy of the estimation of mixing matrix

圖4 稀疏度對混合矩陣估計精度的影響Fig.4 The influence of sparse on the accuracy of the estimation of mixing matrix

4.2.2 實驗二:驗證無約束分離算法的性能

在實驗一得到的混合信號基礎上，假設混合矩陣已知，信噪比為40 dB，α=1的情況下，通過互相關系數(shù)和信干比兩個指標來評價算法的性能，如圖5所示，并從分離信號的信干比、時間頻度、復雜度和運行時間等方面，將本文算法與最短路徑法、最小l1范數(shù)法進行比較，結(jié)果如表1所示。

圖5中信噪比為40 dB時，互相關系數(shù)平均值達到0.991 8，對應的信干比平均值為169.65。整體來說，在信噪比較低的情況下，互相關系數(shù)可達到0.74以上，在信噪比為15 dB以上時，互相關系數(shù)達到0.9以上。在2～40 dB的范圍內(nèi)，信干比的數(shù)值波動在37 dB左右，約占最佳信干比值的21.9%。表明無約束分離算法自身具有分離準確性較高，抗噪聲性能較好的特點。

表1給出了無約束分離算法與其他2個有約束算法的性能比較結(jié)果。其中時間頻度以m個混合信號，n個源信號，采樣值為T的情況進行計算，k表示無約束分離算法的迭代次數(shù)。3種分離算法得到的源信號信干比值均一樣，說明了本文算法的正確性。從時間頻度上看，最短路徑法值最小，最小l1范數(shù)法最大。本文算法與最小l1范數(shù)法相差一個系數(shù)，本文算法的k取值一般為1～3，最小l1范數(shù)法在本實驗中取值為15，運行時間也證明了上述分析。在運算時間上，本文算法不及勢函數(shù)法，但適用范圍要大于勢函數(shù)法，因為勢函數(shù)法只適用于兩個混合信號的情況。因此，綜合運行時間和適用性可以得出:無約束分離算法具有較好的性能。

表1 算法比較結(jié)果(α=1，SNR=40 dB)Table 1 The comparison results of different algorithms(α =1，SNR=40 dB)

圖5 混合分離的分離性能Fig.5 Performance of signal separation

為了驗證“內(nèi)點法”修正了最小夾角法對混合矩陣列向量的選取標準，將無約束分離算法與基于最小夾角法的分離算法進行比較，結(jié)果如圖6所示。

圖6 內(nèi)點法與最小夾角法的性能比較Fig.6 Performance comparison of signal separation between interiorpoint method and minimum angle method

在混合矩陣選取的過程中，不同的混合信號出現(xiàn)列向量散點包含混合信號散點的情況是不確定的，因此，針對不同的混合信號而言，內(nèi)點法較最小夾角法的性能優(yōu)勢表現(xiàn)的也不盡相同。在圖6中，內(nèi)點法的優(yōu)勢較最小夾角法比較突出，且在信噪比較高時，這種優(yōu)勢更明顯。

4 結(jié)論

1)針對現(xiàn)有混合矩陣估計算法易受噪聲干擾，對信號稀疏性要求較高等問題，提出一種Hough加窗算法，利用稀疏信號的特性，用加窗的方法平滑因噪聲或稀疏性差造成的異常值，對得到的聚類區(qū)域分別尋找最大值，達到全局最大值搜索的目的，峰值的個數(shù)即為信源的數(shù)目。

2)針對現(xiàn)有分離算法帶約束的優(yōu)化問題，提出一種無約束的分離方法，并采用“內(nèi)點法”來選擇合適的初始迭代值，修正了現(xiàn)有算法對混合矩陣列向量的選取標準。

仿真結(jié)果不但驗證了算法的有效性，還表明了算法具有較好的抗噪聲性、較低的稀疏敏感性。若信號不滿足稀疏特性，可對其采用適當?shù)淖儞Q使其在變換域中具有稀疏性，該算法在變換域中仍適用。

［1］COMON P.Independent component analysis:a new concept［J］.Signal Processing，Special Issue on Higher Order Statistics，1994，36(3):287-314.

［2］MAMMONE N，La FORESTA F，MORABITO F C.Automatic artifact rejection from multichannel scalp EEG by wavelet ICA［J］.IEEE Sensors Journal，2012，12(3):533-542.

［3］TAKAHATA A K，NADALIN E Z，F(xiàn)ERRARI R，et al.Unsupervised processing of geophysical signals:a review of some key aspects of blind deconvolution and blind source separation［J］.IEEE Signal Processing Magazine，2012，29(4):27-35.

［4］BADAWI W K M，CHIBELUSHI C C，PATWARY M N，et al.Specular-based illumination estimation using blind signal separation techniques［J］.IET Image Processing，2012，6(8):1181-1191.

［5］YANG Xiaomei.The application of blind signal separation in mobile communication［C］//2011 3rd International Conference on Computer Research and Development(ICCRD).Shanghai，China，2011:274-277.

［6］BOFILL P，ZIBULEVSKY M.Underdetermined blind source separation using sparse representations［J］.Signal Processing，2001，81(11):2353-2362.

［7］譚北海，謝勝利.基于源信號數(shù)目估計的欠定盲分離［J］.電子與信息學報，2008，30(4):863-867.TAN Beihai，XIE Shengli.Underdetermined blind separation based on source signals’number estimation［J］.Jour-nal of Electronics＆ Information Technology，2008，30(4):863-867.

［8］張燁，方勇.基于拉普拉斯勢函數(shù)的欠定盲分離中源數(shù)的估計［J］.信號處理，2009，25(11):1719-1725.ZHANG Ye，F(xiàn)ANG Yong.A new method to estimate the number of the sources for underdetermined blind separation based on lapulacial potential function［J］.Signal Processing，2009，25(11):1719-1725.

［9］邱天爽，畢曉輝.稀疏分量分析在欠定盲源分離問題中的研究進展及應用［J］.信號處理，2008，24(6):966-970.QIU Tianshuang，BI Xiaohui.Spares component analysis and application for underdetemined blind source separation［J］.Signal Processing，2008，24(6):966-970.

［10］BOFILL P，ZIBULEVSKY M.Blind separation of more sources than mixtures using sparsity of their short-time fourier transform［C］//Second International Workshop on Independent Component Analysis and Blind Signal Separation.Espoo，F(xiàn)inland，2000:87-92.

［11］CORDERO E，NICOLA F，RODINO L.Sparsity of gabor representation of schrodinger propagators［J］.Applied and Computational Harmonic Analysis，2009，26(3):357-370.

［12］MA Jianwei，PLONKA G，CHAURIS H.A new sparse representation of seismic data using adaptive easy-path wavelet transform［J］.IEEE Geoscience and Remote Sensing Letters，2010，7(3):540-544.

［13］HIROYUKI S，YUZO H.Blind source separation by a geometrical method［C］//Proceedings of the International Joint Conference on Neural Networks.Honolulu，USA，2002:1108-1114.

［14］付寧，彭喜元.K-Hough欠定盲信道估計算法［J］.電子測量與儀器學報，2008，22(5):63-67.FU Ning，PENG Xiyuan.K-Hough underdetermined blind mixing model recovery algorithm［J］.Journal of Electronic Measurement and Instrument，2008，22(5):63-67.

［15］何昭水，謝勝利，傅予力.信號的稀疏性分析［J］.自然科學進展，2006，16(9):1167-1173.HE Shaoshui，XIE Shengli，F(xiàn)U Yuli.Sparse analysis of signal［J］.Progress in Natural Science，2006，16(9):1167-1173.