随机性模型与模拟方法ppt课件

1、随机性模型与模拟方法随机性模型与模拟方法n 随机变量随机变量n 蒙特卡罗方法蒙特卡罗方法n 随机数的生成随机数的生成n 模拟模拟一、随机变量一、随机变量n 何谓随机变量？随机变量是一个其值不可预测的变量。虽然一个随机变量在个别试验中其结果不确定，但在大量重复试验中其结果是具有统计规律的。正是随机变量的这种规律性使我们可以利用它来建模。例如我们可以利用下述的数据：n 得出一个模型。时间时间t秒）秒） 0 1 2 3 4 5 6 7 8 9变量变量X 1 0 2 2 1 2 0 1 0 2n 是一个离散的随机变量并取值于 0，1和2。我们不可能给出 与 的确定的关系式，但是可以通过数 的不同值出现

2、次数来描述这随机型 的规律列表如下：n 这个表给出了随机变量 的变化规律，频率告诉某个特定的事件发生的频繁程度。如果我们需要构造一个含有随机变量的模型，可以假设这个规律总是成立的，模型的假设可以基于这几个数据之上。实际操作时可以把频率分布当作概率函数来处理，但应注意概率是频率的极限值，这两者是有差异的。在处理一个简单的理论模型时，对概率函数X 0 1 2频数频数 3 3 4频率频率 0.3 0.3 0.4tXXXXn必须作出合适的选择。例如，假设在上述问题中的随机变量取三个值时等于可能的，这样其概率函数为n这个例子说明在处理随机变量的模型时有以下两种选择：n （1使用一个理论模型。这在任何一本

3、概率统计的书上都可以找到一些标准的理论模型如二项分布等。每一个都基于一定的假设之下成立的，所以在选用时要特别注意其假设条件。n （2使用基于实际数据的频率表，并不去套用不准理论模型。 P x 0 1 2X131313n使用前者的好处在于能精确地叙述变量的概率，在处理问题时可以充分发挥数理统计的作用。但这一好处把所求模式制约在了处理简单情形。随着复杂性的增加，数学就变的太难。使用后者的好处在于模型时基于观测到的数据而不是基于假设之上。增加复杂性并不成为一大障碍，但我们不再能利用数理统计而得求助于模拟以及模型的统计结果。n在建立随机性模型时，首先要注意，将要处理的是离散还是连续的随机变量。n 1、

4、离散随机变量n 离散随机变量的理论模型是由概率函数 n来刻画的。这个式子说明随机变量 取值 时的概n率。对于离散型的随机变量有下面三种重要的分布 p xP XxXxn（01分布 设随机变量 只可能取0、1两个值，它的分布规律是n n则称 服从01分布。对于一个随机实验，如果它的样本空间只包含两个元素，即 ，我们总能在 上定义一个服从01分布的随机变量n来描述这个随机实验的结果。例如，对新生儿的性别进行登记，检查产品的质量是否合格等都可以用01分布的随机变量来描述。X1(1),0,1(01)kkP Xkppkp12,Se eSX120( )1eeXX eee当当n（2二项分布 设实验 只有两个可

5、能的结果，将 独立地重复地进行 次，则称这一串重复的独立实验为 重贝努利实验。它是一重和重要的数学模型，有着广泛的应用。若用 表示 重贝努利实验中事件 发生的次数， 是一个随机变量，它服从如下的二项分布n n 特别，当 时二项分布就是01分布。 (1),0,1,2,.,nkn kkP xkppkn1n EEnXnAXn（3泊松分布泊松分布 设随机变量设随机变量 所有可能的取值所有可能的取值为为 而取各个值的概率为而取各个值的概率为 n n 其中，其中， 是常数，则称是常数，则称 服从参数为服从参数为 的泊的泊松分布。可以证明当松分布。可以证明当 很小时，以很小时，以 为参数的为参数的二项分布，

