Monte Carlo Simulation Methods(蒙特卡罗模拟方法)主要内容:1.各种随机数的生成方法.2.MCMC方法.12从Buffon 投针问题谈起Buffon 投针问题:平面上画很多平行线，间距为a ．向此平面投掷长 为l (l < a) 的针, 求此针与任一平行线相交的概率p 。

22[0,/2] [0,] sin ,{:sin }.lla X A X随机投针可以理解成针的中心点与最近的平行线的距离X 是均匀地分布在区间 上的r.v.，针与平行线的夹角是均匀地分布在区间 上的r.v.，且X 与相互独立，于是针与平行线相交的充要条件为 即相交3Buffon 投针问题2sin22(sin )2lll p P Xdxdaa 于是有：2lap若我们独立重复地作n 次投针试验，记()nA 为A 发生的次数。

()n f A 为A在n 次中出现的频率。

假如我们取()n f A 作为()p P A 的估计，即ˆ()n pf A 。

然后取2ˆ()n laf A 作为的估计。

根据大数定律，当n 时，..ˆ().a s n pf A p从而有2ˆ()Pn l af A 。

这样可以用随机试验的方法求得的估计。

历史上有如下的试验结果。

43.14159292180834080.831925Lazzarini3.1595148910300.751884Fox 3.15665121832040.601855Smith 3.15956253250000.801850Wolf π的估计值相交次数投针次数针长时间(年)试验者5数值积分问题1()()() ~[0,1] } (1)() n ni i n f x dx Ef X f x X U k n i i d f U nwith probability as n 1k k 计算积分我们可以将此积分看成 的数学期望。

其中(均匀分布)。

