楊 磊 夏亞波 毛欣瑤 廖仙華 方 澄 高 潔

(中國民航大學天津市智能信號與圖像處理重點實驗室 天津 300300)

1 引言

逆合成孔徑雷達(Inverse Synthetic Aperture Radar, ISAR)作為一種主動的微波探測手段，可實現(xiàn)針對非合作運動目標的全天時、全天候、遠距離、高分辨成像探測。ISAR通過發(fā)射并接收大帶寬信號獲得距離向高分辨，通過在相干積累時間(Coherent Processing Interval, CPI)內(nèi)合成一定目標視角變化的多脈沖回波獲得方位向高分辨[1,2]。典型的ISAR系統(tǒng)常放置于地面，用于監(jiān)測空中威脅目標。由于空中目標相對空域稀疏，ISAR接收回波信號具有稀疏特性，可利用壓縮感知(Compressed Sensing, CS)理論，以遠低于奈奎斯特采樣率的回波數(shù)據(jù)恢復ISAR成像結(jié)果[3]。壓縮感知理論由Donoho[4]提出，該理論基于目標稀疏性假設，可實現(xiàn)圖像或信號的稀疏重構(gòu)，從而建立了由低維觀測數(shù)據(jù)精確恢復高維目標信號的理論體系[5]，在雷達信號處理領域應用廣泛。

壓縮感知理論框架下，在目標滿足稀疏性的前提時，可通過有限的低維觀測數(shù)據(jù)恢復高分辨的高維信號[6,7]。其中一類經(jīng)典的信號恢復算法為貪婪算法，該類算法主要利用最小均方誤差的思想，通過人為設置閾值，在每個迭代步驟中選擇局部最優(yōu)解，并在算法收斂后獲得全局最優(yōu)稀疏解。文獻[8]將貪婪算法應用于稀疏信號恢復，從理論方面驗證了此類方法的適用性。盡管該類算法復雜程度較低，運算速度快，但人為設置閾值的方式并不總能保證算法收斂后的解為最優(yōu)稀疏解，有時只是局部最優(yōu)，計算結(jié)果穩(wěn)定性欠佳。

另一類稀疏恢復算法為凸優(yōu)化類算法，該類算法是通過利用優(yōu)化策略對稀疏信號進行有效恢復，主要解決關于目標的 ?0范數(shù)最優(yōu)化問題。已有數(shù)學推導證明?0范數(shù)優(yōu)化問題為非確定性問題(Non-deterministic Problem, NP)[9]，可用最小化?1范數(shù)方法進行近似計算， ?1范數(shù)可認為是最緊致的近似稀疏求解算法。典型的優(yōu)化類算法包括凸優(yōu)化(Con-VeX, CVX)[10]和迭代軟閾值算法[11]等，其中CVX算法基于不斷的優(yōu)化設計最終得到最優(yōu)解，是一種標準的凸優(yōu)化算法，在工程實驗中應用廣泛。但在高維數(shù)據(jù)情況下，典型的CVX算法面臨著運行性能低，耗時過長的缺點。而迭代軟閾值算法通過對目標函數(shù)進行近似替代，并將軟閾值與優(yōu)化策略相結(jié)合，提高了算法計算效率和收斂性能[11]，但由于迭代過程中需要手動調(diào)節(jié)迭代步長，算法自主學習程度不高。文獻[12]將正則化范數(shù)與凸優(yōu)化算法結(jié)合，用于SAR超分辨成像，證明了該類方法在效率和性能方面的優(yōu)勢，但繁瑣的手動調(diào)整正則項系數(shù)步驟嚴重限制了該類算法在超分辨成像的精度及實用性。

