图网络算法—马尔科夫随机场与因子图

在之前的文章中，我们重点介绍了概率图的基本概念与基本定理，感兴趣的读者可以参考我前一篇文章图网络算法——概率图介绍与贝叶斯网络。其中贝叶斯网络是一种比较典型的有向概率图模型。在下面的文章中，我们将来介绍无向概率图的一个代表模型，即马尔科夫随机场。进一步，我们将来介绍因子图的相关概念与基本形式。

1. 马尔科夫随机场(MRF)

1.1 马尔科夫随机场引入

首先，与贝叶斯网络这种有向概率图的一个重要区别是，马尔科夫随机场是一类无向的概率图模型。其基本的组成是 G ( V , E ) G(V,E) G(V,E)，其中V表示的是MRF的节点的集合，而E表示的是MRF中无向边的集合。举一个简单的例子来说，对于一张普通的单通道的图片而言，其是由一个像素点的集合所构成，基进一步，这些像素点可以抽象成一个一般形式的矩阵。对于各个像素点而言，可以为其赋值一些无向边，这样就组成了一个基本的MRF的形式。即如下图所示的形式：
图网络算法—马尔科夫随机场与因子图

1.2 MRF的定义

数学形式上，与贝叶斯网络的定义类似，MRF的定义也是对于联合概率的一种分解形式。对于一般的联合概率 p ( x 1 , x 2 , . . . , x n ) p(x_1,x_2,...,x_n) p(x1,x2,...,xn)，MRF可以将其分解成归一化之后的因子联乘绩的形式。即：
p ( x 1 , x 2 , . . . , x n ) = 1 Z φ ∏ φ i ( D i ) p(x_1,x_2,...,x_n)=\frac{1}{Z_φ}∏φ_i(D_i) p(x1,x2,...,xn)=Zφ1∏φi(Di)
其中， Z φ Z_φ Zφ为联合分布的归一化因子，一般也称为是分配函数。其主要的作用是将这种联乘绩的形式重新归一化到0~1的区间之内。进一步， D i D_i Di是一个随机变量的集合，而 φ i φ_i φi的本质则为一个映射函数，其主要的作用是将随机变量集合 D i D_i Di映射到实数域上的某一个实数之上。其被称为因子或者势函数。

一般情况下，我们可以将概率图中所包含的随机变量分为一系列因子的集合。即 D 1 , D 2 , . . . , D k D_1,D_2,...,D_k D1,D2,...,Dk，进一步，利用势函数，可以将这些因子的集合映射到一个实数之上，即 φ = { φ 1 ( D 2 ) , . . . , φ i ( D i ) , . . . , φ k ( D k ) } φ=\{φ_1(D_2),...,φ_i(D_i),...,φ_k(D_k)\} φ={φ1(D2),...,φi(Di),...,φk(Dk)}，对于归一化因子 Z φ Z_φ Zφ而言，其基本的计算公式如下：
Z φ = ∑ x 1 , x 2 , . . x n p ′ ( x 1 , x 2 , . . . , x n ) Z_φ=∑_{x_1,x_2,..x_n}p'(x_1,x_2,...,x_n) Zφ=x1,x2,..xn∑p′(x1,x2,...,xn)
其中：
p ′ ( x 1 , x 2 , . . . , x n ) = ∏ i = 1 φ i ( D i ) p'(x_1,x_2,...,x_n)=∏_{i=1}φ_i(D_i) p′(x1,x2,...,xn)=i=1∏φi(Di)

下面，我们举一个简单的例子来说明一下MRF的计算流程：

首先，我们先给定一个MRF，即如下的形式：

图网络算法—马尔科夫随机场与因子图
在上述的MRF中，共计包含4个节点和4条无向边。进一步，我们将上述的MRF模型共计分为4个因子，其中 D 1 = { A , B } D_1=\{A,B\} D1={A,B}， D 2 = { B , C } , D 3 = { C , D } , D 4 = { D , A } D_2=\{B,C\},D_3=\{C,D\},D_4=\{D,A\} D2={B,C},D3={C,D},D4={D,A}，则上述的边缘概率可以根据MRF的定义整理成因子连乘积的形式。最后我们给出因变量的取值不同而导致的因子的结果，如下图所示：

图网络算法—马尔科夫随机场与因子图
下面，我们来分析MRF中联合概率的计算流程：