6、当二项分布，当 时趋于以时趋于以 为参数的泊松分为参数的泊松分布，其中布，其中 0,1,2,.,!keP xkk0,1,2,. ,kn0nnpXpX, n p2、连续的随机变量、连续的随机变量n理论模型的连续型随机变量可以由概率密度函数 n 来描述，对所有的 存在n ，且 ，随机变量落在区间 n 的概率可由 来给出，在连续型随机变n 量中下述两种是重要的 。() ( )pdf f xx( )0f x 12( ,x x21( )xxf x dx( )1fx dxn(1)均匀分布 设连续型随机变量 具有概率密度n则称 在区间a，b上服从均匀分布。n 在区间a，b上服从均匀分布的随机变量 ，具有下述

7、意义的等可能性，即它落在区间a，b中任意等长度的子区间内的可能性是相同的，或者说它落在子区间内的概率只依赖于子区间的长度而与子区间的位置无关。n （2正态分布 设连续型随机变量 的概率密度为n其中 为常数，则称 服从参数为 的1,( )0,axbf xba其他X2221( ),2xf xex ,0 , XXXX 正态分布。 连续型随机变量的值如同离散的一样可以用频率表给出，但不同的是离散的随机变量每个频率对应于随机变量的一个值，而对于随机变量每一个频率对应于随机变量的一个取值范围。二、蒙特卡罗方法二、蒙特卡罗方法n蒙特卡罗方法是计算模拟的基础，其名字来源于世界著名的赌城摩纳哥的蒙特卡罗。其思想

8、来源于著名的蒲丰投针问题。n1777年法国科学家蒲丰提出了下述著名问题：平面上画有等距离 的一些平行线，取一根长度为 的针，随机地向有平行线的平面上掷去，求针与平行线相交的概率。n我们用几何概型来解决这一问题。设M为针落下后的中点， 表示中点M到最近一条平行线的距离， 表示针于平行线的交角，如图2.18所示。那么基本时间区域x0a a 0l l ,|0,02axx 图图2.18asin2lx2aoxMxn它为平面上的一个矩形，其面积为 。n 为使针与平行线与 最后的一条平行线相交，其充要条件是n n n n 的面积为 ，这样针与平行线相交的概率为n n 设一共投掷 次（ 是一个事先选好的相当大

9、的自然数），观察到针和直线相交的次数为 。 ( )2aS 020axA01( )sin2S Aldl ()()SApSnmAn/m nn 从上式我们看到，当比值 不变时， 值始终不变。取 为 的近似值，我们可以算出 的近似值。可以想象当投掷次数越来越多时计算的结果就越来越准确。下表时这些实验的有关资料 (此处把 折算为1）：/la 实验者实验者年份年份针长针长投掷次数投掷次数相交次数相交次数 的实验的实验值值pulf18500.8500025323.1596Smith18550.632041218.53.1554De Morggen C 18601.0600382.53.137Fox18840

10、.7510304893.1595Lazzerini19010.83340818083.141592Reina19250.541925208593.1795ppam nn 由此可以看出蒙特卡罗方法的基本步骤：首先，建立一个概率模型，使它的某个参数等于问题的解。然后按照假设的分布，对随机变量选出具体的值这一过程又叫着抽样），从而构造出一个确定性的模型，计算出结果。再通过几次抽样实验的结果，的到参数的统计特性，最终算出解的近似值。n 蒙特卡罗方法主要用再难以定量分析的概率模型，这种模型一般的不到解析的结果，或虽然又解析结果，但计算代价太大以至不可用。也可以用在算不出解析结果的定性模型中。n 用蒙特卡

11、罗方法解题，需要根据随机变量遵循的分布规律选出具体的至，即抽样。随机变量的抽样方法很多，不同的分布采用的方法不尽相同。在计算机上的各种分布的随机数事实上都是按照一定的确定性方法产生的伪随机数。三、随机数的生成三、随机数的生成n我们知道对于丢硬币的随机结果可以用以下的离散随机变量的改里函数来描述 n n 如果我们需要模拟随机变量的以个值或一个集合，可以用丢硬币然后记录其其结果的方法来得到，然而这具又相当的局限性，这里我们用数学程序来产生拟随机变量。即看上去是随机出现的，但并非真正的随家便朗，它们产生于一个梯推公式。不过这些拟随机数并没有明显的规律，当给于适当的伸缩之后，它们非常接近于在 区间的均

12、匀分布。X 0 1P(x) 0.5 0.50,1n这种方法的思想是，设计一个把 和 之间的整数映射到它们自身上的函数 ，然后从 开场，依次计算 例如通过下面的公式可以产生这样的一组随机变量n 给定任意一个初值，如 代入公式得 ，然后用 去除得 ；同样 n 代入公式，可以得 ，重复这一过程可以得到我们所需要的一组随机变量。在程序设计和软件包中通常用 来表示由这样，我们用它来表示从 上的均匀分布所产生的随机变量。0Mf0 x01(),f xx12(),.f xx1973nnXX11(mod1000),1000nnXR071,X1890X 100010.890R 20.333R 0,11890X n

13、 我们可以从它构造出另外的随机变量。例如，可以从 给出区间 上的连续均匀分布的随机变量。如果我们要生成带参数 的指数分布，可以用 。 如果我们要生成平均值未零，标准差为 1 的正态分布，可以用下列公式 n 和n n 来给出 的两个值，令 或 n 可以生成 型的正态分布。()Xab a RND , a b(1/ )ln()XRND1/2212 2ln()cos(2)XRNDRND 1/2112 2ln()cos(2)XRNDRND 1XX2XX( ,) Xn 为了得到离散的随机变量，我们把 分成若干部分。例如设计一个离散的随机变量有下列的概率函数。n取一个RND值：假设 ，那么 ；假设 ，那么

14、；n 假设 ，那么 。n对于连续的随机变量除了取生成的随机变量是每类的中点外，我们可以用同样的思想进行列表分类。如n n n n x 0 1 2 0.3 0.3 0.4 0-10 10-15 15-20频率频率 0.2 0.5 0.3 PxX00.3RND0X 0.30.6RND1X 0.6RND2X n 的一个 值将平移到 。一个更细致的方法是用线性插值而不是取中点，即n n 给出 。n 从已知的 模拟一个连续随机变量的理论分布，可以用以下方法：n 0.36RND12.5X 100.360.250.5X 11.6X pdfn（1逆累积分布函数法 如果随机变量的 n 是 ， 则累积分布函数是

15、n 。如果把它作为一个随机变量， 是 上的均匀分布。从 上的均匀分布取一个 值，解方程 得对应得 的值，例如，设 n n 累积分布函数为n 解 得 。这就是我们所要的由这个分布所生成的 的值pdf( )f x( )( )xF xf t dtF0,1RND( )RNDF x1()xFRND0.5sin0( )0 xxf x否则00( )0.5sin0.5cos0.5(1 cos )xxF xtdttx 0.5(1 cos)RNDXarccos(1 2)XRNDX0,1n(2) 排除法 对于这种方法我们需要用两个 n 值来生成一个 值。设 的值在区间n 外为 ，而 的最大值是 。n 我们可以通过如

16、下的步骤生成 的值。 n 从 上的均匀分布生成 和 ； n 用 计算 ；n 计算 ；n 用 算出 ；n 假设 ，则接受 ，否则排除 n 回到 。 n 对于上面的例子，我们取RNDX , a b0c0,11RND2RND1()xaba RND2ycRND( )yf x).a0,0.5.abc( )f x( )f x1RND( )f x2RNDxx)a)b)c)d) eX四、模拟四、模拟n模拟是现象的模型所产生的再现。所谓数学模拟就是用模型使现象再现。因此，表示现象的部分或总体的基本方程和表示自然规律的数学模型全是数学模拟。然而，狭义地讲主要指的是数字模拟。它是将复杂现象作出可以用数字计算机表达的

