CN111880167A - 一种基于先随机后优化的波达方向估计方法 - Google Patents

Abstract

本发明涉及一种基于先随机后优化的波达方向估计方法，包括：RTO‑MH算法先验转换步骤：构造可逆映射函数T_1D将高斯参考随机变量v与拉普拉斯分布的参数θ相连接，从而使得θ＝T_1D(v)，得到T_1D(v)的转换方程，使得模型中的拉普拉斯参数满足后验形式，进而将Dθ的拉普拉斯先验分布转换为变量为v的标准高斯分布；模型转换步骤：将测量矩阵由复数形式变为实数形式，将实数形式的测量矩阵输入到RTO‑MH算法进行先验转换，在计算结果的归一化功率并选择主瓣作为DOA。本发明的优点在于：通过先验变换算法RTO‑MH提高了DOA估计中BCS的性能，与传统方法相比其准确率得到了很大的提高，且本发明不需要再使用超参数，其处理计算时间也大大地降低。

Description

一种基于先随机后优化的波达方向估计方法

技术领域

本发明涉及通信天线技术领域，尤其涉及一种基于先随机后优化的波达方向估计方法。

背景技术

波达方向(DOA)的估计问题是雷达和天线领域的重要研究任务之一，近年来，压缩感知(CS)理论已经成功应用于DOA估计中，其关键思想是重建稀疏信号；目前有几种稀疏信号的构建方法，例如正交匹配追踪(OMP)和基追踪(BP)等。OMP源自匹配追踪(MP)，为了减轻BP方法的计算负担，有人提出了一种降维方法，称为l₁-SVD；基于CS理论稀疏信号重建方法扩展到了贝叶斯压缩感知(BCS)，贝叶斯压缩感知是基于对信号和噪声的稀疏先验假设从贝叶斯角度提出的；在BCS框架下信号重建主要通过稀疏贝叶斯学习(SBL)实现，在最近的研究中，SBL得到了很大的发展，例如根SBL、变分SBL、联合SBL和网格演化方法，SBL具有多层的假设框架，该框架主要在于迭代学习新的信息并更新超参数。

在DOA估计问题的实际应用中，快拍通常很短甚至有时需要实时处理，但是，有限的数据导致DOA估算的准确性较差；为了降低快拍的需求，现有很多来改进传统空间频谱估计的方法，如伪协方差矩阵方法和空间平滑技术被广泛应用。但是，在稀疏信号重建中快拍对估计精度的影响仍然存在，通常，对基于SBL的DOA方法，估计DOA的准确性很大程度上取决于快拍的数量，如果快拍的数量不足将导致DOA估计的准确度下降。因此，如何在降低快拍需求的同时提高稀疏信号重建的性能是现阶段需要解决的问题。

发明内容

本发明的目的在于克服现有技术的缺点，提供了一种基于先随机后优化的波达方向估计方法，解决了现有DOA估计方法中存在的不足。

本发明的目的通过以下技术方案来实现：一种基于先随机后优化的波达方向估计方法，所述估计方法包括：

RTO-MH算法先验转换步骤：构造可逆映射函数T_1D将高斯参考随机变量v与拉普拉斯先验分布的参数θ相连接，从而使得θ＝T_1D(v)，得到T_1D(v)的转换方程，使得模型中的拉普拉斯参数满足后验形式，进而将Dθ的拉普拉斯先验分布转换为变量为v的标准高斯分布；

模型转换步骤：将测量矩阵由复数形式变为实数形式，将实数形式的测量矩阵输入到RTO-MH算法进行先验转换，在计算结果的归一化功率并选择主瓣作为DOA。

进一步地，所述RTO-MH算法先验转换步骤具体包括：

根据模型中的拉普拉斯单参数的先验分布P(θ)∝exp(-λ|θ|)，得到其后验分布

构造可逆映射函数T_1D将高斯参考随机变量v与拉普拉斯分布的物理参数θ相连接，从而使得θ＝T_1D(v)，得到转换方程

将拉普拉斯单参数的先验变换扩展到拉普拉斯多参数上，得到每个拉普拉斯先验参数的分布元素。

进一步地，所述将拉普拉斯单参数的先验变换扩展到拉普拉斯多参数上，得到每个拉普拉斯先验参数的分布元素包括：

假设拉普拉斯参数θ的先验分布为

进而得到其后验分布

通过所述转换方程定义的一维变换T_1D将Dθ＝T(v)∶＝[T_1D(v₁),…,T_1D(v_n)]^T的每个拉普拉斯先验参数分布元素转换为标准的高斯分布。

