CN115687890A - 一种拖尾噪声背景下的机动目标分布式状态估计方法 - Google Patents

一种拖尾噪声背景下的机动目标分布式状态估计方法 Download PDF

Info

Publication number
CN115687890A
CN115687890A CN202211422756.6A CN202211422756A CN115687890A CN 115687890 A CN115687890 A CN 115687890A CN 202211422756 A CN202211422756 A CN 202211422756A CN 115687890 A CN115687890 A CN 115687890A
Authority
CN
China
Prior art keywords
model
distribution
target
parameter
sensor
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211422756.6A
Other languages
English (en)
Inventor
鲁瑞莲
甘露
高林
李改有
廖红舒
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN202211422756.6A priority Critical patent/CN115687890A/zh
Publication of CN115687890A publication Critical patent/CN115687890A/zh
Pending legal-status Critical Current

Links

Images

Landscapes

  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明属于智能信号处理技术领域，具体的说是一种拖尾噪声背景下的机动目标分布式状态估计方法。本发明提出一种拖尾噪声背景下的机动目标分布式状态估计方法。本发明假设机动目标过程噪声与测量噪声均为拖尾噪声，通过将拖尾噪声建模为高斯‑伽马分布，并利用变分贝叶斯方法估计过程拖尾噪声、量测拖尾噪声相关分布参数，并通过对过程噪声相关分布参数进行一致性融合处理从而提高过程噪声估计精度。接下来对各个机动模型伪似然进行求解，从而解决了拖尾噪声背景下的机动目标多模型状态估计问题，最后基于GCI融合对目标状态后验估计结果进行一致性融合，从而提高目标运动状态估计精度。

Description

