CN107015945A

CN107015945A - 一种基于混合转移分布的高阶交互式多模型滤波方法

Info

Publication number: CN107015945A
Application number: CN201710231055.7A
Authority: CN
Inventors: 周共健; 叶晓平; 许荣庆; 吴立刚
Original assignee: Harbin Institute of Technology
Current assignee: Harbin Institute of Technology
Priority date: 2017-04-10
Filing date: 2017-04-10
Publication date: 2017-08-04
Anticipated expiration: 2037-04-10
Also published as: CN107015945B

Abstract

一种基于混合转移分布的高阶交互式多模型滤波方法，本发明涉及基于混合转移分布的高阶交互式多模型滤波方法。本发明为了解决现有方法中高阶马尔科夫链的设置参数多、设置过程繁琐以及精度较低的问题。本发明包括：一：采用混合转移分布模型得到n阶模型序列转移概率；二：对k时刻进行实时处理；三：对k＝1时的状态进行初始化；四：对k＝2时的状态进行初始化；五：对k进行判断，当k＝n时，则对k时刻的n阶模型序列概率，n阶模型序列状态估计值和相对应的协方差进行初始化；六：对3≤k≤n时的状态进行交互式多模型滤波算法；七：对k>n时的状态进行广义高阶交互式多模型滤波。本发明用于机动目标跟踪领域。

Description

一种基于混合转移分布的高阶交互式多模型滤波方法

技术领域

本发明涉及基于混合转移分布的高阶交互式多模型滤波方法。

背景技术

在目标跟踪的模型不确定问题中，常采用多模型滤波算法解决。其中在H.A.P.Blom,Y.Bar-Shalom.“The interacting multiple model algorithm for systemswith Markovian switching coefficients,”IEEE Transactions on AutomaticControl,vol.33(8),pp.780-783,1988中提出了经典的算法交互式多模型滤波算法(IMM)。该算法虽然能自适应识别当前时刻的模型，但精度却不是很高。

在P.Suchomski,“High-order interacting multiple-model estimation forhybrid systems with Markovian switching parameters,”International Journal ofSystems Science,vol.32(5),pp.669-679,2001中提出了广义高阶交互式多模型滤波方法(IMMn)，利用高阶模型序列更精确地估计目标状态。但这个方法的弊端是其中高阶马尔科夫链的设置比较复杂。其参数的个数随着阶数的增长而呈指数型增长，一来所需设置的参数太多，比较繁琐；二来同时把所有的参数设置得合理而不影响整体的状态估计是一件比较困难的事，因此会影响该算法运用到更高阶上。所以，需要一种更简单的方法来代替高阶马尔科夫链，并应用于高阶交互式多模型滤波方法。

毕欣；杜劲松；王伟；高洁；田星；赵越南；赵乾；丛日刚；仝盼盼；李想；张清石；徐洪庆；高扬，“一种基于自适应转移概率矩阵的交互多模型跟踪方法”，中国，2015-09-02以及韩红；李阳阳；陈兆平；王爽；智建纬；焦李成，“基于模糊推理的交互式多模型方法”中国，2009-07-08这两个专利中都对交互式多模型滤波方法进行了改进，但这些方法的阶数仍处于一阶，并没有利用充分利用高阶模型序列的先验信息，估计精度有待进一步改善。

发明内容

本发明的目的是为了解决现有广义高阶交互式多模型滤波方法中高阶马尔科夫链的设置参数多、设置过程繁琐以及精度较低的问题，而提出一种基于混合转移分布的高阶交互式多模型滤波方法。

一种基于混合转移分布的高阶交互式多模型滤波方法按以下步骤实现：

步骤一：采用混合转移分布模型得到n阶模型序列转移概率ρ(m_k|m_k-n,...,m_k-1)；

步骤二：估计的状态向量为和相对应的协方差为对k时刻进行实时处理；当k＝1时，转至步骤三；当k＝2时，转至步骤四；当3≤k≤n时，转至步骤六；当k>n时，转至步骤七；

步骤三：对k＝1时的状态进行初始化后，转步骤二等待处理下一k＝k+1时刻的雷达观测数据；

步骤四：对k＝2时的状态进行初始化后，转至步骤五；

步骤五：对k进行判断，当k＝n时，则对k时刻的n阶模型序列概率U_k(m_k-n+1,...,m_k)，n阶模型序列状态估计值和与相对应的协方差进行初始化后，转步骤二等待处理下一k＝k+1时刻的雷达观测数据；当k≠n时，直接执行步骤二处理k＝k+1时刻的雷达观测数据；

步骤六：对3≤k≤n时的状态进行交互式多模型滤波后，转至步骤五；

步骤七：对k>n时的状态进行广义高阶交互式多模型滤波后，转步骤二等待处理下一k＝k+1时刻的雷达观测数据。

