第四讲(拓展)_高斯过程(边缘化与条件作用)

Catalogue
  1. 1. 1. 高斯过程(边缘化与条件作用)
    1. 1.1. 1.1. 高斯过程的直观解释
      1. 1.1.1. 1.1.1. 多元高斯分布
        1. 1.1.1.1. 1.1.1.1. 二元高斯分布
      2. 1.1.2. 1.1.2. 边缘化与条件作用
    2. 1.2. 1.2. 相关理论推导
      1. 1.2.1. 1.2.1. 高斯分布的表示
      2. 1.2.2. 1.2.2. 联合高斯分布的分解
    3. 1.3. 1.3. 边缘化操作的作用
      1. 1.3.1. 最小二乘
        1. 1.3.1.1. 高斯牛顿
        2. 1.3.1.2. SLAM问题求解

1. 高斯过程(边缘化与条件作用)

1.1. 高斯过程的直观解释

1.1.1. 多元高斯分布

高斯分布(也叫做正态分布)是高斯过程的基础构件。而我们最感兴趣的是多元高斯分布,其每个随机变量都呈正态分布,联合分布也是高斯的。一般来说,多元高斯分布由均值向量\(\mu\)和协方差矩阵\(\sum\)定义。

  • 均值向量\(\mu\)描述了该分布的期望值,它的每个组件描述了对应维度的均值
  • 协方差矩阵\(\sum\)对每个维度的方差进行建模,并确定不同随机变量之间的关联(协方差矩阵总是对称且半正定的)
  • Σ 的对角线由第\(i\)个随机变量的标准差\(\sigma_i\) 组成,而非对角线的元素则描述了每个元素\(\sigma_{ij}\)之间的相关性
  • \(X\)符合正态分布。协方差矩阵协方差矩阵\(\sum\)描述了该分布的形状

1.1.1.1. 二元高斯分布

从图形上来看,该分布以均值为中心,由协方差矩阵决定其形状。下图展示了这些参数对于一个二维高斯分布的影响。每个随机变量的标准差在协方差矩阵的对角线上,而其它的值则显示了它们之间的协方差。

高斯分布被广泛应用于为真实世界建模,有时在原分布未知的情况下作为替代品,有时用于中心极限定理。接下来我们会进一步讲解如何操纵高斯分布,以及如何从中获得有用的信息。

1.1.2. 边缘化与条件作用

高斯分布有一个很赞的代数性质:它在条件作用和边缘化情况下是封闭的。意思是,经过这些运算后,在结果中得到的分布依旧是高斯分布,这就使得很多统计学和机器学习中的问题变得易解。

  1. 边缘化: 对高斯分布进行边缘化和条件作用进行分解, 即\(P(X,Y)=P(X)P(Y|X)=P(Y)P(X|Y)\): 其中 X 和 Y 代表原始随机变量的子集。

    给定随机变量 X 和 Y 组成的向量的正态概率分布\(P(X,Y)\), 我们可以用以下方法确定每个变量各自的边缘概率分布:

    这个公式所表达的意思很直接了当:X 和 Y 这两个子集各自只依赖于它们 μ 和 Σ 中对应的值。因此,要从高斯分布中边缘化一个随机变量,我们只需把μ 和Σ里那些对应的变量丢掉就行。

    通过边缘化,我们可以获取多元概率分布的一部分信息(这里对y进行边缘化,去掉变量y,得到关于x的边缘分布): \[ \begin{aligned} p_X(x)=\int_y p_{X,Y}(x,y)dy&= \int_y p_{Y|X}(y|x)p_X(x)dy \\ &= p_X(x)\int_y p_{Y|X}(y|x)dy \\ &=p_X(x) \end{aligned} \] 上式子的意思是:

    • 如果我们只想考虑X=x的情况, 即只对X=x的概率感兴趣, 那么我们要考虑变量Y所有可能的值, 需要对变量Y进行积分.(下面有图)
  2. 条件作用 高斯过程的另一个重要运算是条件作用,它可以用于得到一个变量在另一个变量条件下的概率分布。和边缘化类似,这个运算也是封闭的,会得到一个不同的高斯分布。条件运算是高斯过程的基石,它使贝叶斯推断成为可能。条件作用如下定义:

    新的均值只依赖于作为条件的变量,而协方差矩阵则和这个变量无关, 即:

    • \(X \sim \mathcal{N}(\mu_X,\sum_{XX})\)
    • \(P(Y|X) \sim \mathcal{N}(\mu_Y+\sum_{YX}\sum_{XX}^{-1}(X-\mu_{X}),\sum_{YY}-\sum_{YX}\sum_{XX}\sum_{XY})\)
  3. 边缘化与条件作用的直观感受

  • 边缘化可以理解为在高斯分布的一个维度上做累加,这也符合边缘分布的一般定义。
  • 条件作用也有个很好的几何表达——我们可以把它想象成在多元分布上切下一刀,从而获得一个维数更少的高斯分布。

下面的例子是关于Y的边缘分布[MARGINALIZATION(Y)], 实际上就是去掉变量x, 把变量x边缘化掉 (注意哦,与前面的描述是相反的, 这个是去掉变量X)

\[ \begin{aligned} p_Y(y)=\int_x p_{X,Y}(x,y)dx&= \int_x p_{X|Y}(x|y)p_Y(y)dx \\ &= p_Y(y)\int_x p_{X|Y}(x|y)dx \\ &=p_Y(y) \end{aligned} \]

上式的意思就是说, 求出来的边缘分布\(P_Y\)就是给定一个y值, 在原来的联合概率分布\(P(X,Y)\)上, 给定y值之后把变量\(x\)在给定的y值方向上做累加.

  • 图左是 Y 的边缘分布,是边缘化掉变量 X 的结果,类似于沿着 Y 轴把所有的X值做累加
  • 图右是以给定的 X 为条件的关于变量Y的条件分布,类似于在原始分布上切下一刀

图片来自:看得见的高斯过程(原文)


1.2. 相关理论推导

1.2.1. 高斯分布的表示

  • 协方差矩阵+均值
  • 信息矩阵+信息矢量

常见的是"协方差矩阵+均值"形式, \[p(x)=\eta \exp\{-\frac{1}{2}(x-\mu)^T\Sigma^{-1}(x-\mu)\}\]

其中对称正定矩阵Σ为随机变量x的协方差矩阵,μ为x的均值,简记为 \[p(x) = N(\mu, \Sigma)\]

信息矩阵+信息矢量的形式可以由上式推导而来 \[ \begin{aligned} p(x)&=\eta \exp\{-\frac{1}{2}(x-\mu)^T\Sigma^{-1}(x-\mu)\} \\ &=\eta\exp\{-\frac{1}{2}x^T\Sigma^{-1}x+x^T\Sigma^{-1}\mu\} \end{aligned} \] > 其中,运算中产生的常数项都全部吸收到了 η 中.

现在定义信息矩阵\(\Lambda=\Sigma^{-1}\), 信息矢量\(\xi=\Sigma^{-1}\mu=\Lambda\mu\), 则有 \[ p(x)=\eta\exp\{-\frac{1}{2}x^T\Lambda x+x^T\xi\} \]\[p(x) = N^{-1}(\xi, \Lambda)\]

1.2.2. 联合高斯分布的分解

设随机变量\(X, Y\)满足联合高斯分布\(p(X,Y)\)

于是有 \[ p(X,Y)=p(X)p(Y|X)\]

所谓边缘化操作就是求出上面的分布p(X) * 也就是把变量Y边缘化掉(去掉)

用"协方差矩阵+均值"形式给出联合分布的表示: \[ p(X,Y) = N\Bigg(\begin{pmatrix} \mu_X \\ \mu_Y \end{pmatrix}, \begin{pmatrix} \Sigma_{XX} & \Sigma_{XY} \\ \Sigma_{YX} & \Sigma_{YY} \end{pmatrix}\Bigg) \]

补充:舒尔补公式

  1. \[ \begin{aligned} \begin{bmatrix} I & -BD^{-1} \\ 0 & I \end{bmatrix} \begin{bmatrix} A & B \\ C & D \end{bmatrix} = \begin{bmatrix} \Delta_D & 0 \\ C & D \end{bmatrix} \end{aligned} \] \[ \Delta_D=A-BD^{-1}C \]
  2. \[ \begin{aligned} \begin{bmatrix} I & 0 \\ -CA^{-1} & I \end{bmatrix} \begin{bmatrix} A & B \\ C & D \end{bmatrix} = \begin{bmatrix} A & B \\ 0 & \Delta_A \end{bmatrix} \end{aligned} \] \[ \begin{aligned} \begin{bmatrix} A & B \\ C & D \end{bmatrix} \begin{bmatrix} I & -A^{-1}B \\ 0 & I \end{bmatrix} = \begin{bmatrix} A & 0 \\ C & \Delta_A \end{bmatrix} \end{aligned} \] \[\Delta_A=D-CA^{-1}B\]

将上面两式联合起来,就可以把矩阵变成对角型 \[ \begin{aligned} \begin{bmatrix} I & 0 \\ -CA^{-1} & I \end{bmatrix} \begin{bmatrix} A & B \\ C & D \end{bmatrix} \begin{bmatrix} I & -A^{-1}B \\ 0 & I \end{bmatrix} = \begin{bmatrix} A & 0 \\ 0 & \Delta_A \end{bmatrix} \end{aligned} \] 也可以通过对角型恢复出原来的矩阵 \[ \begin{aligned} \begin{bmatrix} I & 0 \\ CA^{-1} & I \end{bmatrix} \begin{bmatrix} A & 0 \\ 0 & \Delta_A \end{bmatrix} \begin{bmatrix} I & A^{-1}B \\ 0 & I \end{bmatrix} = \begin{bmatrix} A & B \\ C & D \end{bmatrix} \end{aligned} \] 就可以快速写出矩阵的逆 \[ \begin{aligned} \begin{bmatrix} A & B \\ C & D \end{bmatrix}^{-1} &= \begin{bmatrix} I & A^{-1}B \\ 0 & I \end{bmatrix}^{-1} \begin{bmatrix} A & 0 \\ 0 & \Delta_A \end{bmatrix}^{-1} \begin{bmatrix} I & 0 \\ CA^{-1} & I \end{bmatrix}^{-1} \\ &= \begin{bmatrix} I & -A^{-1}B \\ 0 & I \end{bmatrix} \begin{bmatrix} A^{-1} & 0 \\ 0 & \Delta_A^{-1} \end{bmatrix} \begin{bmatrix} I & 0 \\ -CA^{-1} & I \end{bmatrix} \end{aligned} \]

联合高斯分布p(X,Y)p(X,Y)的概率密度函数("信息矩阵+信息矢量形式")的表示为: \[ \begin{aligned} p(X, Y)&= \eta \exp\Bigg\{-\frac{1}{2}\begin{pmatrix} X-\mu_X \\ Y-\mu_Y \end{pmatrix}^T\begin{pmatrix} \Sigma_{XX} & \Sigma_{XY} \\ \Sigma_{YX} & \Sigma_{YY} \end{pmatrix}^{-1}\begin{pmatrix} X-\mu_X \\ Y-\mu_Y \end{pmatrix}\Bigg\} \\ &\propto \exp\Bigg\{-\frac{1}{2}\begin{pmatrix} a \\ b \end{pmatrix}^T\begin{pmatrix} A & C^{T} \\ C & D \end{pmatrix}^{-1}\begin{pmatrix} a \\ b \end{pmatrix}\Bigg\} \end{aligned} \]

于是, 就得到了关于X的边缘分布\(p(X)\)和关于Y的条件概率分布\(p(Y|X)\)

  1. 关于X的边缘分布\(p(X)\) \[ \begin{aligned} p(X) &=\eta_1 \exp \Big\{-\frac{1}{2} (X-\mu_X)^{T} \Sigma_{XX}^{-1}(X-\mu_X) \Big\} \sim \mathcal{N}(0,\Sigma_{XX}) \end{aligned} \] > 这说明了,边缘分布\(p(X)\)的协方差就是联合分布\(P(X,Y)\)的协方差矩阵中对应的矩阵块\(\Sigma_{XX}\)

  2. 关于Y的条件概率分布 \[ \begin{aligned} p(Y|X)&= \eta_2\exp\{-\frac{1}{2}[Y-(\mu_Y+\Sigma_{YX}\Sigma_{XX}^{-1}(X-\mu_X))]^T\Theta_{A}^{-1}[Y-(\mu_Y+\Sigma_{YX}\Sigma_{XX}^{-1}(X-\mu_X))]\} \\ & \propto \eta_2\exp\{-\frac{1}{2}[b-CA^{-1}a]^T\Theta_{A}^{-1}[b-CA^{-1}a]\} \\ \end{aligned} \]

其中,\(\Theta_A=\Delta_A=D-CA^{-1}B=D-CA^{-1}C^{T}\) 表示矩阵块A在联合概率分布\(P(X|Y)\)的协方差矩阵中的舒尔补

所以, \[ \Theta_A=\Delta_A=\Sigma_{YY}-\Sigma_{YX}\Sigma_{XX}^{-1}\Sigma_{XY} \]

这说明了,关于Y的条件分布\(p(Y|X)\)的协方差就是联合分布\(P(X,Y)\)的协方差矩阵关于子块\(A=\Sigma_{XX}\)的舒尔补\(\Delta_A=\Sigma_{YY}-\Sigma_{YX}\Sigma_{XX}^{-1}\Sigma_{XY}\)


