1. 从一个实际问题说起:为什么需要因子图
很多人学概率图模型的时候,顺序大概都是这样的:先看贝叶斯网络,觉得有向边挺直观;再看马尔可夫随机场,无向边表示相互影响也能理解;然后突然冒出来一个"因子图",课本上画了一堆小方块和圆圈,说是"二分图",把变量节点和因子节点分开,还讲起了和积算法、最大积算法,一下子就把人劝退了。
我第一次接触因子图是在做传感器融合相关的项目,当时的需求说起来很简单:手上有几个不同来源的观测,想融合成一个统一的状态估计。用卡尔曼滤波硬写也能写,但当观测来源变多、约束关系变复杂之后,代码里到处是矩阵拼接和下标的对应关系,改一处就错一处。后来换用因子图来建模,整个问题一下子清晰了——每个观测写成一条因子,每个未知量是一个变量节点,图一画出来,谁跟谁有关一目了然。
因子图本质上是一种概率图模型的统一表示形式。它做的事情,说穿了就是把一个复杂的联合概率分布拆成一堆局部的函数乘积,这些局部函数就叫"因子",每个因子只跟少数几个变量有关。拆开之后画出二分图,变量是圆圈,因子是方块,边表示这个因子依赖这个变量。有了这张图,很多复杂的推断问题就可以用统一的算法来处理,比如和积算法(sum-product)、最大积算法(max-product)等等。
这篇内容我打算按照我自己的学习路径和实操经验来写,不讲太多泛泛的教科书定义,重点讲清楚三件事:因子图到底在表达什么、怎么从一个实际问题建出因子图、以及在实际工程里怎么用好它。涉及到的实操部分我会尽量给出可以复现的步骤和代码骨架,也会把我在调试过程中踩过的坑和排查经验一起写进去。适合有概率基础、做过一点机器学习或机器人状态估计、想真正把因子图用起来的读者,纯小白也能看懂前面的建模部分,后面偏工程的部分可以按需跳读。
2. 因子图到底在解决什么问题:拆解核心设计思路
2.1 从联合分布到局部函数乘积
先说清楚因子图要表达的东西。任意一个概率模型,核心都是要刻画一个联合分布 $p(x_1, x_2, \dots, x_n)$,这些变量可能是一张图里的像素、一个机器人的位姿序列、或者一组待估的参数。直接写这个联合分布往往很麻烦,因为变量之间关系复杂。但现实中大多数模型都有个很好的性质:联合分布可以分解成若干个只依赖一小部分变量的函数相乘。
举个例子,一个典型的马尔可夫链,联合分布可以写成:
$$p(x_1, x_2, x_3) = \frac{1}{Z}\psi_1(x_1, x_2)\psi_2(x_2, x_3)$$
这里每个 $\psi$ 就是只涉及两个变量的函数,$Z$ 是归一化常数。如果把它画成因子图,就是三个变量圆圈、两个因子方块,第一个方块连 $x_1$ 和 $x_2$,第二个方块连 $x_2$ 和 $x_3$。这就是因子图的核心思想——把一个全局的、高维的联合分布,拆成若干个局部因子的乘积。
为什么这么拆有用?因为局部意味着低维,低维意味着计算量小。一个涉及 100 个变量的联合分布如果直接处理是指数级复杂度,但拆成每两个变量一个因子之后,很多推断算法就变成线性复杂度了。这就是因子图和消息传递算法配合起来这么强大的根本原因。
2.2 因子图相比贝叶斯网络和马尔可夫随机场的优势
有向图(贝叶斯网络)和无向图(马尔可夫随机场)已经能表达很多东西了,为什么还要单独搞一个因子图?我总结是三个原因,都是实操中真正能感受到的。
第一,表达能力强于无向图。马尔可夫随机场的团势函数有个约束:一个团上的势函数必须是该团变量的函数,而且团与团之间不能重复表达同一个关系。因子图没有这个约束,一个因子可以只涉及某几个变量,也可以涉及全部变量,非常灵活。比如一个全局的约束因子,直接连所有变量就完事了,不用纠结团的大小。
第二,把"变量"和"关系"彻底分开。在很多实际系统里,变量就是物理量,因子就是观测或约束,两类东西天然就是分开的。因子图强制你在建模阶段就把它们区分清楚,这个约束看起来是限制,实际上逼着你把模型想明白。我做的状态估计项目里,位姿是变量,里程计读数、GPS 读数、回环检测结果都是因子,分开之后代码结构特别干净。
第三,消息传递算法的推导更统一。在因子图上写和积算法,变量节点和因子节点各自做什么,规则非常清晰,而且对任何结构的图都适用。你不用像在无向图里那样还要区分什么极大团、最小团,推导一遍就能套到所有模型上。
下面这张表是我在实际选型时常用的一个对比,可以直观看出三者的区别:
| 模型类型 | 边/节点结构 | 表达约束 | 推断算法统一性 | 建模直观度 |
|---|---|---|---|---|
| 贝叶斯网络 | 有向边连接变量 | 条件概率 | 需要 moralize 转无向 | 因果性强,直观 |
| 马尔可夫随机场 | 无向边连接变量 | 团势函数 | 依赖团分解 | 团结构复杂时繁琐 |
| 因子图 | 二分图,变量+因子 | 任意因子函数 | 消息传递高度统一 | 变量与约束分离清晰 |
我个人的体会是:如果你的模型偏因果推理、变量关系是单向的,贝叶斯网络更顺手;如果模型里约束关系复杂、需要频繁做推理优化,因子图基本是更好的选择。
2.3 为什么工程上越来越偏爱因子图优化
近些年在机器人、SLAM、传感器融合这些领域,因子图出现的频率明显变高,热搜词里"因子图优化""回环因子加入因子图"就是这种现象的反映。原因很实在。
一是增量式的处理方式天然适合实时系统。系统一边运行一边来新的观测,新的观测就是新的因子,往图里加一个因子就行,不需要把整个问题推倒重来。老的因子可以保留,也可以按一定的策略边缘化掉,这种"增量加因子"的思路在实时估计里特别香。
二是稀疏性好。因子图的风险结构非常稀疏,每个因子只涉及少数变量,这直接决定了求解时的矩阵是稀疏矩阵。稀疏矩阵的求解效率比稠密矩阵高几个数量级,这是因子图能在几百上千个变量规模下实时求解的关键。
三是非线性优化的框架很成熟。用因子图做优化,本质上就是把最大后验估计转化成一个非线性最小二乘问题,然后用高斯-牛顿或者列文伯格-马夸尔特方法迭代求解。这个框架在数值优化领域已经非常成熟,工具库也很多,工程落地的门槛比想象中低不少。
理解了这三点,其实就理解了为什么因子图不只是课本上的一个知识点,而是真正能在工程里提升生产力的工具。接下来的部分,我会从建模到实操,把这个工具怎么用讲清楚。
3. 从零开始建一个因子图:核心细节与实操要点
3.1 识别变量节点和因子节点的实操方法
建模的第一步永远是搞清楚"什么是变量、什么是约束"。我用的方法很简单:列出所有你想估计或推断的量,它们是变量节点;列出所有你手上有的观测、先验、物理规律,它们是因子节点。
拿一个具体场景来说,假设有个移动的小车在一维直线运动,我们想估计它在若干个时刻的位置。手里的数据是:初始位置有先验(大致知道从哪出发)、每个时间步有一个测距传感器读数、相邻时刻之间有速度大致恒定的假设。那变量就是每个时刻的位置 $x_1, x_2, \dots, x_T$;因子就包括:一个先验因子连 $x_1$、若干个测量因子各自连一个 $x_t$、若干个运动约束因子各自连相邻两个位置。
这里有个实操要点:因子不一定要对应观测数据,物理规律、平滑约束、正则项都可以建模成因子。上面那个"相邻位置速度大致恒定"的约束,就是一个没有观测数据的因子,它表达的是一个平滑先验。很多初学者只把观测当因子,结果模型建出来不够约束,求解出来抖动严重,就是这个原因。
再提醒一个细节:变量的粒度要合适。有的系统里是直接估计绝对位姿,有的是估计相对位姿再加一个全局锚点。粒度选粗了约束不够,选细了变量数量爆炸。我一般的做法是先用粗粒度把模型跑通,看结果再决定要不要细化。
3.2 因子的数学形式与选择
因子在数学上就是一个函数,常见的有两类形式,选哪种很关键。
第一类是条件概率形式,因子直接就是某个条件概率 $p(z_t | x_t)$,比如测量模型。这种因子里包含了观测噪声的分布,高斯噪声就是高斯形式。第二类是势函数形式,只保证非负,不要求归一化,用起来更自由。
实际工程中,绝大多数情况会转成负对数形式来做优化。因为我们的目标是求最大后验估计:
$$\hat{x} = \arg\max_x \prod_i f_i(x_i)$$
取负对数之后,乘积变求和,最大化变最小化:
$$\hat{x} = \arg\min_x \sum_i (-\log f_i(x_i))$$
如果每个因子都是高斯形式,$-\log f_i$ 就是一个二次型,整个问题就变成标准的最小二乘。这一步转化是因子图优化能落地的核心,一定要理解透。高斯因子的负对数形式大致长这样:
$$-\log f_i(x) = \frac{1}{2}|h_i(x) - z_i|^2_{\Sigma_i}$$
其中 $h_i(x)$ 是测量预测,$z_i$ 是实际观测,$\Sigma_i$ 是噪声协方差。这个形式意味着:观测和预测差得越多,代价越大;噪声越小,同样的差异代价越大。噪声协方差在这里起到了加权的作用,这是很符合直觉的——不确定的观测权重低,确定的观测权重高。
注意:因子形式选择时,一定要确认噪声协方差矩阵是正定的。我在项目里遇到过因为手工给的协方差矩阵非正定导致求解器报错的情况,排查了大半天,后来发现是某两个观测被错误地设成了完全相关。
3.3 图的连接方式与稀疏性设计要点
因子图的价值很大程度上来自稀疏性,所以连接方式的设计直接影响性能。
基本原则是:因子只连它真正依赖的变量。一个只涉及两个变量的因子,绝不要让它连上第三个变量,哪怕你觉得那个变量"有点关系"。多加一条边,稀疏性就可能被破坏,求解效率可能从线性掉到平方甚至更差。
但这里有个反直觉的点:有时候适当增加约束反而能减少变量数量。比如观测之间如果存在某种解析关系,你可以先把它们合并成一个复合因子,再连到变量上,反而让图更简单。这需要在建模阶段对问题本身有足够理解。
另一个要点是变量的消元顺序。因子图求解通常涉及变量消元,消元顺序对计算的填充量影响巨大。经验做法是优先消元那些连接度低的变量,也就是邻居少的变量先消掉。工具库一般有自动的消元排序策略,但如果你发现求解特别慢,可以手动指定顺序试试。
| 设计维度 | 推荐做法 | 踩坑提醒 |
|---|---|---|
| 因子连接范围 | 只连接真正依赖的变量 | 多连一条边破坏稀疏性 |
| 变量粒度 | 先粗后细,逐步细化 | 一开始就细粒度容易变量爆炸 |
| 消元顺序 | 优先消元低连接度变量 | 顺序不当导致矩阵填充严重 |
| 因子数量 | 够用即可,避免冗余约束 | 冗余因子增加计算量 |
3.4 噪声模型和协方差的实操处理
噪声模型是因子图里最容易被轻视、又最容易出问题的地方。很多人建完图直接给所有因子一个单位协方差,跑出来发现结果不对,然后开始怀疑算法,其实是噪声模型给错了。
正确的思路是:每个因子的协方差应该反映这个观测或约束的真实不确定度。传感器噪声可以从数据手册拿,也可以用静态标定的方法估。实在没有先验信息,就先用经验值跑通,再通过残差分析反推合理的协方差。
我常用的一个技巧是看优化后的残差分布。如果某个因子的残差明显比其他因子大,说明它的协方差可能给得太小(被过度信任),或者这个观测本身有异常。反过来,如果所有残差都很小,可能协方差给得太大了。用这个反馈来调整协方差,通常迭代两三轮就能调到一个合理的水平。
还有一个实用做法是给协方差加一个下限。极端情况下某个方向的协方差接近零,会导致求解时的信息矩阵病态,数值不稳定。给协方差的对角线加一个小量(比如 $10^{-6}$),可以显著提升数值稳定性,代价几乎可以忽略。
4. 消息传递算法:因子图推断的核心实现
4.1 和积算法的推理逻辑与手算示例
和积算法是因子图上最基础的推断算法,目标是精确计算每个变量的边缘分布。它的核心是"消息"这个概念——变量节点和因子节点之间互相传递消息,每个节点把自己收到的东西加工后再传出去。
规矩是这样的:变量节点向因子节点传的消息,是从这个变量能收到的其他所有消息的乘积;因子节点向变量节点传的消息,是这个因子函数乘以它收到的其他所有消息再求和(积掉其他变量)。两个方向交替进行,直到每个变量节点都能收齐所有邻居的消息,然后把这些消息乘起来归一化,就得到该变量的边缘分布。
拿前面的一维小车例子手动算一遍。假设三个位置变量 $x_1, x_2, x_3$,因子是 $f_{12}(x_1,x_2)$ 和 $f_{23}(x_2,x_3)$。要算 $x_2$ 的边缘分布,$x_2$ 只需接收两个消息:一个来自 $f_{12}$,一个来自 $f_{23}$。来自 $f_{12}$ 的消息是 $\sum_{x_1} f_{12}(x_1, x_2) \cdot \mu_{x_1 \to f_{12}}(x_1)$,如果 $x_1$ 有先验,那 $\mu_{x_1}$ 就是先验。同理另一侧。最后 $p(x_2) \propto \mu_{f_{12}\to x_2}(x_2) \cdot \mu_{f_{23}\to x_2}(x_2)$。
手动算一遍最大的收获是理解到:和积算法本质上就是在变量消元,只不过它把所有消元顺序的结果都缓存下来了。所以它比暴力消元高效,代价是消息本身也占内存。
4.2 最大积算法与最大后验估计
很多工程场景里我们不需要整个边缘分布,只需要让后验概率最大的那一组变量取值,也就是最大后验估计(MAP)。这时候用最大积算法更合适。
最大积算法和和积算法的结构几乎一样,唯一区别是把因子到变量的消息里的"求和"换成"求最大值"。为什么?因为我们关心的是最大值而不是积分。从历史角度看,最大积算法就是维特比算法的图模型版本,本质都在做动态规划。
在因子图上,最大积算法的流程是:第一遍从叶子向根方向传消息,每个节点记录下让消息取最大的那个变量取值;第二遍从根向叶子回溯,取出记录下来的取值,就得到全局最优解。这个过程叫"前向-后向"或者"信念传播"的 max 版本。
实际用的时候有个坑要注意:最大积算法得到的是每个变量的最优取值组合,但它不一定给你这个组合的概率值,如果要比较不同模型或者做不确定性度量,还得额外计算。另外如果图里有环,最大积算法不保证收敛到全局最优,只能作为近似。
4.3 有环因子图上的近似推断
真实的因子图往往有环,比如回环检测引入的因子,就把原本链状的结构变成了环形。有环怎么办?
最常用的做法是信念传播在有环图上迭代运行,虽然理论上不保证收敛,但实践里大多数情况都能收敛到一个不错的近似解。收敛判据一般是看消息的变化量,小于某个阈值就停。如果震荡不收敛,可以减小消息的更新步长,或者用阻尼更新。
对于以优化为目标的场景(这也是工程里最常见的),其实可以绕开消息传递,直接用非线性最小二乘求解整个因子图,这就是所谓的"因子图优化"。在这种情况下,环的存在不影响求解框架,只是让信息矩阵变成更复杂的稀疏结构,求解器照样能搞定。
| 算法 | 目标 | 是否要求无环 | 工程常用度 |
|---|---|---|---|
| 和积算法 | 边缘分布 | 树结构精确,有环近似 | 中小规模推理常用 |
| 最大积算法 | 最大后验估计 | 同上 | 决策类问题常用 |
| 因子图优化 | 最大后验估计 | 不要求 | 大规模状态估计首选 |
5. 因子图优化的完整实操流程
5.1 问题建模阶段:从需求到图结构
实操的第一步是把实际问题翻译成因子图。我习惯分四步走:
- 列变量:把所有待估计的量列出来,明确每个变量的维度和物理含义。
- 列因子:把每一个观测、先验、约束列出来,标注它依赖哪些变量。
- 画图:用纸笔或简单的绘图工具把二分图画出来,检查连通性和是否有孤立节点。
- 定初值:给每个变量一个初始猜测,这对接下来的非线性优化很关键。
第三步画图特别重要,能帮你发现很多建模错误。比如某个变量是孤立的,说明它没有约束、无法估计;某个因子连了太多变量,可能是你想偷懒把本应分开的约束合并了。我做过一个姿态估计的模型,画图的时候才发现有两个传感器观测其实依赖于同一组隐含参数,一开始代码里把它们当独立因子处理了,修正之后精度提升很明显。
初值的选择也有讲究。能拿现成的就用现成的,比如上一时刻的估计、里程计积分的预测。如果实在没有,就用一个粗略的均值。非线性优化是局部方法,初值太差有可能收敛到很糟糕的局部极小,这是很常见的坑。
5.2 构建与求解:一个可复现的代码骨架
下面是一个概念性的代码骨架,展示因子图优化的典型结构。这里用类伪代码的方式写,不同语言和库的细节会有差异,但流程是通用的。
# 1. 定义变量 pose = Variable("pose", dimension=3, initial=[0.0, 0.0, 0.0]) landmark = Variable("landmark", dimension=2, initial=estimate_landmark()) # 2. 定义因子,每个因子包含它的误差函数和噪声协方差 prior_factor = Factor(error_fn=pose_prior_error, cov=prior_cov, keys=[pose]) odom_factor = Factor(error_fn=odometry_error, cov=odom_cov, keys=[pose]) # 实际会连两个时刻的位姿 obs_factor = Factor(error_fn=observation_error, cov=obs_cov, keys=[pose, landmark]) # 3. 组装因子图 graph = FactorGraph() graph.add_variables([pose, landmark]) graph.add_factors([prior_factor, odom_factor, obs_factor]) # 4. 求解非线性最小二乘 optimizer = LevenbergMarquardt(graph, max_iter=50, tol=1e-6) result = optimizer.solve()这段骨架里最关键的是每个因子的误差函数定义。以观测因子为例,误差函数大致是:
def observation_error(pose, landmark, measurement): predicted = predict_observation(pose, landmark) return predicted - measurement # 在协方差加权下做最小二乘求解器内部做的事情,是把所有因子的误差函数对变量求导,组装成一个巨大的稀疏线性系统,然后迭代求解。理解这一点,就能明白为什么噪声协方差重要(它决定了信息矩阵的权重)、为什么初值重要(线性化点依赖它)。
5.3 回环因子的加入与图结构变化
热搜词里提到的"回环因子加入因子图",是很典型的一个工程场景,值得单独讲。回环因子是指:系统运行一段时间后,检测到当前位置和很久之前某个位置是同一个地方,于是在两个位置变量之间加一条约束因子。
不加回环因子的图是链状的,每个变量只跟相邻变量有边,结构特别简单。加了回环因子之后,就在链上搭了一座"桥",整个图出现了环。这一步带来的影响是双面的:好处是回环因子提供了一条很强的长距离约束,能修正累积漂移;麻烦是图的结构变复杂了,求解的时候信息矩阵填充增加,计算量上升。
实际操作中加入回环因子的流程是:
- 检测回环:用位置、外观特征等方式判断当前状态和历史上的某个状态是否对应。
- 估计回环约束:计算这两个状态之间的相对变换,以及这个估计的不确定度。
- 加因子:在图里添加一条连接这两个变量节点的因子,协方差取上一步估计的不确定度。
- 重新优化:求解新的因子图,得到修正后的估计。
这里有个重要经验:回环因子的协方差千万不要给得太自信。回环检测本身有误检风险,如果误检且协方差给得小,优化结果会被一条错误的约束拉偏,比不加还糟糕。我一般会给回环因子留一个相对保守(较大)的协方差,让优化过程自己权衡。另外可以做回环检测的一致性校验,加因子之前先做个简单验证,能筛掉不少误检。
6. 常见问题与排查技巧实录
6.1 优化不收敛或震荡的排查思路
这是最常遇到的问题,排查顺序我一般这样走:
先看初值。把每个变量的初值打印出来,和它们的合理范围对比一下。如果某个变量的初值差得离谱,优化基本没戏。解决办法是用更好的初始化方法,比如用里程计或滤波方法先跑一遍给个粗估计。
再看噪声协方差。把所有协方差矩阵打印出来,检查是否有非正定、是否有量级差异过大的情况。量级差异过大是很隐蔽的坑,比如一个因子的协方差是 $10^{-6}$ 量级,另一个是 $10^6$ 量级,信息矩阵条件数极大,数值不稳定,表现出来就是震荡。
再看因子是否有冲突。把优化后的残差按因子打印出来,看看哪些因子残差特别大。持续的大残差往往意味着这个因子和其他因子互相矛盾,可能是观测有误,也可能是模型建错了。定位到具体因子后,逐个隔离验证。
| 现象 | 可能原因 | 排查动作 |
|---|---|---|
| 迭代不下降 | 初值太差或方向错误 | 检查初值与梯度计算 |
| 数值震荡 | 协方差量级失配 | 归一化协方差量级 |
| 收敛到奇怪解 | 局部极小或因子冲突 | 打印各因子残差 |
| 求解报错 | 协方差非正定 | 检查矩阵正定性 |
6.2 残差分布异常的定位方法
残差是因子图调试里最重要的信息源。优化完成后,把每个因子的残差(通常归一化到协方差尺度)统计一下,正常情况下应该大致符合零均值、单位方差的分布。如果某个因子残差系统性偏大,基本可以定位到问题。
我总结的定位逻辑是:如果一个因子残差大且稳定,多半是它的协方差给太小了(系统过度信任它,但它和事实不符,只能靠大残差来抗议);如果残差大且波动,多半是观测本身有噪声或异常值;如果所有因子残差都小得不正常,可能是协方差给太大了,整个系统太"软"。
6.3 实操中的独家避坑技巧
最后分享几个我在实际项目里总结的小技巧,这些在文档里通常不会写。
第一,先在小规模图上验证算法,再上大规模。因子图的很多问题在三个节点的小图上是看不出来的,但大规模一下冒出来。反过来,先用小图把算法逻辑跑通,加节点的时候心里有数,排查效率高很多。
第二,给关键因子加调试开关。在代码里把每个因子做成可开关的,调试时逐个打开,很快就能定位是哪个因子引入的问题。这个习惯帮我省了大量时间。
# 调试时逐个打开因子,观察结果变化 graph.add_factors([ prior_factor if ENABLE_PRIOR else None, odom_factor if ENABLE_ODOM else None, loop_factor if ENABLE_LOOP else None, ])第三,保存每一轮迭代的中间结果。优化不收敛的时候,回看中间过程能看出很多信息,比如某个变量在震荡、某个因子残差在发散。只保存最终结果的调试是盲人摸象。
第四,协方差的对角线加一个微小的下界。前面提过,这里再强调一次,数值稳定性提升非常明显,代价几乎可以忽略。我的经验值是相对量级在 $10^{-6}$ 左右比较合适,具体看问题的尺度。
第五,回环因子保守处理。误检的回环因子对系统的伤害非常大,宁可漏检不要误检。加因子之前做一致性校验,加因子的时候协方差给保守一些,后续如果发现系统表现良好,再逐步收紧。
这些技巧说起来简单,但每一条都是踩过坑之后总结出来的。因子图这个工具,建模阶段想清楚、噪声模型给准确、调试信息留充足,剩下的交给求解器就行。真正的门槛不在算法本身,而在于把实际问题准确翻译成图结构的那一步,这一步做好了,后面基本是工程活。