毛 玲, 趙聯(lián)文, 孟 華, 李雨鍇

(西南交通大學(xué) 數(shù)學(xué)學(xué)院,四川 成都 610031)

盲源分離是在未知源信號(hào)的先驗(yàn)信息及信道參數(shù)的情況下,僅利用傳感器接收到的觀測(cè)信號(hào)對(duì)源信號(hào)進(jìn)行恢復(fù)的過(guò)程[1-3]。由Hyv?rinen等[4-5]提出的快速獨(dú)立成分分析(fast independent component analysis, FastICA)是目前最常用的盲源分離算法。FastICA的使用前提是觀測(cè)信號(hào)個(gè)數(shù)不小于源信號(hào)個(gè)數(shù),并且需要已知源信號(hào)個(gè)數(shù)。而在實(shí)際工程應(yīng)用中,由于受到設(shè)備器材、配置成本等條件的限制,僅能接收到單通道信號(hào)。只有一個(gè)觀測(cè)信號(hào)時(shí)無(wú)法直接使用FastICA進(jìn)行盲源分離,處理此問(wèn)題常用的方法是首先將單通道信號(hào)升維,構(gòu)建虛擬多通道信號(hào),再通過(guò)FastICA對(duì)虛擬多通道信號(hào)進(jìn)行源信號(hào)的恢復(fù)。

由于使用FastICA對(duì)源信號(hào)進(jìn)行恢復(fù)必須正確給出源信號(hào)個(gè)數(shù),在源信號(hào)個(gè)數(shù)未知時(shí)算法失效,因此如何有效估計(jì)盲源信號(hào)個(gè)數(shù)是盲源分離的關(guān)鍵點(diǎn),也是盲源分離效果好壞的關(guān)鍵前提。目前比較常用的源數(shù)估計(jì)算法包括兩種:一種為基于矩陣分析的源數(shù)估計(jì)算法,如主成分分析(principle component analysis, PCA)[6]、奇異值分解(singular value decomposition, SVD)[7]等;另一種為基于信息論的源數(shù)估計(jì)算法,如Akaike[8]提出的赤池信息準(zhǔn)則(Akaike information criterion, AIC)、Rissanen[9]提出的最小描述長(zhǎng)度準(zhǔn)則(minimum description length, MDL)等。雖然這些算法在實(shí)際應(yīng)用中被廣泛使用,但AIC算法得到的估計(jì)結(jié)果不具有一致性,MDL算法在低信噪比條件下誤差率較高[10]。由于在實(shí)際應(yīng)用中,采集到的信號(hào)通常會(huì)受到多種噪音和干擾的影響,因此提出一種有較強(qiáng)抗噪性的源數(shù)估計(jì)算法具有重要意義。

為此,本文提出了一種新的基于信息熵指標(biāo)的源數(shù)估計(jì)算法,無(wú)須知道源信號(hào)的先驗(yàn)信息,直接通過(guò)計(jì)算分離得到的估計(jì)源信號(hào)的信息熵來(lái)度量源信號(hào)的信息量大小,從而確定源信號(hào)的個(gè)數(shù)。本文使用總體經(jīng)驗(yàn)?zāi)B(tài)分解(ensemble empirical mode decomposition, EEMD)[11]對(duì)單通道觀測(cè)數(shù)據(jù)進(jìn)行升維處理,再使用FastICA算法分解得到估計(jì)源信號(hào)。計(jì)算估計(jì)源信號(hào)的信息熵,從而確定源信號(hào)個(gè)數(shù)。

1 基于信息熵的源數(shù)估計(jì)

1.1 信息熵基本理論

1948年Shannon[12]在信息科學(xué)領(lǐng)域提出信息熵(entropy)的概念,來(lái)描述信源的不確定性。不確定性越大,系統(tǒng)包含的混合信息源的信息量越大。

隨機(jī)變量的連續(xù)熵定義如式(1)所示:

(1)

式中:f(x)為隨機(jī)變量X的概率密度。

隨機(jī)變量X和Y的聯(lián)合熵定義如式(2)所示:

(2)

