一种基于L1正则无迹变换的约束多模型滤波方法

摘要

本发明公开一种基于L1正则无迹变换的约束多模型滤波方法，系统测量截断提供了第一个约束，和系统的最新测量序列一起，融入到系统先验的更新过程中；系统观测间隔之间的平滑提供系统的第二个约束，该时空约束的校正和补偿由模糊逻辑算法加权；可行区域由启发式优化方法近似，在其中，通过约束辅助粒子滤波框架对低维流形模型进行建模，并通过基于L1正则化的对模型状态进行无迹卡尔曼平滑滤波。

权利要求

1.一种基于L1正则无迹变换的约束多模型滤波方法，其特征在于，包括如下步骤：步骤1，将多模型随机动力学系统建模为马尔可夫状态空间模型，具体包括系统状态与观测模型方程；所述马尔可夫状态空间模型为由匀速运动的二阶动力学模型以及匀转弯运动的恒定转弯运动模型组成的模型集合，其中，二阶动力学模型的过程噪声为加速度白噪声；

步骤2，通过截断先验约束对系统测量噪声的概率密度函数构建数学模型，以定义系统状态的可行域；

步骤3，对系统状态模型切换之间的平滑反馈构建数学模型，通过补偿约束对状态模型及模型状态向量的动态转移概率进行修正；

步骤4，将系统约束融入重要性密度函数，构建约束辅助粒子滤波，据此选择状态模型粒子，并计算其权重；

步骤5，根据最小二乘准则构建系统状态优化估计的目标函数，通过外围和内在惩罚项实现截断先验约束和平滑补偿约束，以得到系统可信任的观测序列；

步骤6，根据最大相关熵准则，选取模型状态的L1范数作为障碍项，通过数值近似得到可行域内模型状态的最大似然解，即为可行域中心；

步骤7，基于系统截断先验约束和平滑补偿约束，在无迹卡尔曼滤波框架下，对系统可行域的观测序列与模型状态进行预测与更新；

步骤8，根据约束粒子滤波框架下得到的模型权重，对各模型的估计值进行加权融合，计算联合状态估计与协方差矩阵，输出滤波结果。

2.根据权利要求1所述基于L1正则无迹变换的约束多模型滤波方法，其特征在于，步骤

1中，所述将多模型随机动力学系统建模为马尔可夫状态空间模型，具体为：对于机动目标模型，系统状态所处的模式能及时进行切换，将多模型随机动力学系统建模为马尔可夫状态空间模型：

式中， 和 分别代表时间 的 维系统状态向量和 维测量序列； 和 分别是系统状态转移矩阵和测量矩阵； 和 分别是过程噪声和测量噪声； 表示在时间的动力学模型，为 个可能的模式之一，即：历史状态模型的数量随时间呈指数增长，即存在随时间 变化的模型序列 ，系统的结构和/或噪声的统计信息可能因模型而异，假设模型切换过程服从一阶马尔可夫过程，则模型转换概率已知并且建模为：式中，pij表示模型mi转换成模型mj的概率。

3.根据权利要求2所述基于L1正则无迹变换的约束多模型滤波方法，其特征在于，步骤

2中，通过截断先验约束对系统测量噪声的概率密度函数构建数学模型，以定义系统状态的可行域，具体为：

在物理上系统观察噪声具有有限且衰减的能量支持，在数学上对系统测量噪声的概率密度函数在截断的意义上建模，即：式中， 表示系统观测噪声的概率密度函数，ek表示系统观测统计噪声， 表示  维的有界联通区域；

定义模型状态的可行域为：

式中， 表示模型 下状态的标量分量；

给定具有有限噪声支持的测量序列，通过以下递归版本推论模型状态下的后验分布：式中， 表示1到第K时刻的观测， 表示1到第k时刻的状态模型， 表示1到第k时刻的观测周期， 表示概率密度函数， 表示可行域的似然函数，t表示连续观测时间， 是标量分量 的非线性函数， 是归一化常数，是指示函数，定义为：

从上式可知，将截断先验与最新的测量信息一起纳入滤波过程，进而对由测量环境引起的不确定性进行建模和优化控制。

4.根据权利要求3所述基于L1正则无迹变换的约束多模型滤波方法，其特征在于，步骤

3中，对系统状态模型切换之间的平滑反馈构建数学模型，通过补偿约束对状态模型及模型状态向量的动态转移概率进行修正，具体为：对于闭环非线性动力系统，历史更新状态与当前估计状态之间的反馈所产生的平滑补偿是必不可少的环节，代替单次测试统计，衰减记忆的平均值在滑动窗口上表示，滑动窗口定义为：

式中， 和 分别是 和 时刻的采样时长；

根据不变集定理，对于闭环控制系统，在每个采样间隔内，给定适当的初始变量，通过反馈优化来保证预测控制器的均方稳定性，令 其代表目标检测时间，对模型状态向量的动态转移概率修正为如下形式：式中， 表示状态更新概率密度函数， 表示状态模型更新概率密度函数，表示状态模型检测概率密度， 表示状态模型似然更新概率密度函数，是模糊加权因子，表示分别根据过去和当前的观察结果进行的补偿和校正的效果。

说明书

一种基于L1正则无迹变换的约束多模型滤波方法

技术领域