17、数学模型，从数值上进行各种实验。各种方法随着计算机的进步已广泛地应用起来。因此我们所说的模拟主要是指数学模拟。n例 2.18 一列火车大约在下午1点离开 站其规律如下； n 火车从 到 途中所需要的平均时间为 分，由 n 分钟的标准差。如果你要赶的是这趟火车的下一站 ,而你到达 的站的时间分布为n 问你能赶上这列火车的概率是多少？A离站时间离站时间 13.00 13.05 13.10 概率概率 0.7 0.2 0.1时间时间 13.28 13.30 13.32 13.34概率概率 0.3 0.4 0.2 0.1B302ABn为了回答这个问题，我们需要一些随机数。这里我们将采用上面给出的那些随机

18、数，即n 等。而我们所要模拟的是n 火车离站的时间 ；n 火车途中的时间 ；n 你到达车站 的时间 。 n 这样你赶上火车的条件是 。为模拟这个问题只需要生成 ， 和 的值，然后检验这条件。但如何得到 的值是不明显的，因并不知道这个分布。这样，假设一个模型，取平均值为30，标准差为2的正态分布，由所给的条件知 ， n 为离散的，而 为连续的随机变量。0.890,0.333,0.304,0.491,0.630.1t2t3t312ttt1t2t1t2t3t3t)a)b) cBn以分为时间单位，从 的下午以点起算，构造的模型如下n 其中 。n 计算结果为 ， 和 ，这样 。在这种场合你比火车提前到达