进一步地，所述实数形式的模型为t_m＝Φw_m+n_m；其中，t_m＝[Re{y_m},Im{y_m}]^T表示真实的测量矩阵，w_m＝＝[Re{s_i},Im{s_i}]^T表示真实的信号矩阵，n_m＝＝[Re{e_m},Im{e_m}]^T表示真实的噪声矩阵，

表示真实的流型矩阵。

进一步地，所述方法还包括模型建立步骤；所述模型建立步骤包括：。

以K个窄带远场信号s_k(t),k＝1,…,K从方向θ_k,k＝1,…,K撞击M个呈线性阵列分布的传感器，通过相移表示不同传感器的时延，得到先验误差模型y(t)＝A(θ)s(t)+e(t)；

通过线性变换白化先验误差模型，得到具有逆问题的模型表示为：y＝f(θ)+ε,ε～N(0,I_m),θ～(θ₀,I_n)；

RTO通过反复优化随机扰动的代价函数进行后验采样，并通过MH方法校正包含空间谱的采样样本，进而得到参数θ的后验密度

进一步地，所述RTO后验采样包括：

找到一个表示为

的线性化点并进行固定，将该点设置为后验模式；

在线性化点计算F(θ)的雅可比行列式，表示为

对其稀释QR因式分解后进行评估；

根据n维标准高斯分布计算独立样本ξ，并通过

得到提案要点θ_prop。

本发明的有益效果为：本发明通过先验变换算法RTO-MH提高了DOA估计中BCS的性能，与传统方法相比其准确率得到了很大的提高，且本发明不需要再使用超参数，其处理计算时间也大大地降低。

附图说明

图1为本发明的方法流程图；

图2为本发明方法与SBL方法的仿真对比图；

图3为本发明方法与SBL、OMP和l₁-SVD的仿真对比图。

具体实施方式

为使本申请实施例的目的、技术方案和优点更加清楚，下面将结合本申请实施例中附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本申请一部分实施例，而不是全部的实施例。通常在此处附图中描述和示出的本申请实施例的组件可以以各种不同的配置来布置和设计。因此，以下对附图中提供的本申请的实施例的详细描述并非旨在限制要求保护的本申请的保护范围，而是仅仅表示本申请的选定实施例。基于本申请的实施例，本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本申请保护的范围。下面结合附图对本发明做进一步的描述。

如图1所示，本发明涉及一种基于先随机后优化的波达方向估计方法，其包括：

假设有K个窄带远场信号s_k(t),k＝1,…,K从方向θ_k,k＝1,…,K撞击M个呈线性阵列分布的传感器，通过相移表示不同传感器的时延，得到观测模型y(t)＝A(θ)s(t)+e(t)；

其中，y(t)＝[y₁(t),…,y_M(t)]^T、θ＝[θ₁(t),…,θ_K(t)]^T、s(t)＝[s₁(t),…,s_K(t)]^T、e(t)＝[e₁(t),…,e_M(t)]^T·y_m(t)、e_m(t)，m＝1,…,M，表示在时间t时刻第m个传感器的测量值和噪声。A(θ)＝[a(θ₁),…,a(θ_K)]表示阵列流型矩阵，a(θ_K)是第k个信号来源的导向向量，a_m(θ_K)包含第k个信号来源的导向向量到第m个传感器的延迟信息；为了定位稀疏信号的方向，空间域由网格

采样，其中N表示网格编号，N＞＞M>K；采样通常是统一的，因为没有来源的先验信息，测量矢量y(t)和阵列流型矩阵A(θ)是已知的，s(t)需要估算。

公式y(t)＝A(θ)s(t)+e(t)可以看作是贝叶斯逆问题，流型矩阵A(θ)是非线性参数到测量的映射，通常，解决贝叶斯逆问题主要集中在表征后验密度上，对于传统的SBL方法获得后验密度的过程需要更新超参数γ。随着快拍数量的增加，SBL收敛到代价函数的全局最小值的可能性也随之增加；由于全局最小化了SBL超参数，导致最大程度地稀疏，因此可以保证表示形式，并且增加快拍数量可以提高找到这些超参数的可能性，快拍数量太少会导致SBL性能下降。

从贝叶斯压缩感知的另一个角度来看，RTO方法摆脱了快拍数量的依赖，通过有限的快拍解决了贝叶斯逆问题；RTO通过使用提案密度中随机扰动的优化来产生样本，该提案密度将作为Metropolis-Hastings(MH)方法用作Metropolis独立提案。

考虑到具有高斯测量误差的贝叶斯逆问题和高斯先验，可以使用线性变换来白化先验误差模型，从而使逆问题具有以下模型：

