溫江濤，趙倩云，孫潔娣

(1.燕山大學(xué) 河北省測試計量技術(shù)及儀器重點(diǎn)實驗室，河北，秦皇島 066004；2.燕山大學(xué) 信息科學(xué)與工程學(xué)院，河北，秦皇島 066004)

?

基于信號時頻域聚集性的欠定盲分離混合矩陣估計方法

溫江濤1，趙倩云1，孫潔娣2

(1.燕山大學(xué) 河北省測試計量技術(shù)及儀器重點(diǎn)實驗室，河北，秦皇島 066004；2.燕山大學(xué) 信息科學(xué)與工程學(xué)院，河北，秦皇島 066004)

為解決欠定盲分離中混合矩陣估計問題，通過研究觀測信號在時頻域的線性聚集特性，提出一種基于時頻域線性聚集程度差異的混合矩陣估計方法，并著重研究在信號線性聚集程度較弱情況下對混合矩陣的估計. 首先，利用觀測信號或其時頻域中相應(yīng)變換系數(shù)的比值分布衡量信號線性聚集程度；其次，采用優(yōu)化初始中心的K-均值聚類算法估計混合矩陣. 該算法降低了對信號稀疏性的要求，并且可以較高精度地估計出混合矩陣. 仿真實驗結(jié)果表明該方法具有可行性和有效性.

欠定盲源分離；混合矩陣估計；線性聚集性；改進(jìn)K-均值聚類

盲源分離(blind source separation，BSS)是指傳輸信道和源信號都未知時，僅從觀測信號中分離源信號的過程. 在實際問題中，觀測信號的個數(shù)一般少于源信號個數(shù)，即為欠定盲源分離(underdetermined blind source separation，UBSS). 其簡化的數(shù)學(xué)模型為

(1)

式中X(t)=[x1(t)x2(t) …xm(t)]T為由M個傳感器得到的觀測信號；A∈Rm×n未知的混合矩陣；S(t)=[s1(t)s2(t) …sn(t)]T為N個未知的源信號，且M

Lee等[4]提出先估計混合矩陣，然后重構(gòu)源信號的方法，即“兩步法”，它是目前解決稀疏信號分離問題的常用方法. 實現(xiàn)源信號的精確重構(gòu)，需準(zhǔn)確估計混合矩陣A，因此高精度地估計A是兩步法的核心. 近年來在SCA的基礎(chǔ)上，提出了許多解決欠定盲源分離問題的新算法，如退化分離估計技術(shù)(DUET)[5]、時頻混合比算法(TIFROM)[6]、超平面聚類算法[7]、非線性投影列掩蔽算法[8]等，這些算法均假設(shè)信號為稀疏或弱稀疏信號，但并未對信號的稀疏程度做進(jìn)一步的研究. 估計A的前提條件是混合信號須在散點(diǎn)圖中呈現(xiàn)較為明顯的方向性，即信號須有一定的線性聚集性，而實際中很多情況并不滿足此要求.

針對上述問題，本文提出一種基于信號時頻域內(nèi)聚集性估計混合矩陣的算法. 該算法通過統(tǒng)計觀測信號在時頻域內(nèi)變換系數(shù)的比值分布情況衡量觀測信號的線性聚集程度，并根據(jù)不同程度的線性聚集性研究最佳的混合矩陣估計方法. 該方法降低了對信號稀疏性的要求，僅需增強(qiáng)散點(diǎn)圖中線性聚集程度即可較為精確地估計混合矩陣.

1 欠定盲源分離問題

式(1)展開可改寫為如下形式:

(2)

式中：xi(t)是行矢量，表示第i路觀測信號；si(t)也為行矢量，表示第i個源信號. 對式(2)進(jìn)行短時傅里葉變換(STFT)得

(3)

式中xi(t，k)和si(t，k)分別為第i路觀測信號xi(t)和第i個源信號si(t)的STFT系數(shù).

解決欠定盲源分離問題的方法主要基于信號的稀疏性，假設(shè)源信號S為強(qiáng)稀疏信號或其變換系數(shù)具有強(qiáng)稀疏性，即非零值的個數(shù)是唯一的或只有少數(shù)信號具有較大值，大部分?jǐn)?shù)值很小可視為零值. 研究表明多數(shù)信號在變換域中都具有稀疏性[9]，本文主要以混合信號在時頻域中變換系數(shù)的稀疏性為例進(jìn)行理論推導(dǎo).

若在時頻點(diǎn)(tk，kk)僅有sj(tk，kk)為非零值，有

(4)

