在实际的观测环境中经常会存在动力模型不准确和观测模型出现粗差的情况，采用标准Kalman滤波无法抵御粗差对估值的影响，一般会通过引入1个自适应因子来进行调节^[1-3]。本文在抗差Kalman滤波的基础上引入2个自适应因子分别对动态模型和观测模型进行调节，通过这2个自适应因子共同求出状态估值。然而在抗差自适应Kalman滤波中需对每一观测时刻的数据进行迭代计算，从而会使得运算时间变长，对数据的精度也会造成一定的影响。本文通过构建卡方统计量，提出基于卡方检验的抗差自适应Kalman滤波，并通过工程实例对其性能进行验证。

1 抗差自适应Kalman滤波

当观测值存在粗差或者滤波发散等问题出现时，由于Kalman滤波不具备抗差性，易受到影响。针对此问题，杨元喜院士^[4]提出了抗差自适应滤波，其基本思想是在测量值存在粗差时，对测量值使用抗差估计方案处理，在实现抗差估计的前提下，如果动力模型出现非正常误差，通过自适应因子对其调整。即在用Kalman滤波进行运算时，如果测量序列L_k中存在粗差，其滤波结果精度会降低。根据M估计等价权原理对增益矩阵进行分析可知，选择恰当的权函数替换观测噪声的协方差矩阵可以降低粗差对滤波估值精度的影响。利用自适应因子调节增益矩阵可以合理利用观测数据中的信息，减少测量粗差对平差后的结果带来的影响。抗差自适应Kalman滤波的改进部分如下^[5-6]。

1)增益矩阵：

(1)

2) 状态估计：

(2)

3) 估计误差方阵：

(3)

4) 抗差自适应滤波的解：

(4)

式中，观测序列的协方差阵对应的权矩阵,α_k是自适应因子。

2 双自适应因子滤波

当动态模型和观测模型均出现误差时，可以在抗差Kalman滤波中引入1个自适应因子对动态模型和观测模型进行调整，或在抗差自适应Kalman滤波基础上引入2个自适应因子对动态模型和观测模型进行调节，建立双自适应因子滤波算法。

2.1 滤波模型的推导

动态向量的预报误差方程为：

(5)

式中，动态向量预报，状态转移矩阵Φ_k,k－1中包含动力学积分过程，动态向量预报误差V_Xk也反映了动力学的模型误差。

观测向量的误差方程为：

(6)

式中，V_k反映的是观测值的误差。

通过最小二乘原理构造的损失函数为：

(7)

式中，P_k和P_Xk分别为L_k的等价权矩阵和的权矩阵，即,通过对X_k求导，并使其导数为0，可得：

(8)

通过矩阵的恒等变换可知：

(9)

式中，是增益矩阵，P_k为动态估计序列的权矩阵，P_Xk为动态预测向量的权阵，α_k是k时刻的观测模型自适应因子，β_k是k时刻的动力学模型信息的自适应因子^[7]。

2.2 自适应因子选取

自适应因子构造方法是根据指数型的两段法的计算方法，得到观测模型的自适应因子α_k和动力学模型自适应因子β_k：

式中，α_ki∈[0, 1]，β_ki∈[0, 1]。

1) 自适应因子α_k的选取

公式反映了观测值的残差，通过对其第i项分量V_ki进行标准化处理，得到，式中δ_Vki是V_ki的均方差:

(10)

式中，λ_a是一个阈值，取值范围为λ_a∈[1.0, 3.0]。

2) 自适应因子β_k的选取

模型信息的准确性是通过式(5)反映的，将式(9)代入式(5)得：

(11)

式(9)作为动态的预测残差，通过对其第i项分量进行标准化处理，得到是的均方差，或者直接取作为动态的预测残差，也称为简化的动态预测残差。

(12)

式中，λ_β是一个阈值，取值范围为λ_β∈[1.0, 3.0]^[7]。

3 基于卡方检验的抗差自适应Kalman滤波

在本节中，离散系统的动态方程和观测方程与离散系统下的标准Kalman滤波的相同。判断粗差是否存在，主要根据Kalman滤波中的预测残差。预测残差为观测值的预测值与实际观测值之差，即

(13)

式中，状态预测值和观测预测值分别为：

(14)

(15)

如果测量值L_k正常，且预测值也准确，则预报值与观测值的观测误差服从正态分布。同样，残差V_k也服从正态分布，其协方差矩阵是：

(16)

式中，为预测估值的协方差矩阵，D_k为观测向量的协方差矩阵^{[8, 9]}。

