第4讲 随机数的生成及随机变量抽样

y 1
• 定理 设随机变量 U 服从 ( 0 ,1) 上的均匀分布， 则 X F U 的分布函数为 F x 。 • 因此，要产生来自 F x 的随机数，只要先产生 来自 U 0 ,1 的随机数，然后计算 F u 即可。 • 其步骤为 由 U 0 ,1 抽取 u ， •
生成X的随机数
n 1 为常数
并画经验分布函数曲线。
合成法
• 合成法的应用最早见于Butlter 的书中。 • 构思如下： • 如果 X 的密度函数 p x 难于抽样，而 X 关于 Y 的条件密度函数 p x y 以及 Y 的密度函数 g y 均易于抽样，则 X 的随机数可如下产生：
Randnum=exprnd(2,1,10000)+5 cdfplot(Randnum)
二、其他各种分布的随机数的产生
• 基本方法有如下三种： • 逆变换法 • 合成法 • 筛选法
逆变换法
• 设随机变量 X 的分布函数为 F x ，定义
F
1
y inf x : F x y , 0
生成1行1000列501—1000上离散均匀分布的随机 数。 并画经验分布函数曲线。
cdfplot(x)
Randnum=unidrnd(10,1,10000)；cdfplot(Randnum)； pause Randnum=unidrnd(10,1,10000)+10；cdfplot(Randnum)； pause Randnum=unidrnd(500,1,10000)+500； cdfplot(Randnum)
1
1
பைடு நூலகம் 计算 x F
1
u
例3 设密度函数为
1 e f ( x; ) 0
x

x 0, x0
0 为常数
生成Y min{ X 1 , X 2 ,, X n }的随机数
并画经验分布函数曲线。
1 e FX ( x) 0,
随机数的生成及随机变 量抽样 实验目的
学习主要的随机变量抽样方法
实验内容
1、均匀分布U(0,1)的随机数的产生 2、其他各种分布的随机数的产生方法 3、随机数生成实例 4、实验作业
随机数的生成
• 随机数的产生是实现MC计算的先决条件。 而大多数概率分布的随机数的产生都是基于均 匀分布U(0,1)的随机数。 • 首先，介绍服从均匀分布U(0,1)的随机数的产 生方法。 • 其次，介绍服从其他各种分布的随机数的产生 方法。以及服从正态分布的随机数的产生方法。 • 最后，关于随机数的几点注。
1 2
n i 1 i
X 1 2 ln U 1 2 cos 2 U 2 1 X 2 2 ln U 1 2 sin 2 U 2
1
i
例6 生成单位圆上均匀分布的1行10000列随机数， 并画经验分布函数曲线。 Randnum=unifrnd(0,2*pi,1,10000); xRandnum=cos(Randnum) [Y,II] =sort(xRandnum) yRandnum=sin(Randnum) plot(xRandnum(II),yRandnum(II),'.')
离散型随机变量的生成
离散型随机变量X,它的取值是非光滑连续的值，它 只能间断地即离散地取值x1,x2,x3，…，xn，且规定 x1＜x2＜x3＜…<xn。其概率密度函数为
p(xi)=p{X= xi}
概率分布函数为
F ( X ) p{ X x i }

