信息滤波Information Filter这个名字第一次撞见时我以为又是什么新造的概念翻了两页才反应过来——这不就是把卡尔曼滤波整个翻过来写吗。可真动手把推导走完一遍尤其是把它塞进一个多传感器融合的项目里跑起来之后我才明白它换的不只是写法而是整套问题的参照系。标准卡尔曼滤波维护的是状态均值 x̂ 和协方差 P信息滤波维护的则是信息矩阵 Y P⁻¹ 和信息向量 ŷ Y x̂。看上去只是取了个逆代价却完全不对称预测步的公式会变得啰嗦乃至有点丑更新步却简洁到让人怀疑是不是漏了什么——信息形式下一个新观测进来无非是往信息矩阵里加一项 HᵀR⁻¹H往信息向量里加一项 HR⁻¹z没有矩阵求逆没有增益计算。这篇东西我打算从头写。先把信息矩阵和信息向量的定义钉死再说清楚它为什么天然适配高斯分布接着把预测步一步步推完中间用到的矩阵求逆引理会单独拎出来讲来路免得读者看到那个式子直接从椅子上站起来然后是最漂亮的更新步推导以及它与标准卡尔曼滤波同解的证明最后给一份能直接跑的 Python 对拍代码让信息滤波和标准 KF 在同一组数据上跑看两者输出差多少。如果你在做多传感器融合、组合导航、SLAM 后端或者分布式状态估计又或者只是被信息矩阵那套符号绕晕过这篇应该能把中间断掉的推导链条补上。1. 先把问题摆清楚信息形式到底换掉了什么1.1 从一次多传感器融合的返工说起早些年我接过一个多源融合的活儿三个传感器往一个状态估计器里塞数据一个高频低精度的位置源一个低频高精度的位置源还有一个偶尔才吐一次数据的姿态源。方案一开始用的是标准卡尔曼滤波写起来也顺跑起来也对。问题出在扩展上——每接入一个新的传感器就得重排观测矩阵 H、重排噪声阵 R、重排观测向量 z然后重算 S HP⁻Hᵀ R 并且求它的逆。三个传感器还算能忍等到第五个、第六个代码里全是拼接和索引映射改一次错一次。更难受的是采样率不一样低频的那个传感器大部分时间没有数据只能写一堆分支去跳步。后来我把整个估计器改写成信息形式那块拼接代码直接消失了。每个传感器的贡献就是一个独立的加项HᵢᵀRᵢ⁻¹Hᵢ 和 HᵢᵀRᵢ⁻¹zᵢ。传感器有数据就加没数据就跳过不需要为了对齐维度去拼大矩阵。这个改动让我第一次真切感受到信息形式的价值不在数学等价性上而在工程结构上——它把融合这件事从拼装变成了累加。1.2 信息矩阵与信息向量的定义定义本身只有两行但值得写清楚每一项的来路。给定一个高斯分布 x ~ N(x̂, P)并且 P 可逆那么信息矩阵 Y P⁻¹信息向量 Y x̂ P⁻¹x̂把高斯的指数项展开一切就都顺了(x - x̂)P⁻¹(x - x̂) xᵀP⁻¹x - 2x̂ᵀP⁻¹x x̂P⁻¹x̂ xYx - 2ᵀx 常数也就是说一个高斯分布在信息形式下就是 exp(-½xᵀYx ŷᵀx)指数上只剩一个二次型和一个一次项前面的常数不带 x可以不管。这里 Y 必须是半正定乃至正定的对称矩阵这一点后面做数值实现时是个反复出现的坑。Y 的对角元越大表示该分量越确定非对角元则描述两个状态分量之间的耦合强度。默认情况下 Y 的对角占优意味着状态之间耦合不强。这套参数在统计学里叫规范参数或自然参数在估计理论里叫信息形式叫法不同指的是同一件事。1.3 谁适合用这套东西谁不适合先说适合的。第一类是多传感器异步融合这前面已经说过了更新步的加法结构让每个传感器可以本地计算贡献再汇总天然可分布。第二类是初始状态完全未知的场景——标准 KF 要初始化一个庞大的 P₀ 表示我什么都不知道在数值上万般别扭信息形式直接令 Y₀ 0语义干净只要第一个观测进来信息就自然注入。第三类是大维度稀疏问题比如位姿图或者 SLAM 后端信息矩阵的非零结构直接对应变量之间的条件独立关系稀疏求解器用起来非常舒服。再说适合性差一点的地方。如果过程噪声协方差 Q 是奇异的或者系统维度不高但对预测步的求逆特别敏感那信息形式的预测步会成为负担。标准 KF 的预测是纯粹的矩阵乘法 P⁻ FPFᵀ Q不涉及任何求逆信息形式的预测则必须处理一个形如 (Y FᵀMF)⁻¹ 的东西Q 奇异的时候连 M Q⁻¹ 都不存在。所以这不是万能药选之前得看你的瓶颈落在预测还是更新上。2. 推导前必须备好的三件工具2.1 矩阵求逆引理那个让推导能走下去的恒等式整个预测步推导就靠一个恒等式撑着必须先写清楚。对任意合适维度的矩阵 A、B、C、D只要 A 和 C 可逆有(A BCD)⁻¹ A⁻¹ - A⁻¹B(C⁻¹ DA⁻¹B)⁻¹DA⁻¹这个式子在不同教材里名字不太一样叫矩阵求逆引理、Woodbury 恒等式、Sherman-Morrison-Woodbury 公式的都有Sherman-Morrison 是它在一秩更新情形下的特例。它的意义在于把一个先加后求逆的操作换成先求逆再加代价是引入一个新的求逆项但那个求逆的维度通常小得多或者结构好得多。验算一下也不难把右边乘上 (A BCD)逐项展开中间那些 A⁻¹A、DA⁻¹等人为构造出来的项会层层约掉最后剩单位阵。第一次看会觉得这个式子来路不明我的建议是别纠结它怎么被想出来的直接当成工具用等推导做完回头看你会发现它出现的位置恰到好处——它正好把协方差空间里的加法翻译成了信息空间里的减法。2.2 对称正定矩阵与 Cholesky数值实现的地基后面无论推导还是写代码反复要处理形如 Mx b 的线性系统其中 M 对称正定。这时候绝对不要写 np.linalg.inv(M)这是新手最容易踩的坑。正确做法是用 Cholesky 分解或者带主元的求解器理由有三条。其一数值稳定性。显式求逆的误差会随着条件数放大而 Cholesky 分解求解的误差增长慢得多尤其是矩阵接近奇异的时候差别能到几个数量级。其二对称性保持。矩阵求逆的浮点结果往往不再严格对称而 Cholesky 只用到下三角部分输出天然保持结构。其三效率。求逆是 O(n³)Cholesky 也是 O(n³)但常数小一半左右而且分解结果可以复用于多次右端项——在滤波里每个时间步都会遇到同一个 S 配不同右端项复用价值极高。具体到 Python我一般用 scipy.linalg.cho_factor 配合 cho_solve比 numpy 的通用 solve 更快也更稳。矩阵规模小的时候差异不明显一旦状态维度上到几十上百这两者的差距会让你重新审视之前的代码。2.3 稀疏性的来源为什么信息矩阵常常是空的标准 KF 的协方差矩阵 P 通常是稠密的即便初始时稀疏经过几次 FPFᵀ 传播之后非零元会迅速填满。信息矩阵 Y 的行为正好相反它的非零结构直接对应哪些状态之间存在直接约束。举个例子。假设有三个状态变量 a、b、c其中 a 和 b 被一条观测约束在一起b 和 c 被另一条约束联系在一起但 a 和 c 之间没有任何直接观测。在协方差空间里a 和 c 因为共享 b 而变得相关P 里 (a,c) 位置一般非零而在信息空间里Y 的 (a,c) 位置严格为零因为两者之间没有直接的边。换句话说Y 的非零结构就是一张条件独立图图上的边对应非零元。这个性质在 SLAM 后端里被用到了极致。整条轨迹的信息矩阵呈现带状的稀疏结构几千个位姿变量的信息矩阵非零元可能只有百分之几稀疏 Cholesky 分解一顿操作比稠密矩阵快上好几个数量级。这也是为什么大规模位姿图优化几乎清一色用信息形式而非协方差形式。3. 预测步推导从协方差形式硬搬到信息形式3.1 朴素替换会撞上两次求逆线性系统的预测模型是 x_k F x_{k-1} ww ~ N(0, Q)。在标准 KF 里预测一步就是x⁻ F x̂P⁻ F P Fᵀ Q现在我想把 P⁻ 换成信息矩阵 Y⁻ (P⁻)⁻¹。注意 P Y⁻¹代进去Y⁻ (F Y⁻¹ Fᵀ Q)⁻¹这就是最朴素的替换结果。它没错但直接这么算的话每一步要先从 Y 求一次逆得到 P做完两次矩阵乘法再加 Q最后再求一次逆回到 Y。两次求逆而且求逆的矩阵都是稠密的完全丧失了信息形式的任何优势。这显然不是想要的结果。好在有求逆引理可以把这两次求逆压成一次并且把求逆对象换成一个结构更友好的矩阵。3.2 用求逆引理把求逆维度压到状态维度把 Y⁻ (Q F Y⁻¹ Fᵀ)⁻¹ 拿出来对照引理的形式 (A BCD)⁻¹做如下对应A QB FC Y⁻¹D Fᵀ逐项代入引理A⁻¹ Q⁻¹记为 MA⁻¹B MFC⁻¹ YDA⁻¹B FᵀMFDA⁻¹ FᵀM。于是Y⁻ M - MF(Y FᵀMF)⁻¹FM这一步是整个信息滤波预测的核心公式。它值得慢慢看右边的 M 是过程噪声的逆也就是过程噪声能提供的信息量中间那个求逆项里Y 是上一时刻的信息矩阵FᵀMF 可以理解成过程噪声通过状态转移被折射回来的量。两项做减法效果是过程噪声在稀释已有信息——这符合直觉预测永远让不确定性变大信息量必然下降。和朴素形式比好处是三重的。第一只求一次逆而且求逆的矩阵是 Y FᵀMF它的规模就是状态维度 n不是观测维度。第二Y 本身如果是稀疏的FᵀMF 在 Q 为对角阵时也是稀疏的两者相加仍然稀疨稀疏求解器可以全程工作。第三M Q⁻¹ 可以离线算一次反复用如果 Q 是分块对角或者对角M 也是同样结构求起来几乎不花时间。3.3 信息向量预测项的化简有了 Y⁻信息向量是 ŷ⁻ Y⁻x⁻。这里 x⁻ F x̂_{k-1}而 x̂_{k-1} Y⁻¹所以ŷ⁻ Y⁻ F Y⁻¹ ŷ把上一步得到的 Y⁻ 代进去并令 S Y FᵀMFY⁻FY⁻¹ [M - MF S⁻¹ FᵀM] F Y⁻¹ MF Y⁻¹ - MF S⁻¹ (FᵀMF) Y⁻¹关键在第二项。注意 S Y FMF所以 FMF S - Y那么S⁻¹(FᵀMF)Y⁻¹ S⁻¹(S - Y)Y⁻¹ Y⁻¹ - S⁻¹代回Y⁻FY⁻¹ MF Y⁻¹ - MF(Y⁻¹ - S⁻¹) MF S⁻¹所以信息向量的预测式是ŷ⁻ MF(Y FᵀMF)⁻¹ŷ这个结果干净得有点意外。整条式子里只用到一次对 S 的求逆而且这个 S 和上一节算 Y⁻ 时用到的完全是同一个矩阵可以复用分解结果——一次 Cholesky 分解两个右端项成本几乎不增加。从直觉上说预测步在信息形式下是一次线性变换ŷ 被 F 平移了一次又被 M 加权了一次最后被 S⁻¹ 做了个归一化。它不像更新步那样有明确的物理意义这也是信息形式在预测上不如更新上好看的原因。3.4 预测步的完整算法与复杂度把上面两节拼起来信息滤波的预测步是计算 M Q⁻¹若 Q 恒定离线算一次即可组装 S Y_{k-1} FᵀMF对 S 做 Cholesky 分解Y⁻ M - MF S⁻¹ FᵀMŷ⁻ MF S⁻¹ ŷ_{k-1}强制对称化 Y⁻ ← (Y⁻ Y⁻ᵀ)/2复杂度上组装 FᵀMF 是 O(n³)S 的分解也是 O(n³)看起来和标准 KF 的 FPFᵀ Q 差不多。真正的差别在稳定性与稀疏潜力上稀疏条件下上面这些乘法都可以用稀疏矩阵运算而标准 KF 里的 P 很快就稠密了稀疏性用不上。还有一点容易忽略如果 Q 是奇异的M 不存在这套推导直接失效。工程上的处理办法是往 Q 上加一个极小量比如 1e-9 乘以单位阵让它可逆代价是牺牲一点点精度。这个操作在标准 KF 里完全没必要是信息形式特有的功课。4. 更新步推导为什么说信息形式天生适合多源融合4.1 观测更新公式的两行推导更新模型是 z Hx vv ~ N(0, R)。标准 KF 的后验协方差可以写成P P⁻ - P⁻Hᵀ(HP⁻Hᵀ R)⁻¹HP⁻这是卡尔曼滤波的经典结果形式就是先验减掉增益带来的修正。现在对两边取逆用求逆引理处理。对照引理 (A BCD)⁻¹令A P⁻B HC R⁻¹D H那么 A⁻¹ (P⁻)⁻¹ Y⁻C⁻¹ RDA⁻¹B HP⁻HᵀDA⁻¹ HP⁻。直接代入(P⁻ HR⁻¹H)⁻¹ P⁻ - P⁻Hᵀ(R HP⁻Hᵀ)⁻¹HP⁻右边正好就是标准 KF 的后验协方差 P。左边取逆之后是后验信息矩阵也就是 Y。于是立刻得到Y Y⁻ HᵀR⁻¹H一行没有任何附加项。信息向量的推导稍微长一点但同样干净。后验均值是 x̂ x⁻ P⁻Hᵀ(HP⁻Hᵀ R)⁻¹(z - Hx⁻)两边乘 Yŷ Yx̂ (Y⁻ HᵀR⁻¹H)[x⁻ P⁻Hᵀ(HP⁻HᵀR)⁻¹(z - Hx⁻)]把括号展开第一项 (Y⁻ HR⁻¹H)x⁻ ŷ⁻ HᵀR⁻¹Hx⁻。第二项写成 (P⁻)⁻¹P⁻ HᵀR⁻¹HP⁻ 的形式也就是 I HᵀR⁻¹HP⁻再乘上 Hᵀ(HP⁻HᵀR)⁻¹整理得(I HᵀR⁻¹HP⁻)Hᵀ(HP⁻HᵀR)⁻¹(z - Hx⁻) HᵀR⁻¹(R HP⁻Hᵀ)(HP⁻HᵀR)⁻¹(z - Hx⁻) HR⁻¹(z - Hx⁻)两半加起来ŷ ⁻ HR⁻¹Hx⁻ HR⁻¹z - HᵀR⁻¹Hx⁻ ŷ⁻ HR⁻¹z中间那两项精确抵消一点不剩。最终更新式只有两个Y Y⁻ HᵀR⁻¹H ŷ⁻ HᵀR⁻¹z4.2 加法结构带来的工程红利这两个式子没有任何求逆没有增益没有残差乘以增益的步骤。新观测进来就是往信息矩阵和信息向量里各加一个外积项。这个性质在工程上能换来好几样东西。第一样是传感器即插即用。前面提过HᵢᵀR⁻¹H 只和该传感器自己的观测模型与噪声有关跟其他传感器完全解耦。三个传感器融合就是把三个这样的项加起来不需要构造大的拼接矩阵也不需要为了保证维度对齐去写映射代码。低频传感器没数据的时候跳过它的加项就行代码里没有分支。第二样是可分布式计算。假设有若干节点各自持有一部分观测每个节点本地算自己的 HᵢᵀR⁻¹H 和 HᵀR⁻¹z然后把这两块小矩阵发到中心节点汇总中心节点做的只是求和。如果想做无中心的方案求和本身是个交换律成立的运算用一致性协议迭代几轮就能收敛这也是分布式卡尔曼滤波在信息形式下更容易设计的原因。第三样是缺失观测的天然处理。标准 KF 里如果某个时刻没有观测P 更新这一步直接跳过逻辑上没问题但信息形式下即使跳过信息矩阵也会因为预测步的减法而变小物理意义清楚没有观测信息就按模型自然衰减。4.3 与标准卡尔曼滤波的解等价性一个躲不开的问题这两套东西算出来真的一样吗答案是严格一样只要 P 可逆、数值上没有出岔子二者给出的状态估计逐位相等浮点误差范围内。原因是它们描述的是同一个贝叶斯后验只是参数化方式不同。验证的方式很直接跑标准 KF 得到 x̂_kf 和 P_kf跑信息滤波得到 ŷ 和 Y然后比较 x̂_kf 与 Y⁻¹ŷ、P_kf 与 Y⁻¹。如果实现正确差异应该落在 1e-10 量级以下。这个对拍是检验实现正确性的最好方法比盯着公式看半天管用。下面的代码部分我会把完整的对拍脚本给出来跑一遍就知道有没有写错。需要提醒的是等价性只在 P 和 Y 都可逆的严格条件下成立。如果 P₀ 是无穷大表示完全未知标准 KF 根本没法初始化而信息形式可以令 Y₀ 0此时两套东西形式上不等价但信息形式是良定义的标准 KF 不是——这本身就是信息形式的一个优势。5. 代码实测手撸信息滤波并与 KF 对拍5.1 核心实现下面是一份可以直接跑的实现。模型是二维平面的匀速运动模型状态是 [x, vx, y, vy]过程噪声用离散白噪声加速度模型生成观测是两个位置分量。为节省篇幅矩阵构造部分写得紧凑一些。import numpy as np from scipy.linalg import cho_factor, cho_solve def build_model(dt1.0, sigma_a0.05, sigma_z0.3): F np.array([[1, dt, 0, 0], [0, 1, 0, 0], [0, 0, 1, dt], [0, 0, 0, 1]], dtypefloat) G np.array([[0.5 * dt**2, 0], [dt, 0], [0, 0.5 * dt**2], [0, dt]], dtypefloat) Q G np.diag([sigma_a**2, sigma_a**2]) G.T # 给 Q 加一点点抖动保证可逆 Q Q 1e-9 * np.eye(4) H np.array([[1, 0, 0, 0], [0, 0, 1, 0]], dtypefloat) R np.diag([sigma_z**2, sigma_z**2]) return F, Q, H, R def kf_predict(x, P, F, Q): return F x, F P F.T Q def kf_update(x, P, z, H, R): S H P H.T R K np.linalg.solve(S, H P).T # 等价于 P H^T S^{-1}但更稳 x_new x K (z - H x) P_new P - K H P P_new 0.5 * (P_new P_new.T) # 强制对称 return x_new, P_new def if_predict(y_info, Y, F, Q): M np.linalg.inv(Q) S Y F.T M F c cho_factor(S, lowerTrue) Y_pred M - M F cho_solve(c, F.T M) y_pred M F cho_solve(c, y_info) Y_pred 0.5 * (Y_pred Y_pred.T) return y_pred, Y_pred def if_update(y_pred, Y_pred, z, H, R): Rinv np.linalg.inv(R) Y_upd Y_pred H.T Rinv H y_upd y_pred H.T Rinv z Y_upd 0.5 * (Y_upd Y_upd.T) return y_upd, Y_upd有两个细节值得点出来。kf_update 里我用 np.linalg.solve(S, H P).T 而不是显式写 P H.T inv(S)因为前者先解线性系统再转置数值上更稳信息滤波里所有求逆都用 cho_solve 走 Cholesky理由在 2.2 节说过。另外两处都做了强制对称化这一步在长时间递推里几乎是必需的否则浮点误差累积会让矩阵慢慢失去对称性最终导致 Cholesky 分解失败。5.2 对拍脚本与结果解读对拍脚本的思路是先用 KF 跑一遍记录轨迹再用 IF 跑一遍逐时刻比较。初始条件取同一个 P₀ 和 x₀IF 侧用 Y₀ inv(P₀)、₀ Y₀ x₀。def run_compare(T200, seed7): rng np.random.default_rng(seed) F, Q, H, R build_model() # 真值轨迹 x_true np.array([0.0, 1.0, 0.0, 0.5]) xs, zs [], [] for _ in range(T): x_true F x_true rng.multivariate_normal(np.zeros(4), Q) xs.append(x_true.copy()) zs.append(H x_true rng.multivariate_normal(np.zeros(2), R)) xs np.array(xs); zs np.array(zs) # 初值 x0 np.array([0.5, 0.8, -0.3, 0.4]) P0 np.diag([1.0, 1.0, 1.0, 1.0]) # 标准 KF x, P x0.copy(), P0.copy() kf_log [] for z in zs: x, P kf_predict(x, P, F, Q) x, P kf_update(x, P, z, H, R) kf_log.append((x.copy(), P.copy())) # 信息滤波 Y np.linalg.inv(P0) y_info Y x0 if_log [] for z in zs: y_info, Y if_predict(y_info, Y, F, Q) y_info, Y if_update(y_info, Y, z, H, R) x_if np.linalg.solve(Y, y_info) P_if np.linalg.inv(Y) if_log.append((x_if, P_if)) # 对比 dx, dP [], [] for (xk, Pk), (xi, Pi) in zip(kf_log, if_log): dx.append(np.max(np.abs(xk - xi))) dP.append(np.max(np.abs(Pk - Pi))) print(max |dx| , max(dx)) print(max |dP| , max(dP)) if __name__ __main__: run_compare()跑出来 max |dx| 和 max |dP| 通常落在 1e-11 到 1e-9 之间具体取决于随机种子和矩阵条件数。如果结果在 1e-6 以上基本可以断定实现里有 bug最常见的是 Q 的抖动加得太小导致 M 数值病态或者对称化那一步漏掉了。如果直接发散或者出 NaN优先检查预测步里 S 是否真的对称正定以及有没有误把 Fᵀ 写成 F。5.3 条件数、稀疏性与耗时的实测观察对拍通过之后我习惯再观察三个量。第一个是 Y 的迹随时间的走势。在稳定系统里Y 的迹会在几十步之后进入稳态平台——预测步把信息往下压更新步往上抬两者达到平衡。如果迹一直单调增长说明观测噪声设得太小或者过程噪声设得太小信息无限累积数值上迟早出问题如果迹一直往下掉说明观测根本没起作用。第二个是条件数。np.linalg.cond(Y) 在跑了几百步之后通常能到 1e6 到 1e8 量级尤其是速度和位置同时存在、且速度的观测间接依赖位置的时候。cond(P) 和 cond(Y) 数值上互为倒数关系严格说 cond(A⁻¹) cond(A)量级一致所以这不是信息形式独有的问题但信息形式里求逆操作更多受影响更明显。缓解办法有两个用双精度全程别偷懒用 float32以及在信息矩阵上做定期的对称化和数值清洗。第三个是耗时结构。单步耗时上KF 和信息滤波在四维状态这个规模下几乎看不出差别因为矩阵太小了Python 的解释开销盖过了一切。真正能看出差别的是把状态维度拉到几十上百之后KF 更新步要解一个 m×m 的系统m 是观测维度信息滤波更新步只有外积累加完全不含分解。多传感器场景下 m 会随传感器数量增长而信息滤波在这一项上是线性的这就是它的结构性优势所在。6. 踩过的坑与排查清单6.1 对称性与正定性丢失这是信息滤波最常见的故障没有之一。症状是跑了几百步之后 Cholesky 分解报错提示矩阵不是正定或者求解出来的结果莫名其妙地大。原因在于浮点运算不保证对称性Y_pred M - MF S⁻¹ FᵀM 这个式子理论上严格对称但计算机算出来 (1,3) 和 (3,1) 两个位置会差个 1e-15 量级。单次无所谓几千步累积下来不对称部分会被迭代放大最终破坏正定性。解决办法简单粗暴但有效每次构造完 Y 就强制对称化写成 0.5 * (Y Y.T)。这个操作代价几乎为零却能挡住绝大多数发散。更彻底的方案是用平方根信息滤波直接维护 Y 的 Cholesky 因子从结构上保证对称正定代价是公式更复杂代码量翻倍。我个人的取舍是中小规模用对称化就够了几百维以上的长期递推再考虑平方根形式。还有一个小坑是矩阵乘法顺序。np.linalg.solve 和 np.linalg.inv 在某些版本的底层实现里对非对称输入的处理不同导致即使输入只差 1e-16输出也可能差 1e-10。所以条件允许的话先对称化再求解别让求解器替你做决定。6.2 过程噪声协方差奇异Q 奇异意味着某些方向上的过程噪声为零这是很常见的建模选择——比如你认定某个状态例如传感器偏置在短期内完全不变化就会把对应的噪声方差设成 0。标准 KF 里这完全没问题P 照样递推。信息滤波里M Q⁻¹ 直接不存在整套预测公式崩掉。工程上的处理方法有几种按推荐顺序排。最省事的是加抖动Q ← Q εIε 取 1e-9 到 1e-12视状态量的量纲而定量纲差异大的时候按状态分块加不同的 ε别一刀切。稍微讲究一点的做法是对 Q 做特征分解只对零特征值方向注入极小噪声其余方向不动。最正统的做法是改用平方根信息滤波它对奇异的 Q 不敏感因为推导路径完全不同。这里有个经验判断如果你发现为了跑通信息滤波而不得不把 Q 的抖动加到 1e-6 以上那基本说明信息形式不适合这个模型老老实实回去用标准 KF或者在预测步用协方差、更新步用信息形式的混合方案别硬撑。6.3 常见问题速查表现象可能原因排查与处理Cholesky 报非正定长期递推丢失对称性每步强制对称化或改用平方根形式结果与 KF 差 1e-3 以上实现有误或矩阵病态检查 F 转置、Q 抖动、S 是否用同一个分解预测步耗时异常对 S 反复求逆未复用分解一次 cho_factor多个右端项复用状态估计缓慢漂移信息矩阵迹持续增长检查 R 是否过小、观测是否过于频繁初始几拍跳变剧烈Y₀ 过大导致初始信息过强初始 P₀ 设大一些或 Y₀ 直接设零迭代几百步后发散条件数过大全程双精度考虑正则化或换回 KF这张表是我自己踩坑攒出来的每次改代码前扫一眼能省下不少调试时间。还有一条表里没写如果你发现同一个算法在 MATLAB 里跑得好好的移植到 Python 就出问题先查默认数据类型。MATLAB 默认 double而 numpy 里如果不显式指定 dtype某些构造方式会给你 float32精度直接掉一半滤波器的敏感度对这一点非常敏感。7. 从线性到非线性EIF 与稀疏信息滤波7.1 扩展信息滤波的更新形式真实系统大多是非线性的观测写成 z h(x) v。把 h 在当前估计处做一阶泰勒展开得到雅可比矩阵 H ∂h/∂x剩下的推导和线性情形完全一样这就是扩展信息滤波EIF。更新式变成Y Y⁻ HᵀR⁻¹H ŷ⁻ HᵀR⁻¹[z - h(x⁻) Hx⁻]这里括号里的项是观测残差加上线性化修正比标准 EKF 的残差表达式多了一项。这个多出来的 Hx⁻ 是为了让整个式子仍然保持加法结构代价是必须显式地维持一个状态估计 x⁻ 用来做线性化——而纯信息形式本来是不需要维护 x 的这就造成了一点概念上的冗余。实际实现里通常是同时维护 ŷ、Y 和 x̂每次更新后用 x̂ Y⁻¹ŷ 重新同步一次。EIF 在非线性强的场景下不如 EKF 直观因为线性化点必须是一个明确的 x而在信息空间里这个 x 是隐式解出来的。所以业界主流仍是 EKF 或 UKFEIF 主要出现在需要对观测更新做加法式分解的特定场合比如多机器人协同定位。7.2 稀疏化与边缘化信息形式真正的主场信息形式最漂亮的应用不在滤波本身而在批量估计和因子图优化里。一个位姿图的观测约束每一条边贡献的就是一个形如 HᵀR⁻¹H 的加项把所有边加起来得到全局信息矩阵——这和滤波的更新步是一回事只不过这里的状态是整条轨迹的所有位姿维度可能上千。这样得到的信息矩阵呈现稀疏带状结构非零元只出现在有约束关系的变量之间。用稀疏 Cholesky 分解求解复杂度大致和变量数量的线性到一点几次方成正比而不是稠密求解的立方。几百个位姿的轨迹稠密求解就已经卡了稀疏求解还能流畅跑几千个点。更进一步的是边缘化。把旧的状态从问题里消掉比如滑窗优化里丢掉老的位姿在信息形式下就是对信息矩阵做 Schur 补消元后剩下的部分仍然是信息矩阵形式不变只是会引入一些新的非零元fill-in。这个过程在整个估计理论里是自洽的而在协方差空间里做同样的操作会碰到各种数值问题。所以滑窗式的视觉惯性里程计后端几乎清一色用信息形式或者因子图。7.3 什么时候别用信息滤波说几句反过来的话免得有人一头扎进去。如果你只做单传感器、状态维度在十以内、观测频率稳定信息滤波相对标准 KF 没有任何优势反而多了一堆矩阵求逆和对称性维护的麻烦。如果你需要频繁地查询某个子状态的不确定性比如可视化协方差椭圆协方差形式直接读对角块就行信息形式还得求一次逆。如果你的过程噪声协方差接近奇异又不想加抖动那边界会非常难受。我的判断标准很简单多传感器异步接入、状态维度大、信息矩阵天然稀疏、初始状态完全未知——这四条里占两条以上就值得用信息形式一条都不占用标准 KF 省心。别为了形式上的优雅去做工程上的额外功这个亏我吃过。说说我在实际项目里的体会。一开始我迷信两套东西数学等价随便选哪个都行后来发现等价只在纸面上成立落到代码里副作用差异巨大信息形式省掉了观测拼接和增益计算代价是多了一次预测步的求逆和一堆对称性维护。我最后选的是混合方案——滤波器主体用标准 KF只在需要做多源融合的那一层把观测贡献拆成 HᵀR⁻¹H 的累加形式再统一注入。这样既拿到了加法结构带来的扩展性又避开了预测步求逆的坑。另外提一句如果你的状态里有传感器偏置这类量纲和位置速度差好几个数量级的分量记得在组装 Q 和 R 之前做归一化不然信息矩阵的条件数会被最差的那个分量拖着走滤波跑不了多远就开始抖。