19、4分钟。但需要指出，这并不是说我们已经回答了这个问题，要回答这个问题我们要作多次这样的模拟，记下这些结果，算出能赶上火车的频率。通过足够多次的模拟之后我们就可以看出能赶上火车的概率。11100.7,00.70.9,50.91.010RNDtRNDtRNDt333300.3,280.30.7,300.70.9,320.91.0,34RNDtRNDtRNDtRNDt1 212 2ln()cos(2)XRNDRND 0t 15t 229t 330t 1234ttn 一般用在模拟建模时，一次模拟的成功并不能说明什么问题，更不能说我们的主要工作已经完成。你必须多次的进行模拟，然后分析其结果。分析的种类要

20、看模型的对象，而这在模拟的一开始就应该清楚的。在实验的模拟模型的对象是在变化的，但常常包括一下几种：n 对系统的长期性态作出统计；n 比较系统的可选择对象的安排； n 研究参数变化的影响；n 研究模型假设的影响； n 找出系统最优方案；n 上面的例子是相当平凡的，根本不能作为用模拟解决问题的例子。下面我们仅举两个简单的例子以理顺模拟模型的思路。)a)b)c)d) en例例2.19 某个理发店中有两名理发员某个理发店中有两名理发员 和和 ，顾客，顾客随机地来理发，据统计随机地来理发，据统计 的顾客仅需要剪发，化的顾客仅需要剪发，化时时 n 分钟；而有分钟；而有 的顾客即需要剪又需要又需要吹的顾客

21、即需要剪又需要又需要吹风，许花时风，许花时 分。分。n 对任意一个模拟，首先要作的是对任意一个模拟，首先要作的是 n 找出能完全描述任意时刻的系统的状态变找出能完全描述任意时刻的系统的状态变量的集合。量的集合。n 给出能从时刻给出能从时刻 的状态变量算出时刻的状态变量算出时刻 的新的状态变量的程序。的新的状态变量的程序。n 这个例子中有三个状态变量：在等待的顾客这个例子中有三个状态变量：在等待的顾客的人输离散的非负整数）；能否的人输离散的非负整数）；能否 正在工作是或正在工作是或否）；能否否）；能否 正在工作是或否）。正在工作是或否）。)a)bAB40%60%81t 5tABn一次模拟式由始于

22、 ，结束于 的状态变量的值的一系列演算组成的。一个事件是时间中的一点，在这个时刻一个随机变量改变了它的值。在这个例子中的事件有:n 一个顾客到达； 开始服务； 结束服务； 开始服务； 结束服务。n 一个元素是一个离散，或者是系统的长期部分，或者是进入和离开，这里的元素是顾客和两个理发员。n 对研究一个模拟模型来说，有两种程序类型： n （1时间切片 考察状态变量和在时间切片中通常是等时间的切片元素的位置。在每一个时间切片中状态变量可变可不变。0t tENDABAB （2事件序列 考察在每一事件的系统，并不考虑时间之间的时间。 这两种途径我们有时称为“时间传动和“事件传动模型。一般我们用时间传动

23、模型于连续的的确定型系统，事件传动模型于离散的概率模型，但这不时绝对。在这个例子中我们将用“事件传动”。 对于时间切片模型，我们必须决定时间切片的大小，为简单计我们将取 分钟。问题的描述并不包括任何有关顾客到达率的信息。假设在任何一分钟顾客到达的概率是 。实际上有两种不同类型的顾客，取决于是否要吹风。我们通过取服务时间的平均值，即 分，构造一个粗糙的模型。10.6 50.4 86.2 n 为了描述一个顾客是否到来这个随机变量，我们用一个硬币将作为一个随机数的生成器，用n 表示反面， 表示正面。设扔出的序列是 n 。用 表示一个顾客到达，且取初始状态为顾客，运行前 分钟，就有下表的结果：H时间分

24、）时间分） 到达？到达？ A在工作在工作 B在工作在工作 排队排队 0 否否 否否 否否 0 1 是是 是是 否否 0 2 否否 是是 否否 0 3 是是 是是 是是 0 4 是是 是是 是是 0 5 否否 是是 是是 0T, , , , , ,T H T T H T T T H H T10 到 这里，人们将要为我们希望知道什么。通常我们感兴趣的是平均队伍的长度，最长的队伍，顾客等待的平均时间以及两个理发员的忙 闲 程度等，注意到这里有两种不同的平均，即一个是关于时间，而另一个是关于顾客的平均，为回答上述为她我们设 是任意时刻的排队的顾客数。顾客和时间的关系通常可由图 给出。 ( )Q6 是是

25、 是是 是是 0 是是 是是 是是 0 是是 是是 是是 0 否否 是是 是是 010 否否 是是 是是 0WQ2.19 它是一个右连续的阶梯函数是合理的，这是由于只有新顾客到来或有顾客完成服务后离去，函数值才发生变化， 关于时间的平均是 ，其中 图下额面积，设 表示一个时间区间，在其上 保持常数这里 本身是变量）。当我们进行模拟时我们累积其和 。用 记在进行模拟期间到达的顾客数。这样我们所要的两个平均分别为tQQ0ENDt2.19图Q t（）NQQQENDtQtn n 队长平均n 等待时间的平均n n 下面是用来说明累积排队时间 的记录注意这里仅给出 变化的时间）：Q tQEND（）Q tW

26、N（）Q t（）QtQ t时间时间 Q 0 0 0 0 011.584 1 0 0 0 12.935 0 1.351 1.351 1.35117.290 1 4.355 0 1.35117.935 0 0.645 0.645 1.99618.676 1 0.741 0 1.99623.156 0 4.480 4.480 6.47625.217 1 2.061 0 6.47625.327 2 0.110 0.110 6.58625.935 1 0.608 1.216 7.80227.431 2 1.406 1.406 9.208Q t（）n这里所执行的总时间 。在这期间 。队伍的最长长度 。累计

27、排队时间 队伍的平均长度 。平均等待时间n 。n 在我们结束模拟时还有两个顾客，一个是排队的，而另一个是新来的。n 服务的总时间是n 分。因此 忙碌的概率是 n 。 服务的总时间是 分或忙碌时间为 。27.341END 10N max2Q()9.208Q t 9.208/17.3410.34Q 9.208/100.92W 22.185 27.34110081（、）B19.40671%A(50)(13.2858.285)(23.156 15.156)(27.341 23.156)22.185分An 评注1 模拟一个系统的目的不是为了模仿一个现实系统，而是通过解决问题达到优化系统的目的。例如在这个

28、例子中可以分析诸如增雇一个理发员或改变服务时间等对系统的影响。 n 评注2 在我们的模型中，为使问题简单我们已经作了一些假设：n （1假设了在任何一分钟有一个顾客到达的概率是 。n （2默认在同一分钟内的顾客数 。n （3如果两个理发师均空闲，顾客可以任意选。n （4排队的原则是安先后的秩序。如果有预约可以先服务。0.51n（5我们的模型中允许一下情形出现，一个顾客的来到，发现有很多人在等就走啦。也可能是一个顾客在等了一段时间之后等不及了就离开了。意味这允许其中的一个理发员有短暂的休息。 n 例2.20 倒媒台的操作方案 某煤矿公司有一个大型煤台，用于向运媒列车装煤。该倒煤台的容量是 列标准列