一种拖尾噪声背景下的机动目标分布式状态估计方法
技术领域
本发明属于智能信号处理技术领域，具体的说是一种拖尾噪声背景下的分布式机动目标 状态估计方法。
背景技术
跟踪带有敏捷机动目标和异常值损坏的测量的场景常常涉及重尾噪声背景下的状态估计 问题，在这些情况下，即使状态空间模型是线性的，将噪声建模为高斯的卡尔曼滤波器(KF) 的性能恶化，出现了基于重拖尾分布的异常值鲁棒估计量。例如用于建模过程和测量噪声的 Student's-t分布可以被视为广义高斯分布，通过一个可调参数提升分布尾部的自由度(Dof)， 从而保证噪声的拖尾特性。假设状态和过程、测量噪声的联合概率密度函数(PDF)为 Student's-t分布，并基于t分布的线性变换将预测、后验状态PDF也近似为Student's-t， 或者将Student's-t分布扩展为高斯分布与伽马分布，基于此对重拖尾估计问题进行求解。
到目前为止，现有的Student's-t滤波器大多数都是针对非机动目标或者是单个传感器 的基础上实现的，在最小方差(MV)估计意义上寻求Student's-t的最优融合。首先，这些融 合方法更倾向于常规高斯分布，与重拖尾估计的思想缺乏一致性，重拖尾估计的异常值的鲁 棒性源于重尾后验，并且当涉及到通信和计算都受到高度限制的大规模传感器网络时，计算 效率和对节点故障的容忍度也是需要考虑的重要因素。
本发明基于上述问题提出了一种多传感器Student's-t滤波器去解决拖尾噪声背景下的 机动目标分布式状态估计问题，它不像MV那样寻求最优融合，而是与重尾估计一致的次优鲁 棒融合。为此本发明采用GCI融合与变分贝叶斯估计，首先通过变分贝叶斯方法对相关过程 拖尾与量测拖尾分布参数进行估计，并基于估计结果对过程拖尾噪声相关分布参数、目标运 动模型概率等参数进行融合，从而提高拖尾噪声背景下的机动目标分布式状态估计精度与鲁 棒性。
发明内容
针对上述问题，本发明提出一种拖尾噪声背景下的机动目标分布式状态估计方法。本发 明假设机动目标过程噪声与测量噪声均为拖尾噪声，通过将拖尾噪声建模为高斯-伽马分布， 并利用变分贝叶斯方法估计过程拖尾噪声、量测拖尾噪声相关分布参数，并通过对过程噪声 相关分布参数进行一致性融合处理从而提高过程噪声估计精度。接下来对各个机动模型伪似 然进行求解，从而解决了拖尾噪声背景下的机动目标多模型状态估计问题，最后基于GCI融 合对目标状态后验估计结果进行一致性融合，从而提高目标运动状态估计精度。本发明方法 解决实际工程应用中，过程噪声与量测噪声均为拖尾噪声背景下的机动目标分布式状态估计 问题，具有较高的参数估计精度、状态估计精度、较好的对复杂应用背景的适应性和鲁棒性、 较高的计算效率和对节点故障的容忍度等优点，可以满足工程应用需求。本发明大大提高了 算法对复杂场景的适应性和鲁棒性，提高了拖尾噪声背景下的机动目标分布式状态的估计精 度。
本发明的技术方案为：
一种拖尾噪声背景下的机动目标分布式状态估计方法，包括以下步骤：
S1、初始化系统参数，
初始化传感器网络参数：定义传感器网络
Figure BDA0003943165170000021
其中传感器节点S表示接收和 处理数据的节点，通信节点C表示完成数据传输的节点，连接链路
Figure BDA0003943165170000022
表示可通信节 点之间的通信链路，初始化第s个传感器的邻接节点集合为
Figure BDA0003943165170000023
其中p 表示传感器s的邻接节点，初始化传感器一致性加权系数为w；
初始化离散时间状态空间转移模型：定义k时刻目标运动模型集合为
Figure BDA0003943165170000024
其 中M表示目标运动过程中包含运动模型个数，定义k时刻第s个传感器的运动状态为
Figure BDA0003943165170000025
则 有：
Figure BDA0003943165170000026
其中
Figure BDA0003943165170000027
表示第k时刻模型为i，fi(·)表示第i个运动模型的状态转移函数；过程噪声
Figure BDA0003943165170000028
服 从拖尾分布
Figure BDA0003943165170000029
其中p(·)表示概率密度函数，St(x|μ,Σ,v)表示均 值为μ、尺度矩阵为Σ、自由度为v的Student’s-t分布；
初始化量测空间转移模型：
Figure BDA00039431651700000210
其中
Figure BDA0003943165170000031
表示量测值，hs,i(·)表示第s个传感器中第i个运动模型对应的量测转移函数，量测噪 声
Figure BDA0003943165170000032
服从拖尾分布
Figure BDA0003943165170000033
S2、对目标各个模型的运动状态进行模型交互，
结合下式计算第s个传感器目标模型j的模型预测概率
Figure BDA0003943165170000034
与模型交互概率
Figure BDA0003943165170000035
Figure BDA0003943165170000036
其中
Figure BDA0003943165170000037
表示k-1时刻模型i的概率更新值(由上一时刻获取)，
Figure BDA0003943165170000038
表示 k时刻马尔科夫转移概率矩阵(TPM)，其中Pr(·)表示概率值；
结合下式计算目标第i个模型的运动状态交互结果
Figure BDA0003943165170000039
与状态协方差交互结果
Figure BDA00039431651700000310
Figure BDA00039431651700000311
其中，
Figure BDA00039431651700000312
与
Figure BDA00039431651700000313
分别表示k-1时刻模型j的状态估计值与状态估计协方差矩阵(由上一 时刻获取)，∑·表示求和操作。
S3、对目标各个模型的运动状态值、运动状态协方差矩阵进行一步预测计算，
结合下式计算目标第s个传感器第i个模型的运动状态一步预测值
Figure BDA00039431651700000314
与一步预测协方 差矩阵
Figure BDA00039431651700000315
Figure BDA00039431651700000316
其中Fi(·)表示状态转移函数fi(·)的雅克比矩阵，
Figure BDA00039431651700000317
表示k时刻模型i的过程噪声协方差矩阵，(·)T表示求转置操作。
S4、对目标运动状态分布与量测分布进行下述概率密度函数分解，
对目标运动状态分布进行下述分解：
Figure BDA0003943165170000041
其中，
Figure BDA0003943165170000042
表示目标从初始时刻到k-1时刻的各传感器量测集合，N(x|μ,Σ)表示均值为μ、 方差为Σ的高斯分布；其中
Figure BDA0003943165170000043
服从尺度矩阵
Figure BDA0003943165170000044
自由度为
Figure BDA0003943165170000045
的逆-wishart分布，即
Figure BDA0003943165170000046
G(x；a,b)表示形状参数为a、尺度参数为b的伽马分布；∫·表示求 积分操作；
对目标量测分布进行下述分解：
Figure BDA0003943165170000047
其中
Figure BDA0003943165170000048
服从尺度矩阵
Figure BDA0003943165170000049
自由度为
Figure BDA00039431651700000410
的逆-wishart分布，即
Figure BDA00039431651700000411
S5、对目标各个传感器与各运动模型的相关分布先验参数进行模型交互，包括：
结合下式计算步骤S4中目标i个模型的运动状态与量测相关参数逆-wishart分布
Figure BDA00039431651700000412
与
Figure BDA00039431651700000413
参数自由度交互结果
Figure BDA00039431651700000414
与尺度矩阵交互结果
Figure BDA00039431651700000415
Figure BDA00039431651700000416
其中
Figure BDA00039431651700000417
与
Figure BDA00039431651700000418
分别表示k-1时刻逆-wishart分布
Figure BDA00039431651700000419
与
Figure BDA00039431651700000420
的自由度，
Figure BDA00039431651700000421
与
Figure BDA00039431651700000422
分别表示 k-1时刻逆-wishart分布
Figure BDA00039431651700000423
与
Figure BDA00039431651700000424
的尺度矩阵；
结合下式计算第i个运动模型的伽马分布
Figure BDA00039431651700000425
与
Figure BDA00039431651700000426
形状参数交 互结果参数交互结果
Figure BDA00039431651700000427
与尺度参数交互结果
Figure BDA00039431651700000428
Figure BDA0003943165170000051
S6、对目标各个模型潜变量参数进行一步预测计算，包括：
结合下式计算第i个模型的逆-wishart分布
Figure BDA0003943165170000052
与
Figure BDA0003943165170000053
参数自由度预测值
Figure BDA0003943165170000054
与 尺度矩阵预测值
Figure BDA0003943165170000055
Figure BDA0003943165170000056
其中nx表示运动状态向量维度，
Figure BDA0003943165170000057
表示量测状态向量维度，ρ表示遗忘因子。
结合下式计算伽马分布
Figure BDA0003943165170000058
与
Figure BDA0003943165170000059
形状参数预测结果
Figure BDA00039431651700000510
与尺度参数预测结果
Figure BDA00039431651700000511
Figure BDA00039431651700000512
S7、结合变分贝叶斯算法求解各传感器目标状态更新方法与相关参数分布的超参数更新 方法，
设定k时刻各传感器目标运动状态及分布参数集合
Figure BDA00039431651700000513
结合下 式求解上述分布参数集合对应的联合概率密度函数，
Figure BDA00039431651700000514
结合上述联合概率密度函数计算各传感器目标分布参数边缘概率密度(VB-marginal)并 确定对应的超参数更新方法，
结合下式计算第s个传感器目标运动状态
Figure BDA00039431651700000515
边缘对数似然，
Figure BDA0003943165170000061
其中
Figure BDA0003943165170000062
E(·)表示 求期望操作，Eq(x)(x)＝∫xq(x)dx，其中
Figure BDA0003943165170000063
Figure BDA0003943165170000064
(·)-1表示求逆操作，
Figure BDA0003943165170000065
表示与
Figure BDA0003943165170000066
无关的常数项， 结合上式可得目标第i个模型对应的状态更新值
Figure BDA0003943165170000067
与状态更新协方差
Figure BDA0003943165170000068
Figure BDA0003943165170000069
其中，Hs,i表示第s个传感器第i个模型量测转移函数的hs,i(·)的雅克比矩阵，
Figure BDA00039431651700000610
表示新息 协方差矩阵。
结合下式计算
Figure BDA00039431651700000611
边缘似然，
Figure BDA00039431651700000612
其中，det(·)表示求行列式操作，
Figure BDA00039431651700000613
表示联合概率密度函数中与
Figure BDA00039431651700000614
无关的期望常数项；且：
Figure BDA00039431651700000615
综上可得到
Figure BDA00039431651700000616
对应的超参数更新步骤，
Figure BDA00039431651700000617
结合下式计算
Figure BDA0003943165170000071
边缘似然，
Figure BDA0003943165170000072
其中
Figure BDA0003943165170000073
Figure BDA0003943165170000074
表示与
Figure BDA0003943165170000075
无关的常数项；假定
Figure BDA0003943165170000076
结合下式计算
Figure BDA0003943165170000077
超参数更新步骤，
Figure BDA0003943165170000078
结合下式计算
Figure BDA0003943165170000079
边缘似然，
Figure BDA00039431651700000710
确定
Figure BDA00039431651700000711
超参数更新步骤，
Figure BDA00039431651700000712
其中
Figure BDA00039431651700000713
Ψ(·)表示双参数伽马分布。
可获得量测相关分布对应的超参数更新步骤：
确定
Figure BDA00039431651700000714
超参数更新步骤，
Figure BDA00039431651700000715
其中，
Figure BDA00039431651700000716
假定
Figure BDA0003943165170000081
确定
Figure BDA0003943165170000082
超参数更新步骤，
Figure BDA0003943165170000083
其中
Figure BDA0003943165170000084
确定
Figure BDA0003943165170000085
超参数更新步骤:
Figure BDA0003943165170000086
其中
Figure BDA0003943165170000087
S8、对步骤S7中的过程噪声相关分布逆-Wishart分布与伽马分布参数结合下式进行多 传感器一致性融合处理，
Figure BDA0003943165170000088
Figure BDA0003943165170000089
Figure BDA00039431651700000810
Figure BDA00039431651700000811
其中w表示融合权系数，L表示一致性融合次数。
S9、计算目标各传感器与各运动模型对应的伪似然值与k时刻对应的模型概率更新值，
由变分贝叶斯原理结合下式获得伪似然结果，
Figure BDA00039431651700000812
其中DKL(p||q)表示函数p与函数q之间的KL散度，结合S7结果对上式求解获得以下伪似然 值，
Figure BDA0003943165170000091
根据获得的伪似然结果结合下式获得k时刻第i个模型的模型概率更新值
Figure BDA0003943165170000092
Figure BDA0003943165170000093
其中∝表示正比于。
S10、对步骤S9中的模型概率更新值
Figure BDA0003943165170000094
结合下式进行多传感器一致性融合处理，
Figure BDA0003943165170000095
S11、基于目标运动模型概率更新值
Figure BDA0003943165170000096
与状态估计协方差矩阵
Figure BDA0003943165170000097
Figure BDA0003943165170000098
其中
Figure BDA0003943165170000099
表示第L次融合以后的模型概率更新值。
综上所述，得到目标在当前时刻的各传感器状态估计值与状态估计协方差矩阵，完成状 态更新。
本发明的效益是，
1.本发明针对过程与量测噪声拖尾背景下的敏捷机动目标跟踪问题，基于分布式变分贝 叶斯技术对目标的过程噪声与量测噪声协方差分布参数进行实时估计，从而提高了拖尾噪声 背景下的机动目标状态估计精度，弥补了传统卡尔曼滤波针对拖尾噪声估计精度急剧下降的 缺陷；
2.本发明通过在变分贝叶斯估计参数过程中，将各传感器获得的过程噪声分布相关参数 进行一致性融合处理，提高了变分估计的精度，从而提高了目标拖尾噪声背景下的状态估计 精度。
附图说明
图1是本发明的整体流程图；
图2是采用本发明方法的传感器分布图、目标轨迹分布图；
图3是采用最优状态估计、未进行参数估计与采用本发明方法时的分布式目标距离、速 度估计精度对比图，其中蒙特卡洛次数为100次；
图4是采用本发明方法时的单传感器与分布式目标距离估计精度对比图，其中蒙特卡洛 次数为100次；
图5是采用最优状态估计、未进行参数估计与采用本发明方法时的模型概率估计值，其 中蒙特卡洛次数为100次；
图6是采用本发明的方法与单传感器估计结果对比图。
具体实施方式
下面结合附图和实施例，详细描述本发明的技术方案：
实施例
步骤1初始化系统参数，
1.1初始化传感器网络参数：定义传感器网络
Figure BDA0003943165170000101
其中传感器节点S表示接 收和处理数据的节点，通信节点C表示完成数据传输的节点，连接链路
Figure BDA0003943165170000102
表示可通 信节点之间的通信链路，初始化第s个传感器的邻接节点集合为
Figure BDA0003943165170000103
其 中p表示传感器s的邻接节点，初始化传感器一致性加权系数为w；
本实例选用但不限于传感器节点个数为10个，传感器邻接节点判定方式为：当两传感器 之间的距离小于某一判定阈值时，判定两传感器为邻接节点，具体传感器分布及通信关系如 图2所示。
1.2初始化离散时间状态空间转移模型：定义k时刻目标运动模型集合为
Figure BDA0003943165170000111
其中M表示目标运动过程中包含运动模型个数，定义k时刻第s个传感器的运动状态为
Figure BDA0003943165170000112
则有：
Figure BDA0003943165170000113
其中
Figure BDA0003943165170000114
表示第k时刻模型为i，fi(·)表示第i个运动模型的状态转移函数；过程噪声
Figure BDA0003943165170000115
服 从拖尾分布
Figure BDA0003943165170000116
其中p(·)表示概率密度函数，St(x|μ,Σ,v)表示均 值为μ、尺度矩阵为Σ、自由度为v的Student’s-t分布。
本实例选用但不限于目标运动模型个数M＝3，分别为匀速运动、已知转弯率为±0.1rad/s的左右转弯模型，其中模型1表示目标为匀速运动，其状态转移矩阵表示为：
Figure BDA0003943165170000117
T表示采样时间间隔，本实例选用但不限于T＝1s。模型2与3表示已知转弯率的转弯模型， 其状态转移矩阵表示为：
Figure BDA0003943165170000118
本实例选用但不限于F2＝F(w1),w1＝0.1rad/s，F3＝F(w2),w2＝-0.1rad/sw；
本实例选用但不限于三种运动模型对应的运动状态向量为
Figure BDA0003943165170000119
其中x与y分别表 示目标在x轴与y轴的距离，
Figure BDA0003943165170000121
与
Figure BDA0003943165170000122
分别表示目标在x轴与y轴的速度，三种运动状态下的 过程噪声方差均为1m/s2；
1.3初始化量测空间转移模型：
Figure BDA0003943165170000123
其中
Figure BDA0003943165170000124
表示量测值，hs,i(·)表示第s个传感器中第i个运动模型对应的量测转移函数，量测噪 声
Figure BDA0003943165170000125
服从拖尾分布
Figure BDA0003943165170000126
本实例选用但不限于三种运动模型下的量测向量为[R,θ]T，其中R表示目标距离，θ表 示目标方位角，且三种模型下各传感器的量测噪声初始化为σR＝10m,σθ＝0.2°。
步骤2对目标各个模型的运动状态进行模型交互，
2.1结合下式计算第s个传感器目标模型j的模型预测概率
Figure BDA0003943165170000127
与模型交互概率
Figure BDA0003943165170000128
Figure BDA0003943165170000129
其中
Figure BDA00039431651700001210
表示k-1时刻模型i的概率更新值(由上一时刻获取)，
Figure BDA00039431651700001211
表示 k时刻马尔科夫转移概率矩阵(TPM)，其中Pr(·)表示概率值。
本实例中初始时刻模型概率初始化为
Figure BDA00039431651700001212
马尔科夫转移概率矩阵 (TPM)初始化为下值：
Figure BDA00039431651700001213
2.2结合下式计算目标第i个模型的运动状态交互结果
Figure BDA00039431651700001214
与状态协方差交互结果
Figure BDA00039431651700001215
Figure BDA00039431651700001216
其中，
Figure BDA0003943165170000131
与
Figure BDA0003943165170000132
分别表示k-1时刻模型j的状态估计值与状态估计协方差矩阵(由上一 时刻获取)，∑·表示求和操作。
步骤3对目标各个模型的运动状态值、运动状态协方差矩阵进行一步预测计算，
结合下式计算目标第s个传感器第i个模型的运动状态一步预测值
Figure BDA0003943165170000133
与一步预测协方 差矩阵
Figure BDA0003943165170000134
Figure BDA0003943165170000135
其中Fi表示状态转移函数fi(·)的雅克比矩阵，
Figure BDA0003943165170000136
表示k时刻模型i的过程噪声协方差矩阵， (·)T表示求转置操作。
步骤4对目标运动状态分布与量测分布进行下述概率密度函数分解，
4.1对目标运动状态分布进行下述分解：
Figure BDA0003943165170000137
其中，
Figure BDA0003943165170000138
表示目标从初始时刻到k-1时刻的各传感器量测集合，N(x|μ,Σ)表示均值为μ、 方差为Σ的高斯分布；其中
Figure BDA0003943165170000139
服从尺度矩阵
Figure BDA00039431651700001310
自由度为
Figure BDA00039431651700001311
的逆-wishart分布，即
Figure BDA00039431651700001312
G(x；a,b)表示形状参数为a、尺度参数为b的伽马分布；∫·表示求 积分操作；
4.2对目标量测分布进行下述分解：
Figure BDA00039431651700001313
其中
Figure BDA00039431651700001314
服从尺度矩阵
Figure BDA00039431651700001315
自由度为
Figure BDA00039431651700001316
的逆-wishart分布，即
Figure BDA00039431651700001317
步骤5对目标各个传感器与各运动模型的相关分布先验参数进行模型交互，包括：
5.1结合下式计算步骤S4中目标i个模型的运动状态与量测相关参数逆-wishart分布
Figure BDA0003943165170000141
与
Figure BDA0003943165170000142
参数自由度交互结果
Figure BDA0003943165170000143
与尺度矩阵交互结果
Figure BDA0003943165170000144
Figure BDA0003943165170000145
其中
Figure BDA0003943165170000146
与
Figure BDA0003943165170000147
分别表示k-1时刻逆-wishart分布
Figure BDA0003943165170000148
与
Figure BDA0003943165170000149
的自由度，
Figure BDA00039431651700001410
与
Figure BDA00039431651700001411
分别表示 k-1时刻逆-wishart分布
Figure BDA00039431651700001412
与
Figure BDA00039431651700001413
的尺度矩阵。
本实例选用但不限于下述参数初始化结果：
Figure BDA00039431651700001414
其中t表示调谐参数且t＝5，nx＝4，
Figure BDA00039431651700001415
5.2结合下式计算第i个运动模型的伽马分布
Figure BDA00039431651700001416
与
Figure BDA00039431651700001417
形状参 数交互结果参数交互结果
Figure BDA00039431651700001418
与尺度参数交互结果
Figure BDA00039431651700001419
Figure BDA00039431651700001420
本实例选用但不限于下述参数初始化结果：
Figure BDA00039431651700001421
步骤6对目标各个模型潜变量参数进行一步预测计算，包括：
6.1.结合下式计算第i个模型的逆-wishart分布
Figure BDA00039431651700001422
与
Figure BDA00039431651700001423
参数自由度预测值
Figure BDA00039431651700001424
与尺度矩阵预测值
Figure BDA00039431651700001425
Figure BDA0003943165170000151
其中nx表示运动状态向量维度，
Figure BDA0003943165170000152
表示量测状态向量维度，ρ表示遗忘因子，本实例选用但 不限于ρ＝0.98。
6.2.结合下式计算伽马分布
Figure BDA0003943165170000153
与
Figure BDA0003943165170000154
形状参数预测结果
Figure BDA0003943165170000155
Figure BDA0003943165170000156
与尺度参数预测结果
Figure BDA0003943165170000157
Figure BDA0003943165170000158
步骤7结合变分贝叶斯算法求解各传感器目标状态更新方法与相关参数分布的超参数更 新方法，
7.1设定k时刻各传感器目标运动状态及分布参数集合
Figure BDA0003943165170000159
结 合下式求解上述分布参数集合对应的联合概率密度函数，
Figure BDA00039431651700001510
7.2结合上述联合概率密度函数计算各传感器目标分布参数边缘概率密度(VB-marginal) 并确定对应的超参数更新方法，
7.2.1结合下式计算第s个传感器目标运动状态
Figure BDA00039431651700001511
边缘对数似然，
Figure BDA00039431651700001512
其中
Figure BDA00039431651700001513
E(·)表示 求期望操作，Eq(x)(x)＝∫xq(x)dx，其中
Figure BDA0003943165170000161
Figure BDA0003943165170000162
(·)-1表示求逆操作，
Figure BDA0003943165170000163
表示与
Figure BDA0003943165170000164
无关的常数项， 结合上式可得目标第i个模型对应的状态更新值
Figure BDA0003943165170000165
与状态更新协方差
Figure BDA0003943165170000166
Figure BDA0003943165170000167
其中，Hs,i表示第s个传感器第i个模型量测转移函数的hs,i(·)的雅克比矩阵，
Figure BDA0003943165170000168
表示新息 协方差矩阵。
7.2.2结合下式计算
Figure BDA0003943165170000169
边缘似然，
Figure BDA00039431651700001610
其中，det(·)表示求行列式操作，
Figure BDA00039431651700001611
表示联合概率密度函数中与
Figure BDA00039431651700001612
无关的期望常数项；且：
Figure BDA00039431651700001613
综上可得到
Figure BDA00039431651700001614
对应的超参数更新步骤，
Figure BDA00039431651700001615
7.2.3结合下式计算
Figure BDA00039431651700001616
边缘似然，
Figure BDA00039431651700001617
其中
Figure BDA0003943165170000171
Figure BDA0003943165170000172
表示与
Figure BDA0003943165170000173
无关的常数项；假定
Figure BDA0003943165170000174
结合下式计算
Figure BDA0003943165170000175
超参数更新步骤，
Figure BDA0003943165170000176
7.2.4结合下式计算
Figure BDA0003943165170000177
边缘似然，
Figure BDA0003943165170000178
确定
Figure BDA0003943165170000179
超参数更新步骤，
Figure BDA00039431651700001710
其中
Figure BDA00039431651700001711
Ψ(·)表示双参数伽马分布。
7.2.5参考6.2.1～6.2.4可获得量测相关分布对应的超参数更新步骤：
确定
Figure BDA00039431651700001712
超参数更新步骤，
Figure BDA00039431651700001713
其中，
Figure BDA00039431651700001714
假定
Figure BDA00039431651700001715
确定
Figure BDA00039431651700001716
超参数更新步骤，
Figure BDA0003943165170000181
其中
Figure BDA0003943165170000182
确定
Figure BDA0003943165170000183
超参数更新步骤:
Figure BDA0003943165170000184
其中
Figure BDA0003943165170000185
步骤8对步骤7中的过程噪声相关分布逆-Wishart分布与伽马分布参数结合下式进行多 传感器一致性融合处理，
Figure BDA0003943165170000186
Figure BDA0003943165170000187
Figure BDA0003943165170000188
Figure BDA0003943165170000189
其中w表示融合权系数，L表示一致性融合次数，本实例选用但不限于L＝5。
步骤9计算目标各传感器与各运动模型对应的伪似然值与k时刻对应的模型概率更新值，
9.1由变分贝叶斯原理结合下式获得伪似然结果，
Figure BDA00039431651700001810
其中DKL(p||q)表示函数p与函数q之间的KL散度，结合步骤8结果对上式求解获得以下伪 似然值，
Figure BDA0003943165170000191
9.2根据1获得的伪似然结果结合下式获得k时刻第i个模型的模型概率更新值
Figure BDA0003943165170000192
Figure BDA0003943165170000193
其中∝表示正比于。
步骤10对步骤9中的模型概率更新值
Figure BDA0003943165170000194
结合下式进行多传感器一致性融合处理，
Figure BDA0003943165170000195
步骤11基于目标运动模型概率更新值
Figure BDA0003943165170000196
与状态估计协方差矩阵
Figure BDA0003943165170000197
Figure BDA0003943165170000198
其中
Figure BDA0003943165170000199
表示第L次融合以后的模型概率更新值。
综上所述，得到目标在当前时刻的各传感器状态估计值与状态估计协方差矩阵，完成状 态更新。
通过下面的仿真示例证明本发明的实用性：
1.仿真条件及参数
仿真场景为单目标跟踪场景，传感器位置分布图如图2所示，模型编号分别标记为1，2， 3，模型转移过程为：{1,2,3,2,1,3,1}，每个模型持续时间为1→30s,2→40s,3→40s，过 程噪声方差q＝1，量测噪声方差为[σR,σθ]＝[10m,0.2°]，总仿真时间长度为250s，过程噪声 异常值比例为10％，对应的过程噪声方差为100；量测噪声异常值比例为10％，对应的量测噪 声方差为40×[σR,σθ]，融合次数L＝5。
2.仿真内容和结果分析
图1为本发明方法对应的算法处理流程图；
图2是采用本发明方法中的传感器分布图与目标运动轨迹变化图。
图3是不同方法下的距离、速度均方根误差统计图
均方根误差随时间变化曲线图，其中蒙特卡洛次数为100次，从图中可看出本发明方法 中通过参数估计结合多传感器融合后相对于单传感器参数估计方法距离估计误差均有明显的 降低，证明了本发明方法的有效性。
图4本发明方法中的单传感器估计结果与分布式估计结果对比图
蒙特卡洛次数为100次，从图中可看出本发明方法中通过变分估计结合多传感器参数融 合相对于未进行多传感器融合估计误差均有明显的降低，证明了本发明方法的有效性。
图5是采用本发明方法时模型概率估计均方根误差随时间变化曲线图，其中蒙特卡洛次 数为100次，从图中可看出本发明方法具有较好的模型收敛概率。

