一種改進廣義正交匹配追蹤的DOA 估計方法

刁弘揚，胡洲勇，禹永植

1.哈爾濱工程大學 信息與通信工程學院，黑龍江 哈爾濱 150001

2.北京遙感設備研究所，北京 100854

信號的波達方向(direction of arrival,DOA)估計在陣列信號處理領域占據重要地位，廣泛應用于雷達、醫學診療等領域，近年來學者們進行了大量的研究，取得了非常顯著的成果[1?3]。以MUSIC[4]和ESPRIT[5]算法為代表的經典子空間分解類算法，實現簡單且分辨率高，并且在此基礎上發展出了多種改進算法[6?8]，但是這類算法運算量大，需要在信噪比高且快拍數多的情況下才能得到較好的估計性能。壓縮感知[9?10]理論是一個新興起的領域，它能夠突破奈奎斯特采樣定理的限制，壓縮感知理論中的很多重構算法都被用于DOA 估計中，其中平滑l0范數重構算法[11?12]、

l1?SVD算法[13?14]和貪婪算法是廣泛應用的方法。平滑l0范數重構算法采用最速下降法[15]和梯度投影的原理，構造一個函數來近似l0范數，但是此算法不能保證重構的精確度。l1?SVD算法利用l1范數重構算法對信號進行重構，估計出信號的DOA，能夠在低快拍下實現DOA 估計；但是該算法重構過程復雜，會導致較高的計算復雜度，并且在求解過程中，需要人工去設置所需參數，參數選取不好將會對最終的結果產生較大的影響。常用的貪婪算法有匹配追蹤(matching pursuit,MP)、正交匹配追蹤(orthogonal matching pursuit,OMP)[16]、正則OMP[17]等，該類算法需要找到測量矩陣中能夠提供最大能量的某列原子，然后一直進行迭代，直到滿足迭代條件，得到最優解。OMP 算法是一種常用的重構算法，研究者們提出了多種改進算法[18?19]，但是OMP 算法運算時間長，DOA 估計精度有待進一步提高。

