第4讲随机数的生成及随机变量抽样省公共课一等奖全国赛课获奖课件

实验目试验内容学习主要随机变量抽样方法1、均匀分布U(0,1)随机数产生2、其它各种分布随机数产生方法3、随机数生成实例4、试验作业随机数生成及随机变量抽样第1页随机数生成随机数产生是实现MC计算先决条件。而大多数概率分布随机数产生都是基于均匀分布U(0,1)随机数。

首先，介绍服从均匀分布U(0,1)随机数产生方法。其次，介绍服从其它各种分布随机数产生方法。以及服从正态分布随机数产生方法。最终，关于随机数几点注。第2页一、均匀分布U(0,1)随机数产生产生均匀分布标准算法在很多高级计算机语言书都能够看到。算法简单，轻易实现。使用者能够自己手动编程实现。Matlab中也提供给我们用于产生均匀分布各种函数。我们重点是怎样经过均匀分布产生服从其它分布随机数。所以，直接使用Matlab提供可靠安全标准函数，当然不用费事了。

第3页IMSL库中函数使用RNSET:种子设定CALLRNSET(ISEED)RNOPT:

产生器类型设定

CALLRNOPT(IOPT)

RNUN/DRNUN:产生均匀分布随机数CALLRNUN(NR,R)

第4页例1生成1行1000列1—10上离散均匀分布随机数；生成1行1000列21—30上离散均匀分布随机数；生成1行1000列501—1000上离散均匀分布随机数。并画经验分布函数曲线。Randnum=unidrnd(10,1,10000)；cdfplot(Randnum)；pauseRandnum=unidrnd(10,1,10000)+10；cdfplot(Randnum)；pauseRandnum=unidrnd(500,1,10000)+500；cdfplot(Randnum)cdfplot(x)第5页第6页解:由密度函数知

例2设总体X密度函数为其中>0,生成1行10000列随机数.含有均值为指数分布Randnum=exprnd(2,1,10000)+5并画经验分布函数曲线。cdfplot(Randnum)第7页第8页二、其它各种分布随机数产生基本方法有以下三种：逆变换法合成法筛选法第9页逆变换法设随机变量分布函数为

，定义

定理

设随机变量服从上均匀分布，则分布函数为。所以，要产生来自随机数，只要先产生来自

随机数，然后计算即可。其步骤为

第10页为常数例3

设密度函数为并画经验分布函数曲线。第11页例4

设X分布函数为F(X)生成n=201行10000列随机数，并画经验分布函数曲线。第12页n=20Randnum=1-(1-unifrnd(0,1,1,10000)).^(1/n);cdfplot(Randnum)第13页为常数例5

设密度函数为并画经验分布函数曲线。第14页合成法

合成法应用最早见于Butlter书中。构思以下：

假如密度函数难于抽样，而关于条件密度函数以及密度函数均易于抽样，则随机数能够下产生：能够证实由此得到服从。第15页筛选抽样

假设我们要从抽样，假如能够将表示成

，其中是一个密度函数且易于抽样，而，是常数，则抽样能够下进行：

定理

设密度函数，且，其中，，是一个密度函数。令和分别服从和，则在条件下，条件密度为

第16页三、生成标准正态分布随机数随机数产生方法很多。简明介绍三种。

法1、

变换法（Box和Muller1958）设，是独立同分布变量，令则与独立，均服从标准正态分布。法2、结合合成法与筛选法。（略）法3、近似方法（利用中心极限定理）即用个变量产生一个变量。其中是抽自随机数，可近似为一个变量。第17页例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),'.')第18页第19页例7生成单位正方形上均匀分布1行10000列随机数，并画散点图。mm=10000;Randnum=unifrnd(0,4,1,mm);xRandnum=zeros(1,mm);yRandnum=zeros(1,mm);forii=1:mmifRandnum(1,ii)<=1xRandnum(1,ii)=0;yRandnum(1,ii)=Randnum(1,ii);elseifRandnum(1,ii)<=2xRandnum(1,ii)=Randnum(1,ii)-1;yRandnum(1,ii)=1;elseifRandnum(1,ii)<=3xRandnum(1,ii)=1;yRandnum(1,ii)=1-(Randnum(1,ii)-2);elsexRandnum(1,ii)=1-(Randnum(1,ii)-3);yRandnum(1,ii)=0;endendendend[Y,JJ]=sort(xRandnum);plot(xRandnum(JJ),yRandnum(JJ),'.')第20页第21页离散型随机变量生成离散型随机变量X,它取值是非光滑连续值，它只能间断地即离散地取值x1,x2,x3，…，xn，且要求x1＜x2＜x3＜…<xn。其概率密度函数为p(xi)=p{X=xi}概率分布函数为例10