[0001] 本发明涉及非线性滤波技术领域，具体是一种基于L1正则无迹变换的约束多模型滤波方法。

背景技术

[0002] 在目标跟踪的模型不确定问题中，H.A.P.Blom，Y.Bar‑Shalom.“TheinteractingmultiplemodelalgorithmforsystemswithMarkovianswitchingcoefficients，”IEEETran
sactionsonAutomaticControl，vol.33(8)，pp.780‑783，1988提出了经典的交互式多模型
滤波方法。该经典方法使用模型转移概率来自动识别当前使用的模型，进行模型切换，从而
实现在多模型下的自适应滤波估计。但缺点是估计精度不是很高，且依赖于模型转移概率
的设置。当模型不变转移概率的值设置很大时，虽然在模型不变区域误差会减少，但同时却
会增大模型切换区域的误差，并且从模型切换区域到模型不变区域的误差收敛速度会减
慢。

[0003] 例如：

[0004] 专利号CN 110375731 A公开了一种混合交互式多模型滤波方法，其采用多模型参数覆盖噪声协方差阵的方法，解决了单一观测噪声协方差阵在复杂环境下滤波精度下降的
问题，在交互式多模型滤波器主模型中采用基于残差×2检测的鲁棒滤波来抑制了野值干
扰，同时采用残差修正补偿的方法来保证系统在野值点导航不间断。但是其并未考虑截断
先验约束，导致最终的估计精度不是很高；

[0005] 专利号CN 107783944 A公开了一种多模型自校准无迹卡尔曼滤波方法，其将多模型估计理论引入到自校准无迹卡尔曼滤波方法中，推导得到了多模型自校准无迹卡尔曼滤
波方法的完整过程，并自动识别未知输入为零段与非零段，分别发挥无迹卡尔曼滤波方法
与校准无迹卡尔曼滤波方法各自的优势。虽然提高了非线性系统受未知输入影响时的滤波
精度但是在从模型切换区域到模型不变区域的误差收敛速度会相对变慢。

发明内容

[0006] 针对上述现有技术中的不足，本发明提供一种基于L1正则无迹变换的约束多模型滤波方法。

[0007] 为实现上述目的，本发明提供一种基于L1正则无迹变换的约束多模型滤波方法，包括如下步骤：

[0008] 步骤1，将多模型随机动力学系统建模为马尔可夫状态空间模型，具体包括系统状态与观测模型方程；

[0009] 步骤2，通过截断先验约束对系统测量噪声的概率密度函数构建数学模型，以定义系统状态的可行域；

[0010] 步骤3，对系统状态模型切换之间的平滑反馈构建数学模型，通过该补偿约束对状态模型及模型状态向量的动态转移概率进行修正；

[0011] 步骤4，将系统约束融入重要性密度函数，构建约束辅助粒子滤波，据此选择状态模型粒子，并计算其权重；

[0012] 步骤5，根据最小二乘准则构建系统状态优化估计的目标函数，通过外围和内在惩罚项实现截断先验约束和平滑补偿约束，以得到系统可信任的观测序列；

[0013] 步骤6，根据最大相关熵准则，选取模型状态的L1范数作为障碍项，通过数值近似得到可行域内模型状态的最大似然解，即为可行域中心；

[0014] 步骤7，基于系统截断先验约束和平滑补偿约束，在无迹卡尔曼滤波框架下，对系统可行域的观测序列与模型状态进行预测与更新；

[0015] 步骤8，根据约束粒子滤波框架下得到的模型权重，对各模型的估计值进行加权融合，计算联合状态估计与协方差矩阵，输出滤波结果。

[0016] 进一步优选的，步骤1中，所述将多模型随机动力学系统建模为马尔可夫状态空间模型，具体为：

[0017] 对于机动目标模型，系统状态所处的模式可以及时进行切换，将多模型随机动力学系统建模为马尔可夫状态空间模型：

[0018] Xk＝Fk(Mk)Xk‑1+Vk‑1(Mk)

[0019] Zk＝Hk(Mk)Xk+ek(Mk)

[0020] 式中， 和 分别代表时间k的dX维系统状态向量和dZ维测量序列；Fk和Hk分别是系统状态转移矩阵和测量矩阵；Vk和ek分别是过程噪声和测量噪声；Mk表示在
时间的动力学模型，为K个可能的模式之一，即：

[0021] Mi∈{M1，…，MK}，i＝1，…，K

[0022] 历史状态模型的数量随时间呈指数增长，即存在随时间k变化的模型序列Kk，系统的结构和/或噪声的统计信息可能因模型而异，假设模型切换过程服从一阶马尔可夫过程，
则模型转换概率可以已知并且建模为：

[0023]

[0024] 式中，pij表示模型mi转换成模型mj的概率

[0025] 进一步优选的，步骤2中，通过截断先验约束对系统测量噪声的概率密度函数构建数学模型，以定义系统状态的可行域，具体为：

[0026] 在物理上系统观察噪声具有有限且衰减的能量支持，在数学上对系统测量噪声的概率密度函数在截断的意义上建模，即：

[0027]

[0028] 式中， 表示系统观测噪声的概率密度函数，ek表示系统观测统计噪声， 表示dz维的有界联通区域；

[0029] 定义模型条件状态的可行域为：

[0030]

[0031] 式中，am,k表示模型m下状态的标量分量；