式中:f(x,y)為隨機(jī)變量X和Y的聯(lián)合概率密度。

引理1若變量X和Y是相互獨(dú)立的,則聯(lián)合熵等于獨(dú)立熵之和[13]:

H(X,Y)=H(X)+H(Y)。

(3)

推論1聯(lián)合熵不小于其中任意變量的熵:

H(X,Y)≥max(H(X),H(Y))。

(4)

本文使用信息熵來(lái)估計(jì)源數(shù),理論分析如下。

假設(shè)源信號(hào)個(gè)數(shù)為n,第i個(gè)源信號(hào)的概率密度函數(shù)為fi(x),根據(jù)式(1)得到第i個(gè)源信號(hào)的信息熵為

(5)

由于源信號(hào)之間彼此獨(dú)立,根據(jù)引理1得到n個(gè)源信號(hào)的平均信息熵如式(6)所示:

(6)

根據(jù)推論1得

(7)

由FastICA分解得到的各源信號(hào)是相互獨(dú)立的,獨(dú)立信號(hào)混合到一起,其聯(lián)合熵比較大。在引理1和推論1的基礎(chǔ)上如果能分解出相對(duì)純凈的獨(dú)立信號(hào),則平均信息熵會(huì)降低。因此迭代地增加盲源估計(jì)個(gè)數(shù),隨著個(gè)數(shù)的增加,分離的信號(hào)會(huì)越來(lái)越純凈,進(jìn)而信號(hào)平均信息熵越來(lái)越低,當(dāng)過(guò)分解時(shí)(估計(jì)源數(shù)大于真實(shí)源數(shù)),分離結(jié)果會(huì)出現(xiàn)雜亂的噪聲信號(hào),熵值又會(huì)呈現(xiàn)增加的趨勢(shì)。故當(dāng)源數(shù)估計(jì)正確時(shí),分離源信號(hào)不含有其他源的信息,估計(jì)源信號(hào)的平均信息熵會(huì)達(dá)到最小,由此得到最小熵值對(duì)應(yīng)的個(gè)數(shù)即為最佳源數(shù),如式(8)所示:

(8)

基于以上分析,本文提出基于平均信息熵最小的盲源個(gè)數(shù)估計(jì)算法。

1.2 基于GMM的源信號(hào)概率密度擬合

為了計(jì)算源信號(hào)的信息熵,需要已知信號(hào)的概率密度函數(shù)。由于信號(hào)的概率密度未知且較為復(fù)雜,無(wú)法得到其精確的表達(dá)式,根據(jù)統(tǒng)計(jì)理論,任何一個(gè)隨機(jī)變量的分布都可以用高斯混合模型(Gaussian mixed model, GMM)去逼近[14]。因此本文使用GMM來(lái)擬合估計(jì)源信號(hào)的概率密度函數(shù)。

GMM數(shù)學(xué)解析式如式(9)所示:

(9)

式中:K為GMM模型中子模型的個(gè)數(shù);αk為變量x屬于第k個(gè)子模型的概率,

為第k個(gè)子模型的參數(shù);φ(x|θk)為第k個(gè)子模型的概率密度函數(shù),如式(10)所示:

(10)

式中各參數(shù)通常使用期望最大化算法(expectation maximum, EM)[15]迭代計(jì)算得到。

通過(guò)GMM得到源信號(hào)的概率密度,則信息熵進(jìn)一步可以寫(xiě)成數(shù)學(xué)期望如式(11)所示:

=E(-lnf(x|θ))。

(11)

1.3 基于MCMC的信號(hào)熵值計(jì)算

在1.2節(jié)中,通過(guò)GMM可以得到源信號(hào)的概率密度。由于所求信息熵的積分函數(shù)是一個(gè)復(fù)雜的混合模型,無(wú)法求得式(8)的精確值。為了解決復(fù)雜積分求解問(wèn)題,本文使用馬爾可夫鏈蒙特卡羅(Markov chain Monte Carlo, MCMC)算法[16]近似計(jì)算其期望。MCMC的理論基礎(chǔ)為大數(shù)定律:當(dāng)樣本容量足夠多時(shí),樣本均值以概率1收斂于數(shù)學(xué)期望。本文應(yīng)用此定律來(lái)近似計(jì)算期望。根據(jù)MCMC算法,首先生成1組服從目標(biāo)分布的隨機(jī)樣本,然后通過(guò)計(jì)算其樣本均值代入式(11)進(jìn)行估算。

