文章目录
前言
通过模型算法,熟练对Matlab和python的应用。
学习视频链接:
https://www.bilibili.com/video/BV1EK41187QF?p=22&vd_source=67471d3a1b4f517b7a7964093e62f7e6
一、蒙特卡洛
- 蒙特卡洛严格意义上讲并不是一种算法模型,它是一种概率统计的思想方法,没有通用的代码,每个问题对应的代码都是不同的。
- 蒙特卡罗法又称统计模拟法,是一种随机模拟方法,以概率和统计理论方法为基础的一种计算方法,是使用随机数(或更常见伪随机数)来解决很多计算问题的方法。将所求解的问题同一定的概率模型相联系,用电子计算机实现统计模拟或抽样,以获得问题的近似解。为象征性地表明这一方法的概率统计特征,故借用赌城蒙特卡罗命名。
- 只要求解的问题与概率模型有关联,就可以用蒙特卡洛法。
二、经典示例:计算圆周率 π
- 一个经典的蒙特卡洛方法的示例是计算圆周率 π。通过在单位正方形内随机生成点,估计圆的面积来逼近 π 的值。
- 步骤:
1.在正方形[0,2]x[0,2] 内随机生成大量点。
2.计算这些点中有多少点落在单位圆(圆心为 (1,1),半径为 1)内。
3.单位圆的面积是 π,正方形的面积是 4,因此可以通过以下公式估计 π:
π ≈ 圆内点的数量 总点数 × 4 \pi\approx\frac{\text{圆内点的数量}}{\text{总点数}} \times 4 π≈总点数圆内点的数量×4
1.代码实现----Matlab
clear;clc
% 参数初始化:投放10000个点,圆的半径为1,圆心坐标(1,1)
% 初始时还未投放点,有0个点在圆内
p = 100000; r = 1; x0 = 1; y0 = 1; n = 0;
hold on
for i = 1:p
px = rand*2;
py = rand*2;
if (px-1)^2 + (py-1)^2 < 1
plot(px,py,'.','Color',"b")
n = n + 1;
else
plot(px,py,'.','Color',"r")
end
end
axis equal
s = (n/p)*4;
pi0 = s;
运行结果图:
2.代码实现----python
import matplotlib.pyplot as plt
import random
import numpy as np
# 计算圆周率
# 参数初始化:投放10000个点,圆的半径为1,圆心坐标(1,1)
# 初始时还未投放点,有0个点在圆内
p = 100000
r = 1
x0 = 1
y0 = 1
n = 0
# 设置绘图窗口
plt.figure()
plt.title('Monte Carlo Simulation for Esitimation Pi')
plt.xlabel('x')
plt.ylabel('y')
for i in range(p):
px = np.random.rand() * 2
py = np.random.rand() * 2
if (px - x0) ** 2 + (py - y0) ** 2 < r ** 2:
plt.plot(px,py,'.',color = 'b')
n = n + 1
else:
plt.plot(px,py,'.',color = 'r')
plt.axis ('equal')
plt.show()
s = n / p * 4
pi0 = s
运行结果图:
三、示例2:三门问题
- 你参加一档电视节目,节目组提供了A、B、C三扇门,主持人告诉你,其中一扇门后边有辆汽车,其他两扇门后面是一头山羊,你可以选择一扇门打开获得门后的东西。
假如你选择了B门,这时,主持人打开了C门,让你看到C门后是只山羊,然后问你要不要改选A门?(你想要汽车)
这个问题在概率论中是一个条件概率的问题,类似于在排除一个错误选项后,选择到正确选项的概率是多少。求解的问题与概率模型有关,我们使用蒙特卡洛法。
1.代码实现----Matlab
(1)有条件概率
clear;clc
% randi([a,b],m,n) 函数可在指定区间[a,b]内随机取出大小为m*n的整数矩阵
randi([1,5],5,8) % 在区间[1,5]内随机取出大小为5*8的整数矩阵
% 字符串的连接方式:(1)['字符串1','字符串2'] (2)strcat('字符串1','字符串2')
n = 100000;
a = 0; % 不改变注意能赢的次数
b = 0; % 改变注意能赢的次数
for i = 1:n
x = randi([1,3]); % 有车的那扇门
y = randi([1,3]); % 我们选择的那扇门
if x == y
a = a + 1; b = b + 0; % 不改变主意能赢
else
a = a + 0; b = b + 1; % 改变主意能赢
end
end
disp(['使用蒙特卡洛方法得到不改变主意能赢的概率是', num2str(a/n)]);
disp(['使用蒙特卡洛方法得到改变主意才能赢的概率是', num2str(b/n)]);
运行结果:
(2)无条件概率(考虑没有获奖的情况)
clear;clc
% 考虑失败的情况
n = 100000;
a = 0; % 不改变注意能赢的次数
b = 0; % 改变注意能赢的次数
c = 0; % 没有获奖的次数
for i = 1:n
x = randi([1,3]); % 有车的那扇门
y = randi([1,3]); % 我们选择的那扇门
change = randi([0,1]); % 改变主意是1,不改变主意是0
if x == y % 不改变才能赢
if change == 0 % 不改变主意是0
a = a + 1; b = b + 0; c = c + 0;
else
a = a + 0; b = b + 0; c = c + 1;
end
else % 改变才能赢
if change == 0 % 不改变主意是0
a = a + 0; b = b + 0; c = c + 1;
else % 改变主意是1
a = a + 0; b = b + 1; c = c + 0;
end
end
end
disp(['使用蒙特卡洛方法得到不改变主意能赢的概率是', num2str(a/n)]);
disp(['使用蒙特卡洛方法得到改变主意才能赢的概率是', num2str(b/n)]);
disp(['使用蒙特卡洛方法不获奖的概率是', num2str(c/n)]);
运行结果:
2.代码实现----python
(1)有条件概率
import numpy as np
n = 10000
a = 0 # 不改变主意能获胜的次数
b = 0 # 改变主意能获胜的次数
for i in range(n):
x = np.random.randint(1,4) #随机生成1-3的整数
y = np.random.randint(1,4)
if x == y:
a += 1
else:
b += 1
print("使用蒙特卡洛方法得到不改变主意能赢的概率是", a/n)
print("使用蒙特卡洛方法得到改变主意能赢的概率是", b/n)
运行结果:
(2)无条件概率(考虑没有获奖的情况)
import numpy as np
n = 10000
a = 0 # 不改变主意能获胜的次数
b = 0 # 改变主意能获胜的次数
c = 0 # 不获胜的次数
for i in range(n):
x = np.random.randint(1,4) #随机生成1-3的整数
y = np.random.randint(1,4)
change = np.random.randint(0,2) # % 改变主意是1,不改变主意是0
if x == y:
if change == 0:
a += 1
else:
c += 1
else:
if change == 0:
c += 1
else:
b += 1
print("使用蒙特卡洛方法得到不改变主意能赢的概率是", a/n)
print("使用蒙特卡洛方法得到改变主意能赢的概率是", b/n)
print("使用蒙特卡洛方法不能赢的概率是", c/n)
运行结果:
- 蒙特卡洛适用于概率统计的问题,但得出的结果只能用作估计,并不是一个精确的数值。代码中每次生成的随机数不同也会造成结果上的差异。使用蒙特卡洛法,当测试的数量足够大时,得出的结论才会更精确。
总结
本文介绍了使用蒙特卡洛法解决与概率统计有关联的问题。分析了两个示例,并分别使用Matlab和python进行代码编写。
标签:概率,python,else,学习,改变,Matlab,np,蒙特卡洛,主意 From: https://blog.csdn.net/m0_65032457/article/details/140880739