对某车间天天需求某种零件数量历史数据中统计取得表1结果。生成1行1000列零件需求随机数。并画经验分布函数曲线。表1某零件天天需求量X第22页需求量x(件)概率P(x)累积概率F(x)可分配随机数范围X1=10

0.10

F(X1)=0.10（.00--.10]X2=20

0.20

F(X2)=0.30（.10--.30]X3=30

0.40

F(X3)=0.70（.30--.70]X4=40

0.25

F(X4)=0.95

（.70--.95]X5=50

0.05F(X5)=1.00（.95–1)随机变量生成算法为①产生一个u(0，1)，并令i=0；②令i=i+1；③若u＞F(xi)，转回到第②步，不然转至④；④输出得X＝xi。第23页mm=10000;Randnum=unifrnd(0,1,1,mm);xRandnum=zeros(1,mm);forii=1:mmifRandnum(1,ii)<=0.1xRandnum(1,ii)=10;elseifRandnum(1,ii)<=0.3xRandnum(1,ii)=20;elseifRandnum(1,ii)<=0.7xRandnum(1,ii)=30;elseifRandnum(1,ii)<=0.95xRandnum(1,ii)=40;elsexRandnum(1,ii)=50;endendendendendcdfplot(xRandnum)第24页第25页

三角分布(a,m,b)随机变量其密度函数为其分布函数为第26页第27页在用MonteCarlo等方法解应用问题时,随机向量抽样也是经惯用到.若随机向量各分量相互独立,则它等价于多个一元随机变量抽样。随机向量抽样方法

第28页例8生成单位正方形内均匀分布1行10000列随机数，并画散点图。mm=10000xRandnum=unifrnd(0,1,1,mm);yRandnum=unifrnd(0,1,1,mm);plot(xRandnum,yRandnum,'.')第29页第30页第31页mm=100000xRandnum=unifrnd(0,1,1,mm);yRandnum=unifrnd(0,1,1,mm);[Y,JJ]=sort(xRandnum)plot(xRandnum(JJ),yRandnum(JJ),'.')第32页例9生成单位圆内均匀分布1行10000列随机数，并画散点图。mm=10000；Randnum1=unifrnd(-1,1,1,2*mm);Randnum2=unifrnd(-1,1,1,2*mm);xRandnum=zeros(1,mm);yRandnum=zeros(1,mm);s=Randnum1.^2+Randnum2.^2；ii=1;jj=1;whileii<mmifs(1,jj)<=1;xRandnum(1,ii)=Randnum1(1,jj);yRandnum(1,ii)=Randnum2(1,jj);ii=ii+1;endjj=jj+1;endplot(xRandnum,yRandnum,'.')第33页第34页关于随机数几点注注1

因为均匀分布随机数产生总是采取某个确定模型进行，从理论上讲，总会有周期现象出现。初值确定后，全部随机数也随之确定，并不满足真正随机数要求。所以通常把由数学方法产生随机数成为伪随机数。注2

应对所产生伪随机数作各种统计检验，如独立性检验，分布检验，功率谱检验等等。

但其周期又相当长，在实际应用中几乎不可能出现。所以，这种由计算机产生伪随机数能够看成真正随机数来处理。第35页2.设密度函数为1.生成单位球内均匀分布1行10000列随机数，并画散点图。作业:为常数并画经验分布函数曲线。第36页3.生成三角分布(0,1,2)1行10000列随机数，并画散点图。作业:并画经验分布函数曲线。第37页

§5.2随机数与随机变量生成

§5.2.1随机数生成

在系统模拟中只要有随机变量，则在模拟运行每一步中都要对随机变量确定一个详细值。我们将会碰到各种概率分布随机变量，但其中最简单或最基本随机变量是在(0，1)区间上均匀分布随机变量。服从某一分布随机变量都能够经过对(0，1)均匀分布随机变量进行适当转换而得到。(0，1)均匀分布随机变量取值也是在(0，1)区间上均匀分布随机数ui序列(流)独立采样，其密度函数是ui数学期望和方差分别为第38页所以，若能取得(0，1)均匀分布随机数，也就能经过对其适当转换而取得某一要求分布随机变量取值，这就是随机变量生成。为此，首先要掌握(0，1)区间上均匀分布随机数生成方法。均匀分布随机数必须具备均匀性和独立性要求；要生成符合上述要求随机数流，现在多用数学算法来产生，普通是采取递推算法，确定一个初始值(种子数)以后，逐次递推算得随机数流。数学算法取得随机数、常称之为伪随机数(PseudoRandomNumber)序列。数学方法计算产生随机数流必须满足以下要求：(1)尽可能在(0，1)区间均匀分布；(2)含有统计上独立性；(3)产生随机数流能够重复出现，即给以相同初值(种子数)能取得相同随机数流；(4)有足够长周期，即在出现周期性重复之前，能生成足够多个随机数；(5)算法占用计算机内存较少而计算生成速度较快。当前广泛应用算法是线性同余法(LinearcongruentialMethod)，其中又分为：1．混合线性同余法。它是由Lehmer于1951年提出，其算式为第39页xi+1=(axi+c)modmui+1=xi+1/m式中a——乘数(常数)；C——增量(常数)；x0——种子数；m——模数。a,c,m和x0选取对随机数流统计特征和周期长度有极大影响。