對(duì)于用于計(jì)算積分的隨機(jī)樣本,本文使用Metropolis-Hastings(M-H)算法進(jìn)行采樣[17-18],得到服從目標(biāo)分布的樣本點(diǎn)。

M-H采樣算法步驟如下。

輸入:抽樣的目標(biāo)分布的密度函數(shù)f(x|θ)、建議分布函數(shù)q(x,y);

輸出:f(x|θ)的隨機(jī)樣本{xm+1,xm+2,…,xn};

參數(shù):收斂步數(shù)m、迭代步數(shù)n。

① 隨機(jī)初始化x0;

② Fori=1,2,…,n;

從建議分布q(x,xi-1)隨機(jī)生成樣本xi

計(jì)算接受概率:

從均勻分布U(0,1)中生成樣本u

ifα(xi,xi-1)≤u

return (xi-1)

else

return (xi)

得到樣本{xm+1,xm+2,…,xn}。為了保證樣本的收斂性,刪去前m個(gè)樣本,得到目標(biāo)分布的樣本{xm+1,xm+2,…,xn}。

③ end

基于M-H采樣算法得到服從目標(biāo)分布f(x|θ)的隨機(jī)樣本{xm+1,xm+2,…,xn},將這些樣本用來(lái)計(jì)算函數(shù)-lnf(x|θ)的均值,如式(12)所示:

(12)

1.4 基于信源信息熵最小的單通道盲源數(shù)估計(jì)算法

輸入:單通道觀測(cè)信號(hào)數(shù)據(jù)x、分解的源信號(hào)的最大維數(shù)m、GMM模型高斯分布的個(gè)數(shù)K、M-H采樣的樣本個(gè)數(shù)l;

① 使用EEMD將單通道觀測(cè)信號(hào)x分解重構(gòu)成虛擬多通道信號(hào)X;

② Forn=2,3,…,m;

Forp=1,2,…,n

用GMM擬合源信號(hào)Yn的概率密度函數(shù):

使用采樣得到的樣本計(jì)算平均信息熵:

end

end

2 仿真實(shí)驗(yàn)

參考文獻(xiàn)[19]中的仿真數(shù)據(jù)實(shí)驗(yàn),通過(guò)R軟件生成正弦信號(hào)s1(t)、余弦信號(hào)s2(t)、幅值調(diào)制信號(hào)s3(t)這3個(gè)源信號(hào),信號(hào)的表達(dá)式如下式所示:

(13)

式中各信號(hào)的參數(shù)如下:f1=50 Hz、f2=25 Hz、f3=100 Hz、f4=15 Hz,采樣時(shí)間為0.5 s,采樣頻率為1 024 Hz,各源信號(hào)時(shí)域曲線如圖1所示。

圖1 仿真源信號(hào)時(shí)域曲線Figure 1 Time-domain curves of simulated source signals

將源信號(hào)[s1,s2,s3]乘以一個(gè)1×3的服從均勻分布U(0,1)的隨機(jī)矩陣,并添加高斯白噪聲,混合成一個(gè)單通道信號(hào)x(t),時(shí)域曲線如圖2所示。

圖2 仿真觀測(cè)信號(hào)時(shí)域曲線Figure 2 Time-domain curve of simulated observed signal

將單通道觀測(cè)信號(hào)x(t)使用EEMD分解成具有多個(gè)不同瞬時(shí)頻率的IMF分量,并根據(jù)式(14)的Pearson相關(guān)系數(shù)公式,計(jì)算得到IMF分量和x(t)的相關(guān)系數(shù),根據(jù)式(14)計(jì)算的相關(guān)系數(shù)如表1 所示。

表1 IMF分量與觀測(cè)信號(hào)的相關(guān)系數(shù)Table 1 Correlation coefficient between separated signals and observed signals