本发明的有益效果为：

本发明利用接收到的雷达观测数据进行实时处理，从而实现对机动目标的有效跟踪。阶数越高，在模型不变区域的误差就越小，对目标状态有更高的估计精度。与IMM相比，精度高，估计性能更好，提高了10％左右；与IMMn相比，解决了高阶马尔科夫链设置困难的问题，减少了因参数设置不理想而降低滤波效果的可能性。

本发明提出了一种基于混合转移分布的高阶交互式多模型滤波方法。混合转移分布模型由于所需的参数个数少，设置简单的优点，在本发明中被用来代替高阶马尔科夫链来接近高阶模型序列转移概率矩阵，并应用于高阶交互式多模型滤波方法中从而使得高阶算法能够灵活使用。

附图说明

图1为2阶的本发明方法与交互式多模型算法位置均方根误差对比图；

图2为2阶的本发明方法与交互式多模型算法速度均方根误差对比图；

图3为不同阶数的本发明方法位置均方根误差对比图；

图4为不同阶数的本发明方法速度均方根误差对比图。

具体实施方式

具体实施方式一：一种基于混合转移分布的高阶交互式多模型滤波方法包括以下步骤：

步骤四：对k＝2时的状态进行初始化后，转至步骤五；

具体实施方式二：本实施方式与具体实施方式一不同的是：所述步骤一中采用混合转移分布模型得到n阶模型序列转移概率ρ(m_k|m_k-n,...,m_k-1)的具体过程为：

其中m_j为j时刻的模型，j＝k-n,…,k，设模型个数为r个，则m_j的取值范围为1到r；是从模型m_k-g转移到模型m_k的概率，是一阶马尔科夫链中的元素，λ_g是每个步长系数，满足以下条件：

其它步骤及参数与具体实施方式一相同。

具体实施方式三：本实施方式与具体实施方式一或二不同的是：所述模型具体为：

X_k+1＝F_k(m_k)X_k+G_k(m_k)u_k(m_k)+Γ_k(m_k)v_k(m_k)

其中X_k是由k时刻x轴位置x_k，x轴速度y轴位置y_k，y轴速度组成的状态向量。F_k(m_k)表示在k时刻模型m_k下的系统转移矩阵，G_k(m_k)是输入控制矩阵，u_k(m_k)是信号输入，Γ_k(m_k)是噪声系数矩阵，v_k(m_k)是k时刻模型m_k下的零均值白色高斯过程噪声，其协方差为Q_k(m_k)。

其中，T表示采样间隔。

(1)当模型为匀速运动模型时

(2)当模型协同转弯模型时

(3)当模型匀加速运动模型时

其中，a_x，a_y分别是x轴，y轴方向的加速度，x轴为水平方向，y轴为竖直方向，x轴与y轴垂直。

以上列举的三种是在机动目标跟踪中最常用的三种模型，但仅仅作为范例，在实际应用中可以根据需求使用其他的模型。

其它步骤及参数与具体实施方式一或二相同。

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是：所述步骤三中对k＝1时的状态进行初始化的具体过程为

其中z_k＝[x_k y_k]^T表示k时刻的接收到的雷达观测数据，z_k(j)表示z_k的第j个值。r_ij是观测噪声协方差R的第i行第j列元素，即

其它步骤及参数与具体实施方式一至三之一相同。

具体实施方式五：本实施方式与具体实施方式一至四之一不同的是：所述步骤四中对k＝2时的状态进行初始化的具体过程为

然后初始化k时刻模型概率各模型的状态估计相应协方差

其它步骤及参数与具体实施方式一至四之一相同。

具体实施方式六：本实施方式与具体实施方式一至五之一不同的是：所述步骤五中对k进行判断，当k＝n时，则对k时刻的n阶模型序列概率U_k(m_k-n+1,...,m_k)，n阶模型序列状态估计值和与相对应的协方差进行初始化的具体过程为：

其它步骤及参数与具体实施方式一至五之一相同。

具体实施方式七：本实施方式与具体实施方式一至六之一不同的是：所述步骤六中对3≤k≤n时的状态进行交互式多模型滤波的具体过程为：

计算混合概率

其中ρ(m_k|m_k-1)是一阶模型转移概率，可按实际需求预先设置0C₁(m_k)是归一化常数：

计算混合状态估计

计算混合状态估计相对应的协方差

把和作为输入，进行卡尔曼滤波，得到k时刻各模型的状态估计相应协方差和似然函数Λ_k(m_k)。

计算k时刻模型概率

计算k时刻的状态

计算对应的协方差

其它步骤及参数与具体实施方式一至六之一相同。

具体实施方式八：本实施方式与具体实施方式一至七之一不同的是：所述步骤七中对k>n时的状态进行广义高阶交互式多模型滤波的具体过程为：