首先是对于归一化因子的计算过程，回顾我们之前所定义的计算归一化因子的计算公式，有：
Z φ = ∑ x 1 , x 2 , . . x n p ′ ( x 1 , x 2 , . . . , x n ) = ∑ ∏ i = 1 φ i ( D i ) Z_φ=∑_{x_1,x_2,..x_n}p'(x_1,x_2,...,x_n)=∑∏_{i=1}φ_i(D_i) Zφ=x1,x2,..xn∑p′(x1,x2,...,xn)=∑i=1∏φi(Di)
则可以计算出：
Z φ = φ 1 ( A = a 0 , B = b 0 ) ∗ φ 2 ( B = b 0 , C = c 0 ) ∗ φ 3 ( C = c 0 , D = d 0 ) ∗ φ 4 ( D = d 0 , A = a 0 ) + … … + φ 1 ( A = a 1 , B = b 1 ) ∗ φ 2 ( B = b 1 , C = c 1 ) ∗ φ 3 ( C = c 1 , D = d 1 ) ∗ φ 4 ( D = d 1 , A = a 1 ) Z_φ=φ_1(A=a^0,B=b^0)*φ_2(B=b^0,C=c^0)*φ_3(C=c^0,D=d^0)*φ_4(D=d^0,A=a^0) +……\\ +φ_1(A=a^1,B=b^1)*φ_2(B=b^1,C=c^1)*φ_3(C=c^1,D=d^1)*φ_4(D=d^1,A=a^1) Zφ=φ1(A=a0,B=b0)∗φ2(B=b0,C=c0)∗φ3(C=c0,D=d0)∗φ4(D=d0,A=a0)+……+φ1(A=a1,B=b1)∗φ2(B=b1,C=c1)∗φ3(C=c1,D=d1)∗φ4(D=d1,A=a1)
即如取值如下表所示：

图网络算法—马尔科夫随机场与因子图
最终，我们可以计算出 Z φ = 300000 + 300000 + 300000 + 30 + 500 + 500 + , . . . , + 100000 = 7 , 201 , 840 Z_φ=300000+300000+300000+30+500+500+,...,+100000=7,201,840 Zφ=300000+300000+300000+30+500+500+,...,+100000=7,201,840

进一步，我们可以根据归一化因子的结果，计算出各个因子连乘积的归一化概率，即如下的形式：
图网络算法—马尔科夫随机场与因子图

进一步，除了上面的因子之外，我们还可以为MRF增加一些势函数，即如下面的形式：

图网络算法—马尔科夫随机场与因子图
势函数的增加会增加一定的计算量，但是对于最终的概率计算结果不会产生影响。

1.3 MRF的常见应用—图像分割

在CV领域，成对的MRF被广泛的用于对图像进行建模，即如下图所示：

一般情况下，我们会利用成对的MRF来实现最大后验概率推理的问题，即下面的形式：

最后，定义 θ p ( x p ) = − l o g φ p ( x p ) ， θ p q ( x p , x q ) = − l o g φ p q ( x p , x q ) θ_p(x_p)=-logφ_p(x_p)，θ_{pq}(x_p,x_q)=-logφ_{pq}(x_p,x_q) θp(xp)=−logφp(xp)，θpq(xp,xq)=−logφpq(xp,xq)，可以将上述的将乘积的形式准换为求和的形式，并且进一步能够将求解最大后验概率问题转换成成求解最下损失函数的问题，即如下图所示：

进一步，我们举一个两个像素点构成的MRF，求解最大后验概率，即求解最小损失函数的问题，即如下图所示：
进一步，我们可以利用MRF对于图像的建模来实现图像的分割，首先，我们先假设图像中物体的取值为1，而背景的取值为0，即如下的形式：

图网络算法—马尔科夫随机场与因子图
进一步，我们假设图像内部的像素点是具有一定的连续性，我们可以定不同节点之间的边的势函数以及各个像素点的节点势函数为如下形式：

最后，我们根据图的能量函数以及定义好的势函数来对图像中的像素点进行0/1分类，最终所有取值为1的像素点为物体，取值为0的像素点为背景。

1.4 MRF应用—图像去噪

首先，给定一幅带有噪声的图像Y，我们的目标是对噪声进行去噪，从而恢复到原始的图像X，具体一点，我们需要计算最大的后验概率推理来实现图像的去噪过程。即如下的形式：
可以看出，我们可以利用贝叶斯公式来进行对后验概率的推导，进一步，我们可以根据贝叶斯的计算结果获取损失函数E。下面，我们引入MRF模型，对于上述损失函数计算的条件概率，我们使用单个随机变量的势函数来定义。进一步，根据图像是连续变化的特性，即像素点内部是连续变化的特性，对于后面的一般概率，考虑不同节点之间的影响，我们使用两个随机变量的势函数来定义。这样通过多个节点随机变量定义的势函数，可以保留图像的边缘特征，具体的定义形式如下：
最后，我们总结一下上面的建模流程，整个的建模过程分为两个部分，即两个势函数的定义，第一个势函数的目标是令当前节点势函数的取值与噪声图中对应节点的取值相近，即 x p x_p xp与 y p y_p yp相近。而第二个势函数的定义是为了保持图像的边缘特征。即两个节点的之间的距离小于一个阈值。即 θ ( x p , x q ) = m i n ( ∣ x p − x q ∣ , d ) θ(x_p,x_q)=min(|x_p-x_q|,d) θ(xp,xq)=min(∣xp−xq∣,d)。最后，根据loss函数，我们可以计算出每一个像素点的预测值，作为去噪之后的像素点的取值。