Monte Carlo Simulation Methods (蒙特卡罗模拟方法) 主要内容: 1.各种随机数的生成方法. 2.MCMC方法.
从Buffon 投针问题谈起
Buffon 投针问题
试验者 时间(年) 针长 投针次数 相交次数 π的估计值 Wolf 1850 0.80 5000 2532 3.15956 Smith 1855 0.60 3204 1218 3.15665 Fox 1884 0.75 1030 489 3.15951 Lazzarini 1925 0.83 3408 1808 3.14159292
数值积分问题
Monte Carlo数值积分的优点 与一般的数值积分方法比较,Monte Carlo方法 具有以下优点:
随机模拟计算的基本思路 1.针对实际问题建立一个简单且便于实现的概率统计模 型,使所求的量(或解)恰好是该模型某个指标的概率分布或者数字特征。 2.对模型中的随机变量建立抽样方法,在计算机上进行 模拟测试,抽取足够多的随机数,对有关事件进行统计 3.对模拟试验结果加以分析,给出所求解的估计及其精 度(方差)的估计 4.必要时,还应改进模型以降低估计方差和减少试验费 用,提高模拟计算的效率
随机数的生成 1.蒙特卡罗模拟的关键是生成优良的随机数。 2.在计算机实现中,我们是通过确定性的算法生成 随机数,所以这样生成的序列在本质上不是随机 的,只是很好的模仿了随机数的性质(如可以通过 统计检验)。我们通常称之为伪随机数(pseudo-random numbers)。 3.在模拟中,我们需要产生各种概率分布的随机数,而大多数概率分布的随机数产生均基于均匀分布U(0,1)的随机数。
U(0,1)随机数的生成 一个简单的随机数生成器:
一个简单的例子
一个简单的例子(续) 上面的例子中,第一个随机数生成器的周期长度是 10,而后两个生成器的周期长度只有它的一半。我们自然希望生成器的周期越长越好,这样我们得到的分布就更接近于真实的均匀分布。
线性同余生成器 (Linear Congruential Generator )
常用的线性同余生成器 Modulus m Multiplier a Reference 2^31-1 =2147483647 16807 Lewis, Goodman, and Miller 39373 L’Ecuyer 742938285 Fishman and Moore 950706376 1226874159 2147483399 40692 2147483563 40014
复杂一些的生成器(一) 1.Combining Generators:
复杂一些的生成器(二) 2.Multiple recursive generator
算法实现 许多程序语言中都自带生成随机数的方法,如 c 中的 random() 函数,Matlab中的rand()函数等。 但这些生成器生成的随机数效果很不一样,比如 c 中的函数生成的随机数性质就比较差,如果用 c ,最好自己再编一个程序。Matlab 中的 rand() 函数,经过了很多优化。可以产生性质很好的随 机数,可以直接利用。
由rand()函数生成的U[0,1]随机数
由rand函数生成的2维随机点
(Inverse Transform Method) 2.舍取方法 (Acceptance-Rejection Method) 从U(0,1)到其它概率分布的随机数 U(0,1)的均匀分布的随机数,是生成其他概率 分布随机数的基础,下面我们主要介绍两种将 U(0,1)随机数转换为其他分布的随机数的方法。 1.逆变换方法 (Inverse Transform Method) 2.舍取方法 (Acceptance-Rejection Method)
Inverse Transform Method
Inverse Transform Method
几个具体例子(一)
几个具体例子(二)
几个具体例子(三)
标准正态分布随机数的生成 正态分布是概率统计中最重要的分布,在此 我们着重讨论如何生成标准正态分布随机数。 引理:
Box-Muller 算法
逆变换方法(一) 我们无法通过具体的数学表达式计算正态分布函数 的逆函数,我们必须通过数值的方法逼近正态函数 下面我们介绍 Beasley-Springer-Moro 方法。
逆变换方法(二)
逆变换方法(三) 在 matlab 中可以直接通过 norminv() 函数直接 计算标准正态分布函数的逆。 c0=0.3374754822726147 c5=0.0003951896511919 c1=0.9761690190917186 c6=0.0000321767881768 c2=0.1607979714918209 c7=0.0000002888167364 c3=0.0276438810333863 c8=0.0000003960315187 c4=0.0038405729373609 在 matlab 中可以直接通过 norminv() 函数直接 计算标准正态分布函数的逆。
Matlab生成的正态随机数
Acceptance-Rejection Method(一) Acceptance-Rejection 方法最早由 Von Neumann 提出,现在已经广泛应用于各种随机数的生成。 基本思路: 通过一个容易生成的概率分布 g 和一个取舍 准则生成另一个与 g 相近的概率分布 f 。
Acceptance-Rejection Method(二) 具体步骤:
Acceptance-Rejection Method(三) 下面我们验证由上述步骤生成的随机数 Y 确实 具有密度函数 f(x)
Acceptance-Rejection Method(四) 可能的小,也就是使 f 和 g 的分布更为相近。
几个具体例子(一)
几个具体例子(一)
几个具体例子(二)
几个具体例子(二)
随机向量的抽样方法(一)
随机向量的抽样方法(二)
生成多维正态随机数的方法(一)
生成多维正态随机数的方法(二) 生成多维正态随机数的具体步骤:
用Monte Carlo方法求解Laplace方程 参见书上5.8节 P213~P215
马氏链在Monte Carlo随机模拟中的应用 定义 为要模拟服从给定分布的随机变量,用生成一个易于 实现的不可约遍历链 作为随机样本, 使其平稳分布为 的方法,称为马氏链蒙特卡罗方法. 蒙特卡罗方法的一个首要步骤是产生服从给定的概率分布函数 的随机变量(或称为随机样本),由概率论知识,熟知下面的结论.
引理 生成随机变量U,使其分布满足U[0,1],记为U~U[0,1], F(x)是给定的一个分布函数,记 为F(x)的反函数,则X=F-1(U)分布函数为F(x).
米特罗波利斯(Metropolis)等人在1953年最早 给出了通过生成一马氏链实现从分布 中采 样(生成相关的样本)这一重要基本思想.随后, 哈斯汀(Hastings)将其推广到更一般的形式. 下面仅叙述状态空间S为至多可数的情形:
Markov chain Monte Carlo (MCMC) 问题提出:
MCMC方法的基本思路 基本思路: MCMC 是一种简单有效的计算方法,在统计物理, Bayes 统计计算,显著性检验,极大似然估计等领 域都有着广泛的应用。 基本思路:
概率转移核的构造(一) MCMC的方法有很多,在此我们只介绍Metropolis-Hastings方法。 基本思路:
概率转移核的构造(二)
概率转移核的构造(三) Metropolis-Hastings 算法:
(续上页证明)
Metropolis-Hastings 算法的具体步骤
几种常用的 q(x,x’)(一) 1.Metropolis选择:
几种常用的 q(x,x’)(二) 2.独立抽样: