卡尔曼滤波这个词,很多做嵌入式、控制、导航或者信号处理的朋友应该都不陌生。我第一次接触它是在一个姿态解算的项目里,当时拿着一堆陀螺仪和加速度计的原始数据,发现单纯互补滤波出来的角度总是抖得厉害,后来同事甩给我一段卡尔曼滤波的代码,说“你把这个套进去试试”。结果确实稳了不少,但那段代码里五个公式我盯了一下午也没完全搞明白每个矩阵到底在干什么。后来花了很长时间从最小二乘、从贝叶斯估计的角度重新推了一遍,才算是真正把它的来龙去脉理清楚了。
这篇文章我想做的事情很直接:把卡尔曼滤波从“五个公式”还原成一个有血有肉的过程。它到底在解决什么问题,每个公式背后的直觉是什么,为什么是这个形式而不是别的形式,以及在实际写代码的时候哪些地方最容易踩坑。不管你是刚接触状态估计的新手,还是已经用过现成库但一直没搞懂原理的老手,我都尽量用大白话加上必要的推导把这件事讲透。核心关键词卡尔曼滤波、状态估计、协方差矩阵、预测与更新,会贯穿全文反复出现,因为这几个概念就是理解整个算法的钥匙。
1. 卡尔曼滤波到底在解决什么问题
1.1 从一个生活场景说起:你其实一直在做卡尔曼滤波
先别急着看公式。想象你在开车,前方有个路口,你需要估计自己距离路口还有多远。你脑子里有两个信息来源:一是你记得刚才看了一眼路牌,上面写着“前方500米”;二是你凭感觉觉得车速大概60码,开了大概10秒钟。这两个信息都不准——路牌可能立的位置有偏差,你的车速感觉也可能偏快或偏慢。但你不会只信其中一个,你会把两者“综合”一下,得出一个“大概还有300多米”的判断,而且如果路牌很清晰你就更信路牌,如果路牌模糊你就更信自己的速度感。
这就是卡尔曼滤波的核心思想:融合多个不完美的信息源,得到一个比任何单一信息源都更可靠的估计。它做的事情本质上就是“加权平均”,但这个加权不是拍脑袋定的,而是根据每个信息源的不确定程度动态调整的。路牌清晰(测量噪声小),权重就高;车速感觉模糊(过程噪声大),权重就低。卡尔曼滤波把这套直觉用数学严格化了。
在工程上,这个“距离”可以换成任何你想估计的量:无人机的高度、电池的剩余电量、机器人的位置、传感器的真实温度、甚至金融数据的趋势。只要你能建立“状态怎么随时间变化”的模型,并且有“能间接观测到状态”的测量手段,卡尔曼滤波就能派上用场。
1.2 为什么不能直接用测量值:噪声是绕不开的敌人
有人会问,既然有传感器测量,为什么不直接读测量值就完了?答案很简单:测量值本身带着噪声。你用一个温度传感器测水温,读数在25.1、25.3、24.9之间跳,真实温度可能稳定在25.0。你直接取某一次读数,误差可能就有0.3度;但如果你知道温度不会突变,把连续几次读数做平均,误差就能降下来。
问题是,简单的滑动平均有两个致命缺陷。第一,它对所有历史数据一视同仁,但越久远的数据参考价值越低,因为状态可能已经变了。第二,它没有利用“状态如何演化”的先验知识。卡尔曼滤波恰好补上了这两点:它用状态转移模型来预测当前状态(利用演化规律),用协方差来量化不确定性(决定该信预测还是信测量),并且每来一个新测量就更新一次,天然实现了“旧信息权重衰减”。
提示:很多人把卡尔曼滤波当成一个“滤波器”来理解,觉得它是用来去噪的。这个理解不算错,但不够准确。它更本质的身份是状态估计器——去噪只是它估计状态过程中的副产品。
1.3 线性高斯假设:卡尔曼滤波的适用边界
卡尔曼滤波不是万能的,它有两个硬性前提:系统是线性的,噪声是高斯的。线性意味着状态转移和测量关系都能写成矩阵乘法的形式;高斯意味着噪声的分布可以用均值和协方差完全描述。在这两个假设下,卡尔曼滤波被证明是最优估计器(最小均方误差意义下)。
如果系统非线性怎么办?那就得用扩展卡尔曼滤波(EKF)或者无迹卡尔曼滤波(UKF)。EKF的做法是在当前估计点附近做一阶泰勒展开,把非线性函数线性化;UKF则是用一组sigma点来传播分布,不做线性化近似。这两个变种在实际工程中用得非常多,但它们的根基都是标准卡尔曼滤波。所以把标准形式吃透,是理解所有变种的前提。
| 滤波器类型 | 适用系统 | 核心思路 | 典型场景 |
|---|---|---|---|
| 标准KF | 线性系统 | 矩阵运算直接传播 | 线性运动、简单信号 |
| EKF | 弱非线性 | 一阶泰勒展开线性化 | 姿态解算、目标跟踪 |
| UKF | 强非线性 | Sigma点无迹变换 | 高精度导航、复杂模型 |
| 粒子滤波 | 任意非线性非高斯 | 蒙特卡洛采样 | 多峰分布、遮挡跟踪 |
2. 核心思想拆解:预测与更新的双人舞
2.1 两个步骤循环往复:预测靠模型,更新靠测量
卡尔曼滤波的整个流程就是一个不断循环的两步舞:预测(Predict)和更新(Update)。我习惯把它类比成“先猜后改”。
预测阶段,你手里有一个上一时刻的状态估计,还有一个“状态怎么变”的模型。你根据模型往前推一步,得到当前时刻的先验估计。比如你知道上一秒小车在位置10米、速度2米每秒,模型告诉你匀速运动,那这一秒的预测位置就是12米。同时,你的不确定性(协方差)会因为模型本身的不完美而增大——推得越远,越没底。
更新阶段,传感器来了一个测量值。你拿测量值和预测值比较,算出一个残差(也叫新息)。然后根据预测的不确定性和测量的不确定性,算出一个卡尔曼增益,决定这次该信预测多一点还是信测量多一点。最后用增益加权修正预测,得到后验估计,同时把协方差也修正(通常会变小,因为融合了新信息)。
这个循环的美妙之处在于:协方差矩阵会自动调节信任分配。预测不确定性大,增益就偏向测量;测量噪声大,增益就偏向预测。你不需要手动调“信谁多一点”的系数,算法自己会算。
2.2 卡尔曼增益:整个算法的心脏
如果只能记住一个公式,那一定是卡尔曼增益的表达式。它长这样(一维简化形式):
K = P_pred / (P_pred + R)其中 P_pred 是预测的不确定性(先验协方差),R 是测量的不确定性(测量噪声协方差)。你看这个形式,本质上就是一个加权系数:当 P_pred 远大于 R(预测很不准,测量很准),K 趋近于1,几乎完全信测量;当 R 远大于 P_pred(测量很烂,预测很准),K 趋近于0,几乎完全信预测。
多维情况下,K 变成一个矩阵,形式是K = P_pred * H^T * (H * P_pred * H^T + R)^(-1),其中 H 是测量矩阵,负责把状态空间映射到测量空间。那个求逆的括号项叫新息协方差,它衡量的是“残差本身应该有多大波动”。矩阵求逆是卡尔曼滤波计算量最大的地方,测量维度越高越费劲,这也是为什么工程上要尽量精简测量维度。
注意:卡尔曼增益不是常数,它每个时刻都在变。系统刚启动时协方差大,增益也大,收敛快;运行稳定后协方差小,增益也小,输出平滑。这种自适应特性是卡尔曼滤波相比固定系数滤波器的核心优势。
2.3 协方差矩阵:不确定性的“形状”描述
协方差矩阵是很多人觉得抽象的地方。一维情况下它就是一个数(方差),表示估计的离散程度。多维情况下它是一个对称矩阵,对角线是各分量的方差,非对角线是分量之间的相关性。
为什么需要非对角线元素?因为状态分量之间往往不独立。比如你估计小车的位置和速度,位置估计偏大的时候,速度估计往往也偏大(因为你是根据速度积分得到位置的),两者正相关。协方差矩阵的非对角线元素就记录这种相关性。在预测和更新过程中,这些相关性会被正确地传播和修正,这也是卡尔曼滤波比“对每个变量单独滤波”更强大的原因——它利用了变量之间的耦合信息。
实际写代码时,协方差矩阵的初始值 P0 怎么设是个经验活。设太小,滤波器会过度自信,对测量不敏感,收敛慢;设太大,初期波动大,但收敛快。我的习惯是初期给一个偏大的对角阵(比如状态量程的平方),让滤波器快速收敛,跑一段时间后它自己会降下来。
3. 全过程推导:五个公式是怎么来的
3.1 状态空间模型:先把问题写成数学
推导之前,得先把系统用状态空间形式描述出来。离散线性系统的标准写法是:
x_k = F * x_{k-1} + B * u_{k-1} + w_{k-1} (状态方程) z_k = H * x_k + v_k (测量方程)x 是状态向量(我们想估计的东西),F 是状态转移矩阵(描述状态怎么演化),B 和 u 是控制输入(没有控制就省略),w 是过程噪声,服从均值为0、协方差为 Q 的高斯分布。z 是测量向量,H 是测量矩阵,v 是测量噪声,服从均值为0、协方差为 R 的高斯分布。
这里的关键假设是 w 和 v 相互独立,且都是高斯白噪声。有了这个模型,卡尔曼滤波要做的就是:给定 z_1 到 z_k 的测量序列,递推地求出 x_k 的最优估计。所谓“最优”,在最小均方误差意义下,就是求条件期望 E[x_k | z_1...z_k]。
3.2 预测两步:均值往前推,协方差往大长
预测阶段有两个公式。第一个是状态预测:
x_pred = F * x_prev + B * u这个很直观,就是把上一时刻的估计代入状态方程,去掉噪声项(因为噪声均值为0,期望就是0)。第二个是协方差预测:
P_pred = F * P_prev * F^T + Q这个公式值得细说。F * P_prev * F^T 描述的是“上一时刻的不确定性经过状态转移后变成什么样”。矩阵 F 会把不确定性拉伸、旋转、耦合,所以左右各乘一个 F。然后加上 Q,因为过程噪声又引入了新的不确定性。你可以理解为:预测让不确定性只增不减(Q 是半正定的),这符合直觉——没有新信息进来,你只会越来越没底。
实操心得:Q 的设置是卡尔曼滤波调参的重头戏。Q 设太小,滤波器太信任模型,测量来了也不怎么修正,响应迟钝;Q 设太大,滤波器太不信任模型,输出跟着测量抖。我的经验是从“模型误差的物理量级”出发估算,比如速度模型误差大概0.1米每秒,那 Q 里对应项就取0.01量级,然后根据实际效果微调。
3.3 更新三步:算增益、修均值、缩协方差
更新阶段有三个公式。第一步算卡尔曼增益:
K = P_pred * H^T * (H * P_pred * H^T + R)^(-1)前面已经解释过它的加权含义。第二步修正状态估计:
x_new = x_pred + K * (z - H * x_pred)括号里的z - H * x_pred就是残差(新息),表示测量值和预测测量值之间的差。K 乘以残差就是修正量。如果残差为0,说明测量和预测完全一致,不用修正。第三步修正协方差:
P_new = (I - K * H) * P_pred这个公式表示融合测量后不确定性减小。I 是单位阵。有个细节:这个形式在数值计算上不是最优的,因为它可能因为舍入误差导致 P 失去对称正定性。工程上常用Joseph形式:P_new = (I - K*H) * P_pred * (I - K*H)^T + K * R * K^T,虽然计算量大一点,但数值稳定性好得多。我在长时间运行的嵌入式系统里一律用Joseph形式,避免协方差矩阵“跑飞”。
3.4 从贝叶斯角度理解:为什么这五个公式是“最优”的
如果只记公式,很容易忘。从贝叶斯估计的角度看,卡尔曼滤波就是在做递归贝叶斯更新。预测阶段是在算先验分布 p(x_k | z_{1:k-1}),更新阶段是用贝叶斯公式算后验分布 p(x_k | z_{1:k})。在线性高斯假设下,先验和后验都是高斯分布,而高斯分布只需要均值和协方差两个参数就能完全描述,所以整个贝叶斯更新就退化成了均值和协方差的递推——这正是那五个公式在做的事。
从这个角度看,卡尔曼增益的表达式其实就是两个高斯分布相乘后新分布的均值和方差公式的矩阵推广。一维高斯相乘:两个均值按方差倒数加权平均,新方差是两者倒数和再取倒数。你把一维的加权系数写出来,和 K = P/(P+R) 对比一下,会发现完全一致。理解了这一点,那五个公式就不再是“天降公式”,而是高斯分布性质的必然结果。
4. 实操落地:从公式到可运行的代码
4.1 一维温度估计的完整实现
理论说再多不如跑一遍。我用一个最简单的例子:估计一个缓慢变化的温度。状态就是一维温度 x,状态转移矩阵 F=1(假设温度不变),测量矩阵 H=1(温度计直接测温度)。过程噪声 Q 取0.01(温度缓慢漂移),测量噪声 R 取0.25(温度计误差约0.5度,方差0.25)。
import numpy as np # 初始化 x = 24.0 # 初始估计 P = 1.0 # 初始不确定性 F = 1.0 H = 1.0 Q = 0.01 R = 0.25 measurements = [25.1, 24.8, 25.3, 24.9, 25.0, 25.2] for z in measurements: # 预测 x_pred = F * x P_pred = F * P * F + Q # 更新 K = P_pred * H / (H * P_pred * H + R) x = x_pred + K * (z - H * x_pred) P = (1 - K * H) * P_pred print(f"测量={z:.2f} 估计={x:.3f} 增益={K:.3f} 方差={P:.4f}")跑出来的结果会显示:第一次测量时增益约0.8(因为初始P大,更信测量),估计快速靠近测量值;随着P减小,增益逐渐降到0.3左右,估计变得平滑,不再跟着每次测量跳动。这就是卡尔曼滤波“先快后稳”的典型行为。
4.2 多维场景:二维位置速度跟踪
一维太简单,实际项目基本都是多维。假设我们跟踪一个一维运动的小车,状态是 [位置, 速度],测量只有位置。状态转移矩阵 F = [[1, dt], [0, 1]](dt是采样周期),测量矩阵 H = [[1, 0]]。
dt = 0.1 F = np.array([[1, dt], [0, 1]]) H = np.array([[1, 0]]) Q = np.array([[0.01, 0], [0, 0.01]]) # 过程噪声 R = np.array([[0.5]]) # 测量噪声 x = np.array([[0.0], [1.0]]) # 初始位置0,速度1 P = np.eye(2) * 10 # 初始不确定性给大一点 def kalman_step(x, P, z): # 预测 x_pred = F @ x P_pred = F @ P @ F.T + Q # 更新 S = H @ P_pred @ H.T + R K = P_pred @ H.T @ np.linalg.inv(S) x_new = x_pred + K @ (z - H @ x_pred) P_new = (np.eye(2) - K @ H) @ P_pred return x_new, P_new这个例子里有个值得注意的点:虽然我们只测位置,但速度也被间接估计出来了,因为协方差矩阵的非对角线元素把位置和速度关联起来了。测量修正位置的时候,速度估计也跟着被修正。这就是多维卡尔曼滤波“利用耦合信息”的威力,也是它比“对位置和速度分别滤波”强的地方。
4.3 参数整定的实战经验
Q 和 R 的整定是绕不过去的坎。R 相对好办,它直接对应传感器噪声,查数据手册或者做静态测试统计方差就能得到。Q 就麻烦多了,因为它对应的是“模型有多不准”,而这个往往没有现成数据。
我的做法分三步。第一步,先按物理直觉给 Q 一个初值,比如位置模型误差0.1米,那 Q 的位置项就取0.01。第二步,跑实际数据,看残差(新息)的统计特性。如果残差持续偏大且同号,说明模型有系统偏差,Q 给太小了;如果残差剧烈跳动,说明 Q 给太大了。第三步,微调。Q 增大,滤波器响应快但输出抖;Q 减小,输出平滑但可能跟不上真实变化。这个权衡没有标准答案,取决于你对“平滑”和“响应”的优先级。
提示:有个实用技巧是监控“归一化新息平方”(NIS),它等于 残差^T * S^(-1) * 残差,理论上应该服从卡方分布。如果 NIS 长期偏大,说明滤波器过度自信(P 太小或 R 太小),需要调大;长期偏小说明滤波器不够自信。这个方法比肉眼调参靠谱得多。
5. 常见问题与排查实录
5.1 滤波器发散:估计越跑越偏
发散是卡尔曼滤波最经典的故障。表现是估计值逐渐偏离真实值,协方差却不增大,滤波器“死不认错”。根本原因通常是模型和实际系统不匹配,但 Q 给得太小,导致滤波器过度信任错误模型,测量修正量被压得很小。
排查思路:先看残差。如果残差持续偏大且方向一致,基本可以确定是模型偏差。解决办法有两个:一是增大 Q,让滤波器对模型不那么自信;二是检查状态转移矩阵 F 是否写错了(比如符号错误、dt 单位错误)。我遇到过一次 F 里 dt 用了毫秒但实际数据是秒,结果预测速度大了1000倍,滤波器直接飞了。这种低级错误反而最常见,排查时先检查单位。
5.2 协方差矩阵失去正定性
理论上 P 应该始终对称正定,但浮点运算的舍入误差可能破坏这个性质,导致求逆失败或者增益计算出 NaN。长时间运行的嵌入式系统尤其容易出这个问题。
解决办法前面提过:用 Joseph 形式更新协方差。另外可以定期做对称化处理P = (P + P.T) / 2,强制对称。如果条件允许,用平方根形式的卡尔曼滤波(SRKF),它直接传播协方差的平方根,天然保证正定性,数值稳定性最好,代价是计算稍复杂。
5.3 测量延迟与不同步
实际系统里,不同传感器的数据往往不是同一时刻到达的。比如摄像头30帧每秒,IMU 200赫兹,你拿不同时刻的测量去更新同一个状态,就会有时间错位。处理办法有两种:一是做测量缓冲和状态回滚,把状态退回到测量时刻再更新,然后重新预测到当前时刻;二是做延迟补偿,用状态转移矩阵把测量值外推到当前时刻。前者精确但实现复杂,后者简单但有近似误差。小延迟用后者够了,大延迟必须用前者。
| 问题现象 | 可能原因 | 排查方法 | 解决手段 |
|---|---|---|---|
| 估计持续偏离 | 模型偏差、Q过小 | 看残差是否同号偏大 | 增大Q、检查F矩阵 |
| 输出剧烈抖动 | Q过大、R过小 | 看增益是否接近1 | 减小Q、增大R |
| 收敛极慢 | P0过小、R过大 | 看初始增益是否过小 | 增大P0、减小R |
| 数值NaN | 协方差失去正定 | 检查P特征值 | Joseph形式、对称化 |
| 跟不上突变 | Q过小 | 看残差突变后恢复慢 | 增大Q或加自适应机制 |
5.4 非线性系统的坑
标准卡尔曼滤波只适用于线性系统。如果你直接把它套到非线性系统上(比如用角度做状态但测量是sin值),轻则精度下降,重则发散。这时候必须上EKF或UKF。EKF的关键是求雅可比矩阵,手推容易出错,建议用自动微分工具或者数值差分验证。UKF没有雅可比,但sigma点的参数(alpha、beta、kappa)需要调,alpha通常取1e-3到1之间,太小会有数值问题。
我在一个姿态解算项目里踩过的坑是:EKF的线性化点选在预测值还是更新值,结果差别不小。标准EKF在预测值处线性化,但如果预测误差大,线性化误差也大。有些工程实现会在更新后重新线性化一次(迭代EKF),精度更高但计算量翻倍。这个取舍要看你的实时性要求。
6. 写在最后的几点个人体会
卡尔曼滤波这东西,看第一遍觉得公式多,看第二遍觉得道理顺,真正用起来才发现坑都在细节里。我个人最大的体会是:不要把它当黑盒。很多人调包调习惯了,Q 和 R 随便填,跑出来效果不好就换算法,其实问题往往出在对状态模型的理解上。你得清楚每个状态分量的物理意义,清楚每个噪声项对应什么实际误差源,调参才有方向。
另外一个体会是,先在一维上验证,再上多维。一维的卡尔曼滤波可以用纸笔算,每一步的增益、方差变化都能肉眼跟踪,理解透了再推广到矩阵形式,会顺畅很多。我见过不少朋友一上来就啃多维推导,被矩阵求逆和协方差传播绕晕,反而丢了直觉。
最后分享一个我常用的小技巧:在调试阶段,把每一步的预测值、测量值、增益、残差都打印出来,画成曲线。很多时候滤波器行为异常,看一眼增益曲线就知道问题在哪——增益长期接近0说明过度自信,长期接近1说明过度不自信。这比盯着最终输出猜原因高效得多。卡尔曼滤波的五个公式本身不难,难的是建立“不确定性如何流动”的直觉,而这个直觉只能靠一次次实操喂出来。