[0032] 给定具有有限噪声支持的测量序列，可以通过以下递归版本推论模型条件状态下的后验分布：

[0033]

[0034] 式中，Z1∶k表示1到第K时刻的观测，M1∶k表示1到第k时刻的状态模型，T1∶k表示1到第k时刻的观测周期，p(·)表示概率密度函数， 表示可行域的似然函数，t表示连
续观测时间，h(Xm,k)是标量分量Xm,k的非线性函数,εk是归一化常数,pg(·)是指示函数，定
义为：

[0035]

[0036] 从上式可知，将截断先验与最新的测量信息一起纳入滤波过程，进而对由测量环境引起的不确定性进行建模和优化控制。

[0037] 进一步优选的，步骤3中，对系统状态模型切换之间的平滑反馈构建数学模型，通过该补偿约束对状态模型及模型状态向量的动态转移概率进行修正，具体为：

[0038] 对于闭环非线性动力系统，历史更新状态与当前估计状态之间的反馈所产生的平滑补偿是必不可少的环节。代替单次测试统计，衰减记忆的平均值可以在滑动窗口上表示，
滑动窗口可以定义为：

[0039] Tk＝t(k)‑t(k‑1)

[0040] 式中，t(k‑1)和t(k)分别是k和k‑1时刻的采样时长；

[0041] 根据不变集定理，对于闭环控制系统，在每个采样间隔内，给定适当的初始变量，可以通过反馈优化来保证预测控制器的均方稳定性，令 其代表目标检测时间，对
模型条件状态向量的动态转移概率可以修正为如下形式：

[0042]

[0043] 式中，Pkk‑1表示状态更新概率密度函数， 表示状态模型更新概率密度函数，表示状态模型检测概率密度， 表示状态模型似然更新概率密度函数， 是
模糊加权因子，表示分别根据过去和当前的观察结果进行的补偿和校正的效果。

[0044] 进一步优选的，步骤5中，根据最小二乘准则构建系统状态优化估计的目标函数，通过外围和内在惩罚项实现截断先验约束和平滑补偿约束，以得到系统可信任的观测序
列，具体为：

[0045] 对系统测量序列的确定性研究旨在提供状态参数边界，并且可以使用正则化技术来实现约束优化，对应于截断先验和平滑补偿的外围和内在惩罚项构成了系统状态优化估
计的目标函数，给定确定性输入后，根据数学对偶性，可以将动态系统的最优估计问题重写
为如下所示的控制问题，即系统状态优化估计的目标函数为：

[0046]

[0047] 式中，Js是被动传感器数目， 是状态Xk的先验估计， 是误差协方差；其中，过程噪声Qk须满足可逆，以确保最小解存在；第一个外罚项 由外围正则化因子λA控制，反
映了模型逼近的复杂程度；第二个内罚项 由内在结构的正则化因子λI控制，反映了系统
输入空间的几何结构逼近。

[0048] 进一步优选的，步骤6中，根据最大相关熵准则，选取模型状态的L1范数作为障碍项，通过数值近似得到可行域内模型状态的最大似然解，即为可行域中心，具体为：

[0049] 基于系统的过程和测量模型方程，模型条件状态的测量似然函数可以近似为：

[0050]

[0051] 式中，Rk(Mk)和Qk(Mk)分别表示测量噪声和过程噪声的协方差标准差；

[0052] 基于模型条件状态的可行域，模型状态的最大似然解可以近似为：

[0053]

[0054] 进一步优选的，步骤7中，基于系统截断先验约束和平滑补偿约束，在无迹卡尔曼滤波框架下，对系统可行域的观测序列与模型状态进行预测与更新，具体为：

[0055] 步骤7.1，可行域sigma点

[0056] 在k时刻，对于模型m，令 和 分别表示模型状态的第i个sigma点及其权重。修正均值和协方差选择集合：

[0057]

[0058]

[0059]

[0060] 式中 ， 表示模型状态的维 数，λ是可信任区域的半径 ，并 且是均值Cholesky解的第i行，权重因子之和为1，即

[0061] 步骤7.2，约束预测

[0062] k时刻，约束预测均值为：

[0063]

[0064] k时刻，约束预测方差为：

[0065]

[0066] 步骤7.3，模型状态更新

[0067] k时刻，定域观测序列的预测均值 为：

[0068]

[0069] 观测协方差为：

[0070]

[0071] 状态与观测的互协方差为：

[0072]

[0073] 步骤7.4，模型状态估计

[0074] 在上述约束条件下，卡尔曼滤波增益为

[0075]

[0076] 模型状态更新为

[0077]

[0078]

[0079] 进一步优选的，步骤8中，根据约束粒子滤波框架下得到的模型权重，对各模型的估计值进行加权融合，计算联合状态估计与协方差矩阵，输出滤波结果，具体为：

[0080] 将当前测量值的校正和历史更新状态的补偿增强到状态评估中，融合效果分别表示为τ＝0,1，并且可以通过使用模糊高斯函数的隶属度来计算得到权重系数

[0081]

[0082]

[0083] 式中， 是模型状态的新息协方差，从估计精度衡量的意义上讲，当前时刻观测信息更为准确时，αk更接近于1；

[0084] 对于状态向量的最终估计，滤波后的均值和协方差是由模型条件状态估计的加权和组成的融合输出，其计算如下：

[0085]

[0086]

[0087] 式中， 是k时刻模型m的重要性权重。

[0088] 本发明提供的一种基于L1正则无迹变换的约束多模型滤波方法，系统测量截断提供了第一个约束，和系统的最新测量序列一起，融入到系统先验的更新过程中；系统观测间
隔之间的平滑提供系统的第二个约束，该时空约束的校正和补偿由模糊逻辑算法加权；可
行区域由启发式优化方法近似，在其中，通过约束辅助粒子滤波框架对低维流形模型进行
建模，并通过基于L1正则化的对模型状态进行无迹卡尔曼平滑滤波。

附图说明

[0089] 为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本
发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以
根据这些附图示出的结构获得其他的附图。

[0090] 图1为本发明实施例中应用于非线性高斯模型的平滑约束扩展卡尔曼滤波方法的流程示意图；

[0091] 图2为本发明实施例中仿真过程中case1所得RMS位置误差曲线随时间变化的趋势图。

[0092] 图3为本发明实施例中仿真过程中case2所得RMS位置误差曲线随时间变化的趋势图；

[0093] 图4为本发明实施例中仿真过程中case3所得RMS位置误差曲线随时间变化的趋势图。

[0094] 图5为本发明实施例中仿真过程中case4所得RMS位置误差曲线随时间变化的趋势图；

[0095] 图6为本发明实施例中仿真过程中case5所得RMS位置误差曲线随时间变化的趋势图；

[0096] 图7为本发明实施例中仿真过程中RMSE的统计信息图；

[0097] 图8为本发明实施例中仿真过程中RTAMS的统计信息图；

[0098] 图9为本发明实施例中仿真过程中五个滤波器的CV模型切换概率曲线图；

[0099] 图10本发明实施例中仿真过程中2D平面中提出方法的滤波轨迹图；

[0100] 图11本发明实施例中仿真过程中3D空间中提出方法的滤波轨迹图；

[0101] 图12本发明实施例中仿真过程中位置RMSE随时间变化的趋势图；

[0102] 图13本发明实施例中仿真过程中MMUKF的模型切换概率图。

[0103] 本发明目的的实现、功能特点及优点将结合实施例，参照附图做进一步说明。

具体实施方式

[0104] 下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明的一部分实施例，而不是全部的实施例。基
于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其
他实施例，都属于本发明保护的范围。

[0105] 需要说明，本发明实施例中所有方向性指示(诸如上、下、左、右、前、后……)仅用于解释在某一特定姿态(如附图所示)下各部件之间的相对位置关系、运动情况等，如果该
特定姿态发生改变时，则该方向性指示也相应地随之改变。

[0106] 另外，在本发明中如涉及“第一”、“第二”等的描述仅用于描述目的，而不能理解为指示或暗示其相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、
“第二”的特征可以明示或者隐含地包括至少一个该特征。在本发明的描述中，“多个”的含
义是至少两个，例如两个，三个等，除非另有明确具体的限定。

[0107] 在本发明中，除非另有明确的规定和限定，术语“连接”、“固定”等应做广义理解，例如，“固定”可以是固定连接，也可以是可拆卸连接，或成一体；可以是机械连接，也可以是
电连接，还可以是物理连接或无线通信连接；可以是直接相连，也可以通过中间媒介间接相
连，可以是两个元件内部的连通或两个元件的相互作用关系，除非另有明确的限定。对于本
领域的普通技术人员而言，可以根据具体情况理解上述术语在本发明中的具体含义。

[0108] 另外，本发明各个实施例之间的技术方案可以相互结合，但是必须是以本领域普通技术人员能够实现为基础，当技术方案的结合出现相互矛盾或无法实现时应当认为这种
技术方案的结合不存在，也不在本发明要求的保护范围之内。

[0109] 本实施例公开了一种基于L1正则无迹变换的约束多模型滤波方法，参考图1，包括如下步骤：

[0110] 一种基于L1正则无迹变换的约束多模型滤波方法，包括如下步骤：

[0111] 步骤1，将多模型随机动力学系统建模为马尔可夫状态空间模型，具体包括系统状态与观测模型方程；

[0112] 步骤2，通过截断先验约束对系统测量噪声的概率密度函数构建数学模型，以定义系统状态的可行域；

[0113] 步骤3，对系统状态模型切换之间的平滑反馈构建数学模型，通过该补偿约束对状态模型及模型状态向量的动态转移概率进行修正；

[0114] 步骤4，将系统约束融入重要性密度函数，构建约束辅助粒子滤波，据此选择状态模型粒子，并计算其权重；

[0115] 步骤5，根据最小二乘准则构建系统状态优化估计的目标函数，通过外围和内在惩罚项实现截断先验约束和平滑补偿约束，以得到系统可信任的观测序列；

[0116] 步骤6，根据最大相关熵准则，选取模型状态的L1范数作为障碍项，通过数值近似得到可行域内模型状态的最大似然解，即为可行域中心；

[0117] 步骤7，基于系统截断先验约束和平滑补偿约束，在无迹卡尔曼滤波框架下，对系统可行域的观测序列与模型状态进行预测与更新；

[0118] 步骤8，根据约束粒子滤波框架下得到的模型权重，对各模型的估计值进行加权融合，计算联合状态估计与协方差矩阵，输出滤波结果。

