CN114415109B - 一种稀疏贝叶斯学习的直接定位方法 - Google Patents

Abstract

本发明提供了一种稀疏贝叶斯学习的直接定位方法，利用广义双帕累托先验分布比传统分布更强的稀疏特性，保证算法不受相关信号影响的同时，还加强了信号的稀疏性，最终表现为分辨率和定位精度的提高以及低信噪比、小快拍数情况下的稳健性。本发明将稀疏贝叶斯学习引入直接定位系统，从稀疏类方法入手，解决了子空间原理类直接算法无法定位相关辐射源的问题；同时性能上的提升解决了现有最大似然原理类直接定位算法所需采样快拍数较大、分辨率不够、未知辐射源数目时无法准确定位等问题。

Description

一种稀疏贝叶斯学习的直接定位方法

技术领域

本发明涉及无源定位领域，具体涉及一种信号的直接定位方法。

背景技术

无源定位系统指利用目标辐射源自身发射的无线电信号进行定位的系统，在电子侦察、雷达探测等军事领域被广泛运用。直接定位作为无源定位领域一种新的定位技术，无需中间参数估计步骤，在低信噪比下相比传统两步定位方法具有估计精度高、参数自动关联等优点。该技术从2004年被提出之后，先后出现了最大似然原理类算法和子空间类算法：最大似然原理类算法具有定位性能接近性能界克拉美罗下界的优势，子空间原理类算法具有较高的定位分辨率，两类方法均被广泛研究。

然而，随着分布式集群作战概念的提出，无人机蜂群等密集目标辐射源场景对定位系统的精度及分辨力提出了更为严苛的要求。目前的直接定位算法仍存在如下不足：(1)最大似然原理类算法计算复杂度高、所需采样快拍数较大、分辨率低、未知辐射源数目时无法准确定位；(2)子空间原理类算法在辐射源之间存在互相关情况时定位效果不佳，甚至失效。

发明内容

为了克服现有技术的不足，本发明提供一种稀疏贝叶斯学习的直接定位方法。基于稀疏贝叶斯学习的参数估计方法因具备超分辨、高精度、强稳健等优越性能而被广泛应用于空间谱估计中。为了解决现有直接定位技术存在的上述不足，本发明以直接定位技术为核心，提出了一种稀疏贝叶斯学习的直接定位方法，利用广义双帕累托(GDP)先验分布比传统分布更强的稀疏特性，保证算法不受相关信号影响的同时，还加强了信号的稀疏性，最终表现为分辨率和定位精度的提高以及低信噪比、小快拍数情况下的稳健性。

本发明解决其技术问题所采用的技术方案包括以下步骤：

步骤1.构建直接定位接收数据的空间稀疏模型；

(1)阵列接收数据模型建立：二维平面内有N个接收站，在接收站的远场处有Q个固定的信号辐射源，每个信号辐射源的位置为p_q＝[x_q,y_q]^T,q＝1,2,…,Q，p_q为第q个辐射源的位置，x_q为辐射源横坐标，y_q为辐射源纵坐标，各辐射源辐射的信号为互相关的窄带信号，N个接收站接收辐射源辐射的信号并将数据传送到中心处理站实现对辐射源位置的估计；

第n个接收站接收到的K个快拍的数据可建模表示为：

r_n＝A_ns_n+w_n (1)

每个接收站阵元数目为M，采样快拍数为K，为接收站接收到的K快拍数据，/>为阵列流型矩阵，a_n(p_q)为第q个辐射源到达第n个观测站时的导向矢量，主要由目标辐射源与接收站之间的相对位置关系决定，/>为接收站接收到的信号包络，/>为第n个接收站的加性噪声；

在直接定位问题中，辐射源数目Q往往未知，定位目标就是直接利用观测到的数据r_n定位出Q个辐射源的空间位置，n＝1,2,…,N；

(2)构建接收数据的空间稀疏模型：为了将源定位问题转化为一个稀疏表示的问题，受稀疏重构基础理论的启发，本发明考虑将接收数据建立为空域稀疏的模型，将感兴趣的空间范围划分网格点，共G个网格，每个网格表示为潜在的辐射源的位置，划分的网格足够小，真实辐射源的位置在格点上；

将空间谱相应样本排列在G×1的向量中，因为G远大于信号源数Q，空间信号向量是稀疏的，理想情况下，它的大多数元素都接近于0，只有Q个元素与零元素有较大差异，因此采用稀疏恢复的方法来获得辐射源的位置估计；