式中aj為混合矩陣A的第j列，式(4)可改寫為

(5)

從幾何意義來說，式(5)是一個直線方程，故sj(tk，kk)所有采樣時刻將確定一條直線，并且直線的方向為混合矩陣A的第j列向量，此時信號在時頻域中線性聚集性強(qiáng)，可直接采用聚集算法估計混合矩陣. 然而，在實際情況中信號的線性聚集程度是未知的，不確定是否能夠直接采用聚類算法估計混合矩陣. 因此，有必要對信號的線性聚集程度進(jìn)行衡量，并根據(jù)線性聚集程度進(jìn)一步研究精度更高的混合矩陣估計算法.

2 基于信號時頻域聚集性的混合矩陣估計方法

估計混合矩陣A的主要依據(jù)是信號的線性聚集性，信號的線性聚集性強(qiáng)，信號在散點(diǎn)圖中呈現(xiàn)出清晰地方向性，可直接運(yùn)用聚類算法估計混合矩陣；反之，不呈現(xiàn)方向性，直接運(yùn)用聚類算法得不到混合矩陣的估計，需增強(qiáng)線性聚集性，再進(jìn)行混合矩陣的估計. 因此，衡量信號線性聚集程度并增強(qiáng)其程度對精確估計A具有重要影響.

2.1 信號線性聚集程度的衡量

將信號X(t)進(jìn)行STFT，為了更好地衡量信號線性聚集程度，利用式(6)去除能量值較小的變換系數(shù):

(6)

式中ε為閾值. 由此時的變換系數(shù)得到其比值d(t，k)為

(7)

令a=min{d(t，k)}，b=max{d(t，k)}，將區(qū)間[a，b]分為相等的P個子空間，可表示為{a+(i-1)δ，a+iδ}i=1，2，…，P，其中δ=(b-a)/P，且P充分大，統(tǒng)計d(t，k)樣本點(diǎn)落在各子區(qū)間的個數(shù)Di(i=1，2，…，P). 比值d(t，k)相差較小的變換系數(shù)近似匯聚于一條直線，且d(t，k)會落入同一子區(qū)間內(nèi)；反之，落入不同的子區(qū)間. 線性聚集性強(qiáng)時，比值d(t，k)集中落入部分子區(qū)間內(nèi)，對應(yīng)的Di較大，P個子區(qū)間內(nèi)Di呈現(xiàn)極大值和極小值的間隔分布；線性聚集性弱時，比值d(t，k)會散落在相對應(yīng)的子區(qū)間內(nèi)，P個子區(qū)間內(nèi)Di不呈現(xiàn)極大值和極小值的間隔分布. 因此，統(tǒng)計Di分布能夠衡量線性聚集程度.

2.2 增強(qiáng)信號線性聚集程度

通?；旌闲盘栐跁r域中甚至?xí)r頻域中的比值并不會集中落在某些子區(qū)間，即Di在P個子區(qū)間中不呈現(xiàn)峰值點(diǎn)，因此，線性聚集性較弱. 針對此種情況，本文研究了基于時頻域等比系數(shù)增強(qiáng)觀測信號時頻域內(nèi)線性聚集程度的方法.

定義1 觀測信號在時頻域內(nèi)變換系數(shù)的實部之比與虛部之比相等，且為常數(shù)，該常數(shù)是混合矩陣對應(yīng)列向量的比值，稱此類變換系數(shù)為等比系數(shù).

以觀測信號個數(shù)M=2為例，假設(shè)時頻域中某時頻點(diǎn)(tz，kz)上的比值與混合矩陣A第i列的比值相等. 則有

(8)

式(8)可改寫為

(9)

最大概率滿足式(9)的情況是在時頻點(diǎn)(tz，kz)處僅有源信號si(tz，kz)作用，其他源信號對X(tz，kz)不起作用或作用相互抵消. 式(9)可改寫為

(10)

由式(10)可知等比系數(shù)實部的比值與虛部的比值相等，且與混合矩陣對應(yīng)列向量的比值相等. 因此，等比系數(shù)散點(diǎn)圖呈現(xiàn)清晰的方向性，其線性聚集性強(qiáng). 然而實際應(yīng)用中，式(10)中Re[X2(t，k)]/Re[X1(t，k)]和lm[X2(t，k)]/lm[X1(t，k)]的值會存在一定的誤差，并非完全相等，因此設(shè)定一個門限值ω，消除此種誤差，判定等比系數(shù)的公式為

(11)