Claims (1)

1.一种拖尾噪声背景下的机动目标分布式状态估计方法，其特征在于，包括以下步骤：
S1、初始化系统参数，
初始化传感器网络参数：定义传感器网络
Figure FDA0003943165160000011
其中传感器节点S表示接收和处理数据的节点，通信节点C表示完成数据传输的节点，连接链路
Figure FDA0003943165160000012
表示可通信节点之间的通信链路，初始化第s个传感器的邻接节点集合为
Figure FDA0003943165160000013
其中p表示传感器s的邻接节点，初始化传感器一致性加权系数为w；
初始化离散时间状态空间转移模型：定义k时刻目标运动模型集合为
Figure FDA0003943165160000014
其中M表示目标运动过程中包含运动模型个数，定义k时刻第s个传感器的运动状态为
Figure FDA0003943165160000015
则有：
Figure FDA0003943165160000016
其中
Figure FDA0003943165160000017
表示第k时刻模型为i，fi(·)表示第i个运动模型的状态转移函数；过程噪声
Figure FDA0003943165160000018
服从拖尾分布
Figure FDA0003943165160000019
其中p(·)表示概率密度函数，St(x|μ,Σ,v)表示均值为μ、尺度矩阵为Σ、自由度为v的Student’s-t分布；
初始化量测空间转移模型：
Figure FDA00039431651600000110
其中
Figure FDA00039431651600000111
表示量测值，hs,i(·)表示第s个传感器中第i个运动模型对应的量测转移函数，量测噪声
Figure FDA00039431651600000112
服从拖尾分布
Figure FDA00039431651600000113
S2、对目标各个模型的运动状态进行模型交互：
结合下式计算第s个传感器目标模型j的模型预测概率
Figure FDA00039431651600000114
与模型交互概率
Figure FDA00039431651600000115
Figure FDA00039431651600000116
其中
Figure FDA00039431651600000117
表示k-1时刻模型i的概率更新值(由上一时刻获取)，
Figure FDA00039431651600000118
表示k时刻马尔科夫转移概率矩阵(TPM)，其中Pr(·)表示概率值；
结合下式计算目标第i个模型的运动状态交互结果
Figure FDA0003943165160000021
与状态协方差交互结果
Figure FDA0003943165160000022
Figure FDA0003943165160000023
其中，
Figure FDA0003943165160000024
与
Figure FDA0003943165160000025
分别表示k-1时刻模型j的状态估计值与状态估计协方差矩阵(由上一时刻获取)，∑·表示求和操作；
S3、对目标各个模型的运动状态值、运动状态协方差矩阵进行一步预测计算：
结合下式计算目标第s个传感器第i个模型的运动状态一步预测值
Figure FDA0003943165160000026
与一步预测协方差矩阵
Figure FDA0003943165160000027
Figure FDA0003943165160000028
其中Fi(·)表示状态转移函数fi(·)的雅克比矩阵，
Figure FDA0003943165160000029
表示k时刻模型i的过程噪声协方差矩阵，(·)T表示求转置操作；
S4、对目标运动状态分布与量测分布进行下述概率密度函数分解，
对目标运动状态分布进行下述分解：
Figure FDA00039431651600000210
其中，
Figure FDA00039431651600000211
表示目标从初始时刻到k-1时刻的各传感器量测集合，N(x|μ,Σ)表示均值为μ、方差为Σ的高斯分布；其中
Figure FDA00039431651600000212
服从尺度矩阵
Figure FDA00039431651600000213
自由度为
Figure FDA00039431651600000214
的逆-wishart分布，即
Figure FDA00039431651600000215
G(x；a,b)表示形状参数为a、尺度参数为b的伽马分布；∫·表示求积分操作；
对目标量测分布进行下述分解：
Figure FDA0003943165160000031
其中
Figure FDA0003943165160000032
服从尺度矩阵
Figure FDA0003943165160000033
自由度为
Figure FDA0003943165160000034
的逆-wishart分布，即
Figure FDA0003943165160000035
S5、对目标各个传感器与各运动模型的相关分布先验参数进行模型交互，包括：
结合下式计算步骤S4中目标i个模型的运动状态与量测相关参数逆-wishart分布
Figure FDA0003943165160000036
与
Figure FDA0003943165160000037
参数自由度交互结果
Figure FDA0003943165160000038
与尺度矩阵交互结果
Figure FDA0003943165160000039
Figure FDA00039431651600000310
其中
Figure FDA00039431651600000311
与
Figure FDA00039431651600000312
分别表示k-1时刻逆-wishart分布
Figure FDA00039431651600000313
与
Figure FDA00039431651600000314
的自由度，
Figure FDA00039431651600000315
与
Figure FDA00039431651600000316
分别表示k-1时刻逆-wishart分布
Figure FDA00039431651600000317
与
Figure FDA00039431651600000318
的尺度矩阵；
结合下式计算第i个运动模型的伽马分布
Figure FDA00039431651600000319
与
Figure FDA00039431651600000320
形状参数交互结果参数交互结果
Figure FDA00039431651600000321
与尺度参数交互结果
Figure FDA00039431651600000322
Figure FDA00039431651600000323
Figure FDA00039431651600000324
S6、对目标各个模型潜变量参数进行一步预测计算，包括：
结合下式计算第i个模型的逆-wishart分布
Figure FDA00039431651600000325
与
Figure FDA00039431651600000326
参数自由度预测值
Figure FDA00039431651600000327
与尺度矩阵预测值
Figure FDA00039431651600000328
Figure FDA00039431651600000329
其中nx表示运动状态向量维度，
Figure FDA0003943165160000041
表示量测状态向量维度，ρ表示遗忘因子；
结合下式计算伽马分布
Figure FDA0003943165160000042
与
Figure FDA0003943165160000043
形状参数预测结果
Figure FDA0003943165160000044
与尺度参数预测结果
Figure FDA0003943165160000045
Figure FDA0003943165160000046
Figure FDA0003943165160000047
S7、结合变分贝叶斯算法求解各传感器目标状态更新方法与相关参数分布的超参数更新方法，
设定k时刻各传感器目标运动状态及分布参数集合
Figure FDA0003943165160000048
结合下式求解上述分布参数集合对应的联合概率密度函数，
Figure FDA0003943165160000049
结合联合概率密度函数计算各传感器目标分布参数边缘概率密度并确定对应的超参数更新方法，
结合下式计算第s个传感器目标运动状态
Figure FDA00039431651600000410
边缘对数似然，
Figure FDA00039431651600000411
其中
Figure FDA00039431651600000412
E(·)表示求期望操作，Eq(x)(x)＝∫xq(x)dx，其中
Figure FDA00039431651600000413
Figure FDA00039431651600000414
(·)-1表示求逆操作，
Figure FDA00039431651600000415
表示与
Figure FDA00039431651600000416
无关的常数项，结合上式可得目标第i个模型对应的状态更新值
Figure FDA00039431651600000417
与状态更新协方差
Figure FDA00039431651600000418
Figure FDA0003943165160000051
其中，Hs,i表示第s个传感器第i个模型量测转移函数的hs,i(·)的雅克比矩阵，
Figure FDA0003943165160000052
表示新息协方差矩阵；
结合下式计算
Figure FDA0003943165160000053
边缘似然，
Figure FDA0003943165160000054
其中，det(·)表示求行列式操作，
Figure FDA0003943165160000055
表示联合概率密度函数中与
Figure FDA0003943165160000056
无关的期望常数项；且：
Figure FDA0003943165160000057
综上可得到
Figure FDA0003943165160000058
对应的超参数更新步骤，
Figure FDA0003943165160000059
Figure FDA00039431651600000510
结合下式计算
Figure FDA00039431651600000511
边缘似然，
Figure FDA00039431651600000512
其中
Figure FDA00039431651600000513
Figure FDA00039431651600000514
表示与
Figure FDA00039431651600000515
无关的常数项；假定
Figure FDA00039431651600000516
结合下式计算
Figure FDA00039431651600000517
超参数更新步骤，
Figure FDA0003943165160000061
Figure FDA0003943165160000062
结合下式计算
Figure FDA0003943165160000063
边缘似然，
Figure FDA0003943165160000064
确定
Figure FDA0003943165160000065
超参数更新步骤，
Figure FDA0003943165160000066
其中
Figure FDA0003943165160000067
Ψ(·)表示双参数伽马分布；
获得量测相关分布对应的超参数更新步骤：
确定
Figure FDA0003943165160000068
超参数更新步骤，
Figure FDA0003943165160000069
其中，
Figure FDA00039431651600000610
假定
Figure FDA00039431651600000611
确定
Figure FDA00039431651600000612
超参数更新步骤，
Figure FDA00039431651600000613
其中
Figure FDA00039431651600000614
确定
Figure FDA00039431651600000615
超参数更新步骤:
Figure FDA0003943165160000071
其中
Figure FDA0003943165160000072
S8、对步骤S7中的过程噪声相关分布逆-Wishart分布与伽马分布参数结合下式进行多传感器一致性融合处理，
Figure FDA0003943165160000073
Figure FDA0003943165160000074
Figure FDA0003943165160000075
Figure FDA0003943165160000076
其中w表示融合权系数，L表示一致性融合次数；
S9、计算目标各传感器与各运动模型对应的伪似然值与k时刻对应的模型概率更新值，
由变分贝叶斯原理结合下式获得伪似然结果，
Figure FDA0003943165160000077
其中DKL(p||q)表示函数p与函数q之间的KL散度，结合S7结果对上式求解获得以下伪似然值，
Figure FDA0003943165160000081
根据获得的伪似然结果结合下式获得k时刻第i个模型的模型概率更新值
Figure FDA0003943165160000082
Figure FDA0003943165160000083
其中∝表示正比于；
S10、对步骤S9中的模型概率更新值
Figure FDA0003943165160000084
结合下式进行多传感器一致性融合处理，
Figure FDA0003943165160000085
S11、基于目标运动模型概率更新值
Figure FDA0003943165160000086
与状态估计协方差矩阵
Figure FDA0003943165160000087
Figure FDA0003943165160000088
其中
Figure FDA0003943165160000089
表示第L次融合以后的模型概率更新值；
得到目标在当前时刻的各传感器状态估计值与状态估计协方差矩阵，完成状态更新。
CN202211422756.6A 2022-11-15 2022-11-15 一种拖尾噪声背景下的机动目标分布式状态估计方法 Pending CN115687890A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211422756.6A CN115687890A (zh) 2022-11-15 2022-11-15 一种拖尾噪声背景下的机动目标分布式状态估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211422756.6A CN115687890A (zh) 2022-11-15 2022-11-15 一种拖尾噪声背景下的机动目标分布式状态估计方法

