一种非高斯非稳态噪声建模方法

摘要

本发明公开了一种非稳态非高斯观测噪声建模方法，用于目标跟踪系统的滤波器状态估计过程中对观测噪声的处理。本发明的主要特征在于：利用高斯混合模型对非稳态非高斯观测噪声进行近似，混合模型中高斯成员的分布参数融合在滤波器的迭代过程中进行计算和更新，能针对观测噪声统计特性的变化在线调节。在滤波器的每一次迭代处理中，高斯成员的分布参数分为先验参数和后验参数。先验参数由上一次迭代获得的参数估计结果进行计算，后验参数通过极大似然估计求取。利用后验分布参数构建高斯成员分布，并组成观测噪声的高斯混合近似模型，该模型能够保证滤波器精度维持在较高水平。本方法对于观测噪声的建模具有较高精度和较强鲁棒性。

说明书

技术领域

本发明涉及目标跟踪技术领域，特别涉及一种非高斯非稳态噪声建模方法。

背景技术

目标跟踪是将传感器连续收到的目标运动信号数据对应于各种不确定信息源所产生的不 同观测集合进行滤波处理，并相应地估计出机动目标的目标状态参数，如位置、速度和加速 度等。由于传感器量测数据中含有大量的干扰成分，必须对量测信息进行处理，因此目标跟 踪过程也是一个消除误差的处理过程。目标跟踪已经伴随着雷达、红外、激光、声纳等传感 器系统的发展走过了几十年的历程，并由最初的“单目标跟踪”问题逐步发展为多个传感器跟 踪多个目标的“多目标跟踪”问题。

贝叶斯滤波为目标跟踪提供了一个利用概率方法递归进行状态估计的框架，利用所有的 可用信息，以一定信度确定运动目标处于某一确定状态。贝叶斯滤波通过时间更新和状态更 新两个阶段完成目标状态的预测和传感器观测信息的融合，最终实现对目标状态的估计。在 实际应用中，贝叶斯滤波框架有多种实现方式，从概率的表示方式上可将实现方法分成离散 型和连续型。卡尔曼滤波和粒子滤波分别是连续性和离散型方法的典型代表。在过去的几十 年中，研究人员在上述两类方法的基础上进行了深入研究，致力于滤波器估计性能的提高， 并发现，如果能准确获取噪声的统计信息，就能够消除或减弱噪声对状态估计产生的干扰， 提高估计性能。通常情况下，为了简化问题，时变噪声都被近似为零均值高斯白噪声。然而 在很多实际问题中，噪声序列可能是非高斯的，如果近似为高斯噪声将会引入较大的估计误 差。针对这一问题，产生了许多非高斯噪声的建模方法。其中一种典型的方法是将非高斯噪 声的概率密度函数(probability density function,PDF)近似为混合高斯模型，也就是利用一组高 斯概率密度的加权和形式表示非高斯噪声分布。混合高斯模型中高斯成员分布参数的选择决 定了非高斯噪声的概率分布，通常的方法是使用期望最大化(expectation maximization,EM)、 极大似然估计、非高斯自回归处理或相关函数法等技术确定每个高斯成员的分布参数(主要包 括期望、方差、权值)，利用这些高斯成员组成的混合模型对非高斯噪声进行近似，并将其用 于滤波器的处理过程。

上述建模方法能够较好的描述噪声的非高斯特性，但它们有一个共同特点，这些方法都 假设噪声的统计特性是固定不变的，用于近似噪声的高斯混合模型中的每个成员的分布参数 一旦确定，在整个滤波过程中也不再改变，也就是将非高斯噪声作为平稳过程来处理。然而 在实际问题中，传感器量测值常常受到非高斯非稳态噪声的影响。例如在雷达目标跟踪系统 中，雷达量测容易受到闪烁噪声影响，闪烁噪声是一种具有重尾特性的非高斯非稳态噪声， 它由目标表面不同部位产生的两个或以上反射间的干涉形成，使反射中心沿目标长度方向发 生随机游走，导致波面阵传播过程(propogating)产生形状失真，引起测角误差进而干扰跟踪器 的视线(line-of-sight)量测。对于这类非高斯非稳态噪声，以往的高斯混合模型由于使用了确定 不变的分布参数，缺乏对噪声统计特性变化的有效跟踪，最终可能导致滤波器的估计性能下 降，甚至估计发散。

有关文献：

Ibatoulline E A.Iterative method of maximum likelihood for estimation of non-Gaussian  probability density of signals and interferences[C].In Proc.IEEE Int.Symp.EMC,2003.

郝红霞,李红.非高斯非平稳背景的红外弱小目标检测[J].激光与红外，2006.

Liu W,Pokharel P P,Principe J C.Correntropy:Properties and applications in non-Gaussian  signal processing[J].IEEE Trans.Signal Process,2007.

公开号为CN103077303A的专利《电子倍增CCD噪声模型的参数估计方法》公开了一种 电子倍增CCD噪声模型的参数估计方法，首先对噪声分布模型进行初始化设置，以拟作混合 高斯分布模型进行处理；然后对噪声分布模型进行极大似然化迭代，将设定的初始值代入混 合模型中，求样本值来自高斯源的后验概率密度,再将潜在数据代入不完全数据的对数函数中 计算偏导数,求取极值，即可得到参数的迭代估计值；最后进行判断步，将迭代估计值与初始 值比较,依据循环终止条件，判断是否满足终止条件,若满足则迭代停止，若不满足则把迭代值 设为初始值,重新进行极大似然化迭代。该发明以简单的步骤实现了对电子倍增CCD图像噪 声参数的极大似然估计，有效地降低了极大似然法的复杂度,能对电子倍增CCD图像噪声进 行快速、准确的估计。公开号为CN102882820A的专利《认知无线电中非高斯噪声下数字调 制信号识别方法》公开了一种在认知无线电中非高斯噪声背景下基于分数低阶循环谱相关系 数的数字调制识别方法，解决在认知无线电中非高斯噪声背景下识别调制识别性能差、计算 复杂度高的问题。其步骤为：对接收到的信号进行采样；计算采样后信号的分数低阶循环谱 计算信号的分数低阶循环谱在循环频率ε＝0处的截面和在频率f＝0处的截面以及在循环频 率ε面的投影和在频率f面的投影的相关系数ρ1、ρ2、ρ3、ρ4和ρ5；设置信号集的判决门限， 通过采用基于判决树的分类器将不同调制方式的信号识别出来。在非高斯Alpha稳定分布噪 声下，该发明不仅具有较高的识别率和良好的稳健性并且计算复杂度更低，更适合于认知无 线电系统。

综上所述，现有的非高斯噪声模型只能应用于非高斯稳态噪声的近似，而无法表示具有 非稳态特性的噪声，致使滤波器在处理具有非稳态噪声影响的跟踪问题时产生较大的估计误 差。因此如何针对非高斯非稳态噪声设计一种建模方法，是提高滤波器估计性能亟待解决的 问题。

发明内容

本发明的目的是，为具有非高斯非稳态特性的观测噪声提供一种具有自适应在线调节能 力的建模方法。本发明中提供的方法区别于现有方法的显著特征在于：在将非高斯非稳态噪 声建模为高斯混合模型后，各高斯成员的分布参数的取值并不是固定不变的，而是在滤波器 每一次迭代过程中针对实际噪声的统计特性对分布参数进行计算和更新。在滤波器的每一次 迭代处理中，将噪声模型每个高斯成员的分布参数分为先验参数和后验参数。先验参数在上 一次迭代获得的参数估计结果基础上进行初始化，后验参数通过在整个高斯混合模型的观测 似然函数取极大值处获得。利用这样获得的分布参数构建高斯分布成员，并组成对于真实噪 声的高斯混合近似模型，该模型能够有效针对实际噪声的统计特性变化做出调节，保证滤波 器的精度维持在较高水平。

技术方案：

本发明基于高斯混合模型，为具有非高斯非稳态特性的观测噪声提供了一种具有自适应 在线调节能力的建模方法。模型建立的实质是对观测噪声的真实概率密度函数进行估计，其 特征在于，噪声模型的参数确定是在滤波器估计过程的每一次迭代计算中完成的，以迭代过 程中的任意k时刻为例，该方法包含以下步骤：

步骤1、对于任一k时刻的待建模观测噪声，将其近似为高斯混合模型形式，该模型由一 系列服从高斯分布的噪声成员构成；

本发明所涉的混合模型中高斯成员的分布参数主要包括：期望、方差、分布权值，这些 参数在滤波器的每一次迭代分为两种状态：先验参数和后验参数。先验参数由上一次迭代获 得的后验分布参数初始化而来，后验参数根据当前观测噪声的实际分布情况由先验参数更新 而来。以下步骤中，步骤2对高斯成员的先验分布参数进行设定，步骤7将先验分布参数更 新到后验分布参数。

先验分布参数计算：

步骤2、对于滤波过程中的k时刻，若k＝0，则需要进行先验分布参数的初始化，此时 滤波器还未开始迭代；若k>0，则先验分布参数由上一次迭代(k-1时刻)获得的后验分布参 数计算而来；

步骤3、对当前时刻观测噪声进行采样，形成噪声采样集合；

步骤4、利用先验分布参数计算每个高斯成员关于噪声采样的后验条件概率；

步骤5、建立高斯混合模型的似然函数及其对数形式；

步骤6、通过拉格朗日乘子构建对数似然函数的拉格朗日函数；

后验分布参数计算：

步骤7、通过极大似然估计分别求得混合模型中每个高斯成员的后验分布参数，包括以下 步骤：

1)通过令对数似然函数的拉格朗日函数关于高斯成员期望的偏导数为0，求得高斯成员 的期望；

2)通过令拉格朗日函数关于高斯成员标准差的偏导数为0，求得高斯成员的方差；

3)通过令拉格朗日函数关于高斯成员分布权值的偏导数为0，求得高斯成员的分布权值；

步骤8、利用每次迭代的后验高斯成员分布参数构建当前观测噪声的高斯混合近似模型， 这一模型可用于各种滤波算法中对观测噪声的处理，提高滤波器的估计性能；

步骤9、如果滤波器迭代完成，则由滤波器输出状态估计结果，如果迭代未完成，则返回 步骤2。

本发明利用高斯混合模型对非高斯非稳态观测噪声进行建模，模型中高斯成员的分布参 数是在滤波器每次迭代过程中根据噪声的实际统计特性进行在线估计和更新的，最终获得的 噪声模型用于各种滤波器，能够较为精确的反应实际噪声的分布情况，有效提高滤波器估计 性能。

其中主要参数包括：混合高斯模型成员数、观测噪声采样数、拉格朗日乘子、高斯成员 先验分布参数、高斯成员后验分布参数以及滤波器迭代次数。混合高斯模型成员数是指用于 近似非高斯非稳态观测噪声的混合模型包含的高斯分布的数量，观测噪声采样数是指由观测 噪声获得的采样集合包含的采样数量。拉格朗日乘子用于为高斯混合模型的对数似然函数构 建拉格朗日函数。高斯成员先验分布参数是依据滤波器上一次迭代结果为噪声混合模型中每 个高斯成员分配的先验分布参数，包含先验期望、先验方差、先验分布权值。高斯成员后验 分布参数是滤波器每次迭代过程中依据当前观测噪声的实际分布情况对高斯成员先验分布参 数进行更新获得的后验分布参数，包含后验期望、后验方差、后验分布权值。滤波器迭代次 数用于判断滤波器是否终止迭代。

本发明在滤波器每次迭代过程中计算和更新噪声模型中高斯成员的分布参数，并将分布 参数分为先验分布参数和后验分布参数两种状态。这与贝叶斯滤波框架中时间更新和状态更 新两种阶段的划分是一致的，因此本发明提出的噪声模型能广泛适用于贝叶斯滤波的多种实 现方式。

而且，本发明在利用极大似然估计求取高斯成员后验分布参数时，考虑到计算效率和难 度，选取了在同一点处取极值的对数似然函数进行替换。同时通引入拉格朗日乘子，构建拉 格朗日函数，在约束条件下通过求极值的方法获得高斯成员后验分布参数。

有益效果：本发明利用一系列加权高斯分布的混合模型对非高斯非稳态观测噪声进行建 模，将噪声模型的分布参数估计融合在滤波器的迭代过程中进行，结合滤波器迭代中的时间 更新和状态更新，依据实际噪声的统计特性实现噪声分布参数的计算和更新，最终获得的噪 声模型能够更接近真实噪声的分布情况，将此噪声模型用于滤波器状态估计，能提高滤波器 估计精度，加快收敛速度，是滤波器具有更强的鲁棒性。此外，本发明针对非高斯非稳态噪 声进行建模，同时具有广泛的适用范围，可以推广到用于对任意概率分布的噪声进行建模。

附图说明

图1为本发明所述非高斯非稳态方法的系统结构图。

图2为本发明所述非高斯非稳态方法的流程图。

具体实施方式

下面结合附图对本发明的具体实施方式进行详细说明。

本发明公开了一种基于混合高斯模型的非高斯非稳态观测噪声建模方法，该方法的具体 实施包括模型建立、观测噪声采样、拉格朗日函数构建以及高斯成员分布参数估计等关键内 容。本发明所述的非高斯非稳态观测噪声建模方法通过计算机程序实施，图1所示是计算机 实现的系统结构图。下面将按照流程详述本发明提出的技术方案的具体实施方式，流程如图 2所示。本发明所述噪声模型是对观测噪声概率密度函数的近似，该实施方式主要包含以下 关键内容：

对于式(1)描述的目标跟踪系统，为k时刻系统状态变量，服从分布p(x_k|x_k-1)， 表示过程噪声。观测变量是关于状态变量x_k和观测噪声的函数， 服从分布p(z_k|x_k)。f(·)和h(·)分别称为状态转移函数和观测函数。目标跟踪的任务是，在 给定初始状态x₀的情况下，通过不断获得和处理观测值z_k，对每一时刻的状态值x_k进行估计。 本发明所述内容即对非高斯非稳态观测噪声v_k的概率密度函数进行估计。

$\left(\begin{array}{c}{x}_k=f_k\left(x_{k-1}\right)+u_k\\ z_k=h_k\left(x_k\right)+v_k\end{array}\right)---\left(1\right)$

滤波器迭代开始前，执行步骤1建立观测噪声模型，并初始化：

步骤1、利用高斯混合模型对k时刻待建模观测噪声进行建模，将观测噪声的概率密度函 数近似为一系列服从高斯分布的成员密度函数，如式(2)所示。

$p\left(v_k\right)=\underset{i=1}{\overset{M}{\Sigma }}{\alpha }_{i,k}{p}_i\left(v_k\right)=\underset{i=1}{\overset{M}{\Sigma }}{\alpha }_{i,k}N(v_k;{\mu }_{i,k},{\sigma }_{i,k}^2)---\left(2\right)$

其中M为混合模型中高斯成员的数量，通常情况下M＝10可以得到较好的建模效果，并保证 滤波器估计效率。p_i(v_k)表示混合模型中第i个高斯成员的概率密度函数，μ_i,k、和α_i,k分 别为p_i(v_k)对应的期望、方差以及权值，α_i,k满足式

本发明所涉的混合模型中高斯成员的分布参数在每一次滤波器迭代分为两种状态：先验 参数和后验参数。以下步骤中，步骤2对高斯成员的先验分布参数进行设定，步骤7将先验 分布参数更新到后验分布参数。

步骤2、若k＝0，说明滤波器迭代还未开始，为这些高斯成员初始化先验分布参数：设 定每个高斯成员的分布权值为α_i,0＝1/M，期望为μ_i,0＝0，各高斯成员的方差可依据传感器 的先验观测误差特性进行设定，无需考虑实际观测噪声分布情况，造成的偏差可由滤波过程 中的参数确定和更新进行修正。若k>0，每次迭代中高斯成员的先验分布参数由上一次迭代 获得的后验分布参数计算而来，具体方法如下：

先验权值为：

${\alpha }_{i,k^{-}}=\frac{1}{N}---\left(3\right)$

先验期望为：

${\mu }_{i,k^{-}}=\left(\begin{array}{cc}{\mu }_{i,k-1}(1+(i-\frac{M-1}{2})\frac{1}{M}),& \mathrm{if>}\\ {\mu }_{i,k^{-}}={\mu }_{k-1}(1+(i-\frac{M}{2})\frac{1}{M}),& \mathrm{if>}\end{array}\right)---\left(4\right)$

