用蒙特卡洛求pi_什么是马尔可夫链蒙特卡洛（MCMC）？

许多人在学习贝叶斯统计、机器学习的过程中都听说过马尔可夫链蒙特卡洛（Markov Chain Monte Carlo ），但是其概念却很难理解，希望通过这边整理的文章给大家一个直观理解。

蒙特卡洛马尔可夫链是又由两个概念组成的，分别是蒙特卡洛（Monte Carlo）和马尔可夫链（Markov Chains），下面就将分别介绍一下两个概念。

什么是蒙特卡洛？

蒙特卡洛方法来自于摩纳哥的蒙特卡洛赌场，许多纸牌类游戏需要计算其胜利的概率。我们可以将蒙特卡洛理解为简单的模拟，通过模拟的情景来计算其发生的概率。

应用场景：通常情况下，许多概率的计算非常复杂甚至是无法计算的，但是我们可以通过计算机对期待发生的场景进行大量的模拟，从而计算出其发生的概率。简单的公式表达为：

学习过概率的同学应当知道这就是典型的频率学派的观点。

对于蒙特卡洛方法，经典的例子就是计算

值。在（-1，1）之间随机取两个数，如果在单位圆内，则记一次，在圆外则不计入次数。

import random
import numpy as np
from math import sqrt, pi

num1 = 1000 # simulation times 
freq1 = 0
x1 = []
y1 = []
for i in range(1,num1+1):
        x_i, y_i = random.uniform(-1, 1), random.uniform(-1, 1)
        square_distance = x_i**2+y_i**2
        x1.append(x_i)
        y1.append(y_i)
        if square_distance <= 1: # 圆的方程x^2+y^2=1
                freq1 += 1

simulation_pi1 = 4*freq1/num1;

num2 = 10000 # simulation times 
freq2 = 0
x2 = []
y2 = []
for i in range(num2):
        x_i, y_i = random.uniform(-1, 1), random.uniform(-1, 1)
        square_distance = x_i**2+y_i**2
        x2.append(x_i)
        y2.append(y_i)
        if square_distance <= 1: # 圆的方程x^2+y^2=1
                freq2 += 1

simulation_pi2 = 4*freq2/num2;

% matplotlib inline
import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(8,8))

theta = np.linspace(0,2*pi,500)
x,y = np.cos(theta), np.sin(theta)
ax.plot(x, y, color='red', linewidth=1.0)

ax.scatter(x1,y1,color='blue',alpha=0.75)
ax.set_xlim(-1,1);
ax.set_ylim(-1,1);

结果为：

真实的pi值： 3.141592653589793
蒙特卡洛方法下得出的pi值（模拟1000次）： 3.188
相对误差（模拟1000次）： 1.4771917153924754 %
蒙特卡洛方法下得出的pi值（模拟10000次）： 3.1292
相对误差（模拟10000次）： 0.39447041536821975 %

小结：蒙特卡洛就是一种模拟事情发生的方法。

什么是马尔可夫链？

在学习马尔可夫链之前，我们首先了解一下什么是马尔可夫性质（或者叫马尔可夫特性，Markov Property）。

假设你有一个系统，这个系统具有

种可能的状态，并且在状态之间正在移动。真实世界中，我们常常见到这种例子，例如天气从炎热到温和，或者是股票市场从熊市到牛市再到停滞的状态（stagnant states）。

马尔可夫性质是指在给定的一个随机过程（stochastic process）中，在某一个时间点的状态为

，其下一个状态

的概率（

代表

种状态里的任何一种），仅仅取决于现在这个它在给定时刻处于哪一种状态，而与它是如何达到现在的