随机模拟计算随机模拟方法的特点求解确定性问题内容和要求 随机模拟的基本方法；用随机模拟方法求解确定性问题；常用随机过程的模拟；模拟方法在随机服务系统中的应用。1 了解随机模拟方法的特点、基本步骤；2 熟练掌握求解确定性问题的随机模拟方法；3 掌握常见随机过程的模拟方法；4 能够利用模拟方法进行服务系统的随机模拟。Monte-Carlo 方法 一、 MC 的起源和发展Buffon试验排队系统模拟：M/G/1排队系统 二、 MC 的原理 三、随机数的产生原理与 IMSL库均匀分布 U(0,1)的随机数的产生其他各种分布的随机数的产生随机过程模拟 四、 MC的应用举例定积分的 MC计算随机微分方程的数值模拟 五、 EM算法及其 MCMC方法一、 MC 的起源和发展 随机模拟方法，也称为 Monte Carlo方法，是一种基于 “随机数 ”的计算方法。这一方法源于美国在第一次世界大战进行的研制原子弹的 “曼哈顿计划 ”。该计划的主持人之一、数学家冯 诺伊曼用驰名世界的赌城 摩纳哥的 Monte Carlo来命名这种方法，为它蒙上了一层神秘色彩。冯 诺伊曼是公理化方法和计算机体系的领袖人物， Monte Carlo方法也是他的功劳 。 事实上， Monte Carlo方法的基本思想很早以前就被人们所发现和利用。早在17世纪，人们就知道用事件发生的 “频率 ”来决定事件的 “概率 ”。18世纪下半叶的法国学者 Buffon提出用投针试验的方法来确定圆周率的值。这个著名的 Buffon试验是Monte Carlo方法的最早的尝试！ 历史上曾有几位学者相继做过这样的试验。不过呢，他们的试验是费时费力的，同时精度不够高，实施起来也很困难。然而，随着计算机技术的飞速发展 ， 人们不需要具体实施这些试验，而只要在计算机上进行大量的、快速的模拟试验就可以了。 在大众的心目中，科学的代言人是心不在焉的牛顿或者爆炸式发型的爱因斯坦，但这只是传统形象，比他们更了解现代计算技术的冯 诺伊曼是个衣着考究，风度翩翩的人物，他说：纯粹数学和应用数学的许多分支非常需要计算工具，用以打破目前由于纯粹分析的研究方法不能解决非线性问题而形成的停滞状态。 Monte Carlo方法是现代计算技术的最为杰出的成果之一，它在工程领域的作用是不可比拟的。Buffon试验 假设平面上有无数条距离为 1的等距平行线，现向该平面随机投掷一根长度为 的针( )， 则我们可计算该针与任一平行线相交的概率。这里，随机投针指的是：针的中心点与最近的平行线间的距离 均匀的分布在区间 。 上，针与平行线的夹角 （不管相交与否）均匀的分布在区间 上。 因此，针与线相交的充要条件是l1l 21,0x ,021sinxBuffon试验 从而针线相交的概率为 根据上式，若我们做大量的投针试验并记录针与线相交的次数，则由大数定理可以估计出针线相交的概率 ，从而得到 的估计值。 针与线的位置关系：ldxdlXPpl22sin20sin20=xlp排队系统模拟：M/G/1排队系统服务规则： 先到先服务假设： （）顾客到达遵循Poisson分布；（）服务时间服从一般分布；（）到达间隔与服务时间相互独立1排队系统模拟：M/G/1排队系统关心的指标：（）时刻t时，系统中的顾客数；即队长的分布；（）顾客的等待时间；（）服务的忙碌程度；（）用最朴素的 Monte-Carlo方法可以得到这些指标的估计二、 MC 的原理 应用 Monte Carlo方法求解工程技术问题可以分为两类： 确定性问题 随机性问题确定性系统随机性系统模拟自然界Monte-Carlo模拟，即随机模拟（重复 “试验 ”）重复试验计算机模拟思路：1、 针对实际问题建立一个简单且便于实现的概率统计模型，使问题的解对应于该模型中随机变量的概率分布或其某些数字特征，比如，均值和方差等。所构造的模型在主要特征参量方面要与实际问题或系统相一致的。 2、 根据模型中各个随机变量的分布，在计算机上产生随机数，实现一次模拟过程所需的足够数量的随机数。通常先产生均匀分布的随机数，然后生成服从某一分布的随机数，再进行随机模拟试验。思路： 3、 根据概率模型的特点和随机变量的分布特性，设计和选取合适的抽样方法，并对每个随机变量进行抽样（包括直接抽样、分层抽样、相关抽样、重要抽样等）。4、 按照所建立的模型进行仿真试验、计算，求出问题的随机解。 5、 统计分析模拟试验结果，给出问题的估计以及其精度估计。 6、 必要时，还应改进模型以降低估计方差和减少试验费用，提高模拟计算的效率。 收敛性: 由大数定律, Monte-Carlo模拟的收敛是以概率而言的 误差: 用频率估计概率时误差的估计，可由中心极限定理，给定置信水平 的条件下，有： 模拟次数: 由误差公式得NU 2/1|)( XgVar=22/1)(UN三、随机数的产生原理与 IMSL库 随机数的产生是实现 MC计算的先决条件。而大多数概率分布的随机数的产生都是基于均匀分布 U(0,1)的随机数。 首先，介绍服从均匀分布 U(0,1)的随机数的产生方法。 其次，介绍服从其他各种分布的随机数的产生方法。以及服从正态分布的随机数的产生方法。 接下来，介绍随机过程模拟。 最后，关于随机数的几点注。均匀分布 U(0,1)的随机数的产生 产生均匀分布的标准算法在很多高级计算机语言的书都可以看到。算法简单，容易实现。使用者可以自己手动编程实现。 IMSL统计库中也提供给我们用于产生均匀分布的各种函数。我们的重点是怎样通过均匀分布产生服从其他分布的随机数。因此，直接使用 IMSL提供的可靠安全的标准函数，当然不用费事了。IMSL库中的函数使用RNSET: 种子的设定 CALL RNSET (ISEED)RNOPT:产生器的类型的设定 CALL RNOPT (IOPT) RNUN/DRNUN: 产生均匀分布的随机数 CALL RNUN (NR, R) 其他各种分布的随机数的产生 基本方法有如下三种： 逆变换法 合成法 筛选法逆变换法 设随机变量 的分布函数为 ，定义 定理 设随机变量 服从 上的均匀分布，则 的分布函数为 。 因此，要产生来自 的随机数，只要先产生来自 的随机数，然后计算 即可。 其步骤为X ( )xF( ) ( ) 10,:inf1=yyxFxyFU )1,0( )UFX1= ( )xF( )xF( )1,0U( )uF1( )()=uFxuU11,0计算，抽取由合成法 合成法的应用最早见于 Butlter 的书中。 构思如下： 如果 的密度函数 难于抽样，而 关于的条件密度函数 以及 的密度函数均易于抽样，则 的随机数可如下产生： 可以证明由此得到 的服从 。X( )xp X Y( )yxpY( )ygX( )()xyxpyygY抽取由条件分布，抽取的分布由X( )xp筛选抽样 假设我们要从 抽样，如果可以将 表示成 ，其中 是一个密度函数且易于抽样，而 ， 是常数，则 的抽样可如下进行： 定理 设 的密度函数 ，且 ，其中 ， ， 是一个密度函数。令和 分别服从 和 ，则在 的条件 下， 的条件密度为( )xp() ( ) ( )xgxhcxp =( )xp()h( ) 10 =。，回到如果，停止，则如果，抽取，由抽取由1321,01yguyxyguyyhuUX( )xp( ) ( ) ( )xgxhcxp =()h( ) 10 xg1cUY ( )1,0U( )yh( )YgU ( )( ) ( )xpYgUxpY=Y标准正态分布的随机数 的随机数产生方法很多。简要介绍三种。 法 1、 变换法（ Box 和 Muller 1958） 设 ， 是独立同分布的 变量，令 则 与 独立，均服从标准正态分布。 法 2、 结合合成法与筛选法。（略） 法3、 近似方法（利用中心极限定理） 即用 个 变量产生一个 变量。 其中 是抽自 的随机数， 可近似为一个 变量。()1,0N1U2U()1,0U( ) ( )()()=2122112sinln22cosln22121UUXUUX1X2Xn( )1,0U( )1,0N( )2112 = unxiu( )1,0U=niiunu11( )1,0NIMSL库中的函数使用RNSET:种子的设定 CALL RNSET (ISEED)RNOPT:产生器的类型的设定 CALL RNOPT (IOPT)RNNOA: 产生服从标准正态分布的伪随机数 CALL RNNOA (NR, R) 随机过程模拟 高斯白噪声的产生： 一种简单有效的模拟高斯白噪声过程的方法是由独立的单位正态随机序列 ，按照如图所示的方式连接，其中 式中 为白噪声的强度。kktD =kD 由于正态随机数的独立性，当 很小时，按 (a)、 (b)所构成的过程的相关时间很短，从而具有很高的截止频率。 当然， 过小，样本的计算量过大。因此，选取适当小即可。ttt关于随机数的几点注 注1 由于均匀分布的随机数的产生总是采用某个确定的模型进行的，从理论上讲，总会有周期现象出现的。初值确定后，所有随机数也随之确定，并不满足真正随机数的要求。因此通常把由数学方法产生的随机数成为伪随机数。但其周期又相当长，在实际应用中几乎不可能出现。因此，这种由计算机产生的伪随机数可以当作真正的随机数来处理。 注2 应对所产生的伪随机数作各种统计检验，如独立性检验，分布检验，功率谱检验等等。四、 MC的应用举例 例一、定积分的 MC计算随机投点法样本平均值法几种降低估计方差的 MC方法 例二、 随机微分方程的数值模拟定积分的 MC计算 事实上，不少的统计问题，如计算概率、各阶距等，最后都归结为定积分的近似计算问题。 下面考虑一个简单的定积分 为了说明问题，我们首先介绍两种求 的简单的 MC方法，然后给出几种较为复杂而更有效的 MC方法。()dxxfba= 在计算积分上， MC的实用场合是计算重积分 其中 是 维空间的点，当 较大时，用 MC方法比一般的数值方法有优点，主要是它的误差与维数 无关。( )dPPgIk=Pmmm随机投点法方法简述 ： 设 ，有限， ， ，并设 是在 上均匀分布的二维随机变量，其联合密度函数为 。则易见是 中 曲线下方的面积。假设我们向 中进行随机投点，若点落在 下方，(即 称为中的，否则称为不中，则点中的概率为 。若我们进行了 次投点，其中 次中的，则用频率来估计概率 。即 。ab( ) Mxf 0()Mybxayx = 0,:,( )YX,()()MybxaIabM0,1()dxxfba= ( )xfy =( )xfy =( )xfy ()abMp=n0npnnp0= 那么我们可以得到 的一个估计 具体试验步骤为()nnabM01=( )() ()()=+=用上式来估计的个数统计，和，计算，随机数，个独立地产生0,.,11,02nyxfxfMvyabuaxnivuUniiiiiiiii 注1 随机投点法的思想简单明了，且每次投点结果服从二项分布，故 ， 其中