[0119] 下面对上述各步骤的实现方式进行展开说明。

[0120] 对于机动目标模型，系统状态所处的模式可以及时进行切换，将多模型随机动力学系统建模为马尔可夫状态空间模型：

[0121] Xk＝Fk(Mk)Xk‑1+Vk‑1(Mk)，   (1)

[0122] Zk＝Hk(Mk)Xk+ek(Mk)   (2)

[0123] 式中， 和 分别代表时间k的dX维系统状态向量和dZ维测量序列；Fk和Hk分别是系统状态转移矩阵和测量矩阵；Vk和ek分别是过程噪声和测量噪声；Mk表示在时
间的动力学模型，为K个可能的模式之一，即：

[0124] Mi∈{M1，…，MK}，i＝1，…，K   (3)

[0125] 历史状态模型的数量随时间呈指数增长，即存在随时间k变化的模型序列Kk，系统的结构和/或噪声的统计信息可能因模型而异，假设模型切换过程服从一阶马尔可夫过程，
则模型转换概率可以已知并且建模为：

[0126]

[0127] 式中，pij表示模型mi转换成模型mj的概率。

[0128] 通常可以假定水平和垂直运动模型是解耦的,在离散时间变量的3维空间中，为了与静态模型区分开，给模型状态向量增加了一个额外的分量，即转弯速率Ωk，表示为

[0129]T

[0130] 其中,(Xm,k,Ym,k,Zm,k) 和 分别是模型状态的位置和速度分量,上标T表示转置。

[0131] 对于匀速运动，匀速运动的近似恒定速度模型可以建模为具有低方差过程噪声的二阶动力学模型，过程噪声为加速度白噪声，恒速运动的状态转移矩阵为以下模型：

[0132]

[0133] 其中，T是观测间隔，过程噪声矢量vk是一个3×6的独立同分布矩阵，即,其中σe是过程噪声标准差，I是6×6的单位矩阵。过程噪声的平均值可以将操作建模为确定
性输入。

[0134] 对于匀转弯运动，恒定转弯(CT)运动模型(包括顺时针和逆时针)可以分别建模为以下转换坐标系下的状态转移矩阵，由以下式子给出：

[0135]

[0136] 其中Ωk＞0和Ωk＜0分别表示水平方向上的顺时针和逆时针运动。

[0137] 由于外部干扰，在跳跃马尔可夫系统处于变化状态的停留时间内，对于多模型预测控制器，扰动和奇异性是不可避免的且不可预测的。不确定机动的检测是一个具有挑战
性的问题，对于参数估计，可以使用贝叶斯框架，从先验约束开始。

[0138] 通过截断先验约束对系统测量噪声的概率密度函数建模得到模型状态的可行域，具体为：

[0139] 在贝叶斯框架下，状态后验分布的数值近似仅在信任区域内才是准确的，在物理上系统观察噪声具有有限且衰减的能量支持，在数学上对系统测量噪声的概率密度函数在
截断的意义上建模，即：

[0140]

[0141] 式中， 表示系统观测噪声的概率密度函数，ek表示系统观测统计噪声， 表示dz维的有界联通区域；

[0142] 相应地，定义模型状态的可行域为：

[0143]

[0144] 式中，am,k表示模型m下状态的标量分量；

[0145] 给定具有有限噪声支持的测量序列，可以通过以下递归版本推论模型状态下的后验分布：

[0146]

[0147] 式中，Z1∶k表示1到第K时刻的观测，M1∶k表示1到第k时刻的状态模型，T1∶k表示1到第k时刻的观测周期，p(·)表示概率密度函数， 表示可行域的似然函数，t表示连续
观测时间，h(Xm,k)是标量分量Xm,k的非线性函数,εk是归一化常数,pg(·)是指示函数，定义
为：

[0148]

[0149] 从公式(10)可知，将截断先验与最新的测量信息一起纳入滤波过程，进而对由测量环境引起的不确定性进行建模和优化控制。

[0150] 通过平滑补偿约束对模型状态向量的动态转移概率进行修正，具体为：

[0151] 对于闭环非线性动力系统，历史更新状态与当前估计状态之间的反馈所产生的平滑补偿是必不可少的环节。代替单次测试统计，衰减记忆的平均值可以在滑动窗口上表示，
滑动窗口可以定义为：

[0152] Tk＝t(k)‑t(k‑1)   (12)

[0153] 式中，t(k‑1)和t(k)分别是k和k‑1时刻的采样时长；

[0154] 根据不变集定理，对于闭环控制系统，在每个采样间隔内，给定适当的初始变量，可以通过反馈优化来保证预测控制器的均方稳定性，令 其代表目标检测时间，对
模型状态向量的动态转移概率可以修正为如下形式：

[0155]

[0156] 式中，Pk|k‑1表示状态更新概率密度函数， 表示状态模型更新概率密度函数，表示状态模型检测概率密度， 表示状态模型似然更新概率密度函数， 是
模糊加权因子，表示分别根据过去和当前的观察结果进行的补偿和校正的效果。

[0157] 到目前为止，在本实施例中方法的框架中同时考虑了截断和平滑约束，在此过程中，可以通过一系列的启发式优化来实现低方差和滤波性能的稳定性。