计算混合概率

其中，C_n(m_k-n+1,...,m_k)是归一化常数；

计算n阶模型序列的混合状态估计

计算相对应的协方差

把和作为输入，进行卡尔曼滤波，得到k时刻各模型序列的状态估计相应协方差和似然函数Λ_k(m_k-n+1,...,m_k)。

计算k时刻n阶模型序列概率

计算k时刻的模型概率

计算k时刻的状态

计算对应的协方差

其它步骤及参数与具体实施方式一至七之一相同。

具体实施方式九：本实施方式与具体实施方式一至八之一不同的是：以上步骤中用到的卡尔曼滤波通过输入k-1时刻的状态估计和相应的协方差得到k时刻的状态相应的协方差和似然函数Λ_k的具体步骤为：

状态预测：

协方差预测：

计算观测预测值：

其中H_k是观测矩阵；

计算新息：

计算新息协方差：

R是观测噪声协方差，为H_k的转置；

计算似然函数：

Λ_k＝N(z_k；v_k|k-1,S_k|k-1)

其中N(z_k；v_k|k-1,S_k|k-1)表示z_k服从均值为v_k|k-1，协方差为S_k|k-1的高斯分布；

计算增益：

其中是S_k|k-1的逆；

计算状态：

计算协方差：

其它步骤及参数与具体实施方式一至八之一相同。

实施例一：

使用卡尔曼线性滤波模型,观测方程为：

z_k+1＝H_k+1X_k+1+w_k+1

w_k+1是零均值的白色高斯观测噪声，其协方差为R，且与过程噪声互不相关。H_k+1是k+1时刻的观测矩阵。

z_k+1＝[x_k+1,y_k+1]^T

在该仿真场景中，一共使用匀速直线运动和坐标转弯运动这两种模型。具体仿真场景是：机动目标先以X₁＝[1000,50,1000,50]^T的初始状态做匀速运动飞行40s，然后在40-80s的时候以w＝3rad/s转弯运动，最后40s继续做匀速运动。其中，R＝1Im²。采样间隔T＝1s。

基于混合转移分布的高阶交互式多模型滤波方法(阶数n为2)的步骤为：

步骤一：采用混合转移分布模型得到2阶模型序列转移概率ρ(l|i,j)。

ρ(l|i,j)＝λ₁p_jl+λ₂p_il,i,j,l∈{1,2}

其中，λ₁＝0.67,λ₂＝0.23，p_ij i,j∈{1,2}是一阶模型转移概率Π中的元素，

步骤二：估计的状态向量为和相对应的协方差为分以下四种情况对k时刻进行实时处理。

(1)当k＝1时，转步骤三；

(2)当k＝2时，转步骤四；

(3)当3≤k≤n时，转步骤六；

(4)当k>n时，转步骤七；

步骤三：对k＝1时的状态进行初始化

转步骤二等待处理下一k＝k+1时刻的雷达观测数据；

步骤四：对k＝2时的状态进行初始化

然后初始化k时刻模型概率各模型的状态估计相应协方差

然后转步骤五；

步骤五：对k进行判断，当k＝n时，则对U_k(m_k-1,m_k)，和与相对应的进行初始化。

转步骤二等待处理下一k＝k+1时刻的雷达观测数据；

步骤六：对3≤k≤n时的进行交互式多模型滤波。具体步骤如下所示：

计算混合概率：

其中C₁(m_k)是归一化常数；

计算混合状态估计：

计算相对应的协方差：

计算k时刻模型概率：

计算k时刻的状态：

计算对应的协方差：

转步骤五；

步骤七：对k>n时的进行广义高阶交互式多模型滤波。具体步骤如下所示：

计算混合概率：

其中，C₂(m_k-1,m_k)是归一化常数：

计算2阶模型序列的混合状态估计：

计算的协方差：

把和作为输入，进行卡尔曼滤波，得到k时刻各模型序列的状态估计相应协方差和似然函数Λ_k(m_k-1,m_k)。

计算k时刻2阶模型序列概率：

计算k时刻的模型概率：

计算k时刻的状态：

计算对应的协方差：

转步骤二等待处理下一k＝k+1时刻的雷达观测数据；

其中，卡尔曼滤波步骤是通过输入k-1时刻的状态估计和相应的协方差得到k时刻的状态相应的协方差和似然函数Λ_k的具体步骤为：

状态预测：

协方差预测：

计算观测预测值：

计算新息：

计算新息协方差：

计算似然函数：

Λ_k＝N(z_k；v_k|k-1,S_k|k-1)

计算增益：

其中是S_k|k-1的逆；

计算状态：

计算协方差：

500次蒙特卡洛仿真下的交互式多模型算法与基于混合转移分布的2阶交互式多模型滤波方法的均方根误差对比图，如图1和图2所示。

