随机数

#


Monte Carlo方法

定义

Monte Carlo方法亦称为随机模拟(Random simulation)方法,有时也称作随机抽样技术或统计试验方法。此方法对研究的系统进行随机观察抽样,通过对样本值的观察统计,求得所研究系统的某些参数.

基本思想

当试验次数充分多时,某一事件出现的频率近似于该事件发生的概率.

应用举例

例1. 用蒙卡模拟求圆周率π的估计值

设二维随机变量(X,Y)在正方形服从均匀分布,则(X,Y)落在四分之一圆内的概率为 P{X2+Y2<=1}= π/4
现产生n对二维随机点(xi ,yi ), xi ,yi 是区间[0,1]上均匀分布的随机数,若其中有k对满足 xi2+yi2<=1,即相当于做n次投点试验,其中有k点落在1/4圆内。
计算落入1/4圆内的频率为k/n,由大数定律,事件发生的频率依概率收敛于发生的概率,可得圆周率π的估计值为4k/n。次数越多,精度越大。

Matlab代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
clear 
for i=1:5
%仿真5次
k=0; %k用于随机点落在1/4圆内的计数
for j=1:1000000 %样本个数取为N=1000000
a=rand(1,2); %生成区间(0,1)上的均匀分布随机数作取样值
if a(1)^2+a(2)^2<=1 %检查随机数是否满足:
k=k+1;
end
end
PI(i)=4*k/j; %计算π的近似值
end
PI
结果为:PI = 3.1435 3.1402 3.1396 3.1407 3.1400

例9. 报童的策略

*报童每天清晨从报社购进报纸零售,晚上将没有卖掉的报纸退回.每份报纸的购进价为1.3元,零售价为2元,退回价为0.2元.报童售出一份报纸赚0.7元 ,退回一份报纸赔1.1元.报童每天如果购进的报纸太少,不够卖时会少赚钱,如果购得太多卖不完时要赔钱. 试为报童筹划每天应如何确定购进的报纸数使得收益最大.报纸每捆10张,只能整捆购买,报纸可以分为3种类型的新闻日:好、一般、差,它们的概率分别为0.35,0.45和0.2,在这些新闻日中每天对报纸的需求分布的统计结果如下表:

需求量 好新闻的需求概率 一般新闻的需求概率 差新闻的需求概率
40 0.03 0.10 0.44
50 0.05 0.18 0.22
60 0.15 0.40 0.16
70 0.20 0.20 0.12
80 0.35 0.08 0.06
90 0.15 0.04 0.00
100 0.07 0.00 0.00

试确定每天报童应该订购的报纸数量
解:我们通过计算机仿真来解决此问题。最优策略应该是每天的利润最大。
*利润=销售收入-报纸成本+残值

这是一个随机现象的计算机仿真问题,故先确定各种情况的随机数的对应关系。

新闻日和需求量对应的随机数分别如下面两个表格所示:

新闻种类 出现概率 对应的随机数区间
好新闻 0.35 (0,0.35)
一般新闻 0.45 [0.35,0.80)
差新闻 0.20 [0.80,1)

需求量 好新闻的随机数区间 一般新闻的随机数区间 差新闻的随机数区间
40 (0.00,0.03] (0.00,0.10] (0.00,0.44]
50 [0.03,0.08) [0.10,0.28) [0.44,0.66)
60 [0.08,0.23) [0.28,0.68) [0.80,1.00)
70 [0.23,0.43) [0.68,0.88) [0.80,1.00)
80 [0.43,0.78) [0.88,0.96) [0.80,1.00)
90 [0.78,0.93) [0.96,1.00) [100,102]
100 [0.93,1.00] [8,100] [104,125]

计算机仿真的流程:
1)令每天的报纸订购数变化,40—100;
2)让时间从1开始变化(循环)到365;
3)产生新闻种类的随机数,确定当天的新闻类型;
4)产生需求量随机数,确定当天的报纸需求量;
5)计算当天的收入,计算累积利润,
6)比较得出最优定货量。
Matlab编程代码:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
clear all;clc 
x1=rand(365,1);
x2=rand(365,1);
for n=4:10
paper=n*10;
sb(n)=0;
for i=1:365
if x1(i)<0.35
if x2(i)<0.03
news=40;
elseif x2(i)<0.08
news=50;
elseif x2(i)<0.23
news=60;
elseif x2(i)<0.43
news=70;
elseif x2(i)<0.78
news=80;
elseif x2(i)<0.93
news=90;
else news=100;
end
elseif x1(i)<0.8
if x2(i)<0.10
news=40;
elseif x2(i)<0.28
news=50;
elseif x2(i)<0.68
news=60;
elseif x2(i)<0.88
news=70;
elseif x2(i)<0.96
news=80;
else news=90;
end
else
if x2(i)<0.44
news=40;
elseif x2(i)<0.66
news=50;
elseif x2(i)<0.82
news=60
elseif x2(i)<0.94
news=70;
else news=80;
end
end
if paper>=news
sale=news;
remand=paper-news;
else
sale=paper;remand=0;
end
sb(n)=sb(n)+2*sale-
1.3*paper+0.2*remand;
end
end
optnews=40;
optmoney=sb(4);
[40,sb(4)/365]
for n=5:10
if sb(n)>=optmoney
optnews=n*10;
optmoney=sb(n);
end
[n,sb(n)/365]
end
[optnews,optmoney,optmoney/365]

经过计算机仿真后得到最优购货量是每天60份,平均每天利润34.4元,一年的利润是12396元。