Publications (1)

Publication Number Publication Date
CN115687890A true CN115687890A (zh) 2023-02-03

Family

ID=85052930

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211422756.6A Pending CN115687890A (zh) 2022-11-15 2022-11-15 一种拖尾噪声背景下的机动目标分布式状态估计方法

Country Status (1)

Country Link
CN (1) CN115687890A (zh)

Similar Documents

Publication Publication Date Title
CN111985093B (zh) 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法
CN111178385B (zh) 一种鲁棒在线多传感器融合的目标跟踪方法
CN107765347B (zh) 一种高斯过程回归和粒子滤波的短期风速预测方法
Deisenroth et al. Robust filtering and smoothing with Gaussian processes
AU2009289008B2 (en) Estimating a state of at least one target
CN111127523B (zh) 基于量测迭代更新的多传感器gmphd自适应融合方法
CN110225454B (zh) 一种置信度传递的分布式容积卡尔曼滤波协作定位方法
CN111291471B (zh) 一种基于l1正则无迹变换的约束多模型滤波方法
CN114626307B (zh) 一种基于变分贝叶斯的分布式一致性目标状态估计方法
CN111798494A (zh) 广义相关熵准则下的机动目标鲁棒跟踪方法
CN109379702B (zh) 一种三维传感器网络节点定位方法及系统
CN109115228B (zh) 一种基于加权最小二乘容积卡尔曼滤波的目标定位方法
CN115828533A (zh) 一种基于Student’s t分布的交互多模型鲁棒滤波方法
CN114139109A (zh) 一种目标跟踪方法、系统、设备、介质及数据处理终端
CN111340853B (zh) 基于ospa迭代的多传感器gmphd自适应融合方法
CN115687890A (zh) 一种拖尾噪声背景下的机动目标分布式状态估计方法
CN109188424B (zh) 基于量测一致性的分布式多传感器多目标跟踪方法
CN116772867A (zh) 一种基于因子图的节点优选的多auv自适应协同定位方法及系统
CN111262556A (zh) 一种同时估计未知高斯测量噪声统计量的多目标跟踪方法
CN113391285B (zh) 一种量测随机延迟下带闪烁噪声的目标跟踪平滑方法
CN114061592B (zh) 基于多模型的自适应鲁棒auv导航方法
CN115905986A (zh) 一种基于联合策略的稳健卡尔曼滤波方法
CN111695617B (zh) 一种基于改进的协方差交叉算法的分布式火控融合方法
CN114415157A (zh) 一种基于水声传感器网络的水下目标多模型跟踪方法
CN110515069B (zh) 一种用于分布式目标跟踪的自适应一致性信息滤波方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination