23-MCMC-马尔科夫链蒙特卡洛

1.背景

MCMC作为一种随机采样的方法，全名“Markov Chain Monte Carlo”，MCMC有两部分组成，马尔科夫链(Markov Chain简称MC)和蒙特卡洛方法(Monte Carlo Simulation简称MC)组成。
如果我们模型中的P(X)非常的复杂，那么我们是无法直接在P(X)中进行采样的，我们知道的两种采样拒绝采样(Rejection-sampling)和重要性采样(Importance-sampling)是绕过了从P(X)中采样，而采用了提议分布q(X)(propose distribution),借助q(x)进行采样。这两种方法对于提议分布q(x)是有限制条件的，它们要求提议分布q(x)和样本分布P(x)接近，如果两个分布不接近，那么效果就非常的差。但通常情况下P(X)会非常的复杂，那么就会导致我们非常困难的去寻找一个分布q(X)让它跟P(X)接近且q(X)还要简单，这个目标太难了，我们根本没办法去做。为了解决上述找q(x)难的问题，所以我们就用到了MCMC方法进行采样。
MCMC方法是一大种类方法论，我们知道的MH算法和Gibbs算法是基于MCMC方法的两种具体实现。

2.MCMC原理

它借助于构建一套马氏链，用马氏链的性质(即：马氏链在经过若干步骤后，它从初始状态收敛到平稳分布)
我们假设马氏链的状态空间为: { 1 , 2 , . . , K } \{1,2,..,K\} {1,2,..,K},即：我们用下图表示可得：

• 在第m时刻后，$q^{(m)}=q^{(m+1)}=...=C_{常数}$，故满足此性质的分布叫做平稳分布。
• 我们的目标是希望最后平稳分布$q^{(m)}$接近或甚至等于目标分布$P(x)$,这样我们就可以在平稳分布$q^{(m)}$中采集样本 { x 1 , x 2 , . . . , x N } , \{x_1,x_2,...,x_N\}, {x1​,x2​,...,xN​},那么我们就可以用采集的样本近似的代表目标分布$P(X)$.这样我们的采样的任务就完成了。

3.状态转移分布Q

为了更好的表达从 t 时刻到 t+1 时刻的状态变化，我们定义一个状态转移分布Q ,Q矩阵形式如下：
$$Q={\left(\begin{array}{cccc}{Q}_{11}& Q_{12}& ...& Q_{1k}\\ Q_{21}& Q_{22}& ...& Q_{1k}\\ ⋮& ⋮& ⋮& ⋮\\ Q_{k1}& Q_{k2}& ...& Q_{kk}\end{array}\right)}_{k\times k}\text{\hspace{1em}(1)}$$

• $Q_{ij}$表示具体的概率值；一般默认Q分布为随机矩阵，其特征值绝对值小于等于1
• 满足${\sum }_{i=1}^k{Q}_{i1}={\sum }_{i=1}^k{Q}_{i2}=...={\sum }_{i=1}^k{Q}_{ik}=1$;列和为：1
• 满足${\sum }_{j=1}^k{Q}_{1j}={\sum }_{j=1}^k{Q}_{2j}=...={\sum }_{j=1}^k{Q}_{jk}=1$；行和为：1

我们单独分析马氏链从 t 时刻到 t+1 时刻下的分布状态：如下图所示：

• 由下图可得：在 t+1 时刻时的 第 x=j 状态等于在 t 时刻状态下的K种状态的和再乘以转移矩阵；
  $$q^{(t+!)}(x=j)=\sum _{i=1}^K{q}^{(t)}(x=i)\cdot Q_{ij}\text{\hspace{1em}(2)}$$
  我们定义$q^{(t+1)}$为一个行向量：
  $$q^{(t+1)}=(q^{(t+1)}(x=1),q^{(t+1)}(x=2),...,q^{(t+1)}(x=K){)}_{1\times k}\text{\hspace{1em}(3)}$$
  将公式2代入到公式3中可得：
  $$q^{(t+1)}=(\sum _{i=1}^K{q}^{(t)}(x=i)\cdot Q_{i1},\sum _{i=1}^K{q}^{(t)}(x=i)\cdot Q_{i2},...,\sum _{i=1}^K{q}^{(t)}(x=i)\cdot Q_{ik})=q^t\cdot Q\text{\hspace{1em}(4)}$$
  故我们可以得到一个状态转移公式：
  $$q^{(t+1)}=q^t\cdot Q=...=q^{(1)}\cdot Q^t\text{\hspace{1em}(5)}$$