29、车。装满一个空的倒煤台需要一个小组 个小时的时间，费用是 。为提高装煤速度可以以 的代价动用第二个小组。铁道部门每天向这个倒煤台发三列空的标准车。这些列车可在上午 点倒下午 点之间的任何时刻到达。给一列标准车装满煤需要 小时，向倒煤台装煤和从倒煤台向列车装煤不能同时进行。如果列车到达后因等待装煤二停滞，铁道部门将征收每车 的滞费。1.569000元/小时12000元/小时58315000元/小时n此外，每星期四上午 点到下午 点之间还有一列大容量列车到达，其容量为标准的列车的 倍，滞期费为 。请问 n （1安标准规则操作可使装煤费最低？费用使多少？n （2如果标准列车能在指定的时间到达，什么样

30、的调度安排最经济 ？n 这个问题中列车的到达时间使随机因素，适合于建立概率模型，同计算机模拟加以解决。n 首先，模型中需要考虑的费用由两部分组成。一部分使装煤小组向倒煤台装煤的费用，记为 ，另一部分使列车等待装煤的滞期费 。因每天要装的煤数量使固定的， 的大小只受是否使用大二小组影响。25000元/小时LCDCLCn 通过使用第二小组，有可能减少 。模型的主要任务是将总费用 降到最低。故 是模型的目标函数。 n 其次，由于理论上的困难，很难得到最优方案。考虑到这是一个每天重复发生的为她，重要的的是提供一组简单明确的规划，使煤矿公式可以根据规则方便地获得接近最优的解。因此，我们将在方案的优化程度

31、和简明性之间做一个折中。n 设： 为装满列车 所需的煤量； 为倒煤台中剩下煤量； 表示当前时间。其中 和 均以 小时向列车装的煤量为单位。DCLC=CAr0,24)tQ1ADCArQn 根据题意写出下面一些应该遵循的规则：n 有列车等待时，用两个小组装煤节省的滞期费大于增加的装煤费用，此时应使用第二个小组。n 当同时有两列或三列标准列车等待装煤时。应将已装煤量最多的车排在前面先装，已装煤量最少的车排到最后面。可以证明，这样安排滞期费最少。n 当同时有大容量车 和标准车 等待时，先装 后装 的滞期费n n n先装 时的滞期费为 n n当 时，先装 ，否则先装 。AB12225000max(),0

32、15000max(),033DAAABCrQrrrQ22215000max(),025000max(),033DBBABCrQrrrQ12DDCC( )a( )b( )cABABn 设当前待装的车为 ，则用两个小组装倒煤台直到 或 为止，然后装列车。n 周四时，装标准车和大容量列车共需要 小时。即便是倒煤台在周四上午 点以前就已提前装满，当天用两个小组装倒煤台仍需 小时，合计 小时故最快也要到周五早上 点才能完成周四的任务，且此时倒煤台为空。为保证周五正常工作，应马上开始装倒煤台。由以上分析知，周时间最紧张，就始终用两个小组。n 非周四，在此刻 无列车等待，设已知下一列n 车到达时间为 。假设

33、 ，则时间充n 足，可以用一个小组装倒煤台至满或下一列车来。否则用两个小组。AQr4.5Q A151223tt334tQ( )d( ) e( )f5tn 非周四，不知道列车的到达时间。设在时刻 倒煤台中尚有煤量 ，没有列车等待，当天还有 列标准车未到达。假设列车到达时间服从独立的均匀分布，则存在 ，当 时用一个组装煤即可，否则要用两个组。n 的选择应满足最小的原则，因其解析解难以求出，故采用计算机模拟的方法。首先任意取一个 值 ，注意到 ，在上述约束条件下以一定步长 取 n 各种组合，分别用计算机模拟求出平均费用，找出使平均费用最少的一组 ， 和 值，n 作为在该组合给定 下的函数值。选取一系列不同的 的值重复以上过程，就可