還有一類稀疏信號恢復算法是貝葉斯類算法[13,14]，該類算法主要結(jié)合貝葉斯機器學習的思想對信號進行重建。典型的貝葉斯類算法包括變分貝葉斯(Variational Bayes, VB)推斷算法和期望最大化(Expectation-Maximization, EM)算法等。文獻[15]將貝葉斯機器學習結(jié)合應用于SAR地面動目標成像，選取拉普拉斯分布進行先驗建模，利用VB-EM算法進行貝葉斯推理計算，并通過實測數(shù)據(jù)獲得了高精度的恢復結(jié)果。貝葉斯類算法的優(yōu)點在于建模過程應用靈活，而且能提供參數(shù)的區(qū)間估計。但傳統(tǒng)基于貝葉斯機器學習的成像算法在貝葉斯推理過程中面臨復雜的參數(shù)調(diào)整步驟和高維的矩陣逆運算的問題，計算效率較低。

基于以上分析，本文提出一種基于貝葉斯Lasso的稀疏ISAR超分辨成像算法，研究建立凸優(yōu)化Lasso模型與貝葉斯概率模型之間的等價聯(lián)系，解決凸優(yōu)化稀疏算法中面臨的繁瑣手動調(diào)節(jié)正則項系數(shù)的問題。本文通過在分層貝葉斯框架下實現(xiàn)對稀疏Laplace先驗的概率建模，同時考慮對優(yōu)化過程中相應?1范數(shù)正則項稀疏調(diào)節(jié)系數(shù)建立條件概率依賴模型，并利用馬爾科夫鏈蒙特卡洛(Markov Chain Monte Carlo, MCMC)采樣算法從數(shù)據(jù)中自動學習產(chǎn)生正則項系數(shù)，避免繁瑣手動調(diào)節(jié)步驟，最終實現(xiàn)對目標超分辨ISAR成像的精確恢復。針對高維情況下超參數(shù)難以求解的問題，本文采用吉布斯(Gibbs)采樣的方法對建模的參數(shù)進行統(tǒng)計采樣估計。Gibbs采樣方法是一種經(jīng)典的MCMC算法，在統(tǒng)計機器學習方面應用廣泛[16]。本文將貝葉斯Lasso模型與Gibbs采樣方法結(jié)合，避免了復雜的參數(shù)調(diào)整步驟，實現(xiàn)了算法自動化程度提升，增強了稀疏恢復算法的穩(wěn)健性，提高了算法實用性和恢復精度。本文實驗部分應用多種ISAR仿真和實測數(shù)據(jù)，通過與稀疏貝葉斯學習(Sparse Bayesian Learning, SBL)算法和CVX算法進行對比實驗，從抑噪和超分辨性能等不同角度定性的分析了所提算法在性能上的優(yōu)勢，驗證所提算法在ISAR圖像恢復方面的高精度性和實用性。

2 回波信號模型

ISAR成像過程中，雷達保持靜止狀態(tài)，非合作目標通過相對運動形成虛擬孔徑。根據(jù)ISAR成像原理可將目標相對雷達視線方向的運動分為平動和轉(zhuǎn)動，其中平動會造成成像時目標散焦，需要進行平動補償，轉(zhuǎn)動提供了方位向分辨率，是成像過程中需要利用的關鍵部分。經(jīng)過包絡對齊和相位自校準后，ISAR成像可等效為轉(zhuǎn)臺幾何模型如圖1所示，其中目標飛機從A位置運動到B位置。在目標運動過程中，目標觀測角發(fā)生變化，用 θ(t)表示，其中 t為慢時間。坐標系的原點O與目標中心重合，Y軸與雷達觀測目標視線一致，X軸與Y軸垂直，雷達與目標中心點之間的距離以R (t)表示，初始距離設為 R0，則目標散射點與雷達之間的距離可表示為

其中( Xl,Yl)為目標散射點的位置坐標。

