一种基于加权改进平滑l0 范数的DOA 估计方法

王勇，李韬，项建弘

哈尔滨工程大学 信息与通信工程学院，黑龙江 哈尔滨 150001

波达方向（direction of arrival，DOA）估计作为阵列信号处理领域中的重要研究方向，在声呐、雷达和地震勘测等领域具有很大的应用空间[1-3]。其中比较经典的有基于子空间理论的子空间分解类算法，即多重信号分类（multiple signal classification method，MUSIC）算法[4]和信号参数旋转不变（eestimating signal parameter via rotational invariance techniques，ESPRIT）算法[5]。以上2 类算法均需得到信号协方差矩阵，再对其进行分解，在良好条件下，可实现高分辨率DOA 估计，但是需要大量的快拍数据作为估计精度保障，且无法求解相干源信号。尽管以上算法可通过平滑处理的方法来解决相干信号的问题，但会增加计算量，且会损失一定的阵列孔径，使得分辨率下降。另外，最大似然估计（maximum likelihood，ML）方法也具有较好的估计性能，但是需要非线性多维搜索，计算量较大[6]。近年来，压缩感知（compressed sensing，CS）理论的提出和广泛应用[7-9]，为DOA估计提供了新的研究途径，一些学者利用信源在空域的稀疏性，直接将压缩感知的经典算法运用到波达方向估计中。由于计算量小，以正交匹配追踪（orthogonal matching pursuit，OMP）为代表的贪婪类算法[10-11]常被用来解决该问题，该类算法的关键在于经过反复迭代，找到测量矩阵中提供较大能量的原子，从而估计出信号方向。但在低信噪比情况下通常找不到最优解，容易使得DOA估计精度下降。另一类方法是将l0范数最小化问题转化为求解l1范数最小化问题，如L1-SVD 算法[12]估计精度高，需要采用线性规划进行求解，计算量较大，实时处理困难。文献[13]提出平滑l0范数重构算法，其主要思路是通过构造一个带参的高斯平滑函数族来近似信号的l0范数，通过最速下降法求解平滑函数的最优解。由于其计算简单，目前已取得了许多的研究成果，如文献[14]将高斯函数替代为逼近度更高的双曲正切函数去拟合l0范数，并采用修正牛顿法提高了收敛速度。文献[15]和文献[16]中分别提出用一组复合三角函数和修正近似双曲正切函数来逼近l0范数，虽然重构效果有所提高，但函数模型比较复杂，使得迭代计算复杂度较高。

基于平滑l0范数类算法的核心是将求解l0范数问题转化为求解平滑函数的极值问题，但是该类算法通常对噪声敏感，且多用于单快拍下DOA 估计。为了进一步优化算法的适应能力，引入了多观测矢量(multiple measurement vectors，MMV)模型，采用奇异值分解，先对信号降维并提取信号子空间，一定程度上降低了计算量。为了使得算法重构精度更高，提出拟合程度更高的复合优化函数去逼近l0范数，并采用加权机制，加速稀疏解的获取，最后利用最速下降法快速求解。通过实验论证，所提算法在低信噪比、少快拍下具有更优的DOA 估计误差。

1 压缩感知框架下的DOA 估计模型

假设空间中有L个来波信号方向为θi,i∈{1,2,···,L}的远场窄带信号si(t),i∈{1,2,···,L}，入射到阵元数为M的均匀线阵中，其接收范围为 (-0.5π,0.5π)，阵元间距为d。为了获得稀疏信号，将整个空间划分为 2N+1份，则整个空间可表示为 {θ1,θ2,···,θ2N+1}，其中 θk=-0.5π+(k-1)π/2N，k∈{1,2,···,2N+1}。从信号能量角度分析，由于整个空间中只有L(L≪2N+1)个 远场窄带信号，即至多只有L个不同来波方向上才有真实的信号能量分布，因此信号的真实方向相对于整个信号空间则是稀疏的。基于空间等角度均匀划分后所构造的阵列流型矩阵A是一个M×(2N+1)维矩阵，在压缩感知理论中，相当于一个过完备冗余字典，可作为测量矩阵使用。若阵列接收到L个远场窄带信号各自在T个不同时刻的快拍数据，则信号的MMV 接收模型可表示为