y＝f(θ)+ε,ε～N(0,I_m),θ～(θ₀,I_n)

其中，θ₀是先验均值，I_n和I_m分别是大小为n和m的单位矩阵，y是测量向量，f是前向函数，ε是测量误差。

RTO通过反复优化随机扰动的代价函数获得后验样本，为了获得未知参数θ，RTO要求目标分布具有特定形式，这种分布允许在MH过程中使用RTO样本；MH方法会校正包含空间谱的这些后验样本，将目标密度(通常为θ的后验密度)写为：

其中，‖·‖表示标准，上式的右侧定义为参数θ的向量值函数，其形式为：

其中，RTO后验采样包括：首先，找到一个表示为

的线性化点并进行固定，将该点设置为后验模式；其次，在线性化点计算F(θ)的雅可比行列式，表示为

对其稀释QR因式分解后进行评估；最后，根据n维标准高斯分布计算独立样本ξ，并通过

得到提案要点θ_prop。

提案点是根据提案密度进行分配的，表示为：

其中，|·|表示矩阵行列式的绝对值，该提案密度在MH方法中用作独立提案。对于从θ^(i-1)点到

点的更新，MH的接受率为：

其中，

在RTO-MH先验转换中，根据先验误差模型中的拉普拉斯单参数的先验分布P(θ)∝exp(-λ|θ|)，得到其后验分布

其中，σ_obs表示误差标准偏差。

为了满足公式

中的后验形式，构造可逆映射函数T_1D，该函数将高斯参数随机变量v与拉普拉斯分布的物理参数θ相连接，从而θ＝T_1D(v)，转换方程为：

其中

是标准高斯分布的累积分布函数(CDF)，L是拉普拉斯分布的积累分布函数(CDF)。

基于单参数的先验变换，将其扩展到多参数的情况；假设θ的先验分布为：

其中，(Dθ)_i是向量Dθ的第i个元素，D是可逆矩阵，进而得到其后验分布：

其中

表示协方差矩阵。

参考随机Dθ的每个拉普拉斯先验分布元素可以通过转换方程中定义的一维变换T_1D转换为标准高斯分布；所以，Dθ＝T(v)，其中，T(v)∶＝[T_1D(v₁),…,T_1D(v_n)]^T；

其中，v₁是对应的每个参考变量，从而得到参数θ＝D^-1T(v)，变换的雅可比行列式为D^-1J_T，其中，J_T是T的雅可比行列式，记为：

其中,T′_1D表示

的导数，得到

是标准高斯分布的概率密度函数(PDF)；

根据公式θ＝D^-1T(v)的转换，v的后验密度可以表示为：

测量向量是很复杂的，但是，RTO-MH只能用于实值或者实数值，因此，需要将模型转换为实数形式；其表示为：

t_m＝Φw_m+n_m

其中，t_m＝[Re{y_m},Im{y_m}]^T表示真实的测量矩阵，w_m＝＝[Re{s_i},Im{s_i}]^T表示真实的信号矩阵，n_m＝＝[Re{e_m},Im{e_m}]^T表示真实的噪声矩阵，

表示真实的流型矩阵。

通过将测量矩阵由复数形式变为实数形式，将实数形式的测量矩阵输入到RTO-MH算法进行先验转换，在计算结果的归一化功率并选择主瓣作为DOA。

通过仿真实验来验证本发明方法的有效性；将两个互不相关分布来自30°和60°的窄带信号源，撞击在具有16个传感器的均匀线性阵列上；信号的空间范围是通过范围为-90°～90°，网格间隔为1°的均匀网格样本离散化的；首先，将本发明方法与传统SBL方法在快拍数量从1到200增加时进行比较。

如图2所示，垂直虚线表示DOA的真实位置。对于单个快拍，SBL方法只能在30°处估计DOA。而且它的旁瓣很高，这意味着SBL在单个快拍上的性能很差。当快拍达到200张时，SBL的性能将提高。本发明的方法和SBL方法在30°处都有很高的精准度。尽管本发明方法的主瓣在60°处有1°的偏差，主瓣比SBL方法宽，但旁瓣却低得多。

为了说明精度如何表现，将本发明方法的RMSE与OMP、l₁-SVD和SBL进行比较；比较中的信噪比(SNR)范围为[-10，-8，-6，-4，-2、0、2、4、6、8、10]。每个SNR均进行200次独立的蒙特卡洛试验。SBL方法分别模拟了1个快拍和200个快拍。l₁-SVD算法模拟了1个快拍和50个快拍。