式中ω的取值范圍為(0，0.1)，ω是一個經(jīng)驗性取值，滿足式(11)的變換系數(shù)可視為時頻域中的等比系數(shù). 因此，可以進(jìn)一步對等比系數(shù)進(jìn)行分析得到混合矩陣的估計值.

2.3 混合矩陣的估計方法

前文對時頻域中等比系數(shù)進(jìn)行了有效判定，將有效選取的等比系數(shù)采用式(12)進(jìn)行單位化處理，得到集合U為

(12)

式中X(tz，kz)為時頻域上選定的等比系數(shù)，Ωz為所有選定等比系數(shù)的集合. 由于該集合中的向量均為單位向量，進(jìn)一步計算向量之間的比值，即

(13)

運(yùn)用聚類的方法將集合H進(jìn)行分類，若H中任意兩個hm、hn相等或者近似相等，則將其歸為一類.

(14)

式中：sum(Cn)為第n類中所有向量個數(shù)；∑hCn為第n類中所有元素之和，則可通過式(15)估計混合矩陣列向量為

(15)

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

針對欠定盲源分離問題，對本文提出的基于時頻域信號聚集性估計混合矩陣的方法進(jìn)行仿真及分析. 估計誤差采用歸一化均方誤差(NMSE)進(jìn)行評價，并與Bofill的勢函數(shù)算法進(jìn)行對比.

歸一化均方誤差為

(16)

3.1 觀測信號在時頻域內(nèi)線性聚集性較強(qiáng)的情況

取M=2、N=6，即有6路源信號，2路觀測信號. 6路源信號截取自6段長笛聲音信號，取自Bofill[10]的仿真實驗，6路源信號采樣點(diǎn)數(shù)均為32 768點(diǎn). 隨機(jī)選取混合矩陣A，并對其列矢量歸一化為

根據(jù)式(1)混合后得到兩路觀測信號，如圖1所示.

首先衡量觀測信號的線性聚集性，圖2為P=40時觀測信號在時域及時頻域中比值分布情況.

分析圖2中Di的分布可知：觀測信號在時域中Di不具有明顯的極值，因此線性聚集性差；而在時頻域中Di的極值間隔分布，線性聚集性較強(qiáng). 為驗證上述分析的正確性，圖3中給出觀測信號的時域和時頻域散點(diǎn)圖.

由圖3可知，觀測信號在時域散點(diǎn)圖中未呈現(xiàn)方向性，而在時頻域具有明顯的方向性，線性聚集性較強(qiáng).

因此估計混合矩陣A時可以采用改進(jìn)的K-均值聚類算法，得到的估計結(jié)果為

為評估本文算法估計矩陣的精確度，采用Bofill的勢函數(shù)算法[10]估計混合矩陣為

由圖4可知：本文提出的算法估計混合矩陣的歸一化均方誤差明顯小于Bofill勢函數(shù)算法，本文提出算法估計的混合矩陣更為精確.

3.2 觀測信號在時頻域內(nèi)線性聚集性較弱的情況

取M=2、N=4，即4路源信號，2路觀測信號. 采用的4個源信號均為實際語音信號，取自http://www.speech.cs.cmu.edu/cmu_arctic/，每組源信號的采樣點(diǎn)數(shù)為35 920點(diǎn)，混合矩陣為

根據(jù)式(1)混合后得到兩組觀測信號如圖5所示.

首先衡量觀測信號的線性聚集性，圖6是P=40時觀測信號在時域及時頻域中比值分布情況.

分析圖6中Di的分布可知：兩個分布圖均不具有明顯極值，因此，線性聚集性均較弱. 直接運(yùn)用改進(jìn)的K-均值聚類算法估計混合矩陣會有很大的誤差甚至無法得到估計矩陣. 圖7給出了觀測信號在時域、頻域、時頻域的散點(diǎn)圖用于驗證上述分析的正確性.

由圖7可知，觀測信號在時域和頻域散點(diǎn)圖中未呈現(xiàn)方向性，在時頻域散點(diǎn)圖中呈現(xiàn)模糊的方向，線性聚集性依舊較弱. 下面利用本文提出的算法，在時頻域中選取適當(dāng)?shù)缺认禂?shù)增強(qiáng)線性聚集性進(jìn)而估計混合矩陣A. 利用式(11)選出時頻域中的等比系數(shù)，其中選取ω=0.025，如圖8為等比系數(shù)的比值分布情況：所示Di的分布具有明顯的極值，線性聚集性得到了增強(qiáng).

分析圖8中Di的分布可知：等比系數(shù)在時頻域中具有較強(qiáng)的線性聚集性. 圖9為時頻域等比系數(shù)的散點(diǎn)圖.