辐射源位于其中的Q个网格点上，单采样快拍下，接收信号矢量式(1)表示为：

其中，为已知的稀疏表示的过完备阵列流型矩阵，G为划分网格点的格点数，/>为格点上的位置矢量，并假设M<<G、Q<<G，因为目标辐射源仅位于其中的几个格点上，所以/>表示为空间稀疏信号，仅在对应的有信号源的网格点对应的位置上才有数值，其他元素均为0；

由稀疏信号的概念，可知每一列都是是具有稀疏度为Q的稀疏向量，其中仅有Q个非零元，其他元素均为零元素，其中，Q个非零元对应的是真实辐射源的位置；

步骤2.信号的GDP稀疏先验分布假设：

(1)噪声的先验假设；

每个接收站每个阵元接收噪声均为方差的独立的零均值高斯分布，将w_n概率密度函数表示为所有采样快拍下噪声的联合概率密度函数：

噪声方差未知，为简化分析，噪声服从伽马分布：

p(α₀|c,d)＝Gamma(α₀|c,d) (4)

(2)接收数据先验假设；

由噪声功率的假设得到，接收数据的先验分布表示为：

(3)为了对空间稀疏信号建模为GDP先验分布。

建立三层先验分布来实现；首先第一层先验，假设的每一块拍数据服从独立的复高斯分布：

协方差矩阵位Λ＝diag(α)，且表示每个网格点信号的功率；

接下来第二层先验对信号的方差α中的α_g建模，为了得到GDP先验分布，此处对α建立超参数为ξ的伽马分布，具体概率密度函数表示为：

其中，ξ＝[ξ₁,ξ₂,…,ξ_G]^T为超参数，为了得到GDP分布先验，第三层先验假设超参数ξ为伽马分布的随机变量，具备如下的概率分布密度函数：

其中，超参数h为接近于0的正常数；

综合式(5)至(7)三层先验，得到依赖于α和ξ的的概率分布密度函数表达式如下式：

其中，h为可以调节稀疏度的超参数。

步骤3.基于GDP先验分布的稀疏贝叶斯学习参数推导；

基于稀疏贝叶斯学习的直接定位算法主要就是根据接收到的数据矩阵和假设的概率分布依据贝叶斯规则推导出后验概率密度，通过最大化后验概率密度函数得到目标辐射源的位置；

通过前面的先验概率假设，得到后验概率密度函数为：

使得公式(10)中目标函数最大化的稀疏信号中非0元的位置即为目标辐射源的位置。

所述步骤3中，为了最大化式后验概率密度函数，采用EM算法，E步骤主要计算完全似然对数的期望；M步骤使用算法最大化期望值；

(1)E步骤表达为下式：

其中表示计算服从于/>的期望值，并且有

为观测站输出的协方差矩阵，有/>省略掉式(11)中不相关的项，有：

其中，Σ_n,gg表示为Σ_n的第g行g列的元素；

(2)在M步中，分别最大化服从于参数α_g和ξ_g的式(13)，以得到参数α_g和ξ_g的更新迭代公式；

对(13)关于α_g求导得：

采用统计优化中的最大似然准则，令式(14)为0，直接得到更新α_g的表达式；令并带入式(14)有/>令得到α_g的定点更新表达式如下：

其中，g＝1,2,…,G。采用同样的方式令得到ξ_g的更新公式如下：

同时求得噪声精度α₀对数函数：

式(17)中，令上式对α₀求偏导为0，得到α₀的更新公式如下：

通过求解公式，对超参数和稀疏信号的均值、方差进行交替迭代求解，最终获得空间稀疏信号的恢复，通过判断稀疏信号较大功率值对应的位置，即可判断出目标辐射源的空间位置。

本发明的有益效果在于将稀疏贝叶斯学习引入直接定位系统，从稀疏类方法入手，解决了子空间原理类直接算法无法定位相关辐射源的问题；同时性能上的提升解决了现有最大似然原理类直接定位算法所需采样快拍数较大、分辨率不够、未知辐射源数目时无法准确定位等问题。

附图说明

图1为本发明方法的定位实现流程图。

图2为定位场景示意图。

图3为空间稀疏表示示意图。

图4为信号的空间稀疏表示示意图。

图5为不同分布的稀疏性分析示意图。

图6各参数间的关系示意图。

图7为本发明方法定位谱峰图。

图8为本发明方法与传统定位算法随信噪比变化的RMSE曲线图。

图9为本发明方法与传统定位算法随快拍数变化的RMSE曲线图。

具体实施方式

下面结合附图和实施例对本发明进一步说明。

为使本发明的目的、技术方案和优点更加清楚，以下结合附图实施例对本发明作进一步详细描述。

1.构建直接定位接收数据的空间稀疏模型；

2.对信号的分布假设为GDP稀疏先验分布；

3.推导基于GDP先验分布的稀疏贝叶斯学习参数。

实施例：具体步骤如下：

步骤一.构建直接定位接收数据的空间稀疏模型：

(1)阵列接收数据模型建立：如图2所示，二维平面内有N个接收站，在接收站的远场p_q＝[x_q,y_q]^T,q＝1,2,…,Q处有Q个固定的信号辐射源，各辐射源辐射的信号为互相关的窄带信号，N个接收站接收辐射源辐射的信号并将数据传送到中心处理站实现对辐射源位置的估计。

那么第n个接收站接收到的K个快拍的数据可建模表示为：

r_n＝A_ns_n+w_n (19)

每个接收站阵元数目为M，采样快拍数为K，为接收站接收到的K快拍数据，/>为阵列流型矩阵，a_n(p_q)为第q个辐射源到达第n个观测站时的导向矢量，主要由目标辐射源与接收站之间的相对位置关系决定，/>为接收站接收到的信号包络，/>为第n个接收站的加性噪声。

(2)构建接收数据的空间稀疏模型：将空间谱相应样本排列在G×1的向量中，假设辐射源(灰点)刚好位于其中的Q个网格点上，因为G远大于信号源数Q，空间信号向量是稀疏的，理想情况下，它的大多数元素都接近于0，只有Q个元素与零元素有较大差异，因此可以采用稀疏恢复的方法来获得辐射源的位置估计。接收信号矢量式(19)稀疏建模表示为：

其中，为已知的稀疏表示的过完备阵列流型矩阵，G为划分网格点的格点数，/>为格点上的位置矢量，并假设M<<G、Q<<G，因为目标辐射源仅位于其中的几个格点上，所以/>表示为空间稀疏信号，仅在对应的有信号源的网格点对应的位置上才有数值，其他元素均为0，可将其表示为图4的形式。

步骤二.对信号的分布假设为GDP稀疏先验分布：

(1)噪声的先验假设。

每个接收站每个阵元接收噪声均为方差的独立的零均值高斯分布，将w_n概率密度函数表示为所有采样快拍下噪声的联合概率密度函数：

噪声方差未知，为简化分析，假设噪声服从伽马分布：

p(α₀|c,d)＝Gamma(α₀|c,d) (22)

(2)接收数据先验假设。

由噪声功率的假设可以得到，接收数据的先验分布表示为：

(3)为了对空间稀疏信号建模为GDP先验分布。

本发明建立三层先验分布来实现。首先第一层先验，假设的每一块拍数据服从独立的复高斯分布：

协方差矩阵位Λ＝diag(α)，且表示每个网格点信号的功率。

接下来第二层先验对信号的方差α中的α_g建模，为了得到GDP先验分布，此处对α建立超参数为ξ的伽马分布，具体概率密度函数表示为：ξ₁

其中，ξ＝[ξ₁,ξ₂,…,ξ_G]^T为超参数，为了得到GDP分布先验，第三层先验假设超参数ξ为伽马分布的随机变量，具备如下的概率分布密度函数：

其中，超参数h为接近于0的正常数；

综合式(24)至(26)三层先验，得到依赖于α和ξ的的概率分布密度函数表达式如下式：

其中，h为可以调节稀疏度的超参数。

