Monte Carlo模拟
第三章 从概率分布函数的抽样
(Sampling from Probability Distribution
Functions)
舍选抽样法
(acceptance-rejection sampling)
舍选抽样法(舍选抽样法(acceptance-rejection sampling)acceptance-rejection sampling)
直接抽样法的困难:
• 许多随机变量的累积分布函数无法用解析函数给出;
• 有些随机变量的累积分布函数的反函数不存在或难
以求出;
• 即使反函数存在,但计算困难
舍选抽样法(von Neumann):
抽取随机变量x的一个随机序列xi, i=1,2,…, 按一定的
舍选规则从中选出一个子序列,使其满足给定的概率
分布.
Monte Carlo模拟
第三章 从概率分布函数的抽样
(Sampling from Probability Distribution
Functions)
舍选抽样法
(acceptance-rejection sampling)
. 简单舍选抽样法简单舍选抽样法
. 改进的舍选抽样法改进的舍选抽样法
. 典型的例子典型的例子
1. 1. 简单舍选抽样法简单舍选抽样法
舍选法抽样步骤:
1. 产生[a, b]区间内均匀分布的随
机数x: x = (b-a)r1+a, r1 U[0,
1];
2. 产生[0,c]区间内均匀分布的随
机数y: y = cr2, r2 U[0,1];
3. 当y f(x)时,接受x为所需的随
机数,否则,返回到第一步重新
抽取一对(x,y).
Von Neumann rejection method or Hit-and-miss method
设随机变量x的取值区间为x[a,b], 其概率密度函数f(x)
有界,即
抽取r1,r2 U[0,1]
x = a + (b-
a)r1
y = cr2
y f(x)
X = x
>
1. 1. 简单舍选抽样法简单舍选抽样法
a b x
f(x)
c
几何解释:
• 在二维图上,随机选取位于矩形abef内的点[x,y];
• 选取位于曲线f(x)下的那些点,则这些点将服从概率密度为f(x)的分布
e f
1. 1. 简单舍选抽样法简单舍选抽样法
证明:
按舍选抽样法抽出的随机数d的概率:
a b x
f(x)
c
e f
x和y的概率密度函数分别为
联合概率密度函数为
即d的概率函数为f(x)
d
1. 1. 简单舍选抽样法简单舍选抽样法
抽样效率:
对舍选抽样法:欲产生m个随
机变量x的值需产生n对(x,y)
,显然,m n
如果选出某特定分布的一个随机数平均地需要n个随机数r1 U[0, 1]
,则抽样效率定义为
a b x
f(x)
c
e f
d
Monte Carlo模拟
第三章 从概率分布函数的抽样
(Sampling from Probability Distribution
Functions)
舍选抽样法
(acceptance-rejection sampling)
. 简单舍选抽样法简单舍选抽样法
. 改进的舍选抽样法改进的舍选抽样法
. 典型的例子典型的例子
2. 2. 改进的舍选抽样法改进的舍选抽样法
改进的舍选抽样法
简单舍选抽样法的问题:
如果f(x)曲线下的面积占矩形面积的
比例很小,则抽样效率很低,这是因
为随机数x和y是在区间[a, b]和[0,
c]内均匀分布,所产生的大部分投点
不会落在f(x)曲线下
x
c
f(x)
改进方法:
构造一个新的概率密度函数g(x),使它的
形状接近f(x), 且有
式中Cg为常数,而g(x)的抽样相对比较容易。
Cgg(x)
2. 2. 改进的舍选抽样法改进的舍选抽样法
抽样方法:
1. 产生两个随机数
• 产生分布为g(x) 的随机数x ,x[a,b];
• 产生[0, Cgg(x)] 区间上均匀分布的随机数y
,y= Cgg (x) , U[0,1].
2. 接收或舍弃取样值 x.
• 如果 y > f(x),舍弃,返回到1,重复上述过程;
• 否则,接受;
2. 2. 改进的舍选抽样法改进的舍选抽样法
几何解释:
• 在二维图上,随机选取位于曲线Cgg(x)下的点[x,y];
• 选取位于曲线f(x)下的那些点,则这些点将服从概率密度为f(x)的分布
x
c
f(x)
Cgg(x)
2. 2. 改进的舍选抽样法改进的舍选抽样法
证明:
按舍选抽样法抽出的随机数d的概率:
d
x和y的概率密度函数分别为
联合概率密度函数为
即d的概率函数为f(x)
x
c
f(x)
Cgg(x)
2. 2. 改进的舍选抽样法改进的舍选抽样法
抽样效率:
x
c
f(x)
Cgg(x)
常数Cg的选取
• 常数Cg应尽可能地小,因为抽样效率与Cg成反比;
• Cg=max{f(x)/g(x)}, x [a,b]
Monte Carlo模拟
第三章 从概率分布函数的抽样
(Sampling from Probability Distribution
Functions)
舍选抽样法
(acceptance-rejection sampling)
. 简单舍选抽样法简单舍选抽样法
. 改进的舍选抽样法改进的舍选抽样法
. 典型的例子典型的例子
3. 3. 典型的例子典型的例子
例1:标准正态分布的抽样,x[-a,a]
无法用直接抽样法,累积分布函数无解析
表达式
Breit-wigner or Cauchy分布
3. 3. 典型的例子典型的例子
由g(x)抽取x 直接抽样法
抽取u
计算f(x), 如果u<= f(x), 接受x
3. 3. 典型的例子典型的例子
float gaussian_reject(double a)
{
const float c = ;
while(true) {
float eta = randac();
float x = tan(eta * * atan(a)+atan(-a));
float q = c * 1/*
float ksi = randac();
float u = ksi*q;
float p = 1/sqrt(2*)*exp(-x*x/);
if(u <= p) break;
}
return x;
}
3. 3. 典型的例子典型的例子
void test()
{
SetSeed(9,11);
c1 = new TCanvas("c1","Histogram Drawing
Options",200,10,700,900);
c1->Divide(1,2);
TH1F * h1 = new TH1F("h1","h1",100,,);
for(int i=0; i < 5000; i++) {
double x = gaussian_reject();
h1->Fill(x);
}
c1->cd(2);h1->Draw();
}
3. 3. 典型的例子典型的例子
3. 3. 典型的例子典型的例子
3. 3. 典型的例子典型的例子
A
B/2
例2:利用舍选法产生随机数C=cos, S=sin,其中为[0, 2]区
间内均匀分布的随机数
方法1:先产生[0, 2]间均匀分布的随机数: = 2 r,
rU[0,1], 然后直接计算C和S 因需要计算三角函数,
故此方法运算速度慢
方法2:利用舍选法可避免三角函数运算
令A和B为单位圆内直角三角形的两个边,则有
3. 3. 典型的例子典型的例子
因此,只要产生单位圆内的随机坐标A和B,就可用代数运算
求出C和S,算法为
1. 产生两个[0,1]区间上均匀分布的随机数u1和u2;
2. 令v1=2u1-1, v2=u2,则v1U[-1,1],v2U[0,1];
3. 计算r2=v12+v22, 如果r2 > 1,转到1,重新产生;
4. 令A=v1, B=v2, 计算C和S