基于子空間跟蹤的遞推穩健自適應波束形成

摘 要:為了解決實時系統應用中，通用秩信號模型的穩健自適應波束形成算法受計算復雜度高的限制這一問題，基于子空間跟蹤理論，運用梯度下降法提出了一種遞推穩健自適應波束形成算法，有效降低了計算量，提高了系統性能，改善了陣列輸出的信干噪比。仿真結果表明，該算法具有快速收斂性與優越的穩健性。

關鍵詞:穩健自適應波束形成; 通用秩信號模型; 子空間跟蹤; 梯度下降; 計算復雜度

中圖法分類號:TN911.7文獻標志碼:A

文章編號:1001-3695(2010)06-2260-04

doi:10.3969/j.issn.1001-3695.2010.06.075

Recursive robust adaptive beamforming

WANG Zhiqiao， HUANG Hua

(School of Electrical Engineering Information， Sichuan University， Chengdu 610065， China)

Abstract:The robust adaptive beamforming， based on generalrank signal models， is restricted by the great computational complexity in the realtime systems. Based on subspace tracking theory， this paper proposed a recursive robust adaptive beamforming by gradient descent method.The proposed algorithm decreased computational complexity， developed online implementations， made the mean output array SINR. Simulations results demonstrate that the proposed method can converg rapidly and has excellent robust performance.

Key words:robust adaptive beamforming; generalrank signal models; subspace tracking; gradient descent; computational complexity

0 引言

波束形成是陣列信號處理的關鍵技術之一[1]。傳統自適應波束形成算法假定訓練序列中不含期望信號，當系統有誤差時，期望信號會被當做干擾抑制，系統性能嚴重下降[1，2]。為使系統存在誤差時仍能保持良好的性能，出現了大量穩健自適應波束形成算法，主要有線性約束波束形成算法(LCMV)、對角加載波束形成算法(LSMI)、特征空間波束形成算法。目前，基于最差情況下優化和顯式系統誤差建模的穩健自適應波束形成算法得到了廣泛深入的研究[2~12]。文獻[5]中，考慮信道散射效應和信號協方差陣與數據協方差陣存在誤差時，推導出基于通用秩信號模型的穩健自適應波束形成算法。

上述穩健自適應波束形成算法都是批處理算法，即每次數據快拍到來時，需計算波束形成器權向量，計算時要對數據協方差陣求逆或特征值分解，其計算復雜度為O(M3)。文獻[13]在LCMV波束形成器基礎上降秩，得到低計算成本的迭代波束形成器;文獻[14]用約束卡爾曼濾波器降低了文獻[4]中算法的復雜度。這兩種算法只適用于點信源(信號協方差陣秩為1)。為了得到一般信源情況下低計算復雜度的穩健自適應波束形成算法，本文運用信號子空間跟蹤理論，推導出一種梯度下降的遞推穩健自適應波束形成算法。

1 理論基礎

1.1 信號模型

窄帶波束形成器輸出

y(n)=wHx(n) (1)