为了证明GDP先验的稀疏特性，图5比较了不同的对数分布，即文献[D.P.Wipf,B.D.Rao,and S.Nagarajan,Latent variable Bayesian models for promotingsparsity[J].,IEEE Trans.Inf.Theory.,Sep.2011,vol.57,no.9,pp.6236–6255]中的GDP分布、高斯分布，文献[N.L.Pedersen,D.Shutin,C.Manch′on,and B.H.Fleury,Sparseestimation using bayesian hierarchical prior modeling for real and complexmodels[J].,[Online].Available:http://arxiv.org/pdf/1108.4324 2011]中的Laplace1分布。从图中可以看出，GDP分布出现了一个较尖锐的凹峰而其他的分布函数均没有这样的峰值，且随着h的减小，峰值在原点处变得更加集中。虽然高斯分布或拉普拉斯分布的峰值也可以通过改变超参数来更加集中，但它们不会在原点处形成像GDP那样的尖峰。因此，对于复数信号，加入GDP分布比传统分布更能增强信号的稀疏特性。

至此已经建立了直接定位模型中需要的参数概率分布，如图6所示为各参数之间的相互关系。

步骤三.推导基于GDP先验分布的稀疏贝叶斯学习参数：

(1)目标函数表达式。

通过前面的先验概率假设，得到后验概率密度函数为：

使得上述目标函数最大化的稀疏信号中非0元的位置即为目标辐射源的位置。为了最大化式后验概率密度函数，本发明方法采用EM算法：E步骤主要计算完全似然对数的期望；M步使用算法最大化期望值。

(2)E步骤表达为下式：

其中表示计算服从于/>的期望值，并且有

为观测站输出的协方差矩阵，有/>省略掉式(29)中不相关的项，有：

其中，Σ_n,gg表示为Σ_n的第g行g列的元素。

(3)在M步中，分别最大化服从于参数α_g和ξ_g的式(31)，以得到参数α_g和ξ_g的更新迭代公式；

对(31)关于α_g求导得：

令并带入式(32)有/>令/>得到α_g的定点更新表达式如下

其中，g＝1,2,…,G。采用同样的方式令得到ξ_g的更新公式如下：

同时求得噪声精度α₀对数函数

上式中，令上式对α₀求偏导为0，得到α₀的更新公式如下：

通过求解上述公式，对超参数和稀疏信号的均值、方差进行交替迭代求解，最终获得空间稀疏信号的恢复，通过判断稀疏信号较大功率值对应的位置，即可判断出目标辐射源的空间位置。

示例：

采用5个静止的8阵元的均匀线性阵接收站，每个接收站相邻传感器的间距为半波长，目标辐射源分别位于为p₀＝[-1.2,1.2]^T(km)p₁＝[0,0]^T(km)，假设两目标辐射源之间互相关，5个观测站位置分别为u₁＝(-5,-5)^T(km)、u₂＝(-3.5,-5)^T(km)、u₃＝(-2,-5)^T(km)、u₄＝(-0.5,-5)^T(km)、u₅＝(1,-5)^T(km)。

设定SNR＝10dB，每个观测站对接收信号进行快拍数K＝64次采样，得到本发明方法的定位谱峰图如图7所示，由图7可知，谱峰所对应的位置 与设定的辐射源位置p₀和p₁一致，验证了本发明方法的定位准确性，且定位谱峰很尖锐，分辨率高；图8为本发明方法与传统定位算法随信噪比变化的RMSE曲线图，图9为本发明方法与传统定位算法随快拍数变化的RMSE曲线图，由图中可以明显地看出，在定位相关信号源时，传统最大似然类直接定位方法和子空间类直接定位算法均失效，误差水平保持在一个较高的水平，而本发明方法均方误差值要低得多，定位精度高，且在低信噪比和小快拍时定位误差就能达到较低水平，算法稳健。实施例中的图7-图9证明了上述观点。

Claims (2)

1.一种稀疏贝叶斯学习的直接定位方法，其特征在于包括下述步骤：

步骤1.构建直接定位接收数据的空间稀疏模型；

(1)阵列接收数据模型建立：二维平面内有N个接收站，在接收站的远场处有Q个固定的信号辐射源，每个信号辐射源的位置为p_q＝[x_q,y_q]^T,q＝1,2,…,Q，p_q为第q个辐射源的位置，x_q为辐射源横坐标，y_q为辐射源纵坐标，各辐射源辐射的信号为互相关的窄带信号，N个接收站接收辐射源辐射的信号并将数据传送到中心处理站实现对辐射源位置的估计；

第n个接收站接收到的K个快拍的数据可建模表示为：

r_n＝A_ns_n+w_n (1)

每个接收站阵元数目为M，采样快拍数为K，为接收站接收到的K快拍数据，为阵列流型矩阵，a_n(p_q)为第q个辐射源到达第n个观测站时的导向矢量，主要由目标辐射源与接收站之间的相对位置关系决定，/>为接收站接收到的信号包络，/>为第n个接收站的加性噪声；

在直接定位问题中，辐射源数目Q往往未知，定位目标就是直接利用观测到的数据r_n定位出Q个辐射源的空间位置，n＝1,2,…,N；

(2)构建接收数据的空间稀疏模型：将接收数据建立为空域稀疏的模型，将感兴趣的空间范围划分网格点，共G个网格，每个网格表示为潜在的辐射源的位置，划分的网格足够小，真实辐射源的位置在格点上；

将空间谱相应样本排列在G×1的向量中，因为G远大于信号源数Q，空间信号向量是稀疏的，理想情况下，它的大多数元素都接近于0，只有Q个元素与零元素有较大差异，因此采用稀疏恢复的方法来获得辐射源的位置估计；

辐射源位于其中的Q个网格点上，单采样快拍下，接收信号矢量式(1)表示为：

其中，为已知的稀疏表示的过完备阵列流型矩阵，G为划分网格点的格点数，/>为格点上的位置矢量，并假设M<<G、Q<<G，因为目标辐射源仅位于其中的几个格点上，所以/>表示为空间稀疏信号，仅在对应的有信号源的网格点对应的位置上才有数值，其他元素均为0；

由稀疏信号的概念，可知每一列都是是具有稀疏度为Q的稀疏向量，其中仅有Q个非零元，其他元素均为零元素，其中，Q个非零元对应的是真实辐射源的位置；

步骤2.信号的GDP稀疏先验分布假设：

(1)噪声的先验假设；

每个接收站每个阵元接收噪声均为方差的独立的零均值高斯分布，将w_n概率密度函数表示为所有采样快拍下噪声的联合概率密度函数：

噪声方差未知，为简化分析，噪声服从伽马分布：

p(α₀|c,d)＝Gamma(α₀|c,d) (4)

(2)接收数据先验假设；

由噪声功率的假设得到，接收数据的先验分布表示为：

(3)为了对空间稀疏信号建模为GDP先验分布；

建立三层先验分布来实现；首先第一层先验，假设的每一块拍数据服从独立的复高斯分布：

协方差矩阵位Λ＝diag(α)，且表示每个网格点信号的功率；

接下来第二层先验对信号的方差α中的α_g建模，为了得到GDP先验分布，此处对α建立超参数为ξ的伽马分布，具体概率密度函数表示为：

其中，ξ＝[ξ₁,ξ₂,…,ξ_G]^T为超参数，为了得到GDP分布先验，第三层先验假设超参数ξ为伽马分布的随机变量，具备如下的概率分布密度函数：

其中，超参数h为接近于0的正常数；

综合式(5)至(7)三层先验，得到依赖于α和ξ的的概率分布密度函数表达式如下式：

其中，h为可以调节稀疏度的超参数；

步骤3.基于GDP先验分布的稀疏贝叶斯学习参数推导；

基于稀疏贝叶斯学习的直接定位算法主要就是根据接收到的数据矩阵和假设的概率分布依据贝叶斯规则推导出后验概率密度，通过最大化后验概率密度函数得到目标辐射源的位置；

通过前面的先验概率假设，得到后验概率密度函数为：

使得公式(10)中目标函数最大化的稀疏信号中非0元的位置即为目标辐射源的位置。

2.根据权利要求1所述的稀疏贝叶斯学习的直接定位方法，其特征在于：

所述步骤3中，为了最大化式后验概率密度函数，采用EM算法，E步骤主要计算完全似然对数的期望；M步骤使用算法最大化期望值；

(1)E步骤表达为下式：

其中表示计算服从于/>的期望值，并且有

为观测站输出的协方差矩阵，有/>省略掉式(11)中不相关的项，有：

其中，Σ_n,gg表示为Σ_n的第g行g列的元素；

(2)在M步中，分别最大化服从于参数α_g和ξ_g的式(13)，以得到参数α_g和ξ_g的更新迭代公式；

对(13)关于α_g求导得：

采用统计优化中的最大似然准则，令式(14)为0，直接得到更新α_g的表达式；令并带入式(14)有/>令得到α_g的定点更新表达式如下：

其中，g＝1,2,…,G；采用同样的方式令得到ξ_g的更新公式如下：

同时求得噪声精度α₀对数函数：

式(17)中，令上式对α₀求偏导为0，得到α₀的更新公式如下：

通过求解公式，对超参数和稀疏信号的均值、方差进行交替迭代求解，最终获得空间稀疏信号的恢复，通过判断稀疏信号较大功率值对应的位置，即可判断出目标辐射源的空间位置。