随机模拟
合集下载
相关主题
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
• 利用软件产生一组[0,24]上的随机数xi: 10.4 5.3 0.3 9.0 20.1 6.8 7.8 …… 5.4 分别代表甲船到达码头时间;
• 再产生一组[0,24]上的随机数yi: 6.0 3.4 8.1 17.5 0.8 13.4 ……14.0 分别代表乙船到达码头时间;
• 把(xi,yi)(i=1,2…)当作一天中甲船乙船到达的时间 ,统计出能相遇的频数,计算出相遇的频率做为相遇的概率 。
2
Biblioteka Baidu
计算函数:z = norminv(p,mu,sigma)
统计直方图matlab函数hist
2000
• 直方图绘图函数: hist(data,n) 1000
其中,data是需要处理的数据块,
0 1 234 5
绘图原理:利用data中最小数和最大数构成一区间,
将区间等分为n个小区间,统计落入每个小区间的
• 把 [0,1]分成长度为0.34、0.36、0.30的三个区间 [0,0.34]、 (0.34,0.70]、(0.70,1]
• 用rand(1,1)产生1个[0,1]上均匀分布随机数, 如该数在[0,0.34]、 (0.34,0.70]或(0.70,1]内,相当于该天的需 求量相应为500、510和520
三、一些随机模拟的例子
[例1](相遇问题) 甲、乙两船在24小时内独立地随机到达码头. 设两船到达码头时刻都服从[0,24]上的均匀分布,如果甲船到 达码头后停留2小时,乙船到达码头后停留1小时.问两船相遇的 概率有多大?(可用几何概率,此处略)
分析:如果知道甲、乙两船到达的时刻x和y,两船能相遇的条件就
if(x>n) y=y+n*(a-b); else y=y+x*(a-b)-(n-x)*(b-c); end
根据随机数t, 计算需求量x值
根据需求量x, 计算利润并累 加到y中。
End
显示平均利润
零件参数
[例3] 粒子分离器的某关键参数(记作y)由7个零件的 参数(记作x1,x2,…,x7)决定,
二、Matlab中的部分随机数产生命令
注:以下都是产生不同分布m×n 阶随机矩阵
命令
说明
rand(m,n) unifrnd(a,b,m,n)
[0,1]上均匀分布 [a,b]上均匀分布
unifrnd(N,m,n)
1,…,N的等概率分布
randn(m,n)
标准正态分布N(0,1)
exprnd(λ,m,n)
y=0;
for i=1:N
y=y+fun2(x(i));
end
y/N
报童诀窍的简化版
也可完整模拟程序:
售出价a 购进价b 退回价c
a=1;b=0.75;c=0.6;n=510; n=510; N=1000; y=0;
购进数n
总利润y
for i=1:N
模拟天数N
t=rand(1,1); if (t<=0.34) x=500; elseif (t<=0.70) x=510; else x=520; end
计算密度值函数:normpdf(x,mu,sigma)
• 累积分布函数,即积分上限函数
P{ X x} 1 x exp[ (t )2 / 2 2 ]dt
2
计算概率值函数:normcdf(x,mu,sigma)
• 逆累积分布函数值,即已知概率值p,求z使得
P{ X z} 1 z exp[ (t )2 / 2 2 ]dt p
,用平均利润来近似报童的平均收入。
这也是Monte Carlo方法.
问 •如何按分布规律产生随机数据? 题 •随机数据很多时,如何编程?
报童诀窍的简化版 问题1:如何产生以下分布规律的随机数据?
需求量X 500 510 520 注:rand(m,n)可以生成 概率P 0.34 0.36 0.30 [0,1]上均匀分布随机数
else X=520;
end
y(i)=X;
end
文件名为randfun1.m
产生1行N列的全0矩阵,目的 分配好矩阵大小,可以省略
根据随机数t的范围,确定需求量X 值,并保存到数组的相应位置中( 关键部分)
报童诀窍的简化版 问题2:如何从需求量计算利润?售出价a=1, 购进价b=0.75,退回价c=0.6,购进数n=510
是两:船到达的时刻x和y可以随机生x 成 ,y生成x一组2且到达y时刻x,可y以确1定
是否能相遇。 如果重复很多次,统计相遇的比例就可近似为相遇的概率。
这也是Monte Carlo思想.
两船到达码头时刻服从[0,24]上的均匀分布,甲船停留2小时,乙 船停留1小时,相遇概率?
方法一 可以使用软件(excel、C等实现)
x1 x2 x3 x4 x5 x6 x7 标 定 值 0.2 0.3 0.1 0.1 1.5 16 0.75 容差等级 B B B C C B B
容差等级 A=1% B=5% C=15%
均值=标定值 标准差=容差等级*标定值÷3
xi~N (ai , bi2 )
零件损失Q是一个随机变量,平均损失就是期望E(Q)
函数文件名fun2.m function y=fun2(x)
a=1; %售出价 b=0.75; %购进价 c=0.6; %退回价 n=510; %购进数量
if(x>n) y=n*(a-b);
else y=x*(a-b)-(n-x)*(b-c)
end
模拟程n序(a代 b码)
xn
y xNx(==a1r0a0nb0d) ;fu(nn1(xN))(;b c) x n
计算机随机模拟(Monte Carlo)
(建模培训)
吕佳
一、简介
有些问题,由于随机因素很多,用概率论的方法进行求解 可能很难和复杂,这时就需要借助随机模拟方法得到近似解答 。
随机模拟法也叫蒙特卡罗(Monte Carlo)方法,也称为计 算机随机模拟方法。
由于Monte Carlo法计算量大,精度不高,因而需要借助 计算机,并仅适合一些用解析方法或常规数值方法难以解决问 题的低精度求解和验证
• 重复多次就可以若干天的需求量了
生成N=20天的需求量的matlab代码可以为
报童诀窍的简化版
生成N=20天的需求量的matlab代码(定义函数)
function y=randfun1(N)
y=zeros(1,N);
for i=1:N
t=rand(1,1);
if t<=0.34 X=500;
elseif t<=0.70 X=510;
看着是正态随机变量,在后面的标定值及容差等级情况
下,求产品的平均损失?
x1 x2 x3 x4 x5 x6 x7 标 定 值 0.2 0.3 0.1 0.1 1.5 16 0.75 容差等级 B B B C C B B
容差等级 A=1% B=5% C=15%
均值=标定值 标准差=容差等级*标定值÷3
零件参数
模拟:通过产生指定分布的随机数,来代表7个零件的 参数值,计算y值,确定损失大小。多做几次可得均值
零件参数
%编写计算y值的函数:y=funY(x1,x2,x3,x4,x5,x6,x7) %编写计算损失费函数:Q=funQ(y)
D1=[0.1,0.3,0.1,0.1,1.5,16,0.75] %标定值
均值为λ的指数分布
poissrnd(λ,m,n)
均值为λ的泊松分布
normrnd(μ,σ,m,n) 正态分布N(μ,σ2)
有关正态分布matlab计算函数
• 正态分布变量X的数学期望,方差 2 ,密度函数
f ( x) 1 exp[ ( x )2 / 2 2 ) 2
x (, )
两船到达码头时刻服从[0,24]上的均匀分布,甲船停留2小时,乙 船停留1小时,相遇概率?
方法二 %计算机模拟程序
clc; N=1000; c=0; %模拟次数N,相遇次数c清零
for i=1:N
%重复N次到达时间
x=unifrnd(0,24,1,1); %甲船到达时间x(随机数)
y=unifrnd(0,24,1,1); %乙船到达时间y(随机数)
y=174.42( x1 )( x3 )0.85 x5 x2 x1
3
1
2.62
1
0.36(
x4 x2
)0.56
2
(
x4 x2
)1.16
x6 x7
关键参数y的目标值是1.50,当偏离为±0.1时,产品
为次品,损失为1000元;当偏离为±0.3时,产品为废
品,损失为9000元。由于工艺原因,7个零件参数可以
D2=[5,5,5,15,15,5,5,]./100
%容差等级
a=D1 b=D1.*D2./3
零件期望a 零件标准差b
x1=normrnd(a(1),b(1),1,1); 按指定分布产生随机
……
数,作为零件参数
x7=normrnd(a(7),b(7),1,1)
y=funY(x1,x2,x3,x4,x5,x6,x7) Q=funQ(y)
X 500 510 520 P 0.34 0.36 0.30
问:如果报童每天购进报纸为n=510份,每天的平均利润是多少 ?
方法一:概率方法(略) 方法二:如果我们知道每天的需求量,可直接计算利润。 而每天需求量可以按分布生成(随机模拟思想)
售出价a=1,购进价b=0.75,退回价c=0.6,
报童诀窍的简化版
if((x<=y & y<=x+2) |(y<=x & x<=y+1))
c=c+1;
%如果能相遇,则计数器加1
end
end P=c/N %显示相遇的概率近似值
注意:每次运行的结果一般都不一样
报童诀窍的简化版
[例2](报童诀窍的简化版)报童每天清晨从报社购进报纸零
售,晚上退回没有卖掉的报纸.若每份报纸的购进价为b=0.75元 ,售出价为a=1元,退回价为c=0.6元.每天需求量X是离散型随 机变量,其分布为
E(Q) 0 P(Q 0) 1000 P(Q 1000) 9000 P(Q 9000)
1000 P(0.1 | y 1.5 | 0.3) 9000 P(| y 1.5 | 0.3)
由于y的表达式很复杂,要想计算y的分布和上述概 率很困难,我们必须寻找较为有效的近似方法。
购进数量n=510份
需求量X 500 510 520
概率P 0.34 0.36 0.30
按照需求量的分布规律,随机生成N=20个数据:
510 520 500 510 520 500 500 500 500
510 500 500 500 520 510 520 510 510 代表20天的需求量,计算出报童在这20天的总利润和平均利润
三到达B站的时间分布为:
问:到张达三时能赶间上该1趟3:火28车的1概3:率30是多1少3:32 13:34 概率 0.3 0.4 0.2 0.1
计算该产品的关键参 数y值和损失的大小
产生一组零件参数,相当制作一个产品。重复N=1000次
(即生产N个产品),求出损失费用的平均.(代码略)
四、练习
一列火车在下午1点后离开A站前往B站,火车离开A站的规律 如下:
离站时间 13:00 13 : 05 13 :10 火车从概A到率B所需时间平0.均7为30分钟0.,2标准差为02.1分钟。如张
数据量,以数据量为高度绘小矩形,形成直方图。
如果省略参数n,MATLAB将n的默认值取为10。
• 直方图统计计算函数:N=hist(data,n) 计算结果N是n个数的一维数组,分别表示data中各 个小区间的数据量。这种方式只计算而不绘图。
• 其他用法:hist(data,x)或 N= hist(data,x)
• 再产生一组[0,24]上的随机数yi: 6.0 3.4 8.1 17.5 0.8 13.4 ……14.0 分别代表乙船到达码头时间;
• 把(xi,yi)(i=1,2…)当作一天中甲船乙船到达的时间 ,统计出能相遇的频数,计算出相遇的频率做为相遇的概率 。
2
Biblioteka Baidu
计算函数:z = norminv(p,mu,sigma)
统计直方图matlab函数hist
2000
• 直方图绘图函数: hist(data,n) 1000
其中,data是需要处理的数据块,
0 1 234 5
绘图原理:利用data中最小数和最大数构成一区间,
将区间等分为n个小区间,统计落入每个小区间的
• 把 [0,1]分成长度为0.34、0.36、0.30的三个区间 [0,0.34]、 (0.34,0.70]、(0.70,1]
• 用rand(1,1)产生1个[0,1]上均匀分布随机数, 如该数在[0,0.34]、 (0.34,0.70]或(0.70,1]内,相当于该天的需 求量相应为500、510和520
三、一些随机模拟的例子
[例1](相遇问题) 甲、乙两船在24小时内独立地随机到达码头. 设两船到达码头时刻都服从[0,24]上的均匀分布,如果甲船到 达码头后停留2小时,乙船到达码头后停留1小时.问两船相遇的 概率有多大?(可用几何概率,此处略)
分析:如果知道甲、乙两船到达的时刻x和y,两船能相遇的条件就
if(x>n) y=y+n*(a-b); else y=y+x*(a-b)-(n-x)*(b-c); end
根据随机数t, 计算需求量x值
根据需求量x, 计算利润并累 加到y中。
End
显示平均利润
零件参数
[例3] 粒子分离器的某关键参数(记作y)由7个零件的 参数(记作x1,x2,…,x7)决定,
二、Matlab中的部分随机数产生命令
注:以下都是产生不同分布m×n 阶随机矩阵
命令
说明
rand(m,n) unifrnd(a,b,m,n)
[0,1]上均匀分布 [a,b]上均匀分布
unifrnd(N,m,n)
1,…,N的等概率分布
randn(m,n)
标准正态分布N(0,1)
exprnd(λ,m,n)
y=0;
for i=1:N
y=y+fun2(x(i));
end
y/N
报童诀窍的简化版
也可完整模拟程序:
售出价a 购进价b 退回价c
a=1;b=0.75;c=0.6;n=510; n=510; N=1000; y=0;
购进数n
总利润y
for i=1:N
模拟天数N
t=rand(1,1); if (t<=0.34) x=500; elseif (t<=0.70) x=510; else x=520; end
计算密度值函数:normpdf(x,mu,sigma)
• 累积分布函数,即积分上限函数
P{ X x} 1 x exp[ (t )2 / 2 2 ]dt
2
计算概率值函数:normcdf(x,mu,sigma)
• 逆累积分布函数值,即已知概率值p,求z使得
P{ X z} 1 z exp[ (t )2 / 2 2 ]dt p
,用平均利润来近似报童的平均收入。
这也是Monte Carlo方法.
问 •如何按分布规律产生随机数据? 题 •随机数据很多时,如何编程?
报童诀窍的简化版 问题1:如何产生以下分布规律的随机数据?
需求量X 500 510 520 注:rand(m,n)可以生成 概率P 0.34 0.36 0.30 [0,1]上均匀分布随机数
else X=520;
end
y(i)=X;
end
文件名为randfun1.m
产生1行N列的全0矩阵,目的 分配好矩阵大小,可以省略
根据随机数t的范围,确定需求量X 值,并保存到数组的相应位置中( 关键部分)
报童诀窍的简化版 问题2:如何从需求量计算利润?售出价a=1, 购进价b=0.75,退回价c=0.6,购进数n=510
是两:船到达的时刻x和y可以随机生x 成 ,y生成x一组2且到达y时刻x,可y以确1定
是否能相遇。 如果重复很多次,统计相遇的比例就可近似为相遇的概率。
这也是Monte Carlo思想.
两船到达码头时刻服从[0,24]上的均匀分布,甲船停留2小时,乙 船停留1小时,相遇概率?
方法一 可以使用软件(excel、C等实现)
x1 x2 x3 x4 x5 x6 x7 标 定 值 0.2 0.3 0.1 0.1 1.5 16 0.75 容差等级 B B B C C B B
容差等级 A=1% B=5% C=15%
均值=标定值 标准差=容差等级*标定值÷3
xi~N (ai , bi2 )
零件损失Q是一个随机变量,平均损失就是期望E(Q)
函数文件名fun2.m function y=fun2(x)
a=1; %售出价 b=0.75; %购进价 c=0.6; %退回价 n=510; %购进数量
if(x>n) y=n*(a-b);
else y=x*(a-b)-(n-x)*(b-c)
end
模拟程n序(a代 b码)
xn
y xNx(==a1r0a0nb0d) ;fu(nn1(xN))(;b c) x n
计算机随机模拟(Monte Carlo)
(建模培训)
吕佳
一、简介
有些问题,由于随机因素很多,用概率论的方法进行求解 可能很难和复杂,这时就需要借助随机模拟方法得到近似解答 。
随机模拟法也叫蒙特卡罗(Monte Carlo)方法,也称为计 算机随机模拟方法。
由于Monte Carlo法计算量大,精度不高,因而需要借助 计算机,并仅适合一些用解析方法或常规数值方法难以解决问 题的低精度求解和验证
• 重复多次就可以若干天的需求量了
生成N=20天的需求量的matlab代码可以为
报童诀窍的简化版
生成N=20天的需求量的matlab代码(定义函数)
function y=randfun1(N)
y=zeros(1,N);
for i=1:N
t=rand(1,1);
if t<=0.34 X=500;
elseif t<=0.70 X=510;
看着是正态随机变量,在后面的标定值及容差等级情况
下,求产品的平均损失?
x1 x2 x3 x4 x5 x6 x7 标 定 值 0.2 0.3 0.1 0.1 1.5 16 0.75 容差等级 B B B C C B B
容差等级 A=1% B=5% C=15%
均值=标定值 标准差=容差等级*标定值÷3
零件参数
模拟:通过产生指定分布的随机数,来代表7个零件的 参数值,计算y值,确定损失大小。多做几次可得均值
零件参数
%编写计算y值的函数:y=funY(x1,x2,x3,x4,x5,x6,x7) %编写计算损失费函数:Q=funQ(y)
D1=[0.1,0.3,0.1,0.1,1.5,16,0.75] %标定值
均值为λ的指数分布
poissrnd(λ,m,n)
均值为λ的泊松分布
normrnd(μ,σ,m,n) 正态分布N(μ,σ2)
有关正态分布matlab计算函数
• 正态分布变量X的数学期望,方差 2 ,密度函数
f ( x) 1 exp[ ( x )2 / 2 2 ) 2
x (, )
两船到达码头时刻服从[0,24]上的均匀分布,甲船停留2小时,乙 船停留1小时,相遇概率?
方法二 %计算机模拟程序
clc; N=1000; c=0; %模拟次数N,相遇次数c清零
for i=1:N
%重复N次到达时间
x=unifrnd(0,24,1,1); %甲船到达时间x(随机数)
y=unifrnd(0,24,1,1); %乙船到达时间y(随机数)
y=174.42( x1 )( x3 )0.85 x5 x2 x1
3
1
2.62
1
0.36(
x4 x2
)0.56
2
(
x4 x2
)1.16
x6 x7
关键参数y的目标值是1.50,当偏离为±0.1时,产品
为次品,损失为1000元;当偏离为±0.3时,产品为废
品,损失为9000元。由于工艺原因,7个零件参数可以
D2=[5,5,5,15,15,5,5,]./100
%容差等级
a=D1 b=D1.*D2./3
零件期望a 零件标准差b
x1=normrnd(a(1),b(1),1,1); 按指定分布产生随机
……
数,作为零件参数
x7=normrnd(a(7),b(7),1,1)
y=funY(x1,x2,x3,x4,x5,x6,x7) Q=funQ(y)
X 500 510 520 P 0.34 0.36 0.30
问:如果报童每天购进报纸为n=510份,每天的平均利润是多少 ?
方法一:概率方法(略) 方法二:如果我们知道每天的需求量,可直接计算利润。 而每天需求量可以按分布生成(随机模拟思想)
售出价a=1,购进价b=0.75,退回价c=0.6,
报童诀窍的简化版
if((x<=y & y<=x+2) |(y<=x & x<=y+1))
c=c+1;
%如果能相遇,则计数器加1
end
end P=c/N %显示相遇的概率近似值
注意:每次运行的结果一般都不一样
报童诀窍的简化版
[例2](报童诀窍的简化版)报童每天清晨从报社购进报纸零
售,晚上退回没有卖掉的报纸.若每份报纸的购进价为b=0.75元 ,售出价为a=1元,退回价为c=0.6元.每天需求量X是离散型随 机变量,其分布为
E(Q) 0 P(Q 0) 1000 P(Q 1000) 9000 P(Q 9000)
1000 P(0.1 | y 1.5 | 0.3) 9000 P(| y 1.5 | 0.3)
由于y的表达式很复杂,要想计算y的分布和上述概 率很困难,我们必须寻找较为有效的近似方法。
购进数量n=510份
需求量X 500 510 520
概率P 0.34 0.36 0.30
按照需求量的分布规律,随机生成N=20个数据:
510 520 500 510 520 500 500 500 500
510 500 500 500 520 510 520 510 510 代表20天的需求量,计算出报童在这20天的总利润和平均利润
三到达B站的时间分布为:
问:到张达三时能赶间上该1趟3:火28车的1概3:率30是多1少3:32 13:34 概率 0.3 0.4 0.2 0.1
计算该产品的关键参 数y值和损失的大小
产生一组零件参数,相当制作一个产品。重复N=1000次
(即生产N个产品),求出损失费用的平均.(代码略)
四、练习
一列火车在下午1点后离开A站前往B站,火车离开A站的规律 如下:
离站时间 13:00 13 : 05 13 :10 火车从概A到率B所需时间平0.均7为30分钟0.,2标准差为02.1分钟。如张
数据量,以数据量为高度绘小矩形,形成直方图。
如果省略参数n,MATLAB将n的默认值取为10。
• 直方图统计计算函数:N=hist(data,n) 计算结果N是n个数的一维数组,分别表示data中各 个小区间的数据量。这种方式只计算而不绘图。
• 其他用法:hist(data,x)或 N= hist(data,x)