蒙特卡罗随机模拟

Monte Carlo方法也叫随机模拟、随机抽样或者统计实验方法。其主要用途用于模拟一些无法用数值产生的随机系统。比如当系统的各个单元的特征量已知,但系统过于复杂导致无法预测其准确数学模型,这个时候可以用随机模拟法计算系统的相关参数。

蒙特卡罗的基本思想,为了求解数学、物理、生产管理等方面的问题,首先需要建立一个概率模型或者一个抽样试验来计算所求参数的统计特征,最后通过计算机模拟求出近似解。

我们说蒙特卡罗模拟的是随机现象,那么不同的问题的随机现象是不同的。比如掷骰子,各个点数等概率出现,模拟的均匀分布。去某个地方排队,人来的过程是一个随时间陆续变化的,这种一般认为符合指数分布。再有正态分布的例子太多,一般的统计数据没特殊要求都认为是服从正态分布的。每个分布都有对应的函数来表示。

假设一个随机事件在一次实验中发生的概率为P(x),显然0≤p(x)≤1

  • 一次贝努利实验:

P(ξ=x)={p,1?p,x=1x=0

  • . 二项分布

    当把贝努利实验重复多做几次就是二项分布了

    P(X=ak)=Cknpk(1?p)n?k,(p<1)

  • 离散均匀分布

    以相同的概率取所有可能的值

    P(X=ak)=1/n,k=1,2,3...n

  • 泊松分布

    发生概率低,同时次数无限增加的贝努利分布就是泊松分布

    P(X=k)=λkk!e?λ,k=1,2,...,n

再看一些密度函数

  • 均匀分布密度函数

    在区间[a,b]上等概率出现,那么出现的概率:

    P(x)=???1b?a0,x∈[a,b]other

  • 正态分布密度函数

    用的最多

    P(x)=12π??√σexp(?(x?u)22σ2)

  • 指数分布密度函数

    P(x)={λe?λx0x≥0x<0

下面采用matlab进行蒙特卡罗的一些小实验

首先关于matlab产生上述具有特定分布的随机数都有特定的函数,常用的有三种:

  • rand:生成均匀分布的伪随机数。分布在(0~1)之间

    主要语法:rand(m,n)生成m行n列的均匀分布的伪随机数

    rand(m,n,’double’)生成指定精度的均匀分布的伪随机数,参 数还可以是’single’, rand(RandStream,m,n)利用指定的RandStream(我理解为随 机种子)生成伪随机数

  • randn 生成标准正态分布的伪随机数(均值为0,方差为1)

    主要语法:和上面一样

  • randi 生成均匀分布的伪随机整数

    主要语法:randi(iMax)在开区间(0,iMax)生成均匀分布的伪随机整数randi(iMax,m,n)在开区间(0,iMax)生成m*n型随机矩 阵r = randi([iMin,iMax],m,n)在开区间(iMin,iMax)生成 m*n型随机矩阵

还有一些其他的生成函数:

unifrnd():和rand()类似,这个函数生成某个区间内均匀分布的随机数

normrnd():和randn()类似,此函数生成指定均值、标准差的正态分布的随机数

chi2rnd():函数生成服从卡方(Chi-square)分布的随机数。卡方分布只有一个参数:自由度v

frnd():函数生成服从F分布的随机数。F分布有2个参数:v1, v2

trnd():函数生成服从t分布的随机数。t分布有1个参数:自由度v

betarnd():此函数生成服从Beta分布的随机数。Beta分布有两个参数分别是A和B

exprnd():此函数生成服从指数分布的随机数。指数分布只有一个参数: mu

gamrnd():生成服从Gamma分布的随机数。Gamma分布有两个参数:A和B

例如:用蒙塔卡罗方法计算图像的面积,这里为了方便起见计算一下圆的面积,顺便用这种方法预估一下π的值。当然这是规则图像,可以用公式来计算,当遇到的图像不是规则的图形的时候,就非常有用了。

代码如下:

clc
clear
%生成num个(-1,1)的均匀点
num = 100000;
data = 2*rand(2,num)-1;
plot(data(1,:),data(2,:),‘*‘)
hold on;
%
num_in_circle = 0;
radius = 1;%小于1
%plot circle
theta = 0:0.1:2*pi;
plot(radius*cos(theta),radius*sin(theta),‘r‘,‘LineWidth‘,3);
%test samples
for i=1:num
    %随机点在圆内
    if (data(1,i)^2 + data(2,i)^2) <= radius^2
        num_in_circle = num_in_circle + 1;
    end
end
area = 4*num_in_circle/num;%矩形面积为4
PI = area/(radius^2);%
title([‘num=‘,num2str(num),‘ 半径=‘,num2str(radius),‘ 模拟的PI值 = ‘,num2str(PI)])

可以看到num样本越多预测的越准,当再大些时,基本上可以和π相差不多了。

版权声明:本文为博主原创文章,未经博主允许不得转载。

时间: 2024-10-17 01:19:01

蒙特卡罗随机模拟的相关文章

随机模拟(MCMC)

http://cos.name/2013/01/lda-math-mcmc-and-gibbs-sampling/ http://blog.csdn.net/lin360580306/article/details/51240398 随机模拟(或者统计模拟)方法有一个很酷的别名是蒙特卡罗方法(Monte Carlo Simulation).这个方法的发展始于20世纪40年代,和原子弹制造的曼哈顿计划密切相关,当时的几个大牛,包括乌拉姆.冯.诺依曼.费米.费曼.Nicholas Metropoli

【Matlab编程】马氏链随机模拟

本文是利用蒙特卡罗算法对马氏链过程的模拟.假设有10个状态,从每个状态到与之相邻状态的概率是相同的,仿真次数为1000,及进行了1000次状态转移.我们以动画的形式再现了状态转移的过程,并记录了到达每个状态的次数,具体实现如下: close all;clc;clear; figure; s=1; n=1000; r=1; % 圆圈的半径 title('等概率情况的计算机模拟') set(gcf,'doublebuffer','on'); % 设置图形渲染效果 xlabel('Please pre

随机模拟与采样方法

随机模拟方法简介 随机模拟方法又称为蒙特卡罗方法(Monte Carlo Method).蒙特卡洛模拟方法的原理是当问题或对象本身具有概率特征时,可以用计算机模拟的方法产生抽样结果,根据抽样计算统计量或者参数的值:随着模拟次数的增多,可以通过对各次统计量或参数的估计值求平均的方法得到稳定结论.由于涉及到时间序列的反复生成,蒙特卡洛模拟法是以高容量和高速度的计算机为前提条件的,因此只是在近些年才得到广泛推广. 从上面的描述我们不难看出随机模拟方法主要是针对那些确定算法不好解或者解不出来的情况.因此

随机模拟的基本思想和常用采样方法(sampling)

转自:http://blog.csdn.net/xianlingmao/article/details/7768833 引入 我们会遇到很多问题无法用分析的方法来求得精确解,例如由于式子特别,真的解不出来.这时就需要找一种方法求其近似解,并且有手段能测量出这种解的近似程度 (比如渐进性,上下限什么的) 随机模拟的基本思想 现在假设我们有一个矩形的区域R(大小已知),在这个区域中有一个不规则的区域M(即不能通过公式直接计算出来),现在要求取M的面积? 怎么求?近似的方法很多,例如:把这个不规则的区

随机采样和随机模拟:吉布斯采样Gibbs Sampling的具体实现

http://blog.csdn.net/pipisorry/article/details/51525308 吉布斯采样的实现问题 本文主要说明如何通过吉布斯采样进行文档分类(聚类),当然更复杂的实现可以看看吉布斯采样是如何采样LDA主题分布的[主题模型TopicModel:隐含狄利克雷分布LDA]. 关于吉布斯采样的介绍文章都停止在吉布斯采样的详细描述上,如随机采样和随机模拟:吉布斯采样Gibbs Sampling(why)但并没有说明吉布斯采样到底如何实现的(how)? 也就是具体怎么实现

随机采样和随机模拟:吉布斯采样Gibbs Sampling

http://blog.csdn.net/pipisorry/article/details/51373090 马氏链收敛定理 马氏链定理: 如果一个非周期马氏链具有转移概率矩阵P,且它的任何两个状态是连通的,那么limn→∞Pnij 存在且与i无关,记limn→∞Pnij=π(j), 我们有 limn→∞Pn=???????π(1)π(1)?π(1)?π(2)π(2)?π(2)??????π(j)π(j)?π(j)????????????? π(j)=∑i=0∞π(i)Pij π 是方程 πP

python实现简单随机模拟——抛呀抛硬币

还是在上次提到的数据之魅那本书,看到模拟这章,有个python模拟脚本,但书上不全,就自己简单写了下. 流程:在不同的平衡参数p(为0.5时为均匀的)下,模拟60次实验,每次投硬币8次,统计正面朝上的次数,并作图. import random import matplotlib.pyplot as plt repeats, tosses = 60, 8 # p为平衡参数,tosses为每次重复试验中投掷硬币的次数 # 返回当前平衡参数p的情况下,8次实验中正面的次数 def heads(toss

【数据分析 R语言实战】学习笔记 第三章 数据预处理 (下)

3.3缺失值处理 R中缺失值以NA表示,判断数据是否存在缺失值的函数有两个,最基本的函数是is.na()它可以应用于向量.数据框等多种对象,返回逻辑值. > attach(data) The following objects are masked fromdata (pos = 3): city, price, salary > data$salary=replace(salary,salary>5,NA) > is.na(salary) [1] FALSEFALSE TRUE

用蒙特卡罗模拟球π的值

蒙特卡罗(Monte Carlo)方法,又称随机抽样或统计试验方法,属于计算数学的一个分支,它是在本世纪四十年代中期为了适应当时原子能事业的发展而发展起来的.传统的经验方法由于不能逼近真实的物理过程,很难得到满意的结果,而蒙特卡罗方法由于能够真实地模拟实际物理过程,故解决问题与实际非常符合,可以得到很圆满的结果.这也是我们采用该方法的原因. 基本原理及思想编辑 当所要求解的问题是某种事件出现的概率,或者是某个随机变量的期望值时,它们可以通过某种“试验”的方法,得到这种事件出现的频率,或者这个随机