針對稀疏陣列中貪婪算法DOA 估計求解問題，本文引入廣義正交匹配追蹤(generalized orthogonal matching pursuit,GOMP）用于DOA 估計，仿真表明直接將GOMP 算法用于DOA 估計較OMP算法性能有所提高；為了進一步提高DOA 估計的精度，引入最優化理論中的最速下降法，提出基于最速下降法的改進GOMP 算法。

1 基于稀疏表示的DOA 估計模型

假設有一個由M個陣元組成的均勻線陣，陣元間距為d，有K個方向分別為θi(i=1,2,···,K)的遠場窄帶信號xk(t),k∈{1,2,···,K}。通過網格劃分將空間信號稀疏化，整個空間范圍等角度均勻劃分為N份，則所有可能的方向可以表示為θi(i=1,2,···,N)。將空間信號進行稀疏化，只有K個方向上存在真實的目標信號源，其中，N遠遠大于實際信號源的數量K。則空間中所有的方向集合都與一個導向矢量相對應，導向矢量為A=[a(θ1),a(θ2),···,a(θN)]，導向矢量也稱為陣列流型矩陣，則t時刻陣列接收的信號為

式中：x(t)=[x1(t),x2(t),···,xN(t)]T為空間入射信號矢量，其稀疏度為K(K?N)，x(t)中只有K個非零元素，其余N?K個為零元素，非零元素所在的位置即為真實信號源的目標方向；n(t)=[n1(t),n2(t),···,nM(t)]T是均值為0，方差為 σ2的高斯白噪聲；y(t)=[y1(t),y2(t),···yM(t)]T為陣列接收信號矢量。

稀疏表示的DOA 估計就是在已知陣列接收矢量y以及陣列流型矩陣A的條件下，找出x(t)中所有非零元素所在的位置集合，該位置集合即是θi(i=1,2,···,N)中實際信號源所在的角度。

2 基于改進GOMP 算法的DOA 估計

2.1 OMP 算法

由式(1)、(2)可知，x在空間上是K階稀疏的，只有K個非零的元素，而且非零元素對應的位置即為θi(i=1,2,···,N)中實際信號源的角度。

為了找出x中非零元素所在的位置，需要確定陣列流型矩陣A中的有效矢量。OMP 算法通過貪婪搜索的方式來搜索有效列。在每次迭代過程中，找出A中 與陣列接收信號y中剩余部分最相關的一列，剔除該列對于y的貢獻，標記該列并更新殘差。當完成K次迭代后，即可找出參與對y的線性組合中列的標記集合，最后用該標記集合，在集合θi(i=1,2,···,N)中找到信號源所在的真實角度位置。

2.2 GOMP 算法

基于OMP 算法的DOA 估計方法，耗時比較長且估計的精度不高，為了提高DOA 估計精度，采用廣義正交匹配追蹤算法來求解DOA 估計。

OMP 算法和GOMP 算法均是通過貪婪迭代從陣列流型矩陣A中找到能夠提供最大能量的某些列原子，代價函數選用y殘差最小。2 種算法的區別在于進行重構時，OMP 算法每次迭代過程中只選擇與y中剩余部分最相關的一列，而GOMP算法在每次迭代過程中會選擇與y中剩余部分最相關的K列，然后在所有選出的可能的坐標集合中進行搜尋，在索引坐標集合中加入使估計能量值最大的坐標。相對于OMP 算法，GOMP 算法的搜索選擇方式能夠更好地保證從信號能量相近或者不完全和不精確的陣列接收信號中恢復信號，能更準確地找出信號源所在的真實角度位置。

2.3 最速下降法的相關理論

最速下降法基于最優化理論[20]，從梯度的角度來求解最優化問題。函數f(x)的泰勒展開式為

式中：o(x)是高階無窮小量，其在xk的梯度為gk=?f(xk)。記x?xk=akdk，其中，ak是步長，是一個大于零的數值，dk為梯度的方向，則泰勒展開式為f(x)=f(xk)++o(x)。在＜0時f(x)＜f(xk)，此時的函數值是下降的，dk則是下降的方向。

由柯西?施瓦茲不等式可知：

當且僅當dk=gk時，等號成立，最大。所以當且僅當dk=?gk時，最小，以此為搜索方向，此時f(x)是下降最快的。

最速下降法在每次迭代時選擇一個合適的步長，使得目標函數每次減小的程度最大，即每次得到的值和目標函數值之間的差距最小。用f(x)代表目標函數，迭代步長為

梯度的迭代為

因此計算步長的公式就是為找到最好的點xk+1，函數f(x)在這點可以取得最小值，此時f(xk?agk)的極小值就是此次迭代的步長。

最速下降法的過程可以表示為在每次迭代中，從xk出發，按照梯度的負方向去尋找最好的步長值，根據尋找到的步長值來確定下一次迭代過程的出發點xk+1。

將函數的類型擴展至二次型，設函數為

式中：Q是 對稱正定陣；b∈Rn；x∈Rn；令

迭代表示為

步長為

在gk=0時迭代過程停止，在不為零時，設φk(a)=f(xk?agk)，由局部最小點的一階必要條件得到：

2.4 改進GOMP 算法

將最速下降法和GOMP 相結合得到改進GOMP算法，改進算法在信號重構時不使用計算量大的最小二乘法，而是沿著負梯度方向進行搜索而得到最優解。

改進算法中迭代選擇原子及更新索引集的過程和GOMP 算法相一致，在恢復信號過程中，GOMP等貪婪算法均使用最小二乘法進行重構，其矩陣公式為x=(ATA)?1ATy，由于需要進行矩陣的求逆運算，導致運算量很大，會消耗很長的時間；而最速下降法雖然需要進行多次運算，但都是基本運算，相對而言計算量較小，所要耗費的時間也會更少。因此改進算法使用最速下降法代替最小二乘法來恢復信號。

用Φ={φk}表示測量矩陣，此測量矩陣即為DOA 估計模型中的陣列流型矩陣，y表示要恢復的信號，xi表示第i次循環得到的近似信號，ri表示殘差，Λi表示索引集合，Ψ 表示選擇原子的列向量集合，先得到矩陣的負梯度方向為gi=ΨiTri，然后計算每次迭代的步長為

式中ci=Ψigi。從而得到恢復信號：

最后更新殘差：

完成循環，判斷是否滿足迭代停止條件。

改進后的GOMP 算法步驟可總結如下：

1)初始化殘差r0=y，初始化索引集 Λ為空集，初始化迭代次數i=1；

2)計算殘差和測量矩陣 Φ的每一列的內積，找出若干個絕對值最大的元素，并將這些值對應的列序號構成集合P；

3)更新索引集Λi=Λi?1∪P，更新列向量集合Ψi=Ψi?1∪φP；

5)利用式(4)更新殘差；

6)若滿足迭代停止條件，根據索引集 Λ在角度θi(i=1,2,···,N)中找到信號源所在的真實角度位置；若不滿足，則跳轉到步驟2)進行下一次循環。

3 仿真結果分析