假定现在我们只有联合概率分布\(p(X,Y)\)的信息矩阵 现在来讨论上面两个分布的信息矩阵

首先利用舒尔补写出联合概率分布\(p(X,Y)\)的信息矩阵 (也就是协方差矩阵的逆) \[ \begin{aligned} \begin{bmatrix} A & B \\ C & D \end{bmatrix}^{-1} &= \begin{bmatrix} I & -A^{-1}B \\ 0 & I \end{bmatrix} \begin{bmatrix} A^{-1} & 0 \\ 0 & \Delta_A^{-1} \end{bmatrix} \begin{bmatrix} I & 0 \\ -CA^{-1} & I \end{bmatrix} \end{aligned} \]

\[ \begin{aligned} \begin{bmatrix} A & C^{T} \\ C & D \end{bmatrix}^{-1} &= \begin{bmatrix} A^{-1}+A^{-1}C^{T}\Delta_{A}^{-1}CA^{-1} & -A^{-1}C^{T}\Delta_{A}^{-1} \\ -\Delta_{A}^{-1}CA^{-1} & \Delta_{A}^{-1} \end{bmatrix} \\ &= \begin{bmatrix} \Sigma_{XX}^{-1}+\Sigma_{XX}^{-1}\Sigma_{XY}\Theta_{A}^{-1}\Sigma_{YX}\Sigma_{XX}^{-1} & -\Sigma_{XX}^{-1}\Sigma_{XY}\Theta_{A}^{-1} \\ -\Theta_{A}^{-1}\Sigma_{YX}\Sigma_{XX}^{-1} & \Theta_{A}^{-1} \end{bmatrix} \\ &= \begin{bmatrix} \Lambda_{XX} & \Lambda_{XY} \\ \Lambda_{YX} & \Lambda_{YY} \end{bmatrix} \end{aligned} \]

  1. 关于X的边缘分布\(p(X)\)的信息矩阵 > 从前面的讨论, 已经知道了关于X的边缘分布\(p(X)\)的信息矩阵就是联合分布\(P(X,Y)\)的协方差矩阵中对应的矩阵块\(\Sigma_{XX}\), 那么对应的信息矩阵就是\(\Sigma_{XX}^{-1}\)

那么现在的目的就是: 利用原有的联合概率分布信息矩阵, 求出关于X的边缘分布\(p(X)\)的信息矩阵, 也就是边缘化掉变量Y之后的信息矩阵

易得 \[ \begin{aligned} \Sigma_{XX}^{-1}=\Lambda_{XX}-\Lambda_{XY}\Lambda_{YY}^{-1}\Lambda_{YX} \end{aligned} \]

  1. 关于Y的条件概率分布的信息矩阵 > 从前面的讨论, 已经知道了关于Y的条件概率分布的信息矩阵就是联合分布\(P(X,Y)\)的协方差矩阵关于子块\(A=\Sigma_{XX}\)的舒尔补\(\Theta_A=\Delta_A\), 那么对应的信息矩阵就是\(\Theta_A^{-1}\)

易得 \[ \Theta_A^{-1}=\Lambda_{YY} \]


1.3. 边缘化操作的作用

所谓边缘化,就是求某个联合概率分布的边缘分布。比如对于联合概率\(p(X,Y)\),对\(Y\)进行边缘化,就是对\(Y\)在整个空间中积分,即 \[ \begin{aligned} p_X(x)=\int_y p_{X,Y}(x,y)dy&= \int_y p_{Y|X}(y|x)p_X(x)dy \\ &= p_X(x)\int_y p_{Y|X}(y|x)dy \\ &=p_X(x) \end{aligned} \]

因此,\(Y\)边缘化的结果就是X的边缘化分布的\(p(X)\),仍然是一个高斯函数。伴随着边缘化,\(p(Y|X)\)就是\(p(X,Y)\)\(X\)的条件化。

最小二乘

在信息矩阵+信息矢量的表示方式下,边缘化和条件化与最小二乘法有密切关系。在许多基于最小二乘的优化问题中,常有如下形式的优化目标: \[ \min_x\Vert e(x)\Vert_W^2=\min_x~e(x)^TW^{-1}e(x) \] 其中\(W\)\(e(x)\)的协方差矩阵。

高斯牛顿