在粗差混入观测数据的情况下，残差的分布会偏离正态分布。由于各观测值相互独立，则统计量为：

(17)

式中，V_k(i)表示第i个观测残差，D_Vk(i,i)为V_k的协方差矩阵上对角线上的元素^[10]。在L_k中无粗差的情况下，统计量ε_k服从自由度为1的中心χ²分布，反之ε_k服从自由度为1的非中心χ²分布。

检验方法如下。原假设H₀：观测值中无粗差；备选假设H₁：观测值中存在粗差。显著性水平为a，若ε_k≤χ_a²(1)，则接受H₀假设，检验通过；反之，拒绝H₀假设，接受H₁假设。χ_a²(1)为χ²(1)分布上的a分位数。

在观测值中可能存在粗差的情况下，首先采用上述方法进行检验。如果存在粗差，则进行抗差迭代，得到精确状态估值；如果不存在粗差，无需进行抗差迭代，以Kalman滤波得到的估值作为滤波结果。

4 实例分析

在实际观测中异常值不仅仅是含有粗差，还有可能是环境因素使得形变异常，故在进行粗差探测前需对异常值进行处理。我们选取贵广高铁某隧道观测了16期的水准网进行分析，该数据中无异常形变值，其水准网采用二等精密水准测量，在平差后该数据的精度符合相应精度要求。为了检验自适应双因子滤波法的抗差效果，对该水准网中变形监测点某点的第8期加入2 mm(远大于2倍中误差)的粗差，观测周期为7 d。通过使用双自适应因子滤波算法、基于卡方分布的抗差自适应Kalman滤波对加入粗差后的观测数据进行滤波，然后与标准Kalman滤波的滤波估计进行对比，对比结果如图 1~3。

在观测数据中不含粗差时，双自适应因子滤波和标准Kalman滤波的处理结果如图 1。由图 1可知，在观测数据中无粗差时，双自适应因子滤波与标准Kalman滤波的结果十分接近，即使在初始值存在一定误差的情况下也可以快速收敛。因此，这两种方法在不含粗差的情况下是等效的。在观测数据中加入粗差后采用双自适应因子滤波和标准Kalman滤波的处理结果见图 2、图 3。由图可知，当观测数据中存在粗差时，双自适应因子滤波结果与无粗差的观测数据的观测值十分接近，说明双自适应因子能够有效地抑制粗差对滤波估值的影响。

用3种方法(方法1：双自适应因子滤波；方法2：抗差自适应Kalman滤波；方法3：基于卡方分布的抗差自适应Kalman滤波)进行数据处理，其滤波估值的方差见表 1(单位mm²)，各模型的处理结果见图 4~7。

由图 4可知，当观测数据中无粗差时，抗差自适应Kalman滤波估值与标准Kalman滤波估值均随着时间的增加而收敛，但是抗差自适应Kalman滤波的收敛速度明显比标准Kalman滤波慢。这是因为抗差自适应Kalman滤波对每一时刻的观测值均进行了抗差迭代，从而会使得计算时间变长，估计精度也随之降低。

由图 5可知，在观测值中含有粗差时，基于卡方检验的抗差自适应Kalman滤波的估值与无粗差的观测数据十分接近，而且可以快速收敛，其精度也优于抗差自适应Kalman滤波。这是由于该算法综合了统计模型，通过统计检验先对粗差进行判别，使得计算效率和精度均得到了提高。

通过图 6、图 7以及表 1可知，对于小范围内的粗差，双自适应因子滤波模型的滤波精度要优于抗差自适应Kalman滤波。在小范围存在粗差和短时间内，双自适应因子滤波不仅可以抵御粗差对滤波估值的影响，同时可以削弱动态模型不准确的影响。通过均方误差可以看出，双自适应因子滤波与抗差自适应Kalman滤波的抗差性能均比较稳定。通过图 7的残差图和均方误差可知，基于卡方检验的抗差自适应Kalman滤波的效率和稳定性均优于抗差自适应Kalman滤波和双自适应因子滤波。

5 结语

本文将双自适应因子滤波算法、基于卡方检验的抗差自适应Kalman滤波应用于隧道的沉降变形观测。工程实例证明，双自适应因子滤波在小范围内存在粗差的情况下优于抗差自适应Kalman滤波，该算法可以很好地抵制粗差，减弱粗差对滤波估值的影响。而基于卡方检验的抗差自适应Kalman滤波的性能和效率均优于双自适应因子滤波和抗差自适应Kalman滤波。由于这两个自适应因子存在一定的相关性，如何消除该影响仍需进一步研究。