本文通過MATLAB 仿真驗證所提算法的可行性，并與l1?SVD算法、OMP 算法、GOMP 算法進行實驗對比驗證其性能。在實驗仿真中，采用均勻直線陣列，噪聲為高斯白噪聲，均勻線陣的陣元間距為半波長，按照等角度劃分方式將空間信號等角度劃分為181 份，每份角度為1°。均方根誤差表示為

式中：M表示蒙特卡洛實驗次數；K表示信號源的總個數；θkm表示第k個空間信號源的DOA 真實值；表示第k個空間信號源的DOA 估計值。

實驗1殘差二范數對比實驗。設2 個信號源以入射角10°、30°入射到均勻直線陣列上，對結合最速下降法的改進GOMP 算法和原始GOMP算法進行性能對比，在單快拍條件下，比較2 種情況下的殘差二范數，進行500 次蒙特卡洛實驗，仿真結果如圖1 所示。

圖1 改進GOMP 算法與GOMP 算法殘差比較

通過圖1 可以看出，改進GOMP 算法比GOMP算法的殘差值小，隨著信噪比的提升，改進算法的殘差值下降更快，可見對GOMP 算法的優化是有效的。

實驗2不同信噪比下DOA 估計的性能對比實驗。設2 個信號源以入射角10°、30°入射到均勻直線陣列上，均勻線陣陣元數為12，快拍數為8，信噪比以2 dB 為步進從?8～8 dB 變化，在相同的實驗條件下進行500 次蒙特卡洛實驗。當估計偏差小于或等于1°時，認為估計成功，在不同信噪比條件下，統計各個條件下的成功概率和均方根誤差，仿真結果如圖2、3 所示。

圖2 不同信噪比下DOA 估計成功概率

圖3 不同信噪比下DOA 估計均方根誤差

在不同信噪比條件下，由圖2 可以看出，與其他3 種算法相比，改進GOMP 算法具有更高的估計成功概率，在信噪比為?2 dB 時估計成功概率接近100%。由圖3 可以看出，在高信噪比情況下，4 種算法均趨于穩定，改進GOMP 算法具有更小的均方根誤差，特別是在低信噪比情況下，改進GOMP 算法具有明顯的優勢。

實驗3不同快拍數下的DOA 估計的性能對比實驗。設2 個信號源以入射角10°、30°入射到均勻直線陣列上，均勻線陣陣元數為12，信噪比為?4 dB，在相同的實驗條件下進行500 次蒙特卡洛實驗。當估計偏差小于或等于1°時，認為估計成功，在快拍數不同的情況下，統計各個條件下的成功概率和均方根誤差，仿真結果如圖4、5 所示。

圖4 不同快拍數下DOA 估計成功率

圖5 不同快拍數下DOA 估計均方根誤差

在不同快拍數情況下，由圖4 可以看出，改進GOMP 算法具有更高的估計成功概率，且具有很高的穩定性。在快拍數為12 時，已經達到以100%的概率估計出信號的DOA。隨著快拍數的增加，其他3 種算法的成功率緩慢增加至趨于平穩。由圖5 可以看出，4 種算法的均方根誤差均隨著快拍數的增加逐漸減小，改進GOMP 算法一直保持著較小的均方根誤差，特別是在低快拍數情況下，改進GOMP 算法具有明顯的優勢。

實驗4不同快拍數情況下算法運行時間對比實驗。由于l1-SVD 算法運行時間相對于貪婪類算法運行時間長，所以只比較另外3 種算法的運行時間。由圖6 可以看出，隨著快拍數的增加，3 種算法的運行時間均逐步上升，且改進GOMP算法相較于其他2 種算法的運行時間相對較短；在相同的條件下，同其他2 種算法相比，改進GOMP 算法具有更短的運算時間。因此從整體上看改進GOMP 算法的運行時間相對偏低，相比于其他2 種算法具有運行時間上的優勢。

圖6 不同快拍數下算法運行時間比較

4 結論

本文提出了一種新的貪婪DOA 估計方法，該方法可有效提高DOA 估計的精度。GOMP 是一種改進的貪婪恢復算法，其利用更好的原子選擇方式，提高了貪婪算法的穩定性和角度估計精度。

1)利用GOMP 算法求解DOA 問題，相較于原有貪婪DOA 估計方法提高了估計精度和成功率。

2)引入最速下降法對GOMP 算法進行優化，并通過仿真驗證了改進方法具有更高的估計精度和成功率，降低了算法的運行時間，相對于原有貪婪DOA 估計方法的性能有一定的提升。