[0158] 对于动态跳跃多模型系统，基于正则无迹变换的多模型滤波方法旨在提高模型切换概率的准确性，增强非线性滤波的稳定性，有效利用系统统计约束成为关键问题

[0159] 鉴于方程(9)中定义的可行区域，结合测量方程(2)，模型条件下状态的测量似然函数可以递归地近似如下：

[0160]

[0161] 其中，Sk,m是k时刻模型m的新息。

[0162] 在跳变模式系统中，该系统经历了从一种模式到另一种模式的转变，并且机动表现为大的新息，针对此类事件的有效而简单的检测过程基于归一化的平方新息。边缘化处
理是在KF系列滤波器中处理大尺寸测量向量的另一种更为常用的方法，并且测量似然函数
的计算可以分解为

[0163]

[0164] 通过这种方式，状态估计可以分而治之。这种方法的基本原理在于Rao‑Blackwell定理的统计性质，这导致估计量的方差比纯蒙特卡洛抽样小。

[0165] 边缘化过程满足广义正则化理论的基本前提，联合分布函数与边缘分布函数之间存在等价关系。并且可以将主要推理问题简化为对后验分布进行以下类型的期望的计算：

[0166]

[0167] 以这种方式，可以将模型切换沿着低维流形结构进行建模。同时，这种简单的形式使得可以利用系统统计约束选择状态模型，并且可以调用约束辅助粒子滤波方法。

[0168] 概率密度估计是模型切换和信息发现中的基本概念。通常，密度估计方法的执行不受系统测量序列的固有几何结构影响。对于多模型非线性滤波问题，简单的高斯近似不
再合适。在粒子过滤器框架中，每个单目标密度可以由一组带系数的粒子 表示，其中
Ns是采样模型粒子的总数。蒙特卡洛近似的性质保证了约束辅助滤波器(constrained 
auxiliary particle fitlering，CAPF)的可执行性。关键步骤是将截断优先级和平滑度补
偿合并到重要性密度函数的递归推断中。然后，模型粒子的预测密度和更新密度与以下类
型的修正似然函数成比例：

[0169]

[0170] 因此，在空时约束条件下，与标准的顺序重要性重采样(SIR)方法相比，该方法能够有效克服由低噪声引起的粒子衰减和退化问题。修正的次优重要性密度函数应产生适应
性，因此可以消除第二阶段的重采样。

[0171] 相应地，可以通过下式对采样的模型粒子的重要性权重进行更新：

[0172]

[0173] 通过下式进行归一化，

[0174]

[0175] 在有限维可行区域上，可以通过以下公式计算采样模型粒子的概率：

[0176]

[0177] 在上面的CAPF框架中，对于每个单个目标密度的外在本征正则化可以转换为有限维空间上一组加权模型粒子的优化问题。在这一过程中，粒子数量与状态滤波精度可得到
平衡。

[0178] 与标准SIR方法不同，本实施例方法中，状态模型集是在约束辅助粒子滤波框架中估计的。连同最新的测量信息，截断先验和平滑度补偿都被纳入重要性密度函数的构造中。
这样，即使在滤波器存在不稳定性和低方差的情况下，也可以保证模型粒子的多样性和准
确性。如果模式集包含正确的一个且不发生模式跳变，则后验分布的真实概率将收敛到到
统一模式。在上述推导过程中，顺序采样模式粒子被限制在可行区域内，状态模型的概率等
于相应选定模型粒子的归一化重要性权重。因此，所提出的CAPF方法收敛。

[0179] 与闭环控制系统类似，状态估计的一个重要方面是应通过反馈补偿来跟踪系统的有效测量，可以通过正则化方法来确定域测量。

[0180] 对系统测量序列的确定性研究旨在提供状态参数边界，并且可以使用正则化技术来实现约束优化。对应于截断先验和平滑补偿的外围和内在惩罚项构成了系统状态优化估
计的目标函数。给定确定性输入后，根据数学对偶性，可以将动态系统的最优估计问题重写
为如下所示的控制问题：

[0181]

[0182] 其中，Js是被动传感器数目， 是状态Xk的先验估计， 是误差协方差。请注意过程噪声Qk须满足可逆，以确保最小解存在。第一个外罚项 由外围正则化因子λA控制，
反映了模型逼近的复杂程度。第二个内罚项 由内在结构的正则化因子λI控制，反映了系
统输入空间的几何结构逼近。

[0183] 对于某些特定区域，可以通过标量化处理得到凸的目标函数。模型状态变量可分解为目标的位置Xm,k和速度分量 实际上，系统非线性测量函
数可以合理地局部线性化。因此，随时间的变化，历史测量序列可以表示为

[0184]

[0185] 其中 是系统观测噪声。

[0186] 考虑到系统模型公式(1)和(2)，Rk(Mk)和Qk(Mk)分别表示测量噪声和过程噪声的协方差标准差。模型状态的测量似然函数可以近似为

[0187]

[0188] 考虑式(9)中定义的软约束条件，模型状态的最大似然解可以近似为

[0189]

[0190] 对于大规模非线性规划问题，由于传统最小二乘准则对噪声和野值敏感，因此本发明采用基于多变量校准的最大相关熵准则。然后，要求内核数量少于数据点的要求，因
此，降低状态空间复杂度成为确定平滑正则化器的关键因素。为了实现软约束，受接受‑拒
绝方法的思想启发，选择该状态的L1范数的对数形式作为障碍项。从而，增广目标函数可以
转换为以下的L1正则最小二乘的非负对数似然形式