當CPI較大時，ISAR成像分辨率高，但由于高階運動分量的影響，目標運動復雜多變，處于非勻速運動狀態(tài)，使得在長CPI內(nèi)對ISAR的成像分析面臨很大的困難。當CPI較小時，ISAR成像積累時間短，目標可認為處于勻速運動狀態(tài)，目標運動狀態(tài)簡單，因此成像過程易于分析。但由于CPI小時，雷達接收回波數(shù)量過少，會造成成像分辨率低，難以滿足ISAR超分辨的需求。針對此問題，本文利用壓縮感知理論，在有限CPI內(nèi)進行稀疏成像，實現(xiàn)ISAR超分辨，保證最終成像分辨率[17,18]。在CPI較小時雷達觀測角可表示為θ (t)=ωt，其中ω為雷達轉(zhuǎn)動角速度。由于有限CPI內(nèi)雷達回波積累時間相對較短，此時式(1)可進一步表示為

圖1 ISAR成像示意圖

3 貝葉斯Lasso

觀察如式(4)所示的ISAR信號模型，目標成像問題可以解釋為關于散射場 X的線性方程組求解問題。由于ISAR接收回波有時是不完整的，在強噪聲和干擾情況下關于目標 X的線性方程組求解的數(shù)學模型一般是非齊次欠定的，進而造成了解的非唯一性和計算復雜性。針對這個問題，可利用壓縮感知技術提供解決途徑，在稀疏先驗假設條件下，將非齊次欠定線性方程組的求解轉(zhuǎn)變?yōu)閷asso問題進行求解，數(shù)學模型如式(6)所示

為得到散射場 X的后驗概率分布，需要選取符合目標特征的先驗分布。由于ISAR目標具有稀疏性的特點，本文利用拉普拉斯先驗分布對目標場景進行先驗建模。拉普拉斯分布具有良好的重尾特性，可以有效地模擬目標的稀疏特征[7]。但由于拉普拉斯先驗分布與式(7)的高斯似然函數(shù)不共軛，因此無法直接計算后驗分布的閉合解析解。針對這個問題，本文采用層級(hierarchical)貝葉斯模型對先驗進行概率建模，通過引入超參數(shù)構(gòu)建層級之間變量的緊密聯(lián)系，進而實現(xiàn)ISAR成像目標 X的有效估計計算。

由于ISAR成像回波數(shù)據(jù)為復數(shù)形式，本文在分層模型的第1步將 X 建模成依賴于超參數(shù)α 的復高斯分布，建模過程如式(8)所示

其中 α表示方差倒數(shù)矩陣，且α 的第m n個元素αmn>0。由于高斯函數(shù)固有的平滑特性，式(8)不屬于稀疏分布范疇。為促進先驗的稀疏特性并簡化后驗分布的推導，本文選取與高斯分布成對共軛的伽馬(Gamma)分布作為超參數(shù)α 的先驗，并在分層模型的第2步將超參數(shù)α 建模成如式(9)的形式

為提高ISAR成像精度和分辨率，并實現(xiàn)參數(shù)的自動調(diào)整，本文通過對正則化系數(shù) λ設定一個先驗分布并進行分層貝葉斯建模，將懲罰項系數(shù) λ作為超參數(shù)引入分層貝葉斯中。同時求取 λ相應的后驗概率密度分布，進而將懲罰系數(shù)與后驗概率密度分布相聯(lián)系，實現(xiàn)參數(shù)自主學習。為減少外界人為因素的影響，本文將 λ的先驗建模為Gamma分布，概率模型如式(10)所示

圖2 貝葉斯層級模型(有向無環(huán)圖，DAG)

上述貝葉斯分層中參數(shù)與超參數(shù)的關系可用有向無環(huán)圖(DAG)進行表示，對應貝葉斯模型如圖2所示。其中 Y 表示ISAR接收回波數(shù)據(jù)， X為待恢復的成像目標， α和λ 為引入的超參數(shù)。 X和α 的緊密結(jié)合能夠促進先驗的稀疏性，而對 λ的概率建模則將正則項系數(shù)引入層級模型，為參數(shù)自學習提供了條件，虛線方框中的a, b, c, d 和η 均為固定參數(shù)。

常見文獻中對式(12)后驗分布的求解方法可分為兩種，一種是利用VB的方法進行近似求解，另外一種是利用基于統(tǒng)計采樣的MCMC采樣方法進行估計計算。在高維情況下，采用VB的方法進行后驗求解面臨繁雜的參數(shù)調(diào)整和矩陣逆運算，算法自學習能力不強，需要人工手動介入，計算效率較低。而MCMC統(tǒng)計采樣方法依賴馬爾科夫鏈進行遞歸采樣，在進行貝葉斯推理的過程依賴各參數(shù)的條件后驗概率分布，是一個數(shù)據(jù)驅(qū)動的過程，無需過多的參數(shù)調(diào)整步驟，在算法自學習方面具有很大的優(yōu)勢。為提高算法性能和成像精度，減少復雜的人為介入?yún)?shù)調(diào)整，本文在貝葉斯層級模型下建立各參數(shù)依存的概率分布關系，并利用經(jīng)典的Gibbs采樣算法進行后驗概率估計，最終達到對目標散射場 X的高精度恢復。

4 Gibbs采樣算法

Gibbs采樣算法是MCMC方法中一種經(jīng)典的采樣方式，其以變量的后驗概率分布為基礎，構(gòu)造具有馬爾可夫統(tǒng)計特性的轉(zhuǎn)移矩陣，并以轉(zhuǎn)移矩陣為采樣核構(gòu)建采樣鏈，進而對變量及其統(tǒng)計特征進行采樣估計。Gibbs采樣算法通過將高維度計算問題進行分解，轉(zhuǎn)換為低維度問題進行運算，并且能夠?qū)﹄S機采樣樣本充分保留，具有高接受率。因此本文利用Gibbs采樣算法對所提出的貝葉斯Lasso層級模型中的變量和超參數(shù)進行統(tǒng)計采樣計算，并對參數(shù)依賴的后驗分布進行估計，實現(xiàn)對待恢復目標的高精度求解。

在應用Gibbs采樣算法進行采樣計算之前，需要對各參數(shù)服從的后驗概率密度分布進行推導。首先根據(jù)式(11)中噪聲精度 β的先驗分布與式(7)中的似然函數(shù)可得精度β 的后驗分布如式(15)所示

由于采樣算法的馬爾可夫鏈需要一定的迭代次數(shù)才能達到收斂狀態(tài)，收斂后最終得到的樣本就是服從后驗概率分布的樣本估計值，而收斂之前的迭代過程稱為老化期(burn-in period)。處于老化期的采樣結(jié)果具有較大的誤差，在進行估算過程中需對這部分數(shù)據(jù)進行舍棄。表1中步驟(2)設置了老化期門限 T，經(jīng)過一定的迭代后由Gibbs采樣器得到的樣本逐漸收斂到目標分布，此后由采樣得到的樣本可用于參數(shù)和變量的估計計算。為降低誤差的影響，對收斂后的采樣樣本進行均值求解可得出最終X的采樣估計值，即為待恢復的ISAR成像目標散射 場 X。

表1 貝葉斯Lasso算法流程

5 實驗驗證

5.1 仿真數(shù)據(jù)驗證

為進一步驗證貝葉斯Lasso算法的優(yōu)越性和通用性，本文利用仿真數(shù)據(jù)Mig-25對所提模型及算法的成像性能進行驗證。通過在不同降采樣率和信噪比下對目標場景進行恢復重建，并與由距離多普勒(Range-Doppler, RD), CVX和SBL算法得到的成像結(jié)果進行對比，從而驗證所提算法的優(yōu)越性和穩(wěn)健性。圖3給出了信噪比為5 dB降采樣率變化時，通過不同方法所得的實驗結(jié)果。其中圖3第1列和第2列分別表示降采樣率為0.25和未進行降采樣時分別利用RD, SBL, CVX和本文貝葉斯Lasso算法進行處理后的恢復結(jié)果。從圖3中可以看出利用RD方法進行成像時，成像結(jié)果背景噪聲大，目標不明顯。而SBL, CVX和貝葉斯Lasso 3種算法則實現(xiàn)了對噪聲的抑制，以及對成像效果的增強，其中又以貝葉斯Lasso算法的效果最好，從而證明了本文貝葉斯Lasso算法的實用性。

圖3 信噪比5 dB不同降采樣率成像結(jié)果

為驗證本文貝葉斯Lasso算法在抑噪方面的優(yōu)越性能，圖4給出了降采樣率為0.5條件下信噪比變化時通過不同算法恢復的結(jié)果。其中圖4第1列和第2列分別表示信噪比為0 dB和信噪比為5 dB時分別利用RD, SBL ,CVX和該文貝葉斯Lasso算法恢復的成像結(jié)果。從圖4中可以看出RD重建的結(jié)果背景噪聲很大，目標不明顯，利用SBL和CVX算法恢復的圖像雖然有部分抑噪效果，但目標丟失了多數(shù)特征。與這3種算法相比，采用該文貝葉斯Lasso算法恢復的目標散射場飛機輪廓特征突出，背景噪聲消除干凈，性能最優(yōu)。

5.2 實測數(shù)據(jù)驗證

圖4 降采樣率0.5不同信噪比成像結(jié)果

圖5 信噪比5 dB不同降采樣率成像結(jié)果

實測部分采用Yak-42數(shù)據(jù)[12]進行實驗分析，圖5給出了信噪比5 dB條件下降采樣率分別為0.25和未降采樣的恢復結(jié)果。從圖5中可以看出隨著采樣率的降低，RD恢復的目標淹沒在噪聲下，SBL算法恢復的結(jié)果目標丟失部分特征，CVX算法對噪聲的抑制較弱，恢復的目標背景噪聲斑點較多。而本文貝葉斯Lasso算法在低降采樣率下保留目標大部分特征的同時，對噪聲也起到了很好的抑制作用，總體性能優(yōu)于SBL和CVX算法。

為驗證本文貝葉斯Lasso算法在抑噪方面的優(yōu)勢，圖6給出了降采樣為0.5信噪比分別為0 dB和5 dB時不同算法重建的結(jié)果。其中在信噪比為0 dB時，RD算法得到的結(jié)果目標基本識別不清，并且背景噪聲很大，SBL恢復的目標丟失大部分特征，并且機頭部分基本無法識別，而CVX算法抑噪能力不足。從圖6中可以看出，本文采用的貝葉斯Lasso算法在抑制噪聲的同時保留了目標的多數(shù)特征，對目標進行了有效恢復，整體性能明顯優(yōu)于其他3種算法。

圖6 降采樣0.5不同信噪比成像結(jié)果

5.3 地面動目標實驗驗證

本部分通過對地面動目標進行成像分析，驗證貝葉斯Lasso算法在實測數(shù)據(jù)的性能優(yōu)勢。地面動目標成像實驗采用Gotcha-GMTI數(shù)據(jù)集[21]，該數(shù)據(jù)集由美國空軍實驗室發(fā)布，經(jīng)由三通道機載雷達獲取，雷達工作于X波段，作用距離10 km，發(fā)射信號帶寬約640 MHz, PRF約2100 Hz。由于在對地面動目標進行成像過程中目標的運動特性會引起成像的偏移與散焦，本文利用雜波相消并采用呂氏變換方法(LV’s Distribution, LVD)[22]對調(diào)頻發(fā)射信號中的多普勒頻率和線性調(diào)頻率進行協(xié)同估計，以改善動目標成像的聚焦和抑噪性能，提高成像分辨率。

圖7(a)給出了采用RD變換后的參考成像，由圖可見背景噪聲抑制效果不明顯。圖7(b)、圖7(c)分別為利用CVX和SBL算法的成像結(jié)果，從圖中可以看出兩種算法在抑噪方面均起到一定作用，但由于CVX算法難以適應不同距離單元的不同稀疏度，在成像過程中失去了部分目標特征，而SBL算法求解精度有限，目標恢復的效果不佳。圖7(d)為采用本文貝葉斯Lasso算法得到的成像結(jié)果，與其他3種算法相比，由于貝葉斯Lasso算法無需人工干預參數(shù)選擇，能夠自主地適應空變的稀疏度，因此得到的成像結(jié)果信噪比更高，目標散射特征更清晰，從而驗證了本文貝葉斯Lasso算法在成像性能方面的優(yōu)越性。

5.4 相變圖實驗分析

為了定性驗證本文貝葉斯Lasso算法性能的優(yōu)越性，本文利用Donoho等人[23]提出的相位變化圖(Phase Transition Diagram, PTD)進行實驗分析。PTD是對算法精確性和穩(wěn)定性進行評估的一種經(jīng)典方式，其通過多次蒙特卡洛獨立重復實驗，獲得恢復圖像與參考圖像相關量的平均值，并利用均方誤差比較確定算法恢復的準確率，從而對算法的恢復精度進行定性分析。PTD將相變區(qū)域分為不可重建和可重建區(qū)域，并利用不同區(qū)域?qū)λ惴ǖ幕謴托阅苓M行直觀的描述。本文在對算法進行評估分析時，通過100次蒙特卡洛實驗得到重建圖像與目標圖像的平均相關值，并利用歸一化均方誤差(Normalized Mean Square Error, NMSE)進行判別圖像恢復的準確率，NMSE的表達式為

圖7 地面動目標成像

其中 Xk代表第k 次蒙特卡洛實驗所得的恢復結(jié)果，X表示參考目標圖像，閾值設定為0.0001，當NMSE<0.0001時認為圖像重建成功。

圖8 為利用相變圖對S B L,C V X 和貝葉斯Lasso算法進行相變分析所得的結(jié)果。圖8中各子圖的右上方區(qū)域代表可重建區(qū)域，左下方區(qū)域代表不可重建區(qū)域。從圖8中可以看出利用貝葉斯Lasso算法得到的相變圖可重建區(qū)域要明顯大于利用SBL算法得到的結(jié)果，并且與CVX得到的相變圖非常接近，從而可以得出本文貝葉斯Lasso算法對目標的成像恢復能力要明顯優(yōu)于SBL算法，并且能夠達到與CVX標準相近的性能。盡管CVX算法損耗大，運算時間長，并且高維情況下算法實用性不高，但由于CVX是一種標準的凸優(yōu)化恢復算法，其可作為衡量其他算法性能優(yōu)劣的參考方法。從圖8(b)和圖8(c)中可以看出由貝葉斯Lasso算法和CVX算法得到可重建區(qū)域非常接近，表明兩者在算法恢復性能方面具有很高的相似度，從而驗證了貝葉斯Lasso算法的優(yōu)越性。

圖8 相變圖

6 結(jié)束語

本文針對凸優(yōu)化ISAR成像算法繁瑣的參數(shù)調(diào)整過程，研究建立凸優(yōu)化Lasso模型與貝葉斯概率模型之間的聯(lián)系。通過結(jié)合貝葉斯理論，建立貝葉斯Lasso模型，并采用分層方法對參數(shù)建立依賴的概率分布關系，在求取條件后驗分布的同時利用Gibbs采樣算法進行參數(shù)后驗估計，實現(xiàn)了ISAR目標高分辨成像。實驗利用相變圖進行定性比較分析，并利用仿真和實測數(shù)據(jù)進行算法性能驗證，證明了本文采用的貝葉斯Lasso算法的有效性和優(yōu)越性。但是，在研究過程中，筆者發(fā)現(xiàn)標量正則項系數(shù)不能很好地模擬目標場景，有必要研究矩陣化的正則項，并實現(xiàn)矩陣正則項系數(shù)的自動調(diào)整，這將是后續(xù)工作研究的重點。