用计算机怎么做蒙特卡洛模拟,一篇文章讲明白蒙特卡洛模拟

先来个测试题：

用仿真的方法计算圆周率π，请写出计算步骤和简单的公式：

(思考十秒钟，想不出来可以继续往下看，后面有答案)

001 什么是蒙特卡洛模拟

蒙特卡洛模拟，又称为统计实验方法

以概率论和统计理论方法为基础的一种计算方法

通过随机数来解决很多计算问题

主要步骤是：

将实际问题转化为概率模型

通过计算机实现统计模拟，以获得问题的近似解

002 基本原理

蒙特卡洛模拟抓住事物运动的几何数量和几何特征，利用数学方法来模拟，是一种数字模拟实验。它是一个以概率模型为基础，按照这个模型所描绘的过程，通过模拟实验的结果，作为问题的近似解。

可以把蒙特卡洛解法分为三步

1.构造或描绘概率过程

2.实现从已知概率分布抽样

3.建立各种估计量

借助计算机技术，蒙特卡洛模拟方法有两大优点

1.简单，省去了繁复的数学推导和验算过程，使普通人能够理解

2.快速，确定了概率模型，后续运算完全用计算机实现

蒙特卡洛模拟的特点：

随机采样得到的近似解，随着随机采样数值增多，得到正确结果的概率越大

003 怎么用

回到最初的测试题，用蒙特卡洛模拟计算圆周率π

圆周率π的计算推导十分复杂，普通人难以明白，但是我们可以用数值模拟的方式计算出圆周率

在这个正方形内部有个内切圆，他们的面积之比是π/4

现在，在这个正方形内部随机生产10000个点(坐标是(x,y))，计算它们与中心点的距离，从而判断是否落在圆内。

园内的点数/正方形内的点数 = π/4

以此来计算得出π的值

直接上代码(python3.6)

import numpy as np

import pandas as pd

import matplotlib.pyplot as plt

from matplotlib.patches import Circle

%matplotlib inline

n = 10000

r = 1.0

a,b = (0.0,0.0)

xmin, xmax = a-r, a+r

ymin, ymax = b-r, b+r

x = np.random.uniform(xmin,xmax,n)

y = np.random.uniform(ymin,ymax,n)

fig = plt.figure(figsize=(6,6))

axes = fig.add_subplot(1,1,1)#添加子图

#画子图

plt.plot(x,y,'ko',markersize = 1) #plot绘图 markersize表示点的大小；‘ro’r表示red，o表示圆圈

plt.axis('equal')

d = np.sqrt((x-a)**2 + (y-a)**2)

#res 得到圆中的点数

res = sum(np.where(d

pi = res/n*4

print('pi:',pi)

#计算pi的近似值，蒙特卡洛模拟方法，用统计值去近似真实值

#绘制圆形子图

circle = Circle(xy = (a,b), radius = r,alpha = 0.5, color = 'gray')

axes.add_patch(circle)#添加圆形子图

plt.grid(True,linestyle = '--',linewidth = '0.8')

plt.show()

#总结

#蒙特卡洛模拟是用统计值逼近真实值，展示了统计思想

JupyterNotebook运行结果如下

计算圆周率

随着我们设置的点的数量n越多，计算的到的pi值就会越接近真实值3.1415926

好了，蒙特卡洛模拟就讲到这里。

这只是一个简单的例子，关于真实场景下运用蒙特卡洛模拟的例子，等我下一篇文章咯~

下一篇预告：用蒙特卡洛模拟人类财富分配