介绍吉布斯采样前，我们先看一下蒙特卡洛方法。
有很多函数我们无法直接得到他的积分值，但我们可以利用蒙特卡洛方法来进行估计。
比如下面的积分：
∫
a
b
f
(
x
)
d
x
\int_a^b{f(x)}dx
∫abf(x)dx
我们采样一个点作为函数
f
(
x
)
f(x)
f(x)的均值对积分进行估计：
∫
a
b
f
(
x
)
d
x
≈
f
(
x
0
)
(
b
−
a
)
\int_a^b{f(x)}dx\approx f(x_0)(b-a)
∫abf(x)dx≈f(x0)(b−a)
只采样一个点，误差比较大，我们可等间隔的采样n个点：
∫
a
b
f
(
x
)
d
x
≈
b
−
a
n
∑
i
=
0
n
−
1
f
(
x
i
)
\int_a^b{f(x)}dx\approx\frac{b-a}{n}\sum_{i=0}^{n-1} f(x_i)
∫abf(x)dx≈nb−ai=0∑n−1f(xi)
上述我们是等间隔采样，其实隐含着假设：x在(a,b)区间是均匀分布的。如果x在(a,b)区间不是均匀分布，而是按照概率p(x)分布的呢？此时利用采样点进行积分估计的结果如下所示：
∫
a
b
f
(
x
)
d
x
=
∫
a
b
f
(
x
)
p
(
x
)
p
(
x
)
d
x
≈
1
n
∑
i
=
0
n
−
1
f
(
x
i
)
p
(
x
i
)
\int_a^b{f(x)}dx=\int_a^b{\frac{f(x)}{p(x)}p(x)}dx\approx\frac{1}{n}\sum_{i=0}^{n-1} {\frac{f(x_i)}{p(x_i)}}
∫abf(x)dx=∫abp(x)f(x)p(x)dx≈n1i=0∑n−1p(xi)f(xi)
相当于我按照概率p(x)采样n个点，然后根据这n个点对应的
f
(
x
i
)
f(x_i)
f(xi)和
p
(
x
i
)
p(x_i)
p(xi)值来进行估计。均匀分布相当于一种特例。如下：
p
(
x
)
=
1
b
−
a
p(x)=\frac{1}{b-a}
p(x)=b−a1
∫
a
b
f
(
x
)
d
x
≈
1
n
∑
i
=
0
n
−
1
f
(
x
i
)
p
(
x
i
)
=
1
n
∑
i
=
0
n
−
1
f
(
x
i
)
1
/
(
b
−
a
)
=
b
−
a
n
∑
i
=
0
n
−
1
f
(
x
i
)
\int_a^b{f(x)}dx\approx\frac{1}{n}\sum_{i=0}^{n-1} {\frac{f(x_i)}{p(x_i)}}\\=\frac{1}{n}\sum_{i=0}^{n-1} {\frac{f(x_i)}{1/(b-a)}}\\=\frac{b-a}{n}\sum_{i=0}^{n-1} f(x_i)
∫abf(x)dx≈n1i=0∑n−1p(xi)f(xi)=n1i=0∑n−11/(b−a)f(xi)=nb−ai=0∑n−1f(xi)
假设我们知道概率分布p(x)，但是这个概率分布不常见，不知道怎么采样。我们可以借助容易采样的分布q(x)来对p(x)进行采样。方法如下：
如果一个分布，既无法直接采样，也无法借助上文中的方法进行采样呢？如下：
此时，我们就要引入马尔可夫链，借助马尔可夫链的一些特性来解决采样的问题。
简单来说就是一个状态转移过程，且当前状态只跟前一个状态有关。以股市为例，假设有牛市，熊市，横盘三种状态，牛市有0.3的概率保持牛市，有0.2的概率转为熊市，有0.5的概率横盘，这种转移过程就是马尔可夫链。
转移过程会有一个初始状态
π
(
x
o
)
\pi (x_o)
π(xo)，会有一个状态空间，还有一个状态转移矩阵P。
转移过程中，若出现
π
(
x
n
−
1
)
P
=
π
(
x
n
)
=
π
(
x
n
−
1
)
\pi (x_{n-1}) P=\pi (x_{n})=\pi (x_{n-1})
π(xn−1)P=π(xn)=π(xn−1) ，即
π
P
=
π
\pi P=\pi
πP=π，则称概率分布
π
\pi
π是状态转移矩阵P的平稳分布。也称其为平稳的马尔可夫链。
马尔可夫链还有一些比较好的性质：
性质：给定初始状态
π
(
x
0
)
\pi (x_0)
π(x0)，给定转移矩阵P，n轮之后，
π
\pi
π会收敛。且稳定后的
π
\pi
π与初始状态
π
(
x
0
)
\pi (x_0)
π(x0)无关。
性质：对于给定的转移矩阵
P
P
P,
P
n
P^n
Pn在n大于一定值的时候是确定的。
其实马尔可夫链真正对我们有用的是这个公式：
π
P
=
π
\pi P=\pi
πP=π
其中
π
\pi
π是我们想要采样的分布，P是一个状态转移矩阵。
我们只要找到这样的P，就可以把对分布
π
(
x
)
\pi(x)
π(x)进行采样，转化为一个转移过程，因为转移过程中生成的样本点，就是我们想要的符合
π
(
x
)
\pi(x)
π(x)分布的样本点。
有一个细致平稳条件：
π
(
i
)
P
(
i
,
j
)
=
π
(
j
)
P
(
j
,
i
)
\pi (i)P(i,j)=\pi (j)P(j,i)
π(i)P(i,j)=π(j)P(j,i)
由细致平稳条件，我们就能得到
π
P
=
π
\pi P=\pi
πP=π，推导如下：
∑
i
π
(
i
)
P
(
i
,
j
)
=
∑
i
π
(
j
)
P
(
j
,
i
)
=
π
(
j
)
∑
i
P
(
j
,
i
)
=
π
(
j
)
\sum_i\pi (i)P(i,j)=\sum_i\pi (j)P(j,i)=\pi (j)\sum_iP(j,i)=\pi (j)
i∑π(i)P(i,j)=i∑π(j)P(j,i)=π(j)i∑P(j,i)=π(j)
上式就是
π
P
=
π
\pi P=\pi
πP=π 元素展开形式。
所以我们只要找到满足细致平稳条件的P，就能利用马尔可夫链对 π \pi π采样了。
假设 π ( x ) \pi (x) π(x)为我们期望进行采样的目标分布， Q Q Q为状态转移矩阵，其不一定满足细致平稳条件 π ( i ) Q ( i , j ) = π ( j ) Q ( j , i ) \pi (i)Q(i,j)=\pi (j)Q(j,i) π(i)Q(i,j)=π(j)Q(j,i).
为了使其满足细致平稳条件，我们引入一个
α
(
i
,
j
)
\alpha(i,j)
α(i,j),使得下式成立：
π
(
i
)
Q
(
i
,
j
)
α
(
i
,
j
)
=
π
(
j
)
Q
(
j
,
i
)
α
(
j
,
i
)
\pi (i)Q(i,j)\alpha(i,j)=\pi (j)Q(j,i)\alpha(j,i)
π(i)Q(i,j)α(i,j)=π(j)Q(j,i)α(j,i)
此时对应的 ，状态转移矩阵
P
(
i
,
j
)
=
Q
(
i
,
j
)
α
(
i
,
j
)
P(i,j)=Q(i,j)\alpha(i,j)
P(i,j)=Q(i,j)α(i,j)。
怎么才能找到这样的
α
(
i
,
j
)
\alpha(i,j)
α(i,j)呢？其实只要
α
(
i
,
j
)
\alpha(i,j)
α(i,j)满足下面条件即可：
α
(
i
,
j
)
=
π
(
j
)
Q
(
j
,
i
)
α
(
j
,
i
)
=
π
(
i
)
Q
(
i
,
j
)
\alpha(i,j)=\pi (j)Q(j,i)\\ \alpha(j,i)=\pi (i)Q(i,j)
α(i,j)=π(j)Q(j,i)α(j,i)=π(i)Q(i,j)
上面两个式子是等价的，我们只需直接令
α
(
i
,
j
)
=
π
(
j
)
Q
(
j
,
i
)
\alpha(i,j)=\pi (j)Q(j,i)
α(i,j)=π(j)Q(j,i)就行了。
其中， α ( i , j ) \alpha(i,j) α(i,j)我们称之为接受率。有点类似蒙特卡洛采样那里的接受-拒绝方法。
MCMC采样过程如下：
M-H采样是Metropolis Hastings采样的简称。
MCMC已经解决了我们的采样问题，但是其样本接受率可能比较低，就是
α
(
x
t
,
x
∗
)
\alpha(x_t,x_*)
α(xt,x∗)可能为0.1，0.001，这样采样次数就会比较多，浪费了很多时间。怎样提高接受率呢？
观察细致平稳条件
π
(
i
)
Q
(
i
,
j
)
α
(
i
,
j
)
=
π
(
j
)
Q
(
j
,
i
)
α
(
j
,
i
)
,
其中
α
(
i
,
j
)
=
π
(
j
)
Q
(
j
,
i
)
\pi (i)Q(i,j)\alpha(i,j)=\pi (j)Q(j,i)\alpha(j,i) , 其中\alpha(i,j)=\pi (j)Q(j,i)
π(i)Q(i,j)α(i,j)=π(j)Q(j,i)α(j,i),其中α(i,j)=π(j)Q(j,i)
我们知道，在两边同时乘以一个数字，或者除以一个非0数字，细致平稳条件也是满足的。所以对上式两边同乘
min
{
α
(
j
,
i
)
,
α
(
i
,
j
)
}
\min \{\alpha(j,i),\alpha(i,j)\}
min{α(j,i),α(i,j)}，可得下式：
π
(
i
)
Q
(
i
,
j
)
α
(
i
,
j
)
∗
min
{
α
(
j
,
i
)
,
α
(
i
,
j
)
}
=
π
(
j
)
Q
(
j
,
i
)
α
(
j
,
i
)
∗
min
{
α
(
j
,
i
)
,
α
(
i
,
j
)
}
\pi (i)Q(i,j)\alpha(i,j)*\min \{\alpha(j,i),\alpha(i,j)\}=\pi (j)Q(j,i)\alpha(j,i) *\min \{\alpha(j,i),\alpha(i,j)\}
π(i)Q(i,j)α(i,j)∗min{α(j,i),α(i,j)}=π(j)Q(j,i)α(j,i)∗min{α(j,i),α(i,j)}
然后两边再同时除以
α
(
j
,
i
)
∗
α
(
i
,
j
)
\alpha(j,i)*\alpha(i,j)
α(j,i)∗α(i,j),得：
π
(
i
)
Q
(
i
,
j
)
∗
min
{
1
,
α
(
i
,
j
)
α
(
j
,
i
)
}
=
π
(
j
)
Q
(
j
,
i
)
∗
min
{
α
(
j
,
i
)
α
(
i
,
j
)
,
1
}
\pi (i)Q(i,j)*\min \{1,\frac{\alpha(i,j)}{\alpha(j,i)}\}=\pi (j)Q(j,i)*\min \{\frac{\alpha(j,i)}{\alpha(i,j)},1\}
π(i)Q(i,j)∗min{1,α(j,i)α(i,j)}=π(j)Q(j,i)∗min{α(i,j)α(j,i),1}
所以我们得到：
α n e w ( i , j ) = min { α ( i , j ) α ( j , i ) , 1 } = min { π ( j ) Q ( j , i ) π ( i ) Q ( i , j ) , 1 } \alpha_{new}(i,j)=\min \{\frac{\alpha(i,j)}{\alpha(j,i)},1\}=\min \{ \frac{\pi (j)Q(j,i)}{\pi (i)Q(i,j)},1\} αnew(i,j)=min{α(j,i)α(i,j),1}=min{π(i)Q(i,j)π(j)Q(j,i),1}
如果转移矩阵Q是对称的，上式可简化为：
α
n
e
w
(
i
,
j
)
=
min
{
π
(
j
)
π
(
i
)
,
1
}
\alpha_{new}(i,j)=\min \{ \frac{\pi (j)}{\pi (i)},1\}
αnew(i,j)=min{π(i)π(j),1}
M-H采样是对MCMC方法的改进，也属于MCMC采样。
但是M-H采样也有一些问题，比如
- 当特征比较多时， π ( j ) Q ( j , i ) / π ( i ) Q ( i , j ) \pi (j)Q(j,i)/\pi (i)Q(i,j) π(j)Q(j,i)/π(i)Q(i,j)的计算比较耗时。
- 虽然 α n e w ( i , j ) \alpha_{new}(i,j) αnew(i,j)的接受率比较大，但还是小于1的，还是有很多辛苦计算的采样点被拒绝掉。
- 此外，很多时候我们不知道其联合分布，只知道其条件分布，那么这种情况我们就无法对联合分布进行采样。
吉布斯采样主要就是为了解决了上述问题。
在MCMC采样中，我们要寻找一个满足细致平稳条件 π ( i ) P ( i , j ) = π ( j ) P ( j , i ) \pi (i)P(i,j)=\pi (j)P(j,i) π(i)P(i,j)=π(j)P(j,i)的状态转移矩阵P。
下面，我们利用条件概率的变换，来寻找这样的P。首先我们以见到那的二维情况维例：
假设我们有状态
A
=
(
x
1
1
,
x
2
1
)
A=(x_1^1, x_2^1)
A=(x11,x21),
B
=
(
x
1
1
,
x
2
2
)
B=(x_1^1, x_2^2)
B=(x11,x22),
C
=
(
x
1
2
,
x
2
1
)
C=(x_1^2,x_2^1)
C=(x12,x21), 则有：
π
(
A
)
π
(
x
2
2
∣
x
1
1
)
=
π
(
x
1
1
,
x
2
1
)
π
(
x
2
2
∣
x
1
1
)
=
π
(
x
1
1
)
π
(
x
2
1
∣
x
1
1
)
π
(
x
2
2
∣
x
1
1
)
=
π
(
x
1
1
)
π
(
x
2
2
∣
x
1
1
)
π
(
x
2
1
∣
x
1
1
)
=
π
(
x
1
1
,
x
2
2
)
π
(
x
2
1
∣
x
1
1
)
=
π
(
B
)
π
(
x
2
1
∣
x
1
1
)
\pi (A)\pi(x_2^2|x_1^1)=\pi (x_1^1,x_2^1)\pi(x_2^2|x_1^1)\\ =\pi (x_1^1)\pi (x_2^1|x_1^1)\pi(x_2^2|x_1^1) \\ =\pi (x_1^1)\pi(x_2^2|x_1^1)\pi (x_2^1|x_1^1) \\ =\pi (x_1^1,x_2^2)\pi(x_2^1|x_1^1) \\ =\pi (B)\pi(x_2^1|x_1^1)
π(A)π(x22∣x11)=π(x11,x21)π(x22∣x11)=π(x11)π(x21∣x11)π(x22∣x11)=π(x11)π(x22∣x11)π(x21∣x11)=π(x11,x22)π(x21∣x11)=π(B)π(x21∣x11)
π
(
A
)
π
(
x
1
2
∣
x
2
1
)
=
π
(
x
1
1
,
x
2
1
)
π
(
x
1
2
∣
x
2
1
)
=
π
(
x
2
1
)
π
(
x
1
1
∣
x
2
1
)
π
(
x
1
2
∣
x
2
1
)
=
π
(
x
2
1
)
π
(
x
1
2
∣
x
2
1
)
π
(
x
1
1
∣
x
2
1
)
=
π
(
x
1
2
,
x
2
1
)
π
(
x
1
1
∣
x
2
1
)
=
π
(
C
)
π
(
x
1
1
∣
x
2
1
)
\pi (A)\pi(x_1^2|x_2^1)=\pi (x_1^1,x_2^1)\pi(x_1^2|x_2^1)\\ =\pi (x_2^1)\pi (x_1^1|x_2^1)\pi(x_1^2|x_2^1) \\ =\pi (x_2^1)\pi(x_1^2|x_2^1)\pi (x_1^1|x_2^1) \\ =\pi (x_1^2,x_2^1)\pi(x_1^1|x_2^1) \\ =\pi (C)\pi(x_1^1|x_2^1)
π(A)π(x12∣x21)=π(x11,x21)π(x12∣x21)=π(x21)π(x11∣x21)π(x12∣x21)=π(x21)π(x12∣x21)π(x11∣x21)=π(x12,x21)π(x11∣x21)=π(C)π(x11∣x21)
如果我们令状态转移概率矩阵P为下式：
P ( A → B ) = π ( x 2 B ∣ x 1 1 ) ，若 x 1 A = x 1 B = x 1 1 P(A \rightarrow B)=\pi(x_2^B|x_1^{1})，若x_1^A= x_1^B= x_1^1 P(A→B)=π(x2B∣x11)，若x1A=x1B=x11
P ( A → C ) = π ( x 1 C ∣ x 2 1 ) ，若 x 2 A = x 2 C = x 2 1 P(A \rightarrow C)=\pi(x_1^C|x_2^{1})，若x_2^A= x_2^C= x_2^1 P(A→C)=π(x1C∣x21)，若x2A=x2C=x21
P ( A → D ) = 0 ，其他情况 P(A \rightarrow D)=0，其他情况 P(A→D)=0，其他情况
那么对任意两点E,F, 将满足细致平稳条件。
π
(
E
)
P
(
E
→
F
)
=
π
(
F
)
P
(
F
→
E
)
\pi(E)P(E\rightarrow F)=\pi(F)P(F\rightarrow E)
π(E)P(E→F)=π(F)P(F→E)
证明如下：
因为满足
x
1
A
=
x
1
B
=
x
1
1
x_1^A= x_1^B= x_1^1
x1A=x1B=x11，
所以
P
(
A
→
B
)
=
π
(
x
2
B
∣
x
1
1
)
=
π
(
x
2
2
∣
x
1
1
)
P
(
B
→
A
)
=
π
(
x
2
A
∣
x
1
1
)
=
π
(
x
2
1
∣
x
1
1
)
P(A \rightarrow B)=\pi(x_2^B|x_1^{1})=\pi(x_2^2|x_1^{1}) \\ P(B \rightarrow A)=\pi(x_2^A|x_1^{1})=\pi(x_2^1|x_1^{1})
P(A→B)=π(x2B∣x11)=π(x22∣x11)P(B→A)=π(x2A∣x11)=π(x21∣x11)
所以有：
π
(
A
)
P
(
A
→
B
)
=
π
(
A
)
π
(
x
2
2
∣
x
1
1
)
=
π
(
B
)
π
(
x
2
1
∣
x
1
1
)
=
π
(
B
)
P
(
B
→
A
)
\pi(A)P(A\rightarrow B)=\pi (A)\pi(x_2^2|x_1^1) =\pi (B)\pi(x_2^1|x_1^1) =\pi(B)P(B\rightarrow A)
π(A)P(A→B)=π(A)π(x22∣x11)=π(B)π(x21∣x11)=π(B)P(B→A)
因为满足
x
2
A
=
x
2
C
=
x
2
1
x_2^A= x_2^C= x_2^1
x2A=x2C=x21
所以
P
(
A
→
C
)
=
π
(
x
1
C
∣
x
2
1
)
=
π
(
x
1
2
∣
x
2
1
)
P(A \rightarrow C)=\pi(x_1^C|x_2^{1})=\pi(x_1^2|x_2^{1})
P(A→C)=π(x1C∣x21)=π(x12∣x21)
P
(
C
→
A
)
=
π
(
x
1
A
∣
x
2
1
)
=
π
(
x
1
1
∣
x
2
1
)
P(C \rightarrow A)=\pi(x_1^A|x_2^{1})=\pi(x_1^1|x_2^{1})
P(C→A)=π(x1A∣x21)=π(x11∣x21)
所以有：
π
(
A
)
P
(
A
→
C
)
=
π
(
A
)
π
(
x
1
2
∣
x
2
1
)
=
π
(
C
)
π
(
x
1
1
∣
x
2
1
)
=
π
(
C
)
π
(
C
→
A
)
\pi(A)P(A\rightarrow C)=\pi (A)\pi(x_1^2|x_2^1)=\pi (C)\pi(x_1^1|x_2^1)=\pi (C)\pi(C \rightarrow A)
π(A)P(A→C)=π(A)π(x12∣x21)=π(C)π(x11∣x21)=π(C)π(C→A)
因为既不满足
x
1
A
=
x
1
D
=
x
1
1
x_1^A= x_1^D= x_1^1
x1A=x1D=x11
也不满足
x
2
A
=
x
2
D
=
x
2
1
x_2^A= x_2^D= x_2^1
x2A=x2D=x21
所以
P
(
A
→
D
)
=
0
P(A \rightarrow D)=0
P(A→D)=0
P
(
D
→
A
)
=
0
P(D \rightarrow A)=0
P(D→A)=0
所以有：
π
(
A
)
P
(
A
→
D
)
=
0
=
π
(
D
)
π
(
D
→
A
)
\pi(A)P(A\rightarrow D)=0=\pi (D)\pi(D \rightarrow A)
π(A)P(A→D)=0=π(D)π(D→A)
综上，得证任意两点E,F, 有
π
(
E
)
P
(
E
→
F
)
=
π
(
F
)
P
(
F
→
E
)
\pi(E)P(E\rightarrow F)=\pi(F)P(F\rightarrow E)
π(E)P(E→F)=π(F)P(F→E)
假设假设
A
=
(
x
1
1
,
x
2
1
,
.
.
.
,
x
n
1
)
A=(x_1^1,x_2^1,...,x_n^1)
A=(x11,x21,...,xn1),
B
n
=
(
x
1
1
,
x
2
1
,
.
.
.
x
n
−
1
1
,
x
n
2
)
B_n=(x_1^1,x_2^1,...x_{n-1}^1,x_n^2)
Bn=(x11,x21,...xn−11,xn2),我们有：
π
(
A
)
π
(
x
n
2
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
=
π
(
x
1
1
,
x
2
1
,
.
.
.
,
x
n
1
)
π
(
x
n
2
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
=
π
(
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
π
(
x
n
1
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
π
(
x
n
2
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
=
π
(
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
π
(
x
n
2
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
π
(
x
n
1
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
=
π
(
x
1
1
,
x
2
1
,
.
.
.
x
n
−
1
1
,
x
n
2
)
π
(
x
n
1
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
=
π
(
B
)
π
(
x
n
1
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
n
−
1
1
)
\pi (A)\pi(x_n^2|x_1^1,x_2^1,...,x_{n-1}^1)=\pi (x_1^1,x_2^1,...,x_n^1)\pi(x_n^2|x_1^1,x_2^1,...,x_{n-1}^1)\\ =\pi (x_1^1,x_2^1,...,x_{n-1}^1)\pi (x_n^1|x_1^1,x_2^1,...,x_{n-1}^1)\pi(x_n^2|x_1^1,x_2^1,...,x_{n-1}^1) \\ =\pi (x_1^1,x_2^1,...,x_{n-1}^1)\pi(x_n^2|x_1^1,x_2^1,...,x_{n-1}^1) \pi (x_n^1|x_1^1,x_2^1,...,x_{n-1}^1)\\ =\pi (x_1^1,x_2^1,...x_{n-1}^1,x_n^2)\pi(x_n^1|x_1^1,x_2^1,...,x_{n-1}^1) \\ =\pi (B)\pi(x_n^1|x_1^1,x_2^1,...,x_{n-1}^1)
π(A)π(xn2∣x11,x21,...,xn−11)=π(x11,x21,...,xn1)π(xn2∣x11,x21,...,xn−11)=π(x11,x21,...,xn−11)π(xn1∣x11,x21,...,xn−11)π(xn2∣x11,x21,...,xn−11)=π(x11,x21,...,xn−11)π(xn2∣x11,x21,...,xn−11)π(xn1∣x11,x21,...,xn−11)=π(x11,x21,...xn−11,xn2)π(xn1∣x11,x21,...,xn−11)=π(B)π(xn1∣x11,x21,...,xn−11)
因此，如果我们令状态转移概率矩阵为下式：
- P ( A → B n ) = π ( x n B n ∣ x 1 1 , x 2 1 , . . . , x n − 1 1 ) ，若 x i A = x i B n = x i 1 , 1 ≤ i ≤ n 且 i ≠ n . P(A \rightarrow B_n)=\pi(x_n^{B_n}|x_1^1,x_2^1,...,x_{n-1}^1)，若x_i^A= x_i^{B_n}= x_i^1,1\leq i\leq n且i\neq n. P(A→Bn)=π(xnBn∣x11,x21,...,xn−11)，若xiA=xiBn=xi1,1≤i≤n且i=n.
- P ( A → B n − 1 ) = π ( x n − 1 B n − 1 ∣ x 1 1 , x 2 1 , . . . , x n − 2 1 , x n 1 ) ，若 x i A = x i B n − 1 = x i 1 , 1 ≤ i ≤ n 且 i ≠ n − 1. P(A \rightarrow B_{n-1})=\pi(x_{n-1}^{B_{n-1}}|x_1^1,x_2^1,...,x_{n-2}^1,x_{n}^1)，若x_i^A= x_i^{B_{n-1}}= x_i^1,1\leq i\leq n且i\neq n-1. P(A→Bn−1)=π(xn−1Bn−1∣x11,x21,...,xn−21,xn1)，若xiA=xiBn−1=xi1,1≤i≤n且i=n−1.
- …
- P ( A → B j ) = π ( x j B j ∣ x 1 1 , x 2 1 , . . . , x j − 1 1 , x j + 1 1 , . . . , x n 1 ) ，若 x i A = x i B j = x i 1 , 1 ≤ i ≤ n 且 i ≠ j . P(A \rightarrow B_{j})=\pi(x_j^{B_{j}}|x_1^1,x_2^1,...,x_{j-1}^1,x_{j+1}^1,...,x_{n}^1)，若x_i^A= x_i^{B_{j}}= x_i^1,1\leq i\leq n且i\neq j. P(A→Bj)=π(xjBj∣x11,x21,...,xj−11,xj+11,...,xn1)，若xiA=xiBj=xi1,1≤i≤n且i=j.
- …
- P ( A → B 1 ) = π ( x 1 B 1 ∣ x 2 1 , x 3 1 , . . . , x n 1 ) ，若 x i A = x i B 1 = x i 1 , 1 ≤ i ≤ n 且 i ≠ 1. P(A \rightarrow B_{1})=\pi(x_1^{B_{1}}|x_2^1,x_3^1,...,x_{n}^1)，若x_i^A= x_i^{B_{1}}= x_i^1,1\leq i\leq n且i\neq 1. P(A→B1)=π(x1B1∣x21,x31,...,xn1)，若xiA=xiB1=xi1,1≤i≤n且i=1.
- P ( A → D ) = 0 ，其他情况 . P(A \rightarrow D)=0，其他情况. P(A→D)=0，其他情况.
那么对任意两点E，F 将满足细致平稳条件, π ( E ) P ( E → F ) = π ( F ) P ( F → E ) \pi(E)P(E\rightarrow F)=\pi(F)P(F\rightarrow E) π(E)P(E→F)=π(F)P(F→E)
证明过程同二维情况。
综上，状态转移矩阵P就变成了完全条件概率：
A
=
(
x
1
1
,
x
2
1
,
.
.
.
,
x
n
1
)
A=(x_1^1,x_2^1,...,x_n^1)
A=(x11,x21,...,xn1),
B
j
=
(
x
1
1
,
x
2
1
,
.
.
.
x
j
−
1
1
,
x
j
2
,
x
j
+
1
1
,
.
.
.
,
x
n
1
)
B_j=(x_1^1,x_2^1,...x_{j-1}^1,x_{j}^2,x_{j+1}^1,...,x_n^1)
Bj=(x11,x21,...xj−11,xj2,xj+11,...,xn1)
P
(
A
→
B
j
)
=
π
(
x
j
B
j
∣
x
1
1
,
x
2
1
,
.
.
.
,
x
j
−
1
1
,
x
j
+
1
1
,
.
.
.
,
x
n
1
)
P(A \rightarrow B_{j})=\pi(x_j^{B_{j}}|x_1^1,x_2^1,...,x_{j-1}^1,x_{j+1}^1,...,x_{n}^1)
P(A→Bj)=π(xjBj∣x11,x21,...,xj−11,xj+11,...,xn1)
此时，我们的采样过程就变成了按轴转移的过程了。想更新样本点的第j维，使用第j维的完全条件概率分布
π
(
x
j
n
e
w
∣
x
1
o
l
d
,
x
2
o
l
d
,
.
.
.
,
x
j
−
1
o
l
d
,
x
j
+
1
o
l
d
,
.
.
.
,
x
n
o
l
d
)
\pi(x_j^{new}|x_1^{old},x_2^{old},...,x_{j-1}^{old},x_{j+1}^{old},...,x_{n}^{old})
π(xjnew∣x1old,x2old,...,xj−1old,xj+1old,...,xnold),即可采样得到第j维新的值
x
j
n
e
w
x_j^{new}
xjnew。
所有维度的值都更新一遍，就能得到一个新的样本点
(
x
1
n
e
w
,
x
2
n
e
w
,
.
.
.
,
x
n
n
e
w
)
(x_1^{new},x_2^{new},...,x_n^{new})
(x1new,x2new,...,xnnew)
具体采样过程如下：
上面过程是依次轮换坐标轴进行采样。
实际上，按照随机顺序变换坐标轴也是可以的。
该文主要参考下面四篇文章，并加上一些自己的理解和推理细节，感谢作者，写的真好：
MCMC(一)蒙特卡罗方法
MCMC(二)马尔科夫链
MCMC(三)MCMC采样和M-H采样
MCMC(四)Gibbs采样