xi x
p( xi )
例10 对某车间每天需求某种零件的数量历史数据 中统计获得表1的结果。生成1行1000列零件需求的 随机数。并画经验分布函数曲线。
三角分布(a,m,b)的随机变量其密度函数为
0 f ( x) 2(x - a)/[(m - a)(b - a)] 1 当a x或x b时 当a x m时 当m x b时
其分布函数为
当x a时 0 2 (x - a) /[(m - a)(b - a)] 当a x m时 F ( x) 2 1 (b x) /([( b m)(b a )] 当m x b时 1 当x b时
由 Y 的分布 g y 抽取 y ， 由条件分布 p x y 抽取 x
• 可以证明由此得到 X 的服从 p x 。
筛选抽样
• 假设我们要从 p x 抽样，如果可以将 p x 表示 成 p x c h x g x ，其中h 是一个密度函数且易 于抽样，而0 g x 1 ，c 1 是常数，则 X 的抽样 可如下进行： 1 由 U 0 ,1 抽取 u ，由 h y 抽取 y ，
一、均匀分布U(0,1)的随机数的产生
• 产生均匀分布的标准算法在很多高级计算机语 言的书都可以看到。算法简单，容易实现。使用者 可以自己手动编程实现。Matlab 中也提供给我们 用于产生均匀分布的各种函数。我们的重点是怎样 通过均匀分布产生服从其他分布的随机数。因此， 直接使用Matlab提供的可靠安全的标准函数，当然 不用费事了。
例7 生成单位正方形上均匀分布的1行10000列随 机数，并画散点图。
mm=10000;Randnum=unifrnd(0,4,1,mm); xRandnum=zeros(1,mm);yRandnum=zeros(1,mm); for ii=1:mm if Randnum(1,ii)<=1 xRandnum(1,ii)=0; yRandnum(1,ii)=Randnum(1,ii); else if Randnum(1,ii)<=2 xRandnum(1,ii)=Randnum(1,ii)-1; yRandnum(1,ii)=1; else if Randnum(1,ii)<=3 xRandnum(1,ii)=1; yRandnum(1,ii)=1-(Randnum(1,ii)-2); else xRandnum(1,ii)=1-(Randnum(1,ii)-3); yRandnum(1,ii)=0; end end end end [Y,JJ] =sort(xRandnum);plot(xRandnum(JJ),yRandnum(JJ),'.')
三、生成标准正态分布的随机数
• N 0 ,1 的随机数产生方法很多。简要介绍三种。 • 法1、 变换法（Box 和Muller 1958） U • 设 U ， 是独立同分布的 U 0 ,1 变量，令
1 2
• • • • •
则 X 与 X 独立，均服从标准正态分布。 法2、 结合合成法与筛选法。（略） 法3、 近似方法（利用中心极限定理） 1 n 个 U 0 ,1 变量产生一个 N 0 ,1变量。 u u 即用 n x 其中u 是抽自U 0 ,1 的随机数， 12 n u 1 2 可近似为一 个 N 0 ,1 变量。
• 定理 c 其中0 g x 1 ， 1 ，h 是一个密度函数。令 U 和 Y 分别服从U 0 ,1 和 h y ，则在U g Y 的条件 Y • 下， 的条件密度为 p x U g Y p x
Y
2 如果 u g y ，则 x y ，停止， 3 如果 u g y ，回到 1 。 设X 的密度函数 p x ，且 p x c h x g x ，
表1 某零件每天需求量 X
需求量x(件)
X1=10
概率P(x)
0.10
累积概率F(x)
F(X1)=0.10
可分配的随机数范围
（ .00 -- .10]
X2=20
X3=30
0.20
0.40
F(X2)=0.30
F(X3)=0.70
（ .10 -- .30]
（ .30 -- .70]
X4=40
X5=50
例2设总体X的密度函数为
1 ( x ) e , x X ~ f ( x ) 0, 其它
, 为未知参数
其中 >0, 生成 2, 5 1行10000列的随机数. 并画经验分布函数曲线。 解:由密度函数知 X 具有均值为 的指数分布
IMSL库中的函数使用 • RNSET: 种子的设定
• • • CALL RNSET (ISEED) CALL RNOPT (IOPT) CALL RNUN (NR, R)
• RNOPT: 产生器的类型的设定
• RNUN/DRNUN: 产生均匀分布的随机数
例1生成1行1000列的1—10上离散均匀分布的随机 数； 生成1行1000列21—30上离散均匀分布的随机数；
随机向量的抽样方法
在用Monte Carlo等方法解应用问题时,随机 向量的抽样也是经常用到的. 若随机向量各分量相互独立,则它等价于多个 一元随机变量的抽样。
例8 生成单位正方形内均匀分布的1行10000列随机 数，并画散点图。 mm=10000 xRandnum=unifrnd(0,1,1,mm); yRandnum=unifrnd(0,1,1,mm); plot(xRandnum,yRandnum,'.')
) 1 (1 U )
n=20 Randnum=1-(1-unifrnd(0,1,1,10000)).^(1/n); cdfplot(Randnum)
例5 设密度函数为
ci , f X (x) 0, x i x x i 1 其它 x 0 ， i 0 ,...., n 1
0.25
0.05
F(X4)=0.95
F(X5)=1.00
（ .70 -- .95]
（ .95 – 1)
随机变量生成的算法为 ①产生一个u(0，1)，并令i=0； ②令i=i+1； ③若u＞F(xi)，转回到第②步，否则转至④；
mm=10000;Randnum=unifrnd(0,1,1,mm);xRandnum=zeros(1,mm); for ii=1:mm if Randnum(1,ii)<=0.1 xRandnum(1,ii)=10; else if Randnum(1,ii)<=0.3 xRandnum(1,ii)=20; else if Randnum(1,ii)<=0.7 xRandnum(1,ii)=30; else if Randnum(1,ii)<=0.95 xRandnum(1,ii)=40; else xRandnum(1,ii)=50; end end end end end cdfplot(xRandnum)
x

,
ny 1 e , x 0 ， FY ( y ) x0 0,
y 0 y 0
例4 设X分布函数为F(X)
生成Y min{ X 1 , X 2 ,, X n }的随机数
FY ( y ) 1 1 FX ( y ) , Rand Y F
均值为(a+m+b)／3;
2 2 2
方差为(a + b + m – ab – am - bm)／18。 令 u=F(x)，可解得
a ( m a )(b a )u ( 当0 u ( m a ) /(b a )时) x m (b a )u ( m a )(b m) (当(m - a)/(b - a) u 1时)
n 1 X 1 n
(1 (1 U )
)
x, FX ( x) 0,
1 x 0 x0
，
生成n=20的1行10000列随机数，并画经验分布 函数曲线。
FY ( y ) 1 y 1 Rand Y F
1 X n 1 n 1 n
(1 (1 U )