基于稀疏貝葉斯的陣列幅相和互耦誤差聯(lián)合校正方法

王鼎，高衛(wèi)港，吳志東

（信息工程大學(xué)信息系統(tǒng)工程學(xué)院，河南 鄭州 450001）

0 引言

陣列信號處理技術(shù)自面世以來，在雷達(dá)[1]、通信[2]、水聲[3]、醫(yī)學(xué)[4]等領(lǐng)域得到了廣泛的應(yīng)用。同時(shí)，隨著陣列處理理論的不斷完善發(fā)展，空間譜估計(jì)技術(shù)因其在遠(yuǎn)距離探測、低截獲信號處理等方面的優(yōu)勢，受到了國內(nèi)外學(xué)者的關(guān)注。

縱觀空間譜估計(jì)技術(shù)的整個(gè)發(fā)展歷程，從最開始的傳統(tǒng)空間譜估計(jì)，即經(jīng)典傅式譜分析在空域上的直接推廣，這種方法因受限于陣列孔徑無法突破“瑞利限”[5]，到后來的現(xiàn)代空間譜估計(jì)技術(shù)，其中最具代表性的便是以子空間技術(shù)為基礎(chǔ)的多重信號分類（MUSIC,multi signal classification）算法[6]。與傳統(tǒng)空間譜估計(jì)技術(shù)相比，MUSIC 算法的估計(jì)精度和角度分辨率都達(dá)到了新的高度。但是，隨著子空間算法理論的不斷發(fā)展，此類算法的局限性也逐漸暴露出來，對陣列輸出協(xié)方差矩陣以及信號/噪聲子空間的準(zhǔn)確度有高度依賴性導(dǎo)致子空間算法在低信噪比及小快拍情況下效果并不理想。此外，在發(fā)展過程中，還出現(xiàn)了極大似然類測向方法[7]，這種方法通過借助特定維數(shù)模型與觀測數(shù)據(jù)之間的最佳擬合來實(shí)現(xiàn)空間譜估計(jì)，在一定程度上提高了對低信噪比和小樣本情況的適應(yīng)能力，但其又因?yàn)檩^低的計(jì)算效率及對信號個(gè)數(shù)先驗(yàn)的依賴無法很好地發(fā)揮性能優(yōu)勢。

無論是低信噪比還是小樣本情況，陣列信號處理問題本質(zhì)上都屬于從有限的觀測數(shù)據(jù)中獲得相應(yīng)的參數(shù)，而近幾年興起的稀疏重構(gòu)技術(shù)在解決信號缺失條件下信號恢復(fù)問題的過程中展現(xiàn)出了極大優(yōu)勢，目前在壓縮感知[8]、信號分析[9]、圖像處理[10]等領(lǐng)域已有了一定規(guī)模的應(yīng)用，并取得了不錯(cuò)的成果。稀疏重構(gòu)是指利用信號在特定變換域（如空域）上的稀疏性，從包含信號所有可能分量的完備或超完備集合中選取少量幾個(gè)基函數(shù)重構(gòu)原始信號的過程。目前，稀疏重構(gòu)算法主要包括匹配追蹤[11]、Lp范數(shù)[12]和稀疏貝葉斯[13]三類。其中，匹配追蹤類算法對觀測模型的稀疏度有很高的要求，在信號不同分量之間相關(guān)性較強(qiáng)時(shí)性能不好；Lp范數(shù)類算法的局限性則在于正則化因子參數(shù)的選取問題，具體的觀測模型、信號環(huán)境等因素都會對該參數(shù)的選取產(chǎn)生極大的影響；稀疏貝葉斯學(xué)習(xí)方法是由Tipping[14]提出的利用貝葉斯原理與觀測模型先驗(yàn)信息相結(jié)合的一類重構(gòu)方法。目前，稀疏貝葉斯重構(gòu)算法在圖像處理領(lǐng)域的發(fā)展較為完善，在其他領(lǐng)域則相對滯后。