先验方差为：

${\sigma }_{i,k^{-}}^2={\sigma }_{i,k-1}^2---\left(5\right)$

步骤3、对当前时刻噪声进行采样，形成包含N个噪声采样的采样集合，其中第j个噪声 采样表示为N的取值与噪声复杂程度有关，通常情况下可以取N＝10

步骤4、计算每个高斯成员关于噪声采样的后验条件概率，方法如下：

$\left(\begin{array}{c}P\left\{p_i\left(v_k^j\right)\right|v_k^j\}=\frac{P\{p_i\left(v_k^j\right),v_k^j\}}{P\left\{v_k^j\right\}}\\ =\frac{P\left\{v_k^j\right|p_i\left(v_k^j\right)\left\}P\right\{p_i\left(v_k^j\right)\}}{\underset{i=1}{\overset{M}{\Sigma }}p\left\{v_k^j\right|p_i\left(v_k^j\right)\left\}P\right\{p_i\left(v_k^j\right)\}}\\ =\frac{p_i\left\{v_k^j\right|{\mu }_{i,k},{\sigma }_{i,k}^2\}{\alpha }_{i,k}}{\underset{i=1}{\overset{M}{\Sigma }}{p}_i\left\{v_k^j\right|{\mu }_{i,k},{\sigma }_{i,k}^2\}{\alpha }_{i,k}}\end{array}\right)---\left(6\right)$

步骤5、将式(3)、(4)、(5)代入式(6)可以得到高斯成员关于噪声采样的后验条件概率。

建立高斯混合模型的对数似然函数如式(7)所示：

$L_k=\underset{j=1}{\overset{N}{\Sigma }}\ln \underset{i=1}{\overset{M}{\Sigma }}{\alpha }_{i,k}{p}_i\left\{v_k^j\right|{\mu }_{i,k},{\sigma }_{i,k}^2\}---\left(7\right)$

步骤6、引入拉格朗日乘子，并由约束条件得到拉格朗日函数：

$L_k=\underset{j=1}{\overset{N}{\Sigma }}\ln \underset{i=1}{\overset{M}{\Sigma }}{\alpha }_{i,k}{p}_i\left\{v_k^j\right|{\mu }_{i,k},{\sigma }_{i,k}^2\}-\lambda (\underset{i=1}{\overset{M}{\Sigma }}{\alpha }_{i,k}-1)---\left(8\right)$

步骤7、通过极大似然估计分别求得混合模型中每个高斯成员的后验分布参数，包括以下 步骤：

1)求解可以得到k时刻高斯成员的后验期望

${\mu }_{i,k}=\frac{{\Sigma }_{j=1}^{N}P\left\{p_i\left(v_k^j\right)\right|v_k^j\}{v}_k^j}{{\Sigma }_{j=1}^{N}P\left\{p_i\left(v_k^j\right)\right|v_k^j\}}---\left(9\right)$

2)求解可以得到高斯成员的后验方差

${\sigma }_{i,k}^2=\frac{{\Sigma }_{j=1}^{N}P\left\{p_i\left(v_k^j\right)\right|v_k^j\}{(v_k^j-{\mu }_{i,k})}^2}{{\Sigma }_{j=1}^{N}P\left\{p_i\left(v_k^j\right)\right|v_k^j\}}---\left(10\right)$

3)求解并由约束条件得到拉格朗日乘子取值和高斯成员对应的后 验权值

λ＝M  (11)

${\alpha }_{i,k}=\frac{1}{M}\underset{j=1}{\overset{N}{\Sigma }}P\left\{p_i\left(v_k^j\right)\right|v_k^j\}---\left(12\right)$

将式(6)代入分别带入式(9)、(10)和(12)，获得对混合模型中每个高斯成员的后验分布参数。

步骤8、将步骤7获得的后验高斯成员分布参数代入式(2)，构建当前观测噪声的概率密 度函数。这一噪声密度函数可用于各种滤波算法中对观测噪声的处理，提高滤波器的估计性 能。

步骤9、判断滤波器的迭代状态，如果滤波器迭代完成，则由滤波器输出状态估计结果， 如果迭代未完成，则返回步骤2。