實驗2信號直接采用Bofill勢函數(shù)算法，無法估計出混合矩陣；若先進(jìn)行時頻域選取等比系數(shù)的方法增強(qiáng)信號線性聚集性后，再采用Bofill勢函數(shù)算法則可得到估計混合矩陣的估計值為

由圖10可知：本文提出的算法估計混合矩陣的歸一化均方誤差在-62 dB左右，而Bofill勢函數(shù)算法估計混合矩陣的歸一化均方誤差在-42 dB左右，本文提出的算法估計混合矩陣更為精確.

6 結(jié) 論

研究結(jié)果及仿真實驗表明：對于稀疏程度未知的信號，本文算法可以有效鑒別時頻域上的等比系數(shù)，增強(qiáng)信號的線性聚集性，并采用改進(jìn)的聚類算法較為精確地估計出混合矩陣. 本文方法可用于語音信號、機(jī)械振動信號等領(lǐng)域的欠定盲源分離問題；該算法降低了對信號稀疏性的要求且能高精度的估計混合矩陣，為下一步精確重構(gòu)信號奠定了基礎(chǔ).

[1] Li X，Zhang X. Sequential blind extraction adopting second-order statistics[J]. IEEE Signal Processing Letters，2007，14(1):58-61.

[2] Reju V G，Koh Soo Ngee，Soon Ing Yann. An algorithm for mixing matrix estimation in instantaneous blind source separation[J]. Signal Processing，2009，89:1762-1773.

[3] Xie S L，Yang L，Yang J M，et al. Time-frequency approach to underdetermined blind source separation[J]. IEEE Transactions on Neural Networks and Learning Systems，2012，23(2):306-316.

[4] Lee T W，Lewicki M S，Girolami M. Blind source separation of more sources than mixtures using overcomplete representations[J]. IEEE Signal Process，1999，6(4):87-90.

[5] Tengtrairat N，Woo W L. Extension of DUTE to single-channel mixing model and separability analysis[J]. Signal Processing，2014，96(3):261-265.

[6] Abrard F，Deville Y. A time-frequency blind signal separation method applicable to underdetermined mixtures of dependent sources[J]. Signal Process，2005，85(7):211-221.

[7] Jayaraman J，Thiagarajan. Mixing matrix estimation using discriminative clustering for blind source separation[J]. Digital Signal Processing, 2013，23(1):9-18.

[8] Zhou G X，Yang Z Y，Xie S L. Mixing matrix estimation from sparse mixtures with unknown number of source[J]. IEEE Trans on Neural Networks，2011，22(2):211-221.

[9] Bao G Z，Ye Z F，Xu X，et al. A compressed sensing approach to blind separation of speech mixture based on a two-layer sparsity model[J]. IEEE Audio Speech Language Processing，2013，21(5):899-906.

[10] Bofill P，Zibulevesky M. Undetermined blind source separation using sparse representations[J]. Signal Processing，2001，81(11):2353-2362.

(責(zé)任編輯：劉芳)

Mixing Matrix Estimation Based on Cluster Degree of Time-Frequency Signal for Underdetermined Blind Source Separation

WEN Jiang-tao1，ZHAO Qian-yun1，SUN Jie-di2

(1.Key Laboratory of Measurement Technology and Instrumentation of HeBei Province ，Yanshan University，Qinhuangdao，Hebei 066004，China; 2.School of Information Science and Engineering，Yanshan University，Qinhuangdao，Hebei 066004，China)

To solve the mixing matrix estimation problem of underdetermined blind separation，through the study of signal linear aggregation feature in time-frequency domain，an estimation method of mixing matrix based on different signal linear aggregation degree in time-frequency domain was proposed in this paper，and focus on the estimation of mixing matrix under the signal linear aggregation degree in weaker conditions. First，the observed signal or the ratio distribution of the corresponding transformation coefficient in time-frequency domain was used to measure the degree of the signal linear aggregation; second，the improved K-means clustering algorithm was applied to estimate the mixing matrix. The proposed method reduces the requirement for signal sparsity and can estimate the mixing matrix accurately. The simulation results show that the proposed method is feasible and effective.

underdetermined blind sources separation; mixture matrix estimation; linear aggregation; improved K-means clustering

2014-11-15

國家自然科學(xué)基金資助項目(51204145)；河北省自然科學(xué)基金資助項目(E2013203300，E2016203223)

溫江濤(1974—)，男，博士，副教授，E-mail：wens2002@163.com.

TN 911.7

A

1001-0645(2016)07-0733-06

10.15918/j.tbit1001-0645.2016.07.014