> For the complete documentation index, see [llms.txt](https://blessbingo.gitbook.io/garnet/llms.txt). Markdown versions of documentation pages are available by appending `.md` to page URLs; this page is available as [Markdown](https://blessbingo.gitbook.io/garnet/math/mcmc/0x02-cai-yang-shi-li.md).

# 0x02 采样实例

## M-H采样实例

这里的例子是一个**M-H采样**的实例.

* 这里已经得到了目标分布, 即**平稳分布**, 是一个**正态分布**, 均值3, 标准差2.
* 马尔可夫链**状态转移矩阵**$$Q$$的条件转移概率是: 以$$i$$为均值, 标准差为1的正态分布, $$j$$所在位置的概率.

**代码**:

$$\alpha(i,j) = min{ \frac{\pi(j)Q(j,i)}{\pi(i)Q(i,j)},1}$$

由于这里的$$Q$$指的是正态分布, 因此$$Q(i,j)=Q(j,i)$$, 则有$$\alpha(i,j) = min{ \frac{\pi(j)}{\pi(i)},1}$$

```python
import numpy as np
from scipy.stats import norm
import matplotlib.pyplot as plt

def final_dist_prob(value):
    return norm.pdf(value, loc=3, scale=2)


num_samples = 1000
t = 0

pi = 0  # 初始状态值
samples = []

while t < num_samples:
    pi_star = norm.rvs(loc=pi, scale=1, size=1)[0]  # 从条件概率分布中采样
    alpha = min(1, final_dist_prob(pi_star) / (final_dist_prob(pi) + 10e-8))  # 计算接受率

    u = np.random.rand()
    if u < alpha:
        samples.append(pi_star)
        pi = pi_star
        t += 1

x_range = np.linspace(-3, 9, 200)
plt.plot(x_range, final_dist_prob(x_range))
plt.hist(samples, bins=50, normed=1, alpha=0.5)

plt.show()
```

最后采样的结果图:

![image](http://m.qpic.cn/psb?/V12PEvw82GjuY9/sVQ.vVyfzRqWtHz5xJyUh4ZJzok3WgMTwoszIqV1q8E!/b/dIMAAAAAAAAA\&bo=gALZAQAAAAARB2o!\&rf=viewer_4)

可以看到采样得到的样本分布与已知的正态分布十分接近.

## Gibbs二维采样实例

[二维Gibbs采样实例](http://www.cnblogs.com/pinard/p/6645766.html)