• 随机矩阵Q有一个性质：特征值的绝对值小于等于1
  我们对随机矩阵Q进行特征值分解：
  $$Q=A\Lambda A^{-1}\text{\hspace{1em}(6)}$$
  我们假设特征值只有一个为1，其余值得绝对值均小于1：
  $$\Lambda =\left(\begin{array}{cccc}{\lambda }_{11}& & & \\ & ⋮& & \\ & & 1& \\ & & ...& {\lambda }_{kk}\end{array}\right);|{\lambda }_{ii}|\le 1\text{\hspace{1em}(7)}$$
  则可得：
  q ( t + 1 ) = q ( 1 ) ⋅ { A Λ A − 1 } t = q ( 1 ) ⋅ A Λ t A − 1 (8) q^{(t+1)}=q^{(1)}·\{{A\Lambda A^{-1}\}}^t=q^{(1)}·A\Lambda^{t}A^{-1}\tag {8} q(t+1)=q(1)⋅{AΛA−1}t=q(1)⋅AΛtA−1(8)
  $$由随机矩阵性质可得，当t\to \infty 时，{\Lambda }^t=\left(\begin{array}{cccc}0& & & \\ & ⋮& & \\ & & 1& \\ & & ...& 0\end{array}\right);$$
  所以我们代入求得:
  第 m+1 时刻的分布$q^{(m+1)}$
  $$q^{(m+1)}=q^{(1)}\cdot A{\Lambda }^m{A}^{-1}$$
  第 m+1 时刻的分布$q^{(m+1)}$
  $$q^{(m+2)}=q^{(1)}\cdot A{\Lambda }^m{A}^{-1}\cdot A\Lambda A^{-1}=q^{(1)}\cdot A{\Lambda }^{m+1}{A}^{-1}\text{\hspace{1em}(9)}$$

• 因为${\Lambda }^m={\Lambda }^{m+1}$
  可得：
  $$q^{(m+2)}=q^{(m+1)}\text{\hspace{1em}(10)}$$
  $结论：当t\ge m$时，可得：
  $$q^{(m+1)}=q^{(m+2)}=...=q^{(\infty )}\text{\hspace{1em}(10)}$$
  故以上可以证明在第 m 时刻后 ，其分布为同一个概率分布，即趋近于平稳，故为平稳分布。

4. MCMC-采样方法

之前我们讲到的变分推断(Variational Inference),我们知道推断inference的意思是我们看到数据X，我们假设有隐变量Z，我们在隐变量先验P(Z)已知的情况下,求解后验P(Z|X)的过程称之为推断Inference.推断分类如下：

4.1 蒙特卡罗方法(Monte Carlo Method)

蒙特卡罗方法是一种基于采样的随机近似方法。主要用于做数值积分。我们的目标是为了求解后验概率P(Z|X),其中我们定义如下：
$$后验概率：P(Z|X)\text{\hspace{1em}(1)}$$

• Z：隐变量 latent data
• X：观测量 observed data
  当我们知道分布P(z|x)时候，我们常见的做法是为了求在这个分布Z|X的背景下关于函数f(Z)的期望：
  $${\mathbb{E}}_{Z|X}[f(Z)]={\int }_{Z}P(Z|X)f(Z)dZ\approx \frac{1}{N}\sum _{i=1}^{N}f(z_i)\text{\hspace{1em}(2)}$$
• 对于上式，我们问题在于怎样找到所需要的$z^{(1)}，z^{(2)}，...,z^{(N)}\sim P(Z|X)$,那么我们的就集中于如何在复杂分布P(Z|X)中进行采样。