[0191]

[0192] 其中，σj,θ和σj,β是来自第被j个被动传感器的观测序列的高斯核宽。

[0193] 除了截断和平滑约束之外，测量序列也融入到了模型状态的更新中。利用变量的稀疏因子正则化，将高维状态空间的动态规划问题转化为稀疏低维流形结构的最优解问
题，同时当状态向量远离边界时，可以提高系统的抗干扰能力。

[0194] 对于BOT测量系统，由于Hessian矩阵及其正定性。因此，上述基于相关熵最小二乘的目标函数在区域测量区域为凸。公式(19)中的标量处理提供了一种标准技术，用于找到
用于向量优化的帕累托最优点。可以通过使用截断的牛顿内点法来近似最大似然模型状态
变量。即，

[0195]

[0196] 高‑低‑高区间可通过应用回溯搜索方法确定，搜索方向 需满足如下定义

[0197]

[0198] 在上述优化过程中，根据Lasalle不变性原理和Lyapunov稳定性，起始点应在可行范围内。对于纯方位机动目标跟踪，可以通过交叉最小二乘方法获得具有最大渐近性的变
量。

[0199] 最后，约束先验的修正概率密度函数可以近似为如下的高斯分布

[0200]

[0201] 其中， 和 分别是均值和协方差。

[0202] 无迹滤波的前提是：在固定数量的参数的情况下，近似高斯分布要比近似任意非线性函数更容易。矩阵均方根可以通过矩匹配法获得，即

[0203]

[0204]

[0205] 其中下标t和s分别表示截断和平滑。并且

[0206]

[0207] 对于动态系统估计，来自截断先验和平滑度补偿的统计约束将同时完全纳入整个过滤过程。模型状态的约束估计是在无迹卡尔曼滤波框架中实现的，包括以下步骤：

[0208] 可行域sigma点

[0209] 在k时刻，对于模型m，令 和 分别表示模型状态的第i个sigma点及其权重。根据公式(30)中的修正均值和协方差选择集合，

[0210]

[0211]

[0212]

[0213] 其中 表示 模型 状态 的维 数 ，λ是 可信 任区 域的 半径 ，并 且是均值Cholesky解的第i行。权重因子之和为1，即,

[0214] 约束预测

[0215] k时刻，约束预测均值为，

[0216]

[0217] k时刻，约束预测方差为，

[0218]

[0219] 模型状态更新

[0220] k时刻，定域观测序列的预测均值 为

[0221]

[0222] 观测协方差为

[0223]

[0224] 状态与观测的互协方差为

[0225]

[0226] 模型状态估计

[0227] 在上述约束条件下，卡尔曼滤波增益为

[0228]

[0229] 模型状态更新为

[0230]

[0231]

[0232] 将当前测量值的校正和历史更新状态的补偿增强到状态评估中，融合效果分别表示为τ＝0,1。并且可以通过使用模糊高斯函数的隶属度来计算权重系数

[0233]

[0234]

[0235] 其中， 是模型状态的新息协方差，从估计精度衡量的意义上讲，当前时刻观测信息更为准确时，αk更接近于1。

[0236] 对于状态向量的最终估计，滤波后的均值和协方差是由模型状态估计的加权和组成的融合输出，其计算如下

[0237]

[0238]

[0239] 其中， 是k时刻模型m的重要性权重。

[0240] 在边界约束下，通过稀疏正则化确定域的可信任域，并且当系统状态远离边界时，可以适当地提高系统的抗干扰能力。同时，帕累托最优集是原始的多目标问题，并且与经过
测试的现有方法具有竞争力，并且可以有效地平衡系统状态估计的准确性和收敛率。

[0241] 综上，纯方位机动目标跟踪被建模为马尔可夫多模型闭环系统。为了克服模型切换的奇异性问题，本实施例提出了带有L1约束正则化的边缘化多模型无味卡尔曼滤波方法
法。首先，通过边缘化可以有效地降低多模型状态空间的高维。其次，将截断和平滑度先验
考虑到修改后的次优重要性密度函数中，并且模型转换概率矩阵由从约束辅助滤波得出的
粒子权重组成。第三，从可行的感兴趣区域中选择伪域测量序列，并通过测量序列的标量保
证Pareto最优解。最后，可以通过最优反馈稳定地预测控制闭环多模型系统，并通过使用模
糊因子对其进行加权。

[0242] 为了验证本实施例中基于L1正则无迹变换的约束多模型滤波方法，下面结合仿真做出进一步的说明。

[0243] 蒙特卡洛仿真实验在以下场景下运行：S1，一个观测间隔固定的模拟场景，以及S2，观测间隔随时间变化的真实场景。为了比较，还介绍了常规IMMEKF，IMMUKF和MMRBPF的
滤波性能。比较和研究了常规IMMEKF，IMMUKF和MMRBPF的滤波性能。选择两个定量性能参数
矩阵作为参考，它们是(a)均方根误差(RMSE)，测量每次的总平均滤波性能，以及(b)均方根
均方误差(RTAMSE)，测量目标机动后的总平均过滤性能。这些定义分别显示为以下公式
(47)和(48)。

