蒙特卡罗算法算法.pptx
1从Buffon 投针问题谈起 第1页/共60页2Buffon 投针问题第2页/共60页3试验者时间(年)针长投针次数相交次数的估计值Wolf18500.80500025323.15956Smith18550.60320412183.15665Fox18840.7510304893.15951Lazzarini19250.83340818083.14159292第3页/共60页4数值积分问题第4页/共60页5Monte Carlo数值积分的优点与一般的数值积分方法比较,Monte Carlo方法具有以下优点:第5页/共60页6随机模拟计算的基本思路1.针对实际问题建立一个简单且便于实现的概率统计模型,使所求的量(或解)恰好是该模型某个指标的概率分布或者数字特征。2.对模型中的随机变量建立抽样方法,在计算机上进行模拟测试,抽取足够多的随机数,对有关事件进行统计3.对模拟试验结果加以分析,给出所求解的估计及其精度(方差)的估计4.必要时,还应改进模型以降低估计方差和减少试验费用,提高模拟计算的效率第6页/共60页7随机数的生成1.蒙特卡罗模拟的关键是生成优良的随机数。2.在计算机实现中,我们是通过确定性的算法生成 随机数,所以这样生成的序列在本质上不是随机 的,只是很好的模仿了随机数的性质(如可以通过 统计检验)。我们通常称之为伪随机数(pseudo-random numbers)。3.在模拟中,我们需要产生各种概率分布的随机数,而大多数概率分布的随机数产生均基于均匀分布U(0,1)的随机数。第7页/共60页8U(0,1)随机数的生成一个简单的随机数生成器:第8页/共60页9一个简单的例子第9页/共60页10一个简单的例子(续)上面的例子中,第一个随机数生成器的周期长度是 10,而后两个生成器的周期长度只有它的一半。我们自然希望生成器的周期越长越好,这样我们得到的分布就更接近于真实的均匀分布。第10页/共60页11线性同余生成器(Linear Congruential Generator)第11页/共60页12常用的线性同余生成器Modulus mMultiplier aReference231-1=214748364716807Lewis,Goodman,and Miller39373LEcuyer742938285Fishman and Moore950706376Fishman and Moore1226874159Fishman and Moore214748339940692LEcuyer214748356340014LEcuyer第12页/共60页13复杂一些的生成器(一)1.Combining Generators:第13页/共60页14复杂一些的生成器(二)2.Multiple recursive generator第14页/共60页15算法实现许多程序语言中都自带生成随机数的方法,如 c 中的 random()函数,Matlab中的rand()函数等。但这些生成器生成的随机数效果很不一样,比如 c 中的函数生成的随机数性质就比较差,如果用 c,最好自己再编一个程序。Matlab 中的 rand()函数,经过了很多优化。可以产生性质很好的随机数,可以直接利用。第15页/共60页16由rand()函数生成的U0,1随机数第16页/共60页17由rand函数生成的2维随机点第17页/共60页18从U(0,1)到其它概率分布的随机数U(0,1)的均匀分布的随机数,是生成其他概率分布随机数的基础,下面我们主要介绍两种将U(0,1)随机数转换为其他分布的随机数的方法。1.逆变换方法 (Inverse Transform Method)2.舍取方法 (Acceptance-Rejection Method)第18页/共60页19Inverse Transform Method第19页/共60页20Inverse Transform Method第20页/共60页21几个具体例子(一)第21页/共60页22几个具体例子(二)第22页/共60页23几个具体例子(三)第23页/共60页24标准正态分布随机数的生成正态分布是概率统计中最重要的分布,在此我们着重讨论如何生成标准正态分布随机数。引理:第24页/共60页25Box-Muller 算法第25页/共60页26逆变换方法(一)我们无法通过具体的数学表达式计算正态分布函数的逆函数,我们必须通过数值的方法逼近正态函数下面我们介绍 Beasley-Springer-Moro 方法。第26页/共60页27逆变换方法(二)第27页/共60页28逆变换方法(三)a0=2.50662823884b0=-8.47351093090a1=-18.61500062529b1=23.08336743743a2=41.39119773534b2=-21.06224101826a3=-25.44106049637b3=3.13082909833c0=0.3374754822726147c5=0.0003951896511919c1=0.9761690190917186c6=0.0000321767881768c2=0.1607979714918209c7=0.0000002888167364c3=0.0276438810333863c8=0.0000003960315187c4=0.0038405729373609在 matlab 中可以直接通过 norminv()函数直接计算标准正态分布函数的逆。第28页/共60页29Matlab生成的正态随机数第29页/共60页30Acceptance-Rejection Method(一)Acceptance-Rejection 方法最早由 Von Neumann提出,现在已经广泛应用于各种随机数的生成。基本思路:通过一个容易生成的概率分布 g 和一个取舍准则生成另一个与 g 相近的概率分布 f。第30页/共60页31Acceptance-Rejection Method(二)具体步骤:第31页/共60页32Acceptance-Rejection Method(三)下面我们验证由上述步骤生成的随机数 Y 确实具有密度函数 f(x)第32页/共60页33Acceptance-Rejection Method(四)所以为了提高舍取法的效率,我们应该使 c 的取值尽可能的小,也就是使 f 和 g 的分布更为相近。第33页/共60页34几个具体例子(一)第34页/共60页35几个具体例子(一)第35页/共60页36几个具体例子(二)第36页/共60页37几个具体例子(二)第37页/共60页38随机向量的抽样方法(一)第38页/共60页39随机向量的抽样方法(二)第39页/共60页40生成多维正态随机数的方法(一)第40页/共60页41生成多维正态随机数的方法(二)生成多维正态随机数的具体步骤:第41页/共60页42用Monte Carlo方法求解Laplace方程参见书上5.8节 P213P215第42页/共60页43马氏链在Monte Carlo随机模拟中的应用定义定义 为要模拟服从给定分布的随机变量,用生成一个易于实现的不可约遍历链 作为随机样本,使其平稳分布为 的方法,称为马氏链蒙特卡罗方法马氏链蒙特卡罗方法.蒙特卡罗方法的一个首要步骤是产生服从给定的概率分布函数 的随机变量(或称为随机样本),由概率论知识,熟知下面的结论.第43页/共60页44引理引理 生成随机变量U,使其分布满足U0,1,记为UU0,1,F(x)是给定的一个分布函数,记 为F(x)的反函数,则X=F-1(U)分布函数为F(x).第44页/共60页45第45页/共60页46第46页/共60页47米特罗波利斯(Metropolis)等人在1953年最早给出了通过生成一马氏链实现从分布 中采样(生成相关的样本)这一重要基本思想.随后,哈斯汀(Hastings)将其推广到更一般的形式.下面仅叙述状态空间S为至多可数的情形:第47页/共60页48第48页/共60页49第49页/共60页50Markov chain Monte Carlo(MCMC)问题提出:第50页/共60页51MCMC方法的基本思路MCMC 是一种简单有效的计算方法,在统计物理,Bayes 统计计算,显著性检验,极大似然估计等领域都有着广泛的应用。基本思路:第51页/共60页52概率转移核的构造(一)MCMC的方法有很多,在此我们只介绍Metropolis-Hastings方法。基本思路:第52页/共60页53概率转移核的构造(二)第53页/共60页54概率转移核的构造(三)Metropolis-Hastings 算法:第54页/共60页55第55页/共60页56(续上页证明)第56页/共60页57Metropolis-Hastings 算法的具体步骤第57页/共60页58几种常用的 q(x,x)(一)1.Metropolis选择:第58页/共60页59几种常用的 q(x,x)(二)2.独立抽样:第59页/共60页60感谢您的观看!第60页/共60页