式中：Y为M×T维阵列接收信号矩阵；S为(2N+1)×T维远场窄带信号矩阵，其中只有L(L≪2N+1)行不为0，为稀疏信号矩阵。对于S中的列向量s而言，s中只有L(L≪2N+1)个位置不为0，为稀疏向量；G为M×T维的高斯白噪声矩阵；A=[a(θ1)a(θ2) ···a(θ2N+1)]，为M×(2N+1)维空间等角度网格化构造的阵列流型矩阵，其中a(θk)=为导向矢量，k∈{1,2,···,2N+1}，λ 为信号的波长，d为阵元间距。

压缩感知框架下的DOA 估计通常是在获得稀疏网格化构造的阵列流型矩阵A和阵列接收信号矩阵Y之后，利用重构算法对入射信号矩阵S进行恢复，从而找出S中对应的非零位置集合并将其映射到划分的网格空间，最终求解出的θj,j∈{1,2,···,L}即为实际的信号来波方向。

2 改进平滑l0 范数的DOA 估计算法

针对式（1）中信号接收矩阵Y的列向量y，以及稀疏矩阵S中的稀疏列向量s，在空间存在噪声情况时，可转化为下述优化问题进行求解：

图1 和图2 分别在 ρ=0.2 和 ρ=0.1处直观地观察高斯函数fρ(si)与所提复合优化函数gρ(si)对于l0范数的逼近程度，图中H(si)代表函数值。

图1 ρ=0.2时函数对于 l0范数的逼近程度

图2 ρ=0.1时 函数对于 l0范数的逼近程度

可知在 ρ →0这一不断变小的过程中，高斯函数和所提优化函数对于l0范数的逼近程度均有所提高，并且当si=0 时2 个函数的值均为0；在si不为0 的附近，由于所提函数对应的值能更快地趋向1，相对于高斯函数的图形更为“陡峭”，即所提优化函数逼近l0范数程度要高于高斯函数，这也验证了所提优化函数用于逼近l0范数的优越性。因此式（2）可转化为式（3）的最优化问题：

为了加速获得稀疏解，提高算法的收敛性，对所提函数采用加权机制，通常在迭代过程中，对于较大的元素si给与一个较小的权值，而较小的元素si给定一个较大的权值，本文采用的权值函数表形式为

可知此时当si较大时，权值wi较小，符合较大的元素对应较小的权值这一要求。

令

作为权值函数的对角矩阵，在迭代寻优过程中，将权值函数加权到目标稀疏信号

上，可得：

在后续的优化迭代过程中，较小的加权值能够一定程度减缓目标信号中较大元素的下降，而较大的加权值能加速目标信号中较小元素的下降，因此针对稀疏信号s而言，目标信号的大多数元素能更快的趋近于0。由于在噪声环境下，来波信号主要集中在元素值较大的位置，因此引入加权机制在一定程度上提升了算法的收敛性。

通过加权函数的引入，式（3）可转化为式（4）优化问题：

对于式（4），本文采用最速下降法进行求解，稀疏信号向量s的迭代更新方向记为

式中：µ是步长，为一个较小的常数；d为最速下降法的下降方向。

根据目标函数gρ(si)可知其一阶导数为

根据式（4）优化的方向进行优化迭代，并将每次迭代的结果进行梯度投影，其投影方式为

综上所提算法主要包括2 个循环，首先确定一个递减序列 ρ1,ρ2,···,ρj，外循环控制 ρ的取值由ρ1变化到 ρj，内循环则是针对每一个Ysv，计算函数的优化方向d并对信号进行迭代更新和进行梯度投影，从而最终实现信号的重构。

具体的计算步骤：

本文针对的是多快拍模型，由于上述流程只是针对Ysv的某一列向量ysv求得的重构解，因此需要对Ysv所有的列向量分别重复上述步骤求出各自重构解，然后再求平均值，得到多快拍模型下的稀疏重构解，最终再映射到划分网格上，获得DOA 估计值。

3 实验仿真

在Matlab 平台上进行实验仿真验证本文所提算法的可行性，并与SL0、OMP、L1-SVD 算法进行对比，分析所提算法在低信噪比、少快拍下的DOA 估计性能优势。本实验中均采用间距为半波长的均匀线阵，并将空间从 (-90◦,90◦)等角度均匀划分1 81份，每份间隔1°，噪声环境为高斯白噪声。估计指标采用均方根误差和估计成功率，其中DOA 估计的均方根误差定义如下：

