MCMC基础简说-Toy模板网

这篇具有很好参考价值的文章主要介绍了MCMC基础简说。希望对大家有所帮助。如果存在错误或未考虑完全的地方，请大家不吝赐教，您也可以点击"举报违法"按钮提交疑问。

MCMC是一种随机采样方法，由两个MC组成，即蒙特卡罗方法（Monte Carlo Simulation，简称MC）和马尔科夫链（Markov Chain ，也简称MC）。

蒙特卡罗

蒙特卡罗是一种随机模拟方法，这里以求积分为例。

MCMC基础简说对于x在区间[a, b]内均匀分布的积分，我们可以在区间内随机采样n点，用n个点坐标均值来模拟代表所有点的值，可以将积分变为求和式：

MCMC基础简说

可均匀分布只是函数分布的一小部分代表，为了让蒙特卡罗方法适用于更广的范围，我们想引入其他的概率分布模型来使均匀分布对积分式的改变成为特例。假定x在[a, b]的概率分布函数为P(x)，定积分求和就可以表示为：

MCMC基础简说

这样我们就可以将均匀分布变为特例，如P(X)=1/(b-a)就可以和上上图契合。

至此，关键问题来到了函数概率分布函数怎样获得的问题。对于一些常见的概率分布函数的采样还可以通过均匀分布经数式转换得，如高斯分布样本(Z1, Z2)通过均匀分布样本(X1, X2)计算：

MCMC基础简说

但不常见的函数是更多的，我们又进一步提出接受-拒绝采样方法。

接受-拒绝采样方法简单的说就是用已知简便的模型包含要采样的函数模型，再通过过滤掉一部分不符合需采样函数的已知模型的值来逼近需采样函数。 MCMC基础简说

设定一个方便采样的常用概率分布函数 q(x)，以及一个常量 k，使得p(x) 总在 kq(x) 的下方。随机取样一点x0，对y0进行判断，<=p(x0)才会统计为样本。

问题变为对多维分布蒙特卡罗方法无法求得概率分布，这时便是马尔科夫链登上的时候了。

马尔科夫链

马尔可夫链本身是假设某一时刻状态转移的概率只依赖于它的前一个状态，在一个序列模型下更迭。 MCMC基础简说

举一个马尔可夫链模型的具体例子：

MCMC基础简说

牛市（Bull market）, 熊市（Bear market）和横盘（Stagnant market）。各个单元对外指向上的数字代表着该状态转变为另一种状态的可能性，和为1(对内指向没有和为1的限制)。对上图定义状态转移矩阵：矩阵P某一位置P(i,j)的值为P(j|i),即从状态i转化到状态j的概率，并定义牛市为状态0，熊市为状态1, 横盘为状态2.

MCMC基础简说

输入一个初始状态概率分布[0.3, 0.4, 0.3]，令其与P相乘则可以得到下一状态的状态概率分布。为什么可以？具体点说就比如0.3x0.9代表着初始为0下一回合也为0的情况，状态转移矩阵就是转变概率，那么下一回合为0的概率就是0.3x0.9+0.4x0.15(初始为1转变为0) +0.3x0.25，对下回合是1或2同理。

import numpy as np

# 初始概率分布
t = np.array([[0.3, 0.4, 0.3]])

# 状态转移矩阵
p = np.array(
    [[0.9, 0.075, 0.025],
    [0.15, 0.8, 0.05],
    [0.25, 0.25, 0.5]]
)

n = 0
# 迭代，得到稳定概率分布
for i in range(1000):
    t = np.dot(t, p)
    n += 1
    print('t{}时刻的状态分布概率为:\n{}'.format(n, t))

MCMC基础简说