其中:n是時間下標;x(n)是M×1維陣列觀察數據的復向量;w是M×1維波束形成器的復值權向量，最優權向量由最大化信干噪比(SINR)獲得;M是陣元個數;(#8226;)T與(#8226;)H分別代表轉置與共軛轉置。陣列觀察數據向量為

x(n)=s(n)+i(n)+n(n) (2)

其中:s(n)、i(n)與n(n)分別是期望信號、干擾信號與噪聲。假設期望信號、干擾信號與噪聲互不相關，且期望信號與干擾信號為平穩信號，噪聲為高斯白噪聲。

最大化信干噪比(SINR)為

SINR=wHRswwHRi+nw (3)

其中:

Rs=E{s(n)sH(n)} (4)

Ri+n=E{(i(n)+n(n))(i(n)+n(n))H} (5)

分別是M×M維的信號協方差陣與干擾加噪聲協方差陣，E{#8226;}是統計期望。通常，矩陣Rs有任意的秩，即

1≤rank{Rs}≤M

在點信源情況下(Rs秩為1)，有

s(n)=s(n)a(θs)

Rs=σ2sa(θs)a(θs)H

其中:a(θs)是期望信號導向向量;θs是期望信號相對于陣列法線方向的入射角度;s(n)是零均值信號，σ2s=E{|s(n)|2}是s(n)的方差。式(3)可以簡化為

SINR=σ2s|wHa(θs)|2wHRi+nw(6)

在滿秩信源不相干散射情況下，Rs有如下形式

Rs=σ2s∫π/2-π/2ρ(θ)a(θ)aH(θ)dθ (7)

其中:ρ(θ)是歸一化角度功率譜密度[5]。

1.2 傳統自適應波束形成器

通用秩信號模型下，期望信號方向無失真的最小方差波束形成器

minw wH Ri+nws.t.wHRsw=1 (8)

用拉格朗日乘子法解出

wopt=P{R-1i+nRs}(9)

其中:P{#8226;}是生成矩陣主特征向量(最大特征值對應的特征向量)的算子。在點信源的情況下，即為MVDR波束形成器

wopt=R-1i+na(θs)a(θs)HR-1i+na(θs)(10)

當沒有系統誤差時，用數據協方差陣R=E{x(n)xH(n)}=Ri+n+Rs替代式(10)的Ri+n，結果不變[1]。

實際應用中，不能獲得真實數據協方差陣R，只能得到其估計值，常用的無偏估計為

R^=1K∑kn=1x(n)xH(n)(11)

其中:K為快拍數，用估計值R^代替R或Ri+n，即為采樣矩陣求逆波束形成器(SMI)

wSMI=P{R^-1Rs}(12)

當存在期望信號導向向量與角度功率譜密度誤差時，上述自適應波束形成器性能嚴重下降。廣泛使用的穩健技術是對角加載波束形成算法(LSMI)

wLSMI=P{(R^+γI)-1Rs}(13)

其中，最優加載量通常通過數值模擬的方式選擇，沒有更優的選擇方法。

1.3 通用秩信號模型的穩健自適應波束形成算法[5]

假定信號協方差陣與數據協方差陣都存在有限誤差，誤差矩陣都為Hermitian矩陣，即

Rs～=Rs+Δ1 其中，‖Δ1‖≤ξ (14)

R～=R+Δ2 其中，‖Δ2‖≤γ (15)

其中:Rs～與R～分別是實際信號協方差陣與數據協方差陣，Rs與R是對應的假定協方差陣。實際情況下，輸出信干噪比可表示為

SINR=wHR～swwHR～w(16)

最大化SINR優化模型中考慮誤差

minw max‖Δ2‖≤γ wH(R+Δ2)w

s.t. min‖Δ1‖≤ξ wH(Rs+Δ1)w=1 (17)

得到最差情況下存在解析解的最優穩健自適應波束形成器

wSGLW03=P{(R+γI)-1(Rs-ξI)}(18)

其中:wSGLW03記做文獻[5]中算法得到的最優權向量，γ與ξ由環境參數估計。

1.4 信號子空間

假定x是協方差矩陣為C=E{xxH}的復值向量隨機過程，標量函數

J(W)=E{‖x-WWHx‖2}=

tr(C)-2tr(wHCW)+tr(WHCWWHW)(19)

其中:變量W為M×r復值矩陣，不失一般性，假設其秩為r。由文獻[15]中的定理可知，當且僅當W=UrQ時，W是標量函數J(W)的平穩點。其中:Ur是C的任意r個獨立的特征向量，Q是任意的r×r維酉矩陣，而且在每個平穩點處，標量函數J(W)的值都等于其余特征向量對應的特征值之和。因此，當且僅當Ur包含C的r個主特征向量(矩陣特征按從大到小排序后，前r個特征值對應的特征向量)時，標量函數達到最小值，即當標量函數達到最小值時，矩陣變量W收斂為C的r個主特征向量構成子空間的某個正交基上。若已知數據向量是由r個信號加噪聲構成的，最小化標量函數將會得到信號子空間的一個正交基。需要注意的是，由于目標函數的酉不變性，這個正交基不是惟一的，但信號子空間的投影矩陣WWH是惟一的。當然，在r=1時，最小化目標函數將會得到C的歸一化主特征向量。因此，可以通過迭代的方式最小化目標函數，得到矩陣C的歸一化的主特征向量。

2 遞推穩健自適應波束形成算法

2.1 算法推導

存在解析解的穩健自適應波束形成算法 ，計算量集中在計算復雜度都為O(M3)的矩陣求逆與特征值分解，實際應用中并不高效。本文通過矩陣求逆引理與子空間跟蹤理論推導出一種高效算法——遞推穩健自適應波束形成算法。

推導時借助文獻[5]的引理:對于任意的M×M Hermitian矩陣X與任意M×M滿秩Hermitian矩陣Y有如下關系成立

P{XY}=Y-1/2P{Y1/2XY1/2} (20)

顯然，+γI是Hermitian矩陣，雖然在點信源等情況下Hermitian矩陣Rs非滿秩，但對任意ξ>0，矩陣Rs-ξI卻總是滿秩的。需要注意的是，要使矩陣Rs-ξI非負定，ξ須小于Rs的最大特征值，因為當矩陣負定時，式(17)的約束條件不成立。在這種情況下，穩健的自適應波束形成器式(18)(用代替R)為

wSGLW03=(Rs-ξI)-1/2P{(Rs-ξI)1/2(+γI)-1(Rs-ξI)1/2}=

(Rs-ξI)-1/2P{G-1}

(21)

其中，矩陣G定義為

G=(Rs-ξI)-1/2(+γI)(Rs-ξI)-1/2(22)

為了推導遞推算法，用長度為N的固定窗，將第n步對角加載的數據協方差陣(n)展開為

DL(n)=DL(n-1)+1Nx(n)xH(n)-1Nx(n-N)xH(n-N)(23)

其中:DL(n)=(n)+γI。即為秩2更新，γI可以選做DL(n)的初始值，可將G(n)同樣展開表示為

G(n)=G(n-1)+(n)H(n)-(n-N)H(n-N)(24)

其中:變換了的數據快拍表示為

(i)=1N(Rs-ξI)-1/2x(i)(25)

(i)=1N[(Rs-ξI)-1/2]Hx(i)(26)

之所以用(i)與(i)兩個不同的符號表示，是因為一般情況下，矩陣Rs-ξI是不定的，即對于任意的ξ>0且小于Rs的最大特征值，矩陣Rs-ξI同時存在正負特征值，使得其平方根矩陣(Rs-ξI)1/2一定不是Hermitian陣，因此用兩個符號以示區別。

用矩陣求逆引理[16]得到

H-1(n-1)=G-1(n-1)+G-1(n-1)(n-N)H(n-N)G-1(n-1)1-H(n-N)G-1(n-1)(n-N)(27)

G-1(n)=H-1(n-1)-H-1(n-1)(n)H(n)H-1(n-1)1+H(n)H-1(n-1)(n)(28)

其中，中間變量H(n-1)=G(n-1)-(n-N)H(n-N)，G-1(n)的初始值為1γ(Rs-ξI)。

在1.4節中，在當變量W為矩陣C的主特征向量時，無約束目標函數J(W)取得最小值，最直接的迭代求解法是梯度下降類算法。求標量函數J(W)關于W的梯度

J(W)=(-2C+CWWH+WWHC)W(29)

因此，依據梯度下降法的主特征向量更新有

W(n)=W(n-1)-μ[-2C+CW(n-1)WH(n-1)+W(n-1)WH(n-1)C]W(n-1) (30)

用wpro(n)、wpro(n-1)與G-1(n)分別代替上式中的W(n)、W(n-1)與C，得到更新波束形成器權向量的遞推式為

wpro(n)=wpro(n-1)-μ[-2G-1(n)+G-1(n)wpro(n-1)wHpro(n-1)+

wpro(n-1)wHpro(n-1)G

-1(n)]wpro(n-1)

(31)

其中:wpro(n)為本文算法的權向量，μ為步長參數。

2.2 算法描述及計算復雜度分析

將遞推穩健自適應波束形成算法總結如下:

initialize

G-1(0)=1γ(Rs-ξI)

wpro(0)=(Rs-ξI)-1/2P{G-1(0)}

while(new snapshots x(n))

Compute transformed snapshots (2M2-2M)

計算式(25)(26)

if n

G-1(n)=G-1(n-1)-G-1(n-1)(n)H(n)G-1(n-1)1+H(n)G-1(n-1)(n)

else(10M2)

計算式(27)(28)

Compute new weighted vector (7M2-4M)

計算式(31)

end

其中:2M2-2M、5M2、10M2、7M2-4M表示對應計算需要的浮點數操作的次數[16]，每次浮點數操作為一次浮點數加法或浮點數乘法。一般情況下，總的浮點數操作為19M2-6M，計算復雜度為O(M2)，大大降低了文獻[5]算法的計算量。

3 仿真實驗結果

算法推導過程中，并未限定陣列結構，因此遞推穩健自適應波束形成算法適用于任何陣列，但為分析簡便起見，實驗考慮陣元數M=20，陣元間隔為半波長的標準均勻線性陣列[1]， 每一點都是100次的均值。假定只有一個期望信號與一個干擾信號，訓練序列中包含期望信號，干擾噪聲比INR=20 dB。固定窗長度N=20，數據快拍總數K=500。

實驗比較以下算法的輸出信干噪比:采樣矩陣求逆(SMI)，即式(12);對角加載的采樣矩陣求逆(LSMI)，即式(13);文獻[5]的算法，即式(18)，記做SGLW03;本文的遞推算法。

同時在每個圖中顯示最優信干噪比，即式(9)(Ri+n與Rs精確已知的情況)。

LSMI、SGLW03與本文算法的對角加載量為γ=30，按照加載量噪聲比為10 dB的原則選取[1]。

實驗1 考慮點信源且系統沒有誤差。

設參數ξ=16，期望信號與干擾的入射角度為相對于陣列法線方向30°與-30°，圖1顯示了在固定信噪比SNR=0 dB時，本文算法與SGLW03算法的輸出信干噪比?？梢钥闯?，本文算法與SGLW03算法有相似的性能。

實驗2 考慮點信源且期望信號存在DoA估計誤差。

設參數ξ=16，實際的期望信號與干擾的入射角度為30°與-30°，假定的期望信號入射角度為33°，即估計誤差為3°。本文算法使用窗長N=20。圖2給出了SNR=0時平均輸出信干噪比SINR隨快拍數的變化。圖3表示快拍數為200時平均信干噪比SINR隨信噪比SNR的變化。從兩圖中看出，即使DoA失配時，本文的算法和SGLW03算法性能相近，不同快拍和信噪比下都有較高的輸出信干噪比，而SMI隨著快拍增加或信噪比提高性能嚴重下降，高信噪LSMI算法的性能也有所下降。

實驗3 考慮滿秩信源下信號協方差陣下有偏差。

設參數ξ=9，假設實際的期望信號與干擾都存在局部不相干散射[5]，其實際的角度功率譜密度分別服從中心角度為30°與-30°的高斯與均勻分布，角度擴散為4°，相當于高斯分布的標準差為2°，均勻分布的范圍為[-32°，-28°]。而假定的信號協方差陣存在偏差，假設其估計的角度功率譜密度也是服從高斯分布，中心角度為32°，標準差為3°，即估計的中心角度與標準差都有偏差。圖4顯示當SNR=0時，各算法輸出信干噪比隨快拍數的變化。圖5顯示了快拍數為200時輸出信干噪比隨信噪比的變化，其中算法假設窗口長度為N=500，即相當于數據協方差陣采用秩1更新。同實驗2有相似的結果，即使是在信號協方差陣滿秩的情況下，本文的遞推算法和SGLW03算法有同樣優越的穩健性能。

4 結束語

利用信號子空間跟蹤理論，針對通用秩信號模型的穩健自適應波束形成算法計算復雜度高的問題，本文提出了一種遞推穩健自適應波束形成算法，并分析了該算法每次更新權向量只需19M2-6M次浮點數操作，即計算復雜度為O(M2)。從數值仿真中可以看到，本文算法在沒有系統誤差時，能快速收斂于文獻[5]的算法，而且無論是點信源存在DoA估計誤差，還是滿秩信源存在角度功率譜密度參數估計誤差，本文算法都有優越的穩健性能。

參考文獻:

[1]VAN TREES H L. Optimum array processing， detection， estimation， and modulation theory， vol 4[M]. New York: Wiley， 2002.

[2]LI J， STOICA P. Robust adaptive beamforming[M]. Hoboken， NJ: Wiley， 2006.

[3]GERSHMAN A B. Robust adaptive beamforming in sensor arrays[J]. AEUInnt Journal of Electronics and Communications， 1999， 53(6): 305-314.

[4]VOROBYOV S A， GERSHMAN A B， LUO Z Q. Robust adaptive beamforming using worstcase performance optimization: a solutionto the signal mismatch problem[J]. IEEE Trans on Signal Processing， 2003， 51 (2): 313-324.

[5]SHAHBAZPANAHI S， GERSHMAN A B， LUO Z Q. Robust adaptive beamforming for generalrank signal models[J]. IEEE Trans on Signal Processing， 2003， 51(9): 2257-2269.

[6]KIM S J， MAGNANI A， MUTAPCIC A， et al. Robust beamforming via worstcase SINR maximization[J]. IEEE Trans. on Signal Processing， 2008， 56(4): 1539-1547.

[7]YU Z L， SER W， ER M H， et al. Robust response control for adaptive beamformers against arbitrary array imperfections[J]. Signal Processing， 2008， 88(9): 2234-2245.

[8]YU Z L， ER M H， SER W. A novel adaptive beamformer based on semidefinite programming (SDP) with magnitude response constraints[J]. IEEE Trans on Antennas and Propagation， 2008， 56(5): 1297-1307.

[9]ZHANG Y， NG B P， WAN Q. Sidelobe suppression for adaptive beamforming with sparse constraint on beam pattern[J]. IEE Electronics Letters， 2008，44(10): 615-616.

[10] LI J， XIE Y， STOICA P， et al. Beampattern synthesis via a matrix approach for signal power estimation[J]. IEEE Trans on Signal Processing， 2007， 55(12): 5643-5657.

[11]CHEN C， VAIDYANATHAN P P. Quadratically constrained beamforming robust against directionofarrival mismatch[J]. IEEE Trans on Signal Processing， 2007， 55(8): 4139-4150.

[12]SELN Y， ABRAHAMSSON R， STOICA P. Automatic robust adaptive beamforming via ridge regression[J]. Signal Processing， 2008， 88(1): 33-49.

[13]De LAMARE R C. Adaptive reducedrank LCMV beamforming algorithms based on joint iterative optimisation of filters[J]. IEE Electronics Letters， 2008， 44(9): 565-566.

[14]ElKEY A， KIRUBARAJAN T， GERSHMAN A B. Robust adaptive beamforming based on Kalman filter[J]. IEEE Trans on Signal Processing， 2005， 53(8):3032-3041.

[15]YANG B. Projection approximation subspace tracking[J]. IEEE Trans on Signal Processing， 1995， 43(1): 95-107.

[16]GOLUB G H， VAN LOAN C F. Matrix computations[M]. 3rd ed. Baltimore: Johns Hopkins University Press， 1996.