式中：L为空间目标信号源个数，K为蒙特卡洛实验次数，为第k次蒙特卡洛实验中第l个信源的DOA 估计值，θl为第l个信源的真实来波方向。

DOA 估计成功率定义为

式中F(θl)为 入射方向为 θl的信源在满足估计误差RMSE下，能够被估计出来的总次数。本文实验假设估计值与真实值的RMSE小于或等于1°时，则认为估计成功。

实验1验证所提算法在不同快拍条件下的DOA 估计结果有效性。假设阵元数M=16，信噪比 为5 dB，有3 个目标源 信号，分别从-30◦、0◦、20◦入射，其中-30◦和 0◦为相干信号源。分别在单快拍和快拍数为10 的条件下进行仿真，并对空间谱进行了归一化处理，最终DOA 估计结果如图3和图4 所示。

图3 单快拍下DOA 估计结果

图4 快拍数为10 的 DOA 估计结果

由图3 和图4 可知，归一化空间谱在信号的入射方向均出现了峰值，且随着快拍数的增加，归一化空间谱在非信号源方向上的波动更少，在信源方向上的谱峰相对而言更尖锐、宽度更窄，更容易区分信号的来波方向，并且在前2 个为相干信号源的方向上也各自出现了谱峰，因此所提算法在较少快拍，甚至单快拍下均能获取准确的DOA 估计结果，估计算法不受相干源信号影响。

实验2比较在不同信噪比下各算法的DOA估计性能。假设目标信源的入射角度为 0◦、20◦，阵元数量M=16，快拍采样数为10，进行600 次蒙特卡洛独立重复试验。同L1-SVD、OMP、SL0 算法进行比较，与本文所提算法在不同信噪比下的DOA 估计误差和估计成功率的仿真如图5 和图6所示。

图5 不同信噪比下算法的DOA 估计误差

图6 不同信噪比下算法的DOA 估计成功率

在不同信噪比条件下，由图5 可知，高信噪比下，各算法均方根误差均逐渐趋于稳定并收敛到0，当信噪比低于-4 dB 时，所提算法相比其他3 种算法具有明显优势，其均方根误差明显低于其他3 种算法。由图6 可以看出，信噪比为-3 dB 时，各算法估计成功率均高于90%，所提算法趋于100%，在低信噪比下所提算法具有更优的估计性能。

实验3比较在不同快拍下各算法的DOA 估计性能。假设2 个信源来波方向分别是 0◦、20◦，阵元数量16，信噪比RSN=-3 dB，进行600 次蒙特卡洛实验，各算法在不同快拍下的DOA 估计结果如图7 和图8 所示。

图7 不同快拍下算法的DOA 估计误差

图8 不同快拍下算法的DOA 估计成功率

在不同的快拍条件下，由图7 可以看出，即使在-3 dB 低信噪比条件下，随着快拍数的不断增加，各算法的估计误差也能逐渐降低并趋于稳定；在快拍数少于6 时，估计的均方根误差都不同程度的增大，上述所有算法的性能都开始下降，但所提算法相比其他算法的均方根估计误差明显更小。由图8 可以看出，随着快拍数的增加，各算法的估计成功率也在缓慢增加并趋于稳定；当快拍数为10 时，所提算法几乎能以100%的概率估计出信号的来波方向；在快拍数为4 时，其估计成功率也在90%左右。相比其他算法，本文提出的算法在低信噪比、少快拍下具有更好的估计成功率。

4 结论

本文提出了一种加权改进的平滑l0范数最小化算法，并将其应用到多观测矢量DOA 估计中。

1）为了提高算法重构精度，提出了新的复合优化函数来逼近l0范数。经过分析可知，相比原始的SL0 算法，其高斯函数逼近程度更高。

2）为了提高算法的稳定性，本文算法采用了加权机制，使得其在迭代过程中能加速稀疏解的获取，并且在重构前可预先对接收数据进行奇异值分解，提取了信号子空间，使得算法整体的抗噪性能得到提升，也一定程度上降低了计算量。

经过实验对比，本文所提方法相比OMP、原始的SL0 和L1-SVD 算法在低信噪比、少快拍下具有更优的DOA 估计性能。因此在低信噪比、少快拍下的实际应用中，本文所提算法具有更高的可靠性和适用性。但本文只针对一维DOA 估计，后续仍需推广和应用到二维DOA 估计中，以提高算法的适用性。同时在后续研究中也需要考虑阵列流型存在幅相误差这一问题。