上述第一式含义是式中[]表示取整数，a，c，m皆为整常数。2、

乘法线性同余法。若混合线性同余法中c=0，则为乘法线性同余法，其算式为

xi+1=aximodmui+1=xi+1/m(5.7)第40页可参考选取数据有：（1）

a=16807,m=2147483647,x0=123457;（2）

a=655393,m=33554432。

§5.2.2随机数流检验

一、均匀分布性检验1．参数检验。检验ui数字特征，如均值、方差预计值和其理论值差异是否显著。设有u1,u2,…，un随机数流，则它们若ui序列在(0，1)上均匀分布，可假设：u期望和方差分别为第41页2期望和方差分别为则上列假设(8．10)与(8．11)应该成立。据此，可对n个ui计算以下统计量若取显著性水平a=0．05，

当|V1|≤1．96时。则可认为假设(8．10)式成立；当|V2|≤1．96时，则可认为假设(8．11)式成立。因而能够接收此假设，检验经过;不然拒绝接收。2．均匀性检验。它是检验所生成随机数落在(0，1)各子区间频率均匀程度，是否与理论上均匀分布频率有显著性差异。此处介绍惯用方法之一，x2检验方法以下：

将(0，1)区间划分为相等k个子区间，假如落在第i个(i＝l，2，3，…，k)子区间随机数有ni个；而在理论上第i个子区间随机数个数为mi=N／k，其中N为随机数流总个数(拟检验)。由此，可计算x2统计量第42页再按k-1为自由度、显著性水平取0．05，查得

2(a)表值。当算得统计量x2≤

2(a)时，可认为在显著性水平a下能接收为均匀分布假设。

二、独立性检验独立性检验是检验随机数流中前后各数之间是否存在相关性。惯用方法是进行自相关检验。另外，还有PokerTest和Run检验，普通应用较少。

§5.2.3随机变量生成

一、离散型随机变量生成离散型随机变量X,它取值是非光滑连续值，它只能间断地即离散地取值x1,x2,x3，…，xn，且要求x1＜x2＜x3＜…<xn。其概率密度函数为p(xi)=p{X=xi}

概率分布函数为第43页譬如，对某车间天天需求某种零件数量历史数据中统计取得表5．2结果。表5．2某零件天天需求量X需求量x(件)概率P(x)累积概率F(x)可分配随机数范围X1=100.10F(X1)=0.10.00--.09X2=200.20F(X2)=0.30.10--.29X3=300.40F(X3)=0.70.30--.69X4=400.25F(X4)=0.95.70--.94X5=500.05F(X5)=1.00.95--.99随机变量生成算法为①产生一个u(0，1)，并令i=0；②令i=i+1；③若u＞F(xi)，转回到第②步，不然转至④；④输出得X＝xi。第44页二、连续型随机变量生成

以下介绍惯用几个概率分布随机变量生成方法。1．逆变换法。任何概率分布随机变量X，其分布函数F(x)值域是(0,l)；所以，可令

F(x)=u(5．17)假如对上式能解出显式逆函数

x=F-1(u)(5．18)则调用随机数生成算法产生一个随机数u，将其代入显式(8．18)，即可取得随机变量X一次取值x。可惜能求得显式逆函数随机变量，只有以下几个分布。(l)在(a，b)区间均匀分布随机变量X。其概率密度函数为其分布函数为第45页均值为(a十b)／2；方差为(b—a)2／12。两函数图形如图8．1。令u=(x-a)／(b-a)，则可解得x=a+(b–a).u(5．21)(2)三角分布(a，m，b)随机变量其密度函数为其分布函数为

第46页第47页第48页第49页kUP结果00.4430.443P>0.002510.1720.076P>0.002520.4230.032P>0.002530.8190.026P>0.002540.2550.0066P>0.002550.7490.0049P>0.00250.2250.0011P<0.0025,置X=6kUP结果00.4430.443P>0.002510.1720.076P>0.0025