[0244]

[0245]

[0246] 其中te和tmax分别是目标机动的开始和结束时刻。

[0247] 使用控制变量方法研究了来自动态系统本身的固有误差和来自测量环境的外部干扰的影响。

[0248] 固定间隔的仿真场景

[0249] 对于大多数系统，状态误差与过程噪声和测量噪声之间的相关性均为零。引入过程和测量噪声对均值和协方差的影响，其精确度与状态不确定性相同。测量间隔固定为，以
2 2
下是5种情况，case1):σv＝0.03km/s ,σe＝1.5mrad；case2):σv＝0.05km/s ,σe＝1.5mrad；
2 2
case3):σv＝0.1km/s ,σe＝1.5mrad；case4):σv＝0.1km/s ,σe＝3mrad；case5):σv＝0.1km/
2
s ,σe＝5mrad.

[0250] 过程噪声影响

[0251] 在一定范围内，随着过程噪声的增加，图2‑4显示了所得RMS位置误差曲线随时间变化的趋势。

[0252] 当过程噪声为0.005km/s2时，测量噪声相对较大，相比之下，IMMUKF的RMS位置误差曲线显示出最大的总体偏差轨迹差。在平滑的CV运动阶段，MMUKFL1的偏差甚至大于
IMMEKF。这主要是因为当过程噪声足够平稳时，一阶泰勒级数线性化误差相对较小。随着过
程噪声的增加，跟踪性能在平滑CV运动过程中变差，而在机动运动过程中变好。MMUKFL1具
有明显的稳定性和平滑性，这意味着更好的滤波一致性和更强的鲁棒性。

[0253] 观测噪声的影响

[0254] 同样，在一定范围内，随着测量噪声的增加，图3、图5‑6显示了所得RMS位置误差曲线随时间变化的趋势。观测值的不稳定性越大，滤波性能越差。显然，通过约束正则化可以
平滑MMUKFL1的过滤误差。

[0255] 直方图7和8分别报告RMSE和RTAMS的统计信息。重要的是，MMMUKFL1表现出更好的滤波一致性和更强的估计鲁棒性，这主要是由于使用了截断和平滑度测量约束。

[0256] 表1

[0257]

[0258] 此外，表1总结了第二种情况下运行一次MC所需的平均执行时间。与传统的基于IMM的滤波器相比，MMMUKFL1显示出相对较高的计算要求，这主要是由于创新协方差的计算
中增加了测量噪声。显然，与现有技术中分析的基于PF的滤波器相比。本实施例中提出的方
法提供了在线跟踪机动目标的能力。

[0259] 观测间隔的影响

[0260] 接下来研究采样间隔的影响。对于第一种情况，将其增加到T＝3s，其它参数保持不变。图9显示了五个滤波器的CV模型切换概率曲线。

[0261] 尽管正确的转弯模型的模型概率有微小的提高，基于IMM的滤波器无法捕获稀疏样本的目标机动，这主要是因为在交互阶段未考虑当前的测量信息。同时，由于在状态模型
选择阶段考虑当前时刻测量信息，MMRBPF显示出比IMMEKF和IMMUKF更好的性能，但比
MMMUKFL1差。此外，结合着目标动力学特征分析，通过边界优化有效地限制了测量不确定
性，MMMUKFL1的模型切换概率与目标实际机动模型更加一致。

[0262] 真实稀疏观测场景

[0263] 实际的稀疏测量场景由雷达收集的一批飞机方位数据组成。飞行持续时间为107秒，包括40个非周期性采样点。为了研究所提出方法在三维空间中的跟踪特性和质量，假设
高度参数已知为1km。通常，由于传感器的不准确性和测量噪声的影响，无法准确知道状态
的初始值。然而，可以合理推断实际目标初始位置在传感器的观测值附近。它是由高斯分布
给出的 其中

[0264]

[0265] P0|0＝diag[0.15km2 02km2s‑2 0.15km2 02km2s‑2 0.15km2 02km2s‑2]T.

[0266] 系统过程噪声 和观 测噪 声服从高斯分布 ，协方差矩阵分 别为 ：

[0267] 实际上，对于多模型切换情况，基于IMM的滤波器的估计轨迹与真实轨迹之间的偏差大于1.2千米，因此在稀疏观测环境中，选择MMRBPF和CMMPF来与提出的的MMUKFL1在跟踪
性能方面进行比较。

[0268] 图10和11分别显示了2D平面和3D空间中提出方法的滤波轨迹。图12显示了位置RMSE随时间变化的趋势。当在25s发生强机动时，所有滤波器都显示出明显的可见气泡，在
接下来10s的稀疏观测阶段，MMRBPF出现相对较大的偏差，其次是CMMPF和MMUKFL1。之后，在
最后的35s至40s的轻度运动阶段，MMUKF收敛最快。图13显示了上述操作过程中MMUKF的模
型切换概率，仿真结果与真实轨迹变化一致性说明了该方法的优势。鉴于此，认为对于切换
情况的静态方法的修改已经过时。

[0269] 以上所述仅为本发明的优选实施例，并非因此限制本发明的专利范围，凡是在本发明的发明构思下，利用本发明说明书及附图内容所作的等效结构变换，或直接/间接运用
在其他相关的技术领域均包括在本发明的专利保护范围内。