为了寻找上式的最小值,常使用迭代优化的方法,每一次迭代都会寻找一个增量\(\Delta x\)使目标函数减小。为了求增量,往往会将\(e(x)\)在当前\(x\)处展开为一阶近似(这种处理方式即Gauss-Newton Method),即: \[ e(x+\Delta x)\simeq e(x)+J(x)\Delta x \] 其中, \[ J(x)=\frac{\partial e}{\partial x} \] 于是, 可以生成新的优化目标: \[ \min_{\Delta x}~[e(x)+J(x)\Delta x]^TW^{-1}[e(x)+J(x)\Delta x] \] 这是关于\(\Delta x\)的二次函数,对\(\Delta x\)求导,并令导数等于0,有: \[ \begin{aligned} J(x)^TW^{-1}J(x)\Delta x+J(x)^TW^{-1}e(x)=0 \\ \Longrightarrow J(x)^TW^{-1}J(x)\Delta x = -J(x)^TW^{-1}e(x) \end{aligned} \]\(J(x)^TW^{-1}J(x) = H\), 是变量\(\Delta x\)的信息矩阵

\(-J(x)^TW^{-1}e(x)=b\), 则有 \[ H\Delta x=b \] 这就是非线性优化中的增量方程, 通过不断迭代求解\(\Delta x\), 更新变量\(x+=\Delta x\), 使得目标函数不断下降, 直到达到要求.

SLAM问题求解

在很多优化问题中,待优化的变量有明确意义,比如在SLAM或者SfM问题中,要优化的是所有相机的位姿\(c\)以及地图中所有三维点(路标点)的坐标\(p\),设\(\Delta x\)由这两个分量的增量构成,即 \[ \Delta x = \begin{pmatrix} \Delta c \\ \Delta p \end{pmatrix} \] 并且求出了对应的信息矩阵\(\Lambda\) \[ \Lambda=\begin{pmatrix} \Lambda{cc} & \Lambda{cp} \\ \Lambda{pc} & \Lambda{pp} \end{pmatrix}, ~~~ b = \begin{pmatrix} b_c \\ b_p \end{pmatrix} \] 于是有 \[ \begin{pmatrix} \Lambda{cc} & \Lambda{cp} \\ \Lambda{pc} & \Lambda{pp} \end{pmatrix} \begin{pmatrix} \Delta c \\ \Delta p \end{pmatrix}= \begin{pmatrix} b_c \\ b_p \end{pmatrix} \] 直接对信息矩阵求逆可以对上述方程求解, 但是当矩阵维度过大的时候, 直接求逆很难, 通常采用各种分解

可以采用上面的舒尔补公式进行消元, 也就是所说的高斯消元法 >(1) 舒尔补公式 >\[ \begin{aligned} \begin{bmatrix} I & -BD^{-1} \\ 0 & I \end{bmatrix} \begin{bmatrix} A & B \\ C & D \end{bmatrix} = \begin{bmatrix} \Delta_D & 0 \\ C & D \end{bmatrix} \end{aligned} >\] >\[ \Delta_D=A-BD^{-1}C >\] 即对方程两边同时左乘 \[ \begin{bmatrix} I & -\Lambda{cp}\Lambda{pp}^{-1} \\ 0 & I \end{bmatrix} \] 得到一个下三角型的形式 \[ \begin{pmatrix} \Lambda{cc}-\Lambda{cp}\Lambda{pp}^{-1}\Lambda{pc} & 0 \\ \Lambda{pc} & \Lambda{pp} \end{pmatrix} \begin{pmatrix} \Delta c \\ \Delta p \end{pmatrix}= \begin{pmatrix} b_c-\Lambda{cp}\Lambda{pp}^{-1}b_c \\ b_p \end{pmatrix} \] 于是原方程可以转换为两个独立方程 \[ \begin{aligned} (\Lambda{cc}-\Lambda{cp}\Lambda{pp}^{-1}\Lambda{pc})\Delta c &= b_c-\Lambda{cp}\Lambda{pp}^{-1}b_c \\ \Lambda_{pp}\Delta p &= b_p -\Lambda_{pc}\Delta c \end{aligned} \]

可以发现,\(\Delta c\)的系数矩阵,与上文中边缘分布\(p(X)\)的信息矩阵和信息矢量有相同的形式。而\(\Delta p\)的系数矩阵,则与关于Y的条件概率分布\(p(Y|X)\)的信息矩阵有相同形式。

也就是说,这里的高斯消元法,等价于对变量\(\Delta c\)做了边缘化,先将\(\Delta p\)边缘化掉,就是消掉路标点p,单独求相机位姿增量\(\Delta c\),然后再在\(\Delta c\)已知的情况下求\(\Delta p\)


参考: 1. 高斯分布与边缘化 2. 看得见的高斯过程(译文) 3. 看得见的高斯过程(原文)