于是可以将上式积分看作是f (X )的数学期望.若{U,1为U~U[0,1].则可以取作为的估计，由大数定律，可以保证收敛性，即：这表明可以用随机模拟的方法计算积分。

6Monte Carlo 数值积分的优点与一般的数值积分方法比较，Monte Carlo 方法具有以下优点：1. Monte Carlo 一般的数值方法很难推广到高维积分的情形，而方法很容易推广到高维情形2/1/22. ()() dO n O n一般的数值积分方法的收敛阶为 ，而由中心极限定理可以保证 M on te C arl o 方法的收敛阶为 。

此收敛阶与维数无关，且在高维时明显优于一般的数值方法。

随机模拟计算的基本思路1.针对实际问题建立一个简单且便于实现的概率统计模型，使所求的量（或解）恰好是该模型某个指标的概率分布或者数字特征。

2.对模型中的随机变量建立抽样方法，在计算机上进行模拟测试，抽取足够多的随机数，对有关事件进行统计3.对模拟试验结果加以分析，给出所求解的估计及其精度(方差)的估计4.必要时，还应改进模型以降低估计方差和减少试验费用，提高模拟计算的效率7随机数的生成1.蒙特卡罗模拟的关键是生成优良的随机数。

2.在计算机实现中,我们是通过确定性的算法生成随机数,所以这样生成的序列在本质上不是随机的,只是很好的模仿了随机数的性质(如可以通过统计检验)。

我们通常称之为伪随机数(pseudo-random numbers)。

3.在模拟中，我们需要产生各种概率分布的随机数，而大多数概率分布的随机数产生均基于均匀分布U(0,1)的随机数。

89U(0,1)随机数的生成一个简单的随机数生成器：1101 mod ,, /iii i i x ax m u a x x m x m其中 均为整数, 可以任意选取。

111 () , ()i i i i x f x u g x 随机数生成器的一般形式：10一个简单的例子1i+116 mod11, u /11 6,11ii i x x x a m （）1 ,x 1,6,3,7,9,10,5,8,4当时得到序列:,1,6,3.,2 (00)3, 1 ,,1,3,9........3,2,2,1,3,9,5,42,6,7,10,8, 6.......a x ax 如果令 得到序列:如果令 得到序列：11一个简单的例子(续)上面的例子中，第一个随机数生成器的周期长度是 10，而后两个生成器的周期长度只有它的一半。

我们自然希望生成器的周期越长越好，这样我们得到的分布就更接近于真实的均匀分布。

0 (m a x 在给定 的情况下，生成器的周期与和初值种子)选择有关。

12线性同余生成器(Linear Congruential Generator )111 () mod /i ii i x ax c m u x m一般形式：0 c c 1.可以证明 的选择对生成的随机数的均匀性影响不大，所以为了提高计算速度，一般都令 。

02. 1 m m a x m 线性同余器可以达到的最长周期为 ，我们可以通过适当的选择 和，使无论选取怎样的初值 都可以达到最大周期(一般选取 为质数)13常用的线性同余生成器L’Ecuyer 400142147483563L’Ecuyer 406922147483399Fishman and Moore 1226874159Fishman and Moore 950706376Fishman and Moore 742938285L’Ecuyer 39373Lewis, Goodman, and Miller 168072^31-1=2147483647Reference Multiplier a Modulus m14复杂一些的生成器（一）bining Generators:,1,,1,111,11111111111 mod , / 1,2....(1) mod(1)/ 0(1)/ 0( )j i j j i j j i j i j J j ij i j i i i i j J x a x m u x m j J x x m u x m x m m x m m 考虑个简单线性同余生成器：其中为所有中最大的一个15复杂一些的生成器（二）2.Multiple recursive generator1122102(......) mo ,.......)d /ii i k i k i k k i x a x a x a x m u x x x x m需(要选取种子12--1-1 (,......) 1i j i i i k k k x m x x x m m每个有种选择，所以向量 可以取个不同的值，所以这样的随机数生成器的最大周期可以达到 ，大大提高了简单同余生成器的周期。

算法实现许多程序语言中都自带生成随机数的方法，如 c 中的 random() 函数，Matlab中的rand()函数等。

但这些生成器生成的随机数效果很不一样，比如 c 中的函数生成的随机数性质就比较差，如果用c ，最好自己再编一个程序。

Matlab 中的 rand() 函数，经过了很多优化。

可以产生性质很好的随机数，可以直接利用。

16由rand()函数生成的U[0,1]随机数17由rand函数生成的2维随机点18从U(0,1)到其它概率分布的随机数U(0,1)的均匀分布的随机数，是生成其他概率分布随机数的基础，下面我们主要介绍两种将U(0,1)随机数转换为其他分布的随机数的方法。

1.逆变换方法(Inverse Transform Method)2.舍取方法(Acceptance-Rejection Method)1920Inverse Transform Method 11 () ()inf{:()}, 01()()X F x F y x F x y y F y F x 设随机变量 的分布函数为，定义：即是的反函数。

11 (0,1) () () U X F U F x 定理 ：设随机变量服从 上的均匀分布，则的分布函数为 。

11() ()(())(())()F U P X x P F U x P U F x F x 由 的定义和均匀分布的分布函明数可证：得：21Inverse Transform Method 11 () (0,1) () F x U u F u 由定理，要产生来自的随机数，只要先产生来自随机数，然后计算 即可。

具体步骤如下：(1) (0,1)U 上生成均匀分布的随机数。

-1(2) , () () X X U F F x 计算则为来自 分布的随机数.22几个具体例子（一）-11 ~(,) 0 () 1, ()() 01(0,1) (-) (,) X U a b xa x a F x a xb b ax bF y a b a y y U U a b a U U a b 例 ：设，则其分布函数为，生成随机数，则是来自的随机数。

23几个具体例子（二）/12~exp() ()1, 0 ()log(1) log(1)( )1 x X X F x ex F y y X U U U U 例 ：设服从指数分布，则的分布函数为：通过计算得，则：服从指数分布其中服从均匀分布又因为－和有着同样的分布，所以也可以取： log()X U24几个具体例子（三）12013 ()~() ,...... () 0, () n ii i i j j i i X F x c c c P X c p q q p F c q 例 离散分布随机数的生成设为取值的离散随机变量且，令 ，则有，我们按照如下步骤生成随机数：(0,1) ;U (1) 生成上均匀分布随机数1 , 1;k k k q U q k n (2) 寻找满足 ()K X c X F x k k-1k k (3) 令，显然P (X =c )=P (q <U q )=p 则的分布函数是。

25标准正态分布随机数的生成正态分布是概率统计中最重要的分布，在此我们着重讨论如何生成标准正态分布随机数。

引理：121/21121/221212, (0,1) (2ln )cos(2) (2ln )sin(2) U U U Z U U Z U U Z Z 设 是独立同分布的变量，令：则证明：与 独由随机立，且均服从向量的变换公标准式易正态分布。

得，略。

26Box-Muller 算法121. ,U U 生成两个独立的均匀分布随机数 1/21121/221212 (2ln )cos(2) (2ln )sin(2) Z U U Z U U Z Z 2.令 则，为独立的正态随机数27逆变换方法（一）我们无法通过具体的数学表达式计算正态分布函数的逆函数，我们必须通过数值的方法逼近正态函数下面我们介绍 Beasley-Springer-Moro 方法。