三种基本采样方法：

• 概率分布采样
• 拒绝采样
• 重要性采样

4.2 概率分布采样

假设P(z)是我们要采样的概率分布，我们需要有一个均匀分布$U\in [0,1]$,用计算机的随机数来采集N个数据集 { u ( 1 ) , u ( 2 ) , . . . , u ( N ) } ∈ U = [ 0 , 1 ] , x ( i ) \{u^{(1)},u^{(2)},...,u^{(N)}\}\in U=[0,1],x^{(i)} {u(1),u(2),...,u(N)}∈U=[0,1],x(i),我们定义p.d.f为概率密度函数，c.d.f为概率分布函数;

• 概率分布采样步骤如下：
  
  我们对于这种概率采样方法看起来非常的有用，但是我们发现，当我们的概率密度函数十分的复杂的时候，我们就没办法直接求得c.d.f,在用$u^{(i)}$通过c.d.f反向求$x^{(i)}$时出现问题，即：
  $$x^{(i)}=c.d.f^{(-1)}(u^{(i)})\text{\hspace{1em}(3)}$$

• 当我们c.d.f很复杂的时候，我们就没办法直接求得$x^{(i)}$

4.3 拒绝采样 Rejection Sampling

由于对于目标分布p(z)十分的难估计和采样，所以我们假设一个简单的分布q(z)，进行辅助采样，具体如下：

• 我们定义一个叫proposal distribution的分布为q(Z) ,对于任意$\forall z_i$,保证如下：
  $$M\cdot q(z^{(i)})\ge p(z^{(i)})\text{\hspace{1em}(4)}$$
  

• 我们定义接受率为α=L1/L2: 0≤α≤1;
  $$\alpha =\frac{p(z^{(i)})}{M\cdot q(z^{(i)})}\text{\hspace{1em}(5)}$$

采样步骤如下：

• 从q(z)分布中采样得到$z^{(i)}即：z^{(i)}\sim q(z)$
• 从均匀分布U[0,1]从采样得到$u^{(i)}$
• 如果 $u^{(i)}\le \alpha$，那么我们就接受$z^{(i)}$,否则拒绝$z^{(i)}$.

综上所述，我们非常的依赖q(z)的形式，红色部分为我们拒绝区域，绿色部分为我们的拒绝区域。但这个分布又有很多的缺点，当M·q(z)远大于p(z)时，那么拒绝区域就非常的大，导致采样经常失败。就会导致采样率非常的低。所以这个采样方法非常的依赖于M·q(z)的选择，但实际上，我们发现p(z)十分的复杂，我们根本就没有办法来得到那么准确的结果，特别是采样 cost 非常高的话，经常性的采样失败带来的损失是很大的。

4.4 重要性采样

重要性采样并不是对概率进行采样，而是对概率分布的期望进行采样，即：
$${\mathbb{E}}_{p(z)}[f(z)]={\int }_{z}p(z)\cdot f(z)dz={\int }_z\frac{p(z)}{q(z)}\cdot q(z)f(z)dz\text{\hspace{1em}(6)}$$
$$={\int }_{z}f(z)\frac{p(z)}{q(z)}\cdot q(z)dz={\mathbb{E}}_{q(z)}[f(z)\frac{p(z)}{q(z)}]\approx \frac{1}{N}\sum _{i=1}^N[f(z^{(i)})\frac{p(z^{(i)})}{q(z^{(i)})}]\text{\hspace{1em}(7)}$$

• 我们假设样本$z^{(i)}$是在分布q(z)中进行采样的。即：$z^{(i)}\sim q(z),i=1,2,...,N$
• 我们定义权重$weight=\frac{p(z^{(i)})}{q(z^{(i)})}$
• 重要性采样也是严重的依赖于p(z)与q(z)的相似性程度。用来平衡不同的概率密度值之间的差距。同样重要性采样也可能
  会出现一些问题，就是两个分布之间的差距太大了话，总是采样采不到重要的样本，采的可能都是实际分布概率值小的部分。也就是采样效率不均匀的问题。在这个基础上，我们进一步提出了 SamplingImportance Resampling。