如图3所示，本发明方法在1个快拍的RMSE级别明显比SBL和l₁-SVD低，在200个快拍所有方法的精度几乎相同，这说明本发明的方法现对于传统的方法，在低快拍数量的处理方面具有很大的优势，在快拍数量为50时，本发明方法的准确性超过了传统的方法。

本发明的方法与传统的BCS方法(例如SBL)的准确性高度取决于快拍数量相比，该方法不需要很多样本即可更新超参数。仿真结果表明，所提方法具有比现有方法更好的估计精度。SBL，OMP和l₁-SVD当快拍数量为单个时，有效地减少了处理时间。它的计算负担接近OMP算法，低于l₁-SVD和SBL方法。

以上所述仅是本发明的优选实施方式，应当理解本发明并非局限于本文所披露的形式，不应看作是对其他实施例的排除，而可用于各种其他组合、修改和环境，并能够在本文所述构想范围内，通过上述教导或相关领域的技术或知识进行改动。而本领域人员所进行的改动和变化不脱离本发明的精神和范围，则都应在本发明所附权利要求的保护范围内。

Claims (6)

1.一种基于先随机后优化的波达方向估计方法，其特征在于：所述估计方法包括：

RTO-MH算法先验转换步骤：构造可逆映射函数T_1D将高斯参考随机变量v与拉普拉斯先验的参数θ相连接，从而使得θ＝T_1D(v)，得到T_1D(v)的转换方程，使得模型中的拉普拉斯参数满足后验形式，进而将Dθ的拉普拉斯先验分布转换为变量为v的标准高斯分布；

模型转换步骤：将测量矩阵由复数形式变为实数形式，将实数形式的测量矩阵输入到RTO-MH算法进行先验转换，在计算结果的归一化功率并选择主瓣作为DOA。

2.根据权利要求1所述的一种基于先随机后优化的波达方向估计方法，其特征在于：所述RTO-MH算法先验转换步骤具体包括：

根据模型中的拉普拉斯单参数的先验分布P(θ)∝exp(-λ|θ|)，得到其后验分布

构造可逆映射函数T_1D将高斯参考随机变量v与拉普拉斯先验的物理参数θ相连接，从而使得θ＝T_1D(v)，得到转换方程

将拉普拉斯单参数的先验变换扩展到拉普拉斯多参数上，得到每个拉普拉斯先验参数的分布元素。

3.根据权利要求2所述的一种基于先随机后优化的波达方向估计方法，其特征在于：所述将拉普拉斯单参数的先验变换扩展到拉普拉斯多参数上，得到每个拉普拉斯先验参数的分布元素包括：

假设拉普拉斯参数θ的先验分布为

进而得到其后验分布

通过所述转换方程定义的一维变换T_1D将Dθ＝T(v)：＝[T_1D(v₁)，…，T_1D(v_n)]^T的每个拉普拉斯先验参数分布转换为标准高斯分布。

4.根据权利要求1所述的一种基于先随机后优化的波达方向估计方法，其特征在于：所述实数形式的模型为t_m＝Φw_m+n_m；其中，t_m＝[Re{y_m}，Im{y_m}]^T表示真实的测量矩阵，w_m＝＝[Re{s_i}，Im{s_i}]^T表示真实的信号矩阵，n_m＝＝[Re{e_m}，Im{e_m}]^T表示真实的噪声矩阵，

表示真实的流型矩阵。

5.根据权利要求1-4中任意一项所述的一种基于先随机后优化的波达方向估计方法，其特征在于：所述方法还包括模型建立步骤；所述模型建立步骤包括：。

以K个窄带远场信号s_k(t)，k＝1，…，K从方向θ_k，k＝1，…，K撞击M个呈线性阵列分布的传感器，通过相移表示不同传感器的时延，得到先验误差模型y(t)＝A(θ)s(t)+e(t)；

通过线性变换白化先验误差模型，得到具有逆问题的模型表示为：y＝f(θ)+ε，ε～N(0，I_m)，θ～(θ₀，I_n)；

RTO通过反复优化随机扰动的代价函数进行后验采样，并通过MH方法校正包含空间谱的采样样本，进而得到参数θ的后验密度

6.根据权利要求5所述的一种基于先随机后优化的波达方向估计方法，其特征在于：所述RTO后验采样包括：

找到一个表示为

的线性化点并进行固定，将该点设置为后验模式；

在线性化点计算F(θ)的雅可比行列式，表示为

对其稀释QR因式分解后进行评估；

根据n维标准高斯分布计算独立样本ξ，并通过

得到提案要点θ_prop。

Images