實(shí)際中，無論是具有超分辨率的子空間算法，還是稀疏重構(gòu)算法，都無法達(dá)到理論估計(jì)精度，這是因?yàn)樵趯?shí)際中無法避免由于各種原因?qū)е碌年嚵姓`差，包括各個(gè)陣元和通道的幅度和相位不具有一致性[15]、天線之間相互影響產(chǎn)生互耦效應(yīng)[16]、陣元安裝位置出現(xiàn)偏差引入的位置誤差[17]等，都會使陣列流型和理想的陣列流型之間存在較大的偏差，從而使各種算法的性能出現(xiàn)一定程度的惡化，甚至完全失效。針對幅相誤差和互耦誤差共存的情況，研究人員[18-20]提出了一系列的聯(lián)合校正方法，其中，文獻(xiàn)[18]利用均勻線陣互耦矩陣的特殊結(jié)構(gòu)，通過交替迭代的方式對陣列誤差參數(shù)進(jìn)行優(yōu)化校正。文獻(xiàn)[19]利用陣列協(xié)方差矩陣的關(guān)系式構(gòu)造了一種二次代價(jià)函數(shù)，并通過交替迭代進(jìn)行參數(shù)估計(jì)。文獻(xiàn)[20]提出了基于協(xié)方差匹配技術(shù)的誤差聯(lián)合校正算法。文獻(xiàn)[21]利用MIMO 系統(tǒng)，引入若干個(gè)經(jīng)過精確校準(zhǔn)的輔助陣元對誤差及波達(dá)方向進(jìn)行聯(lián)合估計(jì)。

基于對當(dāng)前研究現(xiàn)狀的分析，本文針對幅相誤差和互耦誤差共存的情況，提出了基于信號空域稀疏性的貝葉斯測向與誤差校正方法，并分析了信號測向精度及相應(yīng)誤差校準(zhǔn)的理論下限，仿真結(jié)果驗(yàn)證了方法的有效性。

本文用到的一些數(shù)字符號說明如下。向量和矩陣均用粗斜體表示。xT和xH分別表示向量的轉(zhuǎn)置和共軛轉(zhuǎn)置分別表示L1 范數(shù)、L2 范數(shù)、Frobenous 范數(shù)。表示矩陣A的行列式。tr(A) 表示矩陣A的跡。(A)j,:表示矩陣A的第j行，(A):,j表示矩陣A的第j列，(A)i,j表示矩陣A的第(i,j) 個(gè)元素。xj表示向量x的第j個(gè)元素。diag(A) 表示由A的對角線元素構(gòu)成的向量；diag(x) 表示對角矩陣，其對角向量對應(yīng)于向量x；toeplitz(x) 表示由x形成的托普利茲矩陣。x'(?)表示x(?)關(guān)于?的導(dǎo)數(shù)。Re(·) 和Im(·) 分別表示對復(fù)數(shù)進(jìn)行取實(shí)和取虛操作。x⊙y表示x和y的Hadamard 積。表示給定分布概率條件下對應(yīng)的期望運(yùn)算符。

1 信號模型

假設(shè)在空域中有K個(gè)遠(yuǎn)場窄帶信號入射到M元陣列，入射方向?yàn)椋敲碝×1維的陣列接收數(shù)據(jù)為

圖1 陣列信號接收模型

在實(shí)際中，式(1)所對應(yīng)的理想情況是不存在的，更普遍的是幅相誤差與互耦誤差同時(shí)存在，在這種情況下，陣列輸出的計(jì)算式為

采用文獻(xiàn)[22]中使用的離格模型，式(2)可以寫為

為了進(jìn)一步明確陣列輸出與互耦誤差系數(shù)和幅相誤差參數(shù)之間的關(guān)系，對式(3)做如下變換

通過對空間角度域進(jìn)行均勻采樣，將式(3)和式(4)進(jìn)行空域擴(kuò)展，得到陣列輸出的超完備表示形式為

2 模型求解

2.1 參數(shù)求解

引入中間變量γ=[γ1,…,γ L]T表示方 向?1,…,? L處對應(yīng)的信號分量幅度的方差，即其中R=diag(γ)。

陣列輸出的概率密度函數(shù)為

結(jié)合式(6)和式(7)，得到給定γ,c,δ2,β條件下的概率為

令式(11)的導(dǎo)數(shù)為零，可得c第q次迭代時(shí)的值為

令式(13)的導(dǎo)數(shù)為零，可得f第q次迭代時(shí)的值為

令式(15)的導(dǎo)數(shù)為零，可得第q次迭代時(shí)δ2的值為

令式(17)的導(dǎo)數(shù)為零，可得第q次迭代時(shí)的值為

其中，constant 是一個(gè)常數(shù)余項(xiàng)，P是一個(gè)半正定矩陣，計(jì)算式為

其中，μt由式(25)得到，實(shí)際上，β與是聯(lián)合稀疏的，對應(yīng)的K個(gè)源位置處有K個(gè)非零值，只需計(jì)算對應(yīng)γl的最大K個(gè)位置處的β值，其余值設(shè)為零。因此，β、ν和P的大小可以分別減小為K、K和K×K維。

由式(19)可知

對于式(12)、式(13)、式(16)和式(18)中給出的期望值取決于條件概率密度函數(shù)并且的形式不是直觀給出的，因此有必要對上述參數(shù)更新方法進(jìn)行簡化。

由于在t時(shí)刻，信號幅度僅與該時(shí)刻的陣列觀測值x(t) 有關(guān)，根據(jù)式(6)和式(7)可以得到關(guān)于x(t) 的后驗(yàn)概率密度函數(shù)為

將式(26)～式(30)代入?yún)?shù)的更新式便可實(shí)現(xiàn)對參數(shù)c,δ2,γ,β的更新。而后估計(jì)出信號幅度的均值與方差，進(jìn)入下一次的EM 迭代。由文獻(xiàn)[23]可知，EM 算法具有嚴(yán)格的收斂性，該迭代過程逐漸逼近穩(wěn)定的估計(jì)值，當(dāng)?shù)諗繒r(shí)，由各變量的估計(jì)值便可以獲得信號的空域分布、陣列誤差模型和噪聲功率等。

3 參數(shù)估計(jì)方差的克拉美羅界

N個(gè)陣列觀測樣本x1,…,xN的概率密度函數(shù)為

對概率密度函數(shù)（式(31)）取對數(shù)可以得到如下對數(shù)似然函數(shù)

其中，?和δ2是實(shí)數(shù)域上的參數(shù)，s(t) 和c是復(fù)向量。因此，將s(t) 和c分解為實(shí)部和虛部之和的形式，分別為其中，分別代表復(fù)向量的實(shí)部與虛部。經(jīng)過上述分解之后，得到的參數(shù)集為

根據(jù)式(32)的似然函數(shù)，可以得到關(guān)于噪聲方差δ2的克拉美-羅下界（CRLB,Cramer-Rao lower bound）為

在存在陣列互耦或者陣列幅相誤差時(shí)，式(35)給出的CRLB 只是簡單地將誤差向量分解為實(shí)部與虛部，這并不能很好地反映出相關(guān)參數(shù)的物理屬性，因此，結(jié)合文獻(xiàn)[24]中的結(jié)論，可以將的CRLB 轉(zhuǎn)化為的CRLB。這樣就可以更加清晰地得出復(fù)向量c的估計(jì)精度。

4 仿真實(shí)驗(yàn)與結(jié)果分析

本節(jié)將通過模擬仿真實(shí)驗(yàn)對本文提出的多種誤差同時(shí)存在時(shí)的校正方法進(jìn)行驗(yàn)證，下面對一些固定參數(shù)進(jìn)行說明，陣列為8 元均勻直線陣，半徑波長比為0.5，信號入射角度為(30.3° 36.8°)，算法最大迭代次數(shù)為1 500 次，停止閾值為10-3。采用均方根誤差作為評判指標(biāo)，其中，N為蒙特卡羅仿真次數(shù)，在本文中設(shè)置為為第i次迭代的估計(jì)值，?為真實(shí)值。互耦系數(shù)為c=(0.6 +j0.4,0.2 +j0.1)T。值得一提的是，在這里假設(shè)第一個(gè)陣元完全校準(zhǔn)，因此，幅度誤差分別為(0 0.15 0.30 0.20 0.25 -0.20 -0.20 -0.30)，相位誤差分別為(0° -20°-30° -40° 35° 25° 20° 40°) 。

4.1 校正前后的空間譜圖

本節(jié)驗(yàn)證幅相誤差和互耦誤差對信號空間譜的影響，并探究本文方法對信號的分辨能力。假設(shè)2 個(gè)不相關(guān)信號以(30.3° 36.8°)的方向入射至8 元均勻直線陣，信噪比為10 dB，快拍數(shù)為100，其余參數(shù)不變。圖2 給出利用本文方法誤差校正前后的空間譜。其中，虛線輔助線為信號真實(shí)方向。

圖2 誤差校正前后的空間譜

從圖2 中可以看出，在多種誤差同時(shí)存在的情況下，空間譜的峰值被削弱，譜峰也發(fā)生相應(yīng)的偏移，導(dǎo)致分辨率下降。根據(jù)文獻(xiàn)[25]，在統(tǒng)計(jì)意義下，幅相誤差對靠近來波方向的空間譜峰影響較大，其中，幅度誤差只影響空間譜的峰值尖銳程度，不會改變譜峰的位置，相位誤差則會導(dǎo)致譜峰發(fā)生偏移，而互耦誤差則對兩者皆有影響。利用本文方法對幅相誤差和互耦誤差進(jìn)行聯(lián)合估計(jì)，可以很好地對來波信號方向進(jìn)行估計(jì)，提高空間譜的分辨能力。

4.2 性能分析

本節(jié)比較的方法包括離格稀疏貝葉斯推斷（OGSBI,off-grid sparse Bayesian inference）、求根稀疏貝葉斯學(xué)習(xí)（rootSBL,root sparse Bayesian learning）[26]和L1SVD 算法，分別從信噪比和快拍數(shù)2 個(gè)角度對本文方法的性能進(jìn)行探究，結(jié)果分別如圖3～圖6 所示。其中，當(dāng)探究信噪比對算法的影響時(shí)，信噪比的范圍為-10～10 dB，快拍數(shù)固定為100；當(dāng)探究快拍數(shù)對算法的影響時(shí)，快拍數(shù)范圍為10～100，信噪比固定為0，其余參數(shù)不變。

圖3 信號方位估計(jì)均方根誤差隨信噪比變化曲線

圖4 信號方位估計(jì)隨信噪比變化箱線

圖5 信號方位估計(jì)均方根誤差隨快拍數(shù)變化曲線

圖6 信號方位估計(jì)隨快拍數(shù)變化箱線

從圖3 和圖5 中可以看出，當(dāng)分別采用信噪比和快拍數(shù)為變量時(shí)，本文方法均能達(dá)到很好的效果，且逼近克拉美羅界。在信噪比較小和快拍數(shù)較低的情況下，稀疏貝葉斯類算法中的rootSBL 和OGSBI 與本文方法表現(xiàn)大致相同，但隨著信噪比和快拍數(shù)的增加，rootSBL 和OGSBI 算法雖然考慮了離散網(wǎng)格對信號方位估計(jì)造成的影響，但互耦誤差的存在使估計(jì)的空間譜偽峰較高，一定程度上限制了這2 種算法的性能。L1SVD 算法未考慮網(wǎng)格誤差的影響，效果最差。

從圖4 和圖6 中可以看出，隨著信噪比和快拍數(shù)的提升，箱線圖中的中位線逐漸逼近真實(shí)的來波方向，整體的波動程度迅速下降。2 次實(shí)驗(yàn)的蒙特卡羅次數(shù)均為100 次，直觀來看，當(dāng)信噪比提升時(shí)，箱線圖的波動程度相比于快拍數(shù)增大時(shí)收斂得更快，箱體更小。

將快拍數(shù)設(shè)置為 10～100，信噪比設(shè)置為-10～10 dB，其余條件保持不變，繪制互耦誤差和幅相誤差隨信噪比快拍數(shù)的變化曲面，分別如圖7 和圖8 所示。

圖7 互耦誤差隨信噪比快拍數(shù)變化曲面

圖8 幅相誤差隨信噪比快拍變化曲面

從圖7 和圖8 中可以看出，隨著信噪比和快拍數(shù)的提高，互耦誤差和幅相誤差曲面均能迅速下降，其中，互耦誤差曲面比幅相誤差曲面整體偏低的原因有2 個(gè)：1) 在觀測模型中，幅相誤差與信號方向之間的相互耦合效果要比互耦誤差與信號方向之間的耦合效果強(qiáng)；2) 采用實(shí)驗(yàn)值與真實(shí)值之間差值的二范數(shù)作為評價(jià)指標(biāo)，本文中互耦誤差系數(shù)只有2 個(gè)，幅相誤差系數(shù)則有8 個(gè)。

表1 和表2 是在信噪比為10 dB、快拍數(shù)為100、其余條件與上文保持一致條件下的互耦系數(shù)和幅相系數(shù)估計(jì)結(jié)果。

表1 互耦系數(shù)估計(jì)結(jié)果

表2 幅相系數(shù)估計(jì)結(jié)果

5 結(jié)束語

本文提出了幅相誤差和互耦誤差同時(shí)存在時(shí)基于稀疏貝葉斯算法的求解方法。該方法共包含2 個(gè)階段，第一階段，在空域稀疏的前提下構(gòu)建了誤差存在時(shí)接收信號的超完備模型，得到接收信號的后驗(yàn)概率密度函數(shù)。由于密度函數(shù)具有極強(qiáng)的非線性，直接優(yōu)化十分困難，第二階段采用EM 算法對該函數(shù)進(jìn)行迭代優(yōu)化求解相應(yīng)參數(shù)。此外，本文還推導(dǎo)了相應(yīng)參數(shù)的克拉美羅界。仿真實(shí)驗(yàn)結(jié)果表明，本文方法測向性能優(yōu)于OGSBI、rootSB 和L1SVD 中的方法，同時(shí)能夠?qū)Ψ嗾`差和互耦誤差進(jìn)行有效的校正。