4.5 采样重要性重采样

经过重要性采样后，我们得到了 N 个样本点，以及对应的权重。那么我用权重来作为采样的概率，重新测采样出 N 个样本。也就是如下图所示

通过二次采样可以降低采样不平衡的问题。至于为什么呢？大家想一想，我在这里表达一下自己的看法。p(zi)/q(zi) 是 Weight，如果 Weight 比较大的话，说明 p(zi) 比较大而 q(zi) 比较的小，也就是我们通过 q(zi) 采出来的数量比较少。那么我们按权重再来采一次，就可以增加采到重要性样本的概率，成功的弥补了重要性采样带来的缺陷，有效的弥补采样不均衡的问题。

5. MCMC-马尔科夫链

我们知道蒙特卡罗方法[Monte Carlo Method]是一个基于随机采样的近似方法，我们定义马尔科夫链[Markov Chain]是时间和状态都是离散的链。随机过程就是研究对象是一个随机对象序列 { X 1 , X 2 , . . . , X N } ， 而 不 是 一 个 单 独 的 随 机 变 量 X . \{X_1,X_2,...,X_N\}，而不是一个单独的随机变量X. {X1​,X2​,...,XN​}，而不是一个单独的随机变量X.

5.1 马尔科夫性质

马氏链定义是它的随机变量序列满足马尔科夫性质；马尔科夫性质定义：当一个随机过程在给定现在状态及所有过去状态情况下，其未来状态的条件概率分布仅依赖于当前状态；换句话说，在给定现在状态时，它与过去状态（即该过程的历史路径）是条件独立的，那么此随机过程即具有马尔可夫性质。具有马尔可夫性质的过程通常称之为马尔可夫过程；简单来说：在给定当前状态下，未来只依赖于当前，未来跟过去无关，就是英雄不问出处，未来只取决于现在，跟过去没任何关系。

• 齐次(一阶)马尔科夫链的数学表达：
  $$P(X_{t+1}=x|X_1,X_2,...,X_t)=P(X_{t+1}=x|X_t)\text{\hspace{1em}(8)}$$

• 转移矩阵$P,每一个转移矩阵为[P_{ij}]$
  $$P_{ij}=P(X_{t+1}=j|X_t=i)\text{\hspace{1em}(9)}$$

• 马尔科夫链图形化表达：
  
  我们来局部分析概率从$X到X^{\ast }$
  

• 假如存在一个向量$\pi =[{\pi }^{(1)},{\pi }^{(2)},...,{\pi }^{(N)},...],假如有一个概率分布\pi 满足：{\sum }_{i=1}^{\infty }{\pi }^{(i)}=1$，且满足如下公式：那么我们就称这个序列 { π ( k ) } \{\pi(k)\} {π(k)}为马氏链$X_t$的平稳分布。
  $${\pi }_{t+1}(X^{\ast })=\int {\pi }_t(X)\cdot P(X⟼X^{\ast })dx\text{\hspace{1em}(10)}$$

• 重点：我们引入平稳分布的目的是当我们要通过采样的方式来求出后验分布P(Z)时，如果我们把P(Z)看作平稳分布的 { π ( k ) } \{\pi(k)\} {π(k)}，我们就可以构造一系列的 { X 1 , X 2 , . . . . , X N } \{X_1,X_2,....,X_N\} {X1​,X2​,....,XN​}的马氏链，让这个马氏链的 { π ( k ) } \{\pi(k)\} {π(k)}来逼近于P(Z)。

• 细致平稳概念(detailed balance)
  $$\pi (X)\cdot P(X⟼X^{\ast })=\pi (X^{\ast })\cdot P(X^{\ast }⟼X)\text{\hspace{1em}(11)}$$

• detailed-balance是平稳分布的充分非必要条件。

• $P(X⟼X^{\ast })：表示条件概率，相当于P(X|X^{\ast })$
  

5.2 细致分布证明平稳分布

$\int \pi (X)\cdot P(X⟼X^{\ast })dx$

我们知道细致分布为：$\pi (X)\cdot P(X⟼X^{\ast })=\pi (X^{\ast })\cdot P(X^{\ast }⟼X)$

$=\int \pi (X^{\ast })\cdot P(X^{\ast }⟼X)dx$

$=\pi (X^{\ast })\cdot \underset{=1}{\underbrace{\int P(X^{\ast }⟼X)dx}}$

$=\pi (X^{\ast })$
$$结论：\pi (X^{\ast })=\int \pi (X)\cdot P(X⟼X^{\ast })dx；[平稳分布公式]$$

• 注：为了找到平稳分布，我们不方便用平稳分布公式来找，我们曲线救国的方式用细致分布公式来找到平稳分布。

6.MH算法(Metropolis-Hastings)

我们关心的是我们的后验概率P(Z)，为了求关于后验概率相关函数的期望即：
$${\mathbb{E}}_{P(Z)}[f(z)]={\int }_{z}f(z)\cdot p(z)dz\approx \frac{1}{N}\sum _{i=1}^{N}f(z^{(i)})\text{\hspace{1em}(12)}$$

• 我们在概率分布P(Z)中不断的采样，得到N个$z^{(i)}\sim p(z)$.但实际上我们从P(Z)中采样是十分的困难，所以我们就用MCMC法进行采样。

我们由细致分布公式可知：$\pi (X)\cdot P(X⟼X^{\ast })=\pi (X^{\ast })\cdot P(X^{\ast }⟼X)$,我们所求最终的平稳分布为$\pi (X)$，所以我们把目标后验分布P(Z)看作平稳分布的$\pi (X)$，通过构造一个马氏链，让这个马氏链最终趋向于平稳分布P(Z);那么我们的问题转换为怎么找到detailed balance里的转移矩阵，使得这个马氏链趋向于平稳P(Z).

• 构造转移矩阵Q使得满足detailed-balance：
  $$P(Z)\cdot \underset{P(Z⟼Z^{\ast })}{\underbrace{Q(Z⟼Z^{\ast })\cdot \alpha (z,z^{\ast })}}=P(Z^{\ast })\cdot \underset{P(Z^{\ast }⟼Z)}{\underbrace{Q(Z^{\ast }⟼Z)\cdot \alpha (z^{\ast },z)}}\text{\hspace{1em}(13)}$$

• $\alpha (z,z^{\ast }):$表示接受率。
  $$\alpha (z,z^{\ast })=min(1，\frac{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}{p(Z)\cdot Q(Z⟼Z^{\ast })})\text{\hspace{1em}(14)}$$

• $Q(Z^{\ast }⟼Z)=Q(Z|Z^{\ast });Q(Z⟼Z^{\ast })=Q(Z^{\ast }|Z)$

$$P(Z)\cdot P(Z⟼Z^{\ast })=P(Z^{\ast })\cdot P(Z^{\ast }⟼Z)\text{\hspace{1em}(15)}$$

证：
$$P(Z)Q(Z⟼Z^{\ast })\cdot \alpha (z,z^{\ast })=P(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)\cdot \alpha (z^{\ast },z)\text{\hspace{1em}(16)}$$
证明如下：

$P(Z)Q(Z⟼Z^{\ast })\cdot \alpha (z,z^{\ast })$

$=P(Z)Q(Z⟼Z^{\ast })\cdot min(1，\frac{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}{p(Z)\cdot Q(Z⟼Z^{\ast })})$

$=\cdot min(P(Z)Q(Z⟼Z^{\ast })，p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z))$

$=p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z))\cdot min(\frac{p(Z)\cdot Q(Z⟼Z^{\ast })}{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)},1)$

由于可得：
$$\alpha (z^{\ast },z)=min(\frac{p(Z)\cdot Q(Z⟼Z^{\ast })}{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}，1)\text{\hspace{1em}(17)}$$
所以原式可得：
$$P(Z)Q(Z⟼Z^{\ast })\cdot \alpha (z,z^{\ast })=p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z))\cdot \alpha (z^{\ast },z)\text{\hspace{1em}(18)}$$
证毕，公式16成立！
结论：只要我们取接受率$\alpha (z^{\ast },z)$,，那么马氏链就满足detailed-balance 的条件,从而得到平稳分布。
$$\alpha (z,z^{\ast })=min(1，\frac{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}{p(Z)\cdot Q(Z⟼Z^{\ast })})\text{\hspace{1em}(19)}$$

• P(Z) 在转移矩阵 Q(Z|Z∗)α(Z∗, Z) 下是一个平稳分布，也就是一个马尔可夫链，通过在这个马尔可夫链中采样就可以得到我们的相应的数据样本点了。实际上这就是大名鼎鼎的 MetropolisHastings 采样法。

Metropolis-Hastings Sampling (M.H算法具体采样步骤)：

• 1 > 我们用计算机随机数从均匀分布U(0,1)中随机采样得到$u^{(i)},即可得：u^{(i)}\sim U(0,1)$
• 2 > 我们在分布$Q(Z|Z^{(i-1)})$中采样得到样本$Z^{\ast }$
• 3 > 我们定义$\alpha (z,z^{\ast })=min(1，\frac{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}{p(Z)\cdot Q(Z⟼Z^{\ast })})$
• 4 > 如果$u^{(i)}\le \alpha 时，可得Z^{(i)}=Z^{\ast }，否则Z^{(i)}=Z^{(i-1)}$，这样执行了 N 次，就可以得到 {Z(1), Z(2)
  , · · · , Z(N)} 个样本点

注意：对于后验概率P(Z)来说，我们隐含了一个重要的条件，即：
$$P(Z)=\frac{\hat{p}(Z)}{Z_p}\text{\hspace{1em}(20)}$$

• $Z_p:$表示归一化因子，是求不出来的，但在接受率上正好分子与分布约掉了，所以没影响接受率。Pˆ(Z) = likelihood × prior。所以这里的 P(Z) 和 P(Z∗) 就是 Pˆ(Z) 和 Pˆ(Z∗)。由于归一化因子被抵消了，干脆就直接写成了$P(Z)和P(Z^{\ast })$。

7.吉布斯采样算法Gibbs

吉布斯算法采样本质上也是属于MH采样，不过是一种特殊的MH采样。Gibbs采样假设要采样的分布P(Z)维度非常的高。即：
$$P(Z)=P(Z_1,Z_2,...,Z_N)\text{\hspace{1em}(21)}$$

• 采样第 i 维时，我们固定其他N-1个维度；即：
  $$Z_i\sim P(Z_i|Z_{1,2,...,i-1,i+1,...N})\text{\hspace{1em}(22)}$$
  举例说明：假如$P(Z)=P(Z_1,Z_2,Z_3)$，三维数据分布
• 首先我们给定初始值 $Z_1^{(0)},Z_2^{(0)},Z_3^{(0)}$
• 在 t + 1 时刻采样
  

我们知道MH算法中的接受率α表示如下，我们现在可以用Gibbs算法来类比MH算法的接受率：

• MH 算法接受率α：
  $$\alpha (z,z^{\ast })=min(1，\frac{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}{p(Z)\cdot Q(Z⟼Z^{\ast })})\text{\hspace{1em}(23)}$$

$$\frac{p(Z^{\ast })\cdot Q(Z^{\ast }⟼Z)}{p(Z)\cdot Q(Z⟼Z^{\ast })})$$

• $p(Z^{\ast })=p(Z^{\ast }|Z_{-i}^{\ast })P(Z_{-i}^{\ast })$概率拆分
• $Q(Z^{\ast }⟼Z)=P(Z|Z^{\ast })$条件概率
• $P(Z)=P(Z_i|Z_{-i})\cdot P(Z_{-i})$概率拆分
• $Q(Z⟼Z^{\ast })=P(Z^{\ast }|Z)$条件概率
  $$=\frac{p(Z^{\ast }|Z_{-i}^{\ast })P(Z_{-i}^{\ast })\cdot P(Z|Z^{\ast })}{P(Z_i|Z_{-i})\cdot P(Z_{-i})\cdot P(Z^{\ast }|Z)}\text{\hspace{1em}(24)}$$
• $P(Z|Z^{\ast })=P(Z_i|Z_{-i}^{\ast })$,对于Gibbs算法来说成立
• $P(Z^{\ast }|Z)=P(Z_i^{\ast }|Z_{-i})$,对于Gibbs算法来说成立
• $P(Z_{-i}^{\ast })=P(Z_{-i}):因为每次采样都是固定其他分量，所以概率相等，且两项等价：Z_{-i}^{\ast }=Z_{-i}$
  $$=\frac{p(Z^{\ast }|Z_{-i}^{\ast })P(Z_{-i}^{\ast })\cdot P(Z_i|Z_{-i}^{\ast })}{P(Z_i|Z_{-i}^{\ast })\cdot P(Z_{-i})\cdot P(Z_i^{\ast }|Z_{-i}^{\ast })}=1\text{\hspace{1em}(25)}$$
  所以，对于吉布斯采样算法，其接受率为1，所以吉布斯采样就是特殊的MH算法，且接受率为1，每次从状态转移矩阵$Q(Z|Z^{(t)})=P(Z_t^{\ast }|Z_{-i})$中采样得到$Z_i^{\ast }$即可。这里我们需要额外提醒一下，使用 Gibbs Sampling 是有使用前提的，也就是固定其他维度后的一
  维分布时方便进行采样的，如果固定其他维度的时候得到的一维分布仍然是非常难进行采样的，那么
  使用 Gibbs Sampling 也是没有用的。

8.MCMC-回顾

8.1采样的动机

蒙特卡洛方法主要是为了做采样的。那我们采样的动机是什么？在机器学习中我们为什么要去采样。

• 采样的本身就是常见的任务
  在我们的机器学习中，我们经常出现我们已经学出了一个模型P(X)，这个模型是一堆数据，有一种任务是让你从概率分布模型中产生一些样本，这个就是采样的过程。例如，P(X)是一个生成模型，表示一堆图片，让你生成一堆新的图片，那这个其实就是当你把模型学好后就可以从这个分布中采样新的数据样本即可。
• 求和或求积分的简化
  $${\int }_{x}P(X)f(x)dx={\mathbb{E}}_{p(x)}[f(x)]\text{\hspace{1em}(26)}$$
  从概率分布P(X)中采样出N个样本 { X ( 1 ) , X ( 2 ) , . . . X ( N ) , } ∈ P ( X ) \{X^{(1)},X^{(2)},...X^{(N)},\}\in P(X) {X(1),X(2),...X(N),}∈P(X),那么我们可以用如下公式代替期望：
  $${\mathbb{E}}_{p(x)}[f(x)]\approx \frac{1}{N}\sum _{i=1}^{N}f(X^{(i)})\text{\hspace{1em}(27)}$$

8.2 好的样本定义

• 样本趋向于高概率区域
  我们希望在分布P(X)中采样的分布都是在高概率区域，这样采样才有意义，比如，我们希望选取的人大代表应该始终代表最广大人民的根本利益，所以选择的应该是在大部人群中选取，而不是在总在那些富豪之间选取，毕竟富豪在整个人口分布中占据小部分。
• 样本与样本之间最好是相互独立的。
  我们最好是在不同种类中选取相互独立的样本，这样才能更好的表达整体数据，而不是总选择相关的样本。

但在实际过程中，我们采样是十分困难的，比如无向图模型中的P(X)中，往往有一个归一化因子Z,满足如下：
$$P(X)=\frac{1\cdot \hat{P}(X)}{Z}\text{\hspace{1em}(28)}$$

• Z:划分函数：partition-function;归一化常数Z,当维度很大的时候，往往难计算出来$Z={\int }_x\hat{P}(X)dx$,那么我们就很难直接在P(X)中进行采样【详见概率分布采样】。

所以我们遇到的困难如下：

• 划分函数Z是十分难求出来的。
• 即使划分函数Z求出来了，也会因为高维数据遇到高维的维度诅咒问题。
  比如： X ∈ { 1 , 2 , . . . , K } p ， 那 么 总 的 状 态 空 间 ： K P X\in\{1,2,...,K\}^p，那么总的状态空间：K^P X∈{1,2,...,K}p，那么总的状态空间：KP,这个维度就呈现指数级增长，那么我们就没办法遍历所有数据。

正因为我们发现直接采样P(X)十分的困难，所以我们才另辟蹊径采取了拒绝采样和重要性采样。借助简单的辅助分布q(X)来代替复杂的P(X)来进行采样，来绕过直接采样P(X);
MCMC的思想是试图构建一个马氏链去逼近最终的目标分布，比如MH算法采样和Gibbs算法采样。

9. MCMC的采样思路以及面临的困难

蒙特卡洛方法的马尔科夫链经过若干周期后会趋向于平稳分布，我们假设状态空间为离散 { 1 , 2 , . . . , K } \{1,2,...,K\} {1,2,...,K},定义状态转移矩阵为$Q=[Q_{ij}{]}_{k\times k}$,即假定到达m步后马氏链趋向于平稳分布。因此启发我们和采样方法结合起来产生MCMC。MCMC其实是一套方法论，我们常用的是MH算法采样和吉布斯算法采样。

• MCMC思路：
  利用马尔科夫链收敛于平稳分布，我们需要设计转移矩阵Q，使得平稳分布q(X)近似等于目标分布P(X)
• 对于MCMC来说，我们不需要显示的去寻找一个分布q(x),而是通过构造一个转移矩阵$Q=[Q_{ij}{]}_{k\times k}$和一个初始的样本$q^{(1)}(X)$,经过若干步骤后就达到平稳分布q(X),但我们不需要知道q(X)的形式.我们只需构造马氏链即可，马氏链必备转移矩阵Q和状态空间 { 1 , 2 , . . . , K } \{1,2,...,K\} {1,2,...,K},那么这个马氏链就可以了。
• 具体采样步骤如下：假定我们已知状态转移矩阵Q，它的状态空间为K个，那么初始的$q^{(1)}(X)$是一个离散型的概率分布，这个分布我们可以任意给，$q^{(1)}(X)$可以是均匀分布，也可以从1~K中任意的选择一个状态，这是马氏链的本身的性质,
• 马氏链的性质：无论初始状态选择1~K中任意的一个状态，经过有限步骤后，分布会收敛到同一个平稳分布。平稳分布q(X)与状态转移矩阵相关，与初始状态无关

9.1 MCMC问题

• 理论只保证收敛性，但无法知道何时收敛
• 混合时间过长，难以收敛，导致发散；根本原因是P(X)太复杂，维度太高且维度之间有相关性。
• 在平稳分布以后的样本$X^{(m)}$与$X^{(m+1)}$是相关的，即样本之间有一定的相关性，为了解决上述问题，我们应该每隔几步采样。
  举例说明：
  对于混合高斯模型来说，会有一个什么问题呢？就是样本都趋向于一个峰值附近，很有可能会过不了低谷，导致样本都聚
  集在一个峰值附近。这个问题出现的原因我们可以从能量的角度来解释这个问题。在无向图中，我们
  常用下列公式来进行表示：
  $$P(X)=\frac{1}{Z}\hat{P}(X)=\frac{1}{Z}ex{p}^{-E(X)}\text{\hspace{1em}(29)}$$
  
  实际上这里的 E(X) 指的就是能量函数，能量和概率是成反比的，概率越大意味着能量越低，能量越低，越难发生跳跃的现象。所以，采样很容易陷入到一个峰值附近。并且，多峰还可以分为均匀和陡峭，陡峭的情况中，能量差实在是太大了，就很难发生跳跃。就像孙悟空翻出如来佛祖的五指山一样，佛祖的维度很高，孙悟空在翻跟头的时候，一直在一个低维里面不同的打转，根本就跳不出来，就是来自佛祖的降维打击。所以，在高维情况下，很容易发生在一个峰值附近不停的采样，根本就跳不出来，导致采到的样本的多样性低，样本之间的关联性大，独立性低。