(14)

由式(15)計(jì)算IMF分量的相關(guān)系數(shù)累計(jì)貢獻(xiàn)率R,選擇累計(jì)貢獻(xiàn)率超過(guò)閾值η(本文選擇閾值為η=90%)的各分量和x(t)組成新的觀測(cè)信號(hào)X(t)。

(15)

式中:ri為按降序排列的第i個(gè)相關(guān)系數(shù);K為分量總個(gè)數(shù)。根據(jù)累計(jì)貢獻(xiàn)率,最終選擇了5個(gè)IMF分量和x(t)重構(gòu)成新的觀測(cè)信號(hào)X(t)=[x(t),imf2,imf3,imf4,imf1,imf5]T。

由于FastICA要求源信號(hào)個(gè)數(shù)不大于觀測(cè)信號(hào)個(gè)數(shù),因此選擇源數(shù)從2到6。為了降低隨機(jī)誤差的影響,共進(jìn)行仿真實(shí)驗(yàn)50次,分別計(jì)算出這50次實(shí)驗(yàn)中不同可能源數(shù)下的估計(jì)源信號(hào)的平均信息熵值,結(jié)果如圖3所示。

圖3 50次實(shí)驗(yàn)下不同源數(shù)的源信號(hào)平均信息熵值Figure 3 Average information entropy of the source signals with different source numbers in 50 experiments

圖3中當(dāng)源信號(hào)個(gè)數(shù)為3時(shí),計(jì)算得到的估計(jì)源信號(hào)信息熵值最小,因此認(rèn)為源信號(hào)的最佳個(gè)數(shù)為3。結(jié)果表明,本文提出的基于信息熵的源數(shù)估計(jì)算法估計(jì)得到的最佳源數(shù)和設(shè)定的仿真源信號(hào)個(gè)數(shù)相等。在50次實(shí)驗(yàn)中成功次數(shù)為46,得到本文源數(shù)估計(jì)算法的估計(jì)準(zhǔn)確率為92%。

根據(jù)以上實(shí)驗(yàn)結(jié)果選擇源信號(hào)個(gè)數(shù)為3,并使用FastICA算法對(duì)盲源信號(hào)進(jìn)行恢復(fù),得到估計(jì)源信號(hào)如圖4所示。

圖4 仿真信號(hào)的估計(jì)信號(hào)時(shí)域曲線Figure 4 Estimated signals′ time-domain curves of simulated signals

為考察本文算法對(duì)源信號(hào)的估計(jì)效果,根據(jù)式(14)計(jì)算估計(jì)信號(hào)與源信號(hào)之間的Pearson相關(guān)系數(shù),結(jié)果如表2所示。

表2 估計(jì)信號(hào)和源信號(hào)相關(guān)系數(shù)Table 2 Correlation coefficient between estimated signals and source signals

由圖4和表2結(jié)果可知,分離得到的估計(jì)源信號(hào)和仿真源信號(hào)的波形非常接近,且估計(jì)源信號(hào)和仿真源信號(hào)之間的相關(guān)系數(shù)絕對(duì)值均大于0.9。結(jié)果表明,本文提出的基于信息熵最小的源數(shù)估計(jì)算法估計(jì)能夠精確地估計(jì)出盲源個(gè)數(shù),為EEMD和FastICA算法能夠有效分離混合信號(hào)得到目標(biāo)信號(hào)奠定了基礎(chǔ)。

為了進(jìn)一步驗(yàn)證算法的可靠性,將仿真實(shí)驗(yàn)設(shè)置為3個(gè)源、4個(gè)源、5個(gè)源和6個(gè)源混合的4種情況,并使用上述算法對(duì)分離結(jié)果的源信號(hào)平均信息熵值進(jìn)行計(jì)算,結(jié)果如圖5所示。

圖5結(jié)果表明,當(dāng)?shù)磾?shù)等于真實(shí)源數(shù)時(shí),估計(jì)源信號(hào)的平均信息熵會(huì)達(dá)到最小,進(jìn)一步說(shuō)明了本文提出的基于信息熵最小的源數(shù)估計(jì)算法的有效性。