500次蒙特卡洛仿真下不同阶数的基于混合转移分布的高阶交互式多模型滤波方法的均方根误差对比图，如图3和图4所示。

通过混合转移分布模型和广义高阶交互式多模型滤波方法相结合的方法，使得高阶模型转移概率能比较容易设置，从而能灵活运用高阶滤波算法，而不受阶数的影响。

本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，本领域技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。

Claims

1.一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述基于混合转移分布的高阶交互式多模型滤波方法包括以下步骤：

步骤三：对k＝1时的状态进行初始化后，转步骤二处理k＝k+1时刻的雷达观测数据；

步骤四：对k＝2时的状态进行初始化后，转至步骤五；

步骤五：对k进行判断，当k＝n时，则对k时刻的n阶模型序列概率U_k(m_k-n+1,...,m_k)，n阶模型序列状态估计值和与相对应的协方差进行初始化后，转步骤二处理k＝k+1时刻的雷达观测数据；当k≠n时，直接执行步骤二处理k＝k+1时刻的雷达观测数据；

步骤六：对3≤k≤n时的状态进行交互式多模型滤波算法后，转至步骤五；

步骤七：对k>n时的状态进行广义高阶交互式多模型滤波后，转步骤二处理k＝k+1时刻的雷达观测数据。

2.根据权利要求1所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤一中采用混合转移分布模型得到n阶模型序列转移概率ρ(m_k|m_k-n,...,m_k-1)的具体过程为：

3.根据权利要求2所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述模型具体为：

X_k+1＝F_k(m_k)X_k+G_k(m_k)u_k(m_k)+Γ_k(m_k)v_k(m_k)

其中X_k是由k时刻x轴位置x_k，x轴速度y轴位置y_k，y轴速度组成的状态向量；F_k(m_k)表示在k时刻模型m_k下的系统转移矩阵，G_k(m_k)是输入控制矩阵，u_k(m_k)是信号输入，Γ_k(m_k)是噪声系数矩阵，v_k(m_k)是k时刻模型m_k下的零均值白色高斯过程噪声，其协方差为Q_k(m_k)；

其中，T表示采样间隔；

(1)当模型为匀速运动模型时：

(2)当模型协同转弯模型时：

(3)当模型匀加速运动模型时：

其中，a_x，a_y分别是x轴，y轴方向的加速度。

4.根据权利要求3所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤三中对k＝1时的状态进行初始化的具体过程为：

其中z_k＝[x_k y_k]^T表示k时刻接收到的雷达观测数据，z_k(j)表示z_k的第j个值；r_ij是观测噪声协方差R的第i行第j列元素，即

5.根据权利要求4所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤四中对k＝2时的状态进行初始化的具体过程为：

初始化k时刻模型概率各模型的状态估计相应协方差

6.根据权利要求5所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤五中对k进行判断，当k＝n时，则对k时刻的状态n阶模型序列概率U_k(m_k-n+1,...,m_k)，n阶模型序列状态估计值和与相对应的协方差进行初始化的具体过程为：

7.根据权利要求6所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤六中对3≤k≤n时的状态进行交互式多模型滤波算法的具体过程为：

计算混合概率：

其中ρ(m_k|m_k-1)是一阶模型转移概率，C₁(m_k)是归一化常数：

计算混合状态估计：

计算混合状态估计相对应的协方差：

把和作为输入，进行卡尔曼滤波，得到k时刻各模型的状态估计相应协方差和似然函数Λ_k(m_k)；

计算k时刻模型概率：

计算k时刻的状态

计算对应的协方差

8.根据权利要求7所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤七中对k>n时的进行广义高阶交互式多模型滤波的具体过程为：

计算混合概率：

其中，C_n(m_k-n+1,...,m_k)是归一化常数；

计算n阶模型序列的混合状态估计：

计算相对应的协方差：

把和作为输入，进行卡尔曼滤波，得到k时刻各模型序列的状态估计相应协方差和似然函数Λ_k(m_k-n+1,...,m_k)；

计算k时刻n阶模型序列概率：

计算k时刻的模型概率

计算k时刻的状态：

计算对应的协方差：

9.根据权利要求8所述的一种基于混合转移分布的高阶交互式多模型滤波方法，其特征在于：所述步骤六和步骤七中进行卡尔曼滤波的具体步骤为：

状态预测：

其中是卡尔曼滤波的输入；

协方差预测：

其中是卡尔曼滤波的输入；

计算观测预测值：

其中H_k是观测矩阵；

计算新息：

计算新息协方差：

R是观测噪声协方差，为H_k的转置；

计算似然函数：

其中表示z_k服从均值为协方差为S_k|k-1的高斯分布；

计算增益：

计算状态：

计算协方差：