孫延鵬，王藝霖，屈樂樂

(沈陽航空航天大學(xué) 電子信息工程學(xué)院，沈陽 110136)

?

基于貝葉斯壓縮感知的頻率步進(jìn)探地雷達(dá)成像算法

孫延鵬，王藝霖，屈樂樂

(沈陽航空航天大學(xué) 電子信息工程學(xué)院，沈陽 110136)

針對(duì)傳統(tǒng)壓縮感知頻率步進(jìn)探地雷達(dá)成像算法存在計(jì)算量大和對(duì)噪聲和重建正則化參數(shù)敏感的問題，提出一種基于稀疏貝葉斯學(xué)習(xí)的貝葉斯壓縮感知成像算法。該成像算法的核心通過稀疏貝葉斯線性回歸模型中相關(guān)向量機(jī)的學(xué)習(xí)來實(shí)現(xiàn)對(duì)探測(cè)場景反射系數(shù)的重構(gòu)。仿真結(jié)果表明，相比其他的經(jīng)典算法，所提成像算法能夠更好地利用了探測(cè)場景的統(tǒng)計(jì)先驗(yàn)信息，能夠更好地兼顧重構(gòu)精度和計(jì)算效率。

頻率步進(jìn)探地雷達(dá)；貝葉斯壓縮感知；成像算法

頻率步進(jìn)體制是探地雷達(dá)技術(shù)的工作體制之一。頻率步進(jìn)探地雷達(dá)是通過發(fā)射工作頻率按固定間隔逐步上升的連續(xù)波，并對(duì)不同頻率點(diǎn)回波幅度與相位信息進(jìn)行數(shù)據(jù)采集的連續(xù)波雷達(dá)。隨著探地雷達(dá)技術(shù)應(yīng)用的進(jìn)一步推進(jìn)，探地雷達(dá)不斷朝著多通道、多波段和高分辨率等方向發(fā)展[1]，這就使得探地雷達(dá)系統(tǒng)的數(shù)據(jù)采集量大大增加，成像復(fù)雜度大大提高。隨著壓縮感知(Compressed sensing,即CS)理論[2]的提出，這一問題正得到解決，該理論證明[3]:當(dāng)信號(hào)在某個(gè)變換域是稀疏或可壓縮的，稀疏信號(hào)采樣速率不再取決于信號(hào)的帶寬，而主要取決于兩個(gè)基本準(zhǔn)則[4]:稀疏性和非相干性，突破了Nyquist理論對(duì)信號(hào)采樣的限制。

目前CS 理論依舊在不斷地發(fā)展,其重建算法主要包括兩類[5]：一類是以基追蹤算法BP(basis pursuit)[6]為代表l1的范數(shù)線性優(yōu)化算法，該算法是在所有測(cè)量向量的不同組合中尋求最優(yōu)化的解，因此計(jì)算復(fù)雜性較高；另一類是以正交匹配追蹤算法OMP(Orthogonal Matching Pursuit)[7]為代表的貪婪算法，該算法是通過遞歸的方法，對(duì)已選擇的原子集合進(jìn)行正交化從而保證了迭代的最優(yōu)性。

傳統(tǒng)的壓縮感知頻率步進(jìn)探地雷達(dá)存在數(shù)據(jù)計(jì)算量大和對(duì)噪聲敏感的問題。針對(duì)這一問題，本文提出一種基于稀疏貝葉斯學(xué)習(xí)的貝葉斯壓縮感知(Bayesian compressed sensing,BCS)頻率步進(jìn)探地雷達(dá)成像算法[8-9]。該算法基于貝葉斯統(tǒng)計(jì)學(xué)的理論[10]，依據(jù)先驗(yàn)知識(shí)對(duì)系數(shù)向量的每一個(gè)元素賦予一個(gè)先驗(yàn)概率分布用以限制模型的復(fù)雜度，并引入超參數(shù)；通過最大化超參數(shù)的邊緣對(duì)數(shù)似然函數(shù)，獲得每個(gè)系數(shù)向量對(duì)應(yīng)的先驗(yàn)概率分布的參數(shù)值。與傳統(tǒng)CS成像重建算法相比，基于BCS的頻率步進(jìn)探地雷達(dá)成像重建算法能夠較好地兼顧計(jì)算效率和重構(gòu)精度，因此具有較好的成像重建性能。

1 信號(hào)模型

頻率步進(jìn)探地雷達(dá)成像模型[11]如圖1所示，假定發(fā)射天線、接收天線與成像區(qū)域位于同一平面內(nèi)。成像區(qū)域沿水平方向(x軸)和距離方向(z軸)劃分為K×L個(gè)均勻的網(wǎng)格，每一個(gè)網(wǎng)格對(duì)應(yīng)一個(gè)像素點(diǎn)。

假設(shè)測(cè)量過程中共有M個(gè)天線測(cè)量位置，發(fā)射天線在一個(gè)掃描周期中，頻率從f0到fN-l共N個(gè)掃描頻點(diǎn)，探測(cè)區(qū)域總共有P個(gè)實(shí)際散射目標(biāo)數(shù)，則在第m=1,…,M個(gè)位置對(duì)復(fù)基帶信號(hào)進(jìn)行采樣后得到復(fù)信號(hào)如下：

圖1 頻率步進(jìn)探地雷達(dá)成像模型

(1)

ρp表示第p個(gè)點(diǎn)目標(biāo)的復(fù)反射系數(shù)，τp,m表示為電磁波在第m個(gè)天線測(cè)量位置與第p個(gè)點(diǎn)目標(biāo)的雙程傳輸延時(shí)。

為了能夠用壓縮感知的模型求解出地下待成像區(qū)域的2維反射率分布r(K,L)，將r(K,L)通過列堆疊轉(zhuǎn)換成一個(gè)KL×1維的向量，記為wx，如圖2所示。

圖2 網(wǎng)格劃分示意圖

系統(tǒng)在第m個(gè)探測(cè)位置接收到的回波信號(hào)可以表示為一個(gè)基矩陣Bm和反射系數(shù)向量wx相乘的形式，如下：

Sm=Bmwx

(2)

Sm=[s0(f0),…,sm(fN-1)]T為N×1維頻域測(cè)量數(shù)據(jù)向量，[·]T表示轉(zhuǎn)置，sm(fn)為當(dāng)工作頻率為fn時(shí)，天線在第m個(gè)測(cè)量位置采集到的回波數(shù)據(jù)。基矩陣Bm是大小為N×KL的基矩陣，其第j列為：

[Bm]j=[e-j2πf0τm,j,…,e-j2πfN-1τm,j]T

(3)

Sx=Bxwx

(4)

一般來說，反射系數(shù)向量wx是稀疏的，滿足CS理論對(duì)信號(hào)的要求。為了更好地滿足約束等容性條件RIP(Restricted Isometric Property)[18]，我們從M個(gè)天線測(cè)量位置隨機(jī)選擇Q1個(gè)位置，再從選中的測(cè)量位置隨機(jī)選取Q2個(gè)數(shù)據(jù)，這樣構(gòu)成了一個(gè)Q1Q2×NM大小的測(cè)量矩陣Ψ，并且Q1Q2?NM，用測(cè)量矩陣Ψ對(duì)Sx進(jìn)行投影可以得到：

tx=ΨSx=ΨBxwx=Φxwx

(5)

式中，Φx=ΨBx是一個(gè)大小為Q1Q2×KL的投影矩陣。

2 貝葉斯的壓縮感知(BCS)

由于模型中投影矩陣的元素與回波數(shù)據(jù)均為復(fù)數(shù)，無法直接對(duì)其使用BCS算法。所以根據(jù)復(fù)數(shù)的運(yùn)算法則，對(duì)式(5)中的矩陣和向量做如下變化：

(6)

(7)

其中Re(·)與Im(·)分別代表復(fù)數(shù)的實(shí)部和虛部?？紤]到測(cè)量噪聲，由式(5)(6)(7)可以得到：

t=Φw+n

(8)

其中測(cè)量噪聲n的元素可以近似為均值是0，方差是σ2的正態(tài)分布。因此可以寫出t的高斯似然函數(shù)[12]：

(9)

‖·‖表示l2范數(shù)。這樣就將w的求解轉(zhuǎn)換為w是稀疏的這一先驗(yàn)條件下的線性回歸問題。通過賦予參數(shù)w={ωi}i=1,2KL一個(gè)先驗(yàn)分布，可以避免過擬合情況的發(fā)生。設(shè)參數(shù)ωi的先驗(yàn)分布是均值為零,逆方差為?i的正態(tài)分布，所以p(w/a)的密度函數(shù)為：

(10)

其中a=[?1,?2,…,?2KL]T，根據(jù)先驗(yàn)分布與后驗(yàn)分布的共軛關(guān)系，可以得到w的后驗(yàn)分布的形式為多元正態(tài)分布，它的均值向量μ和協(xié)方差矩陣∑如公式(11)、(12)(設(shè)?0=σ-2)：

μ=?0∑ΦTt

(11)

∑=(?0ΦTΦ+A)-1

(12)

其中A=diag(a)，其中diag(·)表示以其中的向量作為對(duì)角元素的對(duì)角矩陣。通過觀察可以發(fā)現(xiàn)對(duì)均值向量μ和協(xié)方差矩陣∑的求解也就是對(duì)超參數(shù)?0和a的求解，在相關(guān)向量機(jī)[13]框架下，通過第二型最大似然估計(jì)方法(ML)[14]或者最大期望算法(EM)可以求出超參數(shù)的極值。本文采用ML的方法得到：

(13)

(14)

3 快速貝葉斯壓縮感知算法

針對(duì)BCS算法計(jì)算復(fù)雜度大的問題，文獻(xiàn)[14]提出了一種快速BCS算法，該算法對(duì)w進(jìn)行邊緣化處理，得到關(guān)于超參數(shù)的邊緣似然函數(shù)：

L(a,σ2)=logp(t/a,σ2)

(15)

其中C=σ2I+ΦA(chǔ)-1ΦT，并以此作為著手點(diǎn)，把邊緣似然函數(shù)分解成兩部分：

(16)

通過極值定理可知，當(dāng)?i取如下的值時(shí)，可以讓l(?i)達(dá)到最大值：

(17)

(18)

當(dāng)?i使每一個(gè)l(?i)都取最大值的時(shí)候，那么對(duì)應(yīng)的a=[?1,?2,…,?2KL]T就是使邊緣似然函數(shù)取最大值的超參數(shù)。

通過以上分析，快速貝葉斯壓縮感知算法的本質(zhì)就是根據(jù)公式(17)(18)，使邊緣似然函L(a,σ2)最大化。操作步驟如下：

(2)計(jì)算均值向量μ和協(xié)方差矩陣∑和所有投影向量對(duì)應(yīng)hi和qi。

(4)根據(jù)公式(11)(12)(14)重新計(jì)算μ，∑，σ2，hi和qi的值，并返回第三步，直到兩次的迭代結(jié)果小于門限值或迭代次數(shù)超過一定范圍，退出循環(huán)。

貝葉斯壓縮感知與其他算法的不同之處在于所求得的稀疏系數(shù)最終表現(xiàn)為某概率密度函數(shù)的形式。例如本文的稀疏系數(shù)為均值向量是μ和協(xié)方差矩陣是∑的多元正態(tài)分布，因此可以認(rèn)為均值向量就是對(duì)應(yīng)的稀疏系數(shù)，而且可以將協(xié)方差矩陣作為誤差棒來刻畫稀疏系數(shù)的準(zhǔn)確程度。

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

在這部分，我們對(duì)BCS成像算法進(jìn)行了仿真，并分別采用基追蹤算法(BP)、正交匹配追蹤算法(OMP)及BCS算法重構(gòu)稀疏目標(biāo)。以驗(yàn)證BCS算法的性能表現(xiàn)。BP算法是將凸優(yōu)化問題轉(zhuǎn)化為線性規(guī)劃(Linear Programming)問題求解[15]。信號(hào)重建誤差小、精度高，但相對(duì)比較復(fù)雜，計(jì)算復(fù)雜度較高，重建速度比較慢。OMP算法為基于l0范數(shù)的貪婪算法中的匹配追蹤系列算法，它是基于信號(hào)稀疏分解而構(gòu)造的一類算法[16]。這類算法的計(jì)算簡單、速度較快，但是信號(hào)重建效果較差。

在仿真中，假設(shè)系統(tǒng)發(fā)射電磁波的起始頻率f0=1 GHz，最高頻率fN-1=3 GHz，頻率步進(jìn)點(diǎn)數(shù)N=101，天線測(cè)量位置總數(shù)M=30，信噪比為10 dB。三個(gè)目標(biāo)點(diǎn)的深度和水平位置坐標(biāo)分別為(40,10)、(40,20)、(40,30)(單位cm)。不考慮地面回波的影響下，成像區(qū)域沿水平方向從1 cm到40 cm，沿深度方向從30 cm到50 cm。成像區(qū)域被均勻地劃分為21×40個(gè)大小為1 cm的網(wǎng)格目標(biāo)，成像結(jié)果如圖3所示。

通過表1對(duì)比可以看出OMP算法的速度最快，BP 算法的速度最慢。這主要是由算法本身的性能造成的：BCS 算法和 BP 算法在每次迭代時(shí)都需要計(jì)算矩陣的逆，這使得它們需要更多的運(yùn)算時(shí)間，但是BCS 方法能夠更快地收斂，因此速度明顯比 BP算法更快；OMP 算法每次迭代的主要步驟是計(jì)算殘差與采樣矩陣的內(nèi)積，且迭代次數(shù)主要取決于重構(gòu)目標(biāo)的稀疏度，因此其速度最快。對(duì)于各種算法的重構(gòu)精度，從圖中可以看出，在信噪比為10 dB時(shí)，BP和BCS的重構(gòu)精度相近，但是BP所消耗的時(shí)間最長；OMP的恢復(fù)精度較差。

圖4為在信噪比為10 dB的條件下，測(cè)量個(gè)數(shù)與算法運(yùn)行時(shí)間的關(guān)系，可以看出隨著測(cè)量點(diǎn)的增加，算法運(yùn)行時(shí)間也不斷增加。

從圖5可以看出隨著測(cè)量點(diǎn)的增加，三種算法的均方誤差都有變小的趨勢(shì)，其中BP算法的均方誤差值要比BCS算法的小一些，但相差不多，BP算法和BCS算法的均方誤差要明顯小于OMP算法。

圖3 當(dāng)信噪比為10 dB時(shí)，各算法的重構(gòu)結(jié)果

對(duì)比算法OMPBPBCS時(shí)間(s)0 03135 110 13

圖4 測(cè)量個(gè)數(shù)對(duì)重構(gòu)時(shí)間的影響

圖6所示為當(dāng)采樣點(diǎn)為450時(shí)，噪聲對(duì)各種算法性能的影響。可以看出，當(dāng)噪聲較大對(duì)，BCS算法的重構(gòu)效果較差，不過當(dāng)信噪比大于10dB之后，BCS算法與OMP算法的重構(gòu)效果相近，且要好于BP算法。

圖5 測(cè)量個(gè)數(shù)對(duì)重建均方誤差的影響

圖6 噪聲對(duì)均方誤差的影響

5 結(jié)論

本文主要介紹了一種基于稀疏貝葉斯學(xué)習(xí)的貝葉斯壓縮感知頻率步進(jìn)探地雷達(dá)成像重建算法，所提成像算法采用分層貝葉斯模型，能夠從較少的壓縮采樣數(shù)據(jù)中重建目標(biāo)的位置和反射系數(shù)。仿真結(jié)果表明，相比其他的經(jīng)典算法，所提成像算法能夠更好地利用探測(cè)場景的統(tǒng)計(jì)先驗(yàn)信息，更好地兼顧重構(gòu)精度和計(jì)算效率。本文雖然將貝葉斯壓縮感知的方法應(yīng)用到復(fù)數(shù)域，不過模型中復(fù)數(shù)的實(shí)部和虛部并不是相互獨(dú)立的，采用何種建模能更精確的體現(xiàn)目標(biāo)的稀疏特性，將是今后研究的方向。

[1]T Counts,C Gurbuz.W R Scott.Multistatic ground penetra-ting radar experiments[J].IEEE Transactions on Geoscience and Remote Sensing,2007,45(8):2544-2553.

[2]D L Donoho.Compressed sensing[J].IEEE Transactionson Information Theory,2006,52(4):1289-1306.

[3]S Ji,Y Xue,L Carin.Bayesian compressive sensing[J].IEEE Transactions on Signal Processing,2008,56(6):2346-2356.