實(shí)驗(yàn)結(jié)果證明,當(dāng)估計(jì)源數(shù)小于真正源數(shù)時(shí),隨著迭代源數(shù)增加,源信號(hào)的平均信息熵值減小;而當(dāng)?shù)磾?shù)大于真正源數(shù)時(shí),源信號(hào)的平均信息熵值又在增大;只有當(dāng)平均信息熵達(dá)到最小時(shí),對(duì)應(yīng)迭代源數(shù)則為正確源數(shù)。

3 與傳統(tǒng)源數(shù)估計(jì)算法的對(duì)比實(shí)驗(yàn)

為驗(yàn)證本文算法的魯棒性和高估計(jì)準(zhǔn)確率,將本文算法與傳統(tǒng)源數(shù)估計(jì)算法PCA和AIC進(jìn)行對(duì)比。由于噪聲能量會(huì)影響算法的性能,故在仿真實(shí)驗(yàn)中,通過(guò)考察不同估計(jì)算法隨信噪比(SNR)的變化來(lái)分析其估計(jì)準(zhǔn)確性能。為便于分析,實(shí)驗(yàn)以算法的估計(jì)準(zhǔn)確率為性能的評(píng)價(jià)指標(biāo),以估計(jì)正確頻率表示各算法的估計(jì)準(zhǔn)確率A如式(16)所示:

(16)

式中:n為估計(jì)正確次數(shù);N為模擬總次數(shù)。

數(shù)據(jù)采用第2節(jié)的仿真信號(hào),設(shè)置SNR為-10～10 dB,步長(zhǎng)為1 dB,進(jìn)行1 000次模擬實(shí)驗(yàn),計(jì)算不同信噪比下各算法的估計(jì)準(zhǔn)確率,得到實(shí)驗(yàn)結(jié)果如圖6所示。

圖6 3種算法的估計(jì)準(zhǔn)確率隨信噪比變化的情況Figure 6 Variation of 3 algorithms′ accuracy with different SNR

由圖6知,在高斯噪聲環(huán)境中,當(dāng)SNR在-10～-7 dB之間時(shí),本文算法和PCA算法的估計(jì)準(zhǔn)確率均高于AIC算法,且本文算法比另外兩種算法的準(zhǔn)確率都要高;當(dāng)SNR大于-7 dB時(shí),本文算法和AIC算法的估計(jì)準(zhǔn)確率均高于PCA;但當(dāng)SNR逐漸增大時(shí),AIC算法的估計(jì)準(zhǔn)確率會(huì)呈現(xiàn)下降趨勢(shì),估計(jì)結(jié)果不具有一致性。本文算法整體表現(xiàn)性能均高于另外兩種算法,且最終估計(jì)準(zhǔn)確率穩(wěn)定在94%左右,因此本文算法更具有魯棒性,且估計(jì)準(zhǔn)確率更高。

4 算法在TETRA通信信號(hào)中的驗(yàn)證

通信信號(hào)在日常傳輸過(guò)程中,經(jīng)常會(huì)受到不同類型干擾信號(hào)的干擾[20]。為了保證通信正常,需要對(duì)混有干擾的數(shù)據(jù)進(jìn)行目標(biāo)信號(hào)提取,從而達(dá)到將目標(biāo)信號(hào)與干擾分開(kāi)的目的。使用TETRA手持電臺(tái)發(fā)射信號(hào)進(jìn)行實(shí)時(shí)采集,并對(duì)其進(jìn)行頻譜搬移變?yōu)榛鶐盘?hào)。信號(hào)參數(shù)為T(mén)ETRA通信信號(hào)采樣率為9.21 MHz,中心頻率為0 Hz,信號(hào)帶寬為25 kHz,幅值為1 V,采樣時(shí)長(zhǎng)為9 999/9 210 s。疊加兩類通信干擾信號(hào):單音干擾和寬帶線性掃頻干擾信號(hào)。干擾信號(hào)參數(shù)為單音干擾信號(hào)中心頻率為1 kHz、幅值為1 V;寬帶線性掃頻干擾信號(hào)的起始頻率為5 Hz,截止頻率為100 Hz,掃頻斜率為95 Hz/s,信號(hào)幅度為1 V。

首先基于本文提出的源數(shù)估計(jì)算法估計(jì)出源信號(hào)個(gè)數(shù),然后結(jié)合EEMD和FastICA將源信號(hào)分離出來(lái)?；煊懈蓴_的TETRA信號(hào)如圖7所示。

圖7 混有干擾信號(hào)的TETRA觀測(cè)信號(hào)時(shí)域曲線Figure 7 Time-domain curve of TETRA observation signal with interference signals

將觀測(cè)信號(hào)進(jìn)行EEMD分解,計(jì)算IMF分量和TETRA信號(hào)的Pearson相關(guān)系數(shù)如表3所示。

表3 IMF分量與觀測(cè)信號(hào)的相關(guān)系數(shù)Table 3 Correlation coefficients between IMF components and observed signal

計(jì)算IMF分量的相關(guān)系數(shù)累計(jì)貢獻(xiàn)率,得到imf1～imf7分量的累計(jì)貢獻(xiàn)率超過(guò)90%,故選擇這7個(gè)IMF分量和x(t)組成新的觀測(cè)信號(hào),即X(t)=[x(t),imf2,imf1,imf5,imf6,imf4,imf7,imf3]T,使用FastICA對(duì)X(t)進(jìn)行盲源分離,通過(guò)本文算法計(jì)算得到的源信號(hào)平均信息熵值如表4所示。

表4 估計(jì)源信號(hào)平均信息熵值Table 4 Average information entropy of the estimated source signals

根據(jù)表4結(jié)果,當(dāng)源信號(hào)個(gè)數(shù)為3時(shí),估計(jì)源信號(hào)平均信息熵值最小,因此得到源信號(hào)的最佳個(gè)數(shù)為3。確定源信號(hào)個(gè)數(shù)后,再使用FastICA算法,得到TETRA信號(hào)盲分離結(jié)果如圖8所示。

圖8 TETRA信號(hào)盲分離結(jié)果Figure 8 BSS result of TETRA signal

計(jì)算分離得到的估計(jì)信號(hào)與源信號(hào)之間的Pearson相關(guān)系數(shù),如表5所示。

表5 估計(jì)信號(hào)和源信號(hào)相關(guān)系數(shù)Table 5 Correlation coefficient between estimated signals and source signals

根據(jù)表5結(jié)果可知,估計(jì)信號(hào)y1對(duì)應(yīng)了單音干擾信號(hào),y2對(duì)應(yīng)了TETRA通信信號(hào),y3對(duì)應(yīng)了掃頻干擾信號(hào),相關(guān)系數(shù)均達(dá)到0.9以上。結(jié)果驗(yàn)證了本文提出的基于信源信息熵最小的單通道盲源數(shù)估計(jì)算法的有效性,能準(zhǔn)確獲取盲源信號(hào)個(gè)數(shù)。

5 結(jié)論

本文通過(guò)對(duì)估計(jì)源信號(hào)的分布信息進(jìn)行估計(jì),提出了一種將信息熵作為評(píng)價(jià)指標(biāo)的盲源信號(hào)源數(shù)估計(jì)算法。該算法運(yùn)用遞歸的方式檢驗(yàn)不同可能盲源數(shù)下通過(guò)FastICA分離得到的信號(hào)的信息量,通過(guò)分析和比較信息熵,篩選出最佳的盲源數(shù)作為真實(shí)盲源數(shù)的估計(jì)值。文中的算法充分利用了源信號(hào)的信息特征,進(jìn)而給出了源信號(hào)個(gè)數(shù)估計(jì)的算法。基于仿真數(shù)據(jù)和真實(shí)數(shù)據(jù)的數(shù)字實(shí)驗(yàn)結(jié)果也驗(yàn)證了本文算法的有效性和較強(qiáng)的抗噪性,為源數(shù)未知條件下的盲源分離提供了有效技術(shù)支持。