[4]楊海蓉,張成,丁大為,等.壓縮傳感理論與重構(gòu)算法[J].電子學(xué)報(bào),2011,39(1):142-148.

[5]方紅,王年,章權(quán)兵.基于稀疏貝葉斯學(xué)習(xí)的圖像重建方法[J].中國圖象圖形學(xué)報(bào),2009,14(6):1064-1069.

[6]D L Donoho,M A Saunders.Atomic decomposition by basis pursuit[J].SIAM Journal on Scientific Computing,1998,20(1):33-61.

[7]白凌云,梁志毅,徐志軍.基于壓縮感知信號(hào)重建的自適應(yīng)正交多匹配追蹤算法[J].計(jì)算機(jī)應(yīng)用研究,2011,28(11):4060-4063.

[8]V N Vapnik.Statistical learning theory[M].Wiley,1998.

[9]S D Babacan,R Monila.Bayesian compressive sensing using Laplace priors[J].IEEE Transactions on Image Processing,2010,19(1):53-63.

[10]茆詩松.貝葉斯統(tǒng)計(jì)[M].北京：中國統(tǒng)計(jì)出版社，1999.

[11]屈樂樂,方廣有,楊天虹.壓縮感知理論在頻率步進(jìn)探地雷達(dá)偏移成像中的應(yīng)用[J].電子與信息學(xué)報(bào),2011,33(1):21-26.

[12]徐建平,皮亦鳴,曹宗杰.基于貝葉斯壓縮感知的合成孔徑雷達(dá)高分辨成像[J].電子與信息學(xué)報(bào),2011,33(12):2864-2868.

[13]C Cortes,N Vapnik.Support vector networks[J].Machine Learning,1995,20(3):273-297.

[14]M E Tipping,A C Faul.Fast marginal likelihood maximization for sparse Bayesian models[C].Proceedings of the Ninth International Work-shop onArtificial Intelligence and Statistics.Key West,FL,2003.

[15]李少東,楊軍,胡國旗.一種改進(jìn)的壓縮感知信號(hào)重構(gòu)算法[J].信號(hào)處理,2012,28(5):744-749.

[16]J A Tropp,A C Gilbert.Signal recovery from randommeasurements via orthogonal matching pursuit[J].IEEE Transactions on Information Theory,2007,53(12):4655-4666.

[17]郭鵬.基于貝葉斯壓縮感知的自適應(yīng)測(cè)量算法[J].計(jì)算機(jī)工程與應(yīng)用,2013,49(9):200-203.

[18]Candes E J,Romberg J,Tao T.Robust uncer-tainty principles:exact signal reconstruction from highly incomplete frequency information[J].IEEE Transactions on Information Theory,2006,52(2):489-509.

(責(zé)任編輯：劉劃 英文審校：林嘉)

Stepped-frequency GPR imaging algorithm based on bayesian compressive sensing

SUN Yan-peng，WANG Yi-lin，QU Le-le

(College of Electronic and Information Engineering,Shenyang Aerospace University,Shenyang 110136,China)

The stepped-frequency ground penetrating radar(GPR) based on the traditional compressive sensing is rather computationally intensive and sensitive to the regularization parameter.To solve the above problem,animaging algorithm based on Bayesian compressive sensing(BCS) is proposed in the paper.Within the sparse Bayesian linear regression model,the proposed BCS-based imaging algorithm uses the relevance vector machine to reconstruct the reflectivity of the investigation scene.The numerical simulation results show that the proposed imaging algorithm can take advantage of the prior statistical information of the scene reflectivies and achieve both reconstruction accuracy and computation efficiency compared with the traditional CS-based imaging methods.

stepped-frequency ground penetrating radar;Bayesian compressive sensing;imaging algorithm

2015-07-16

國家自然科學(xué)基金(項(xiàng)目編號(hào)：61302172)，遼寧省自然科學(xué)基金(項(xiàng)目編號(hào)：2014024002)，遼寧省博士啟動(dòng)基金(項(xiàng)目編號(hào)：20121035)

孫延鵬(1973-)，男，山東茌平人，教授，主要研究方向：航空電子信息系統(tǒng)，E-mail：syp@syiae.edu.cn。

2095-1248(2015)05-0068-06

TN958

A

10.3969/j.issn.2095-1248.2015.05.010

信息科學(xué)與工程