实验数据处理方法
第二部分,Monte Carlo模拟
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
第八章 从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
? Monte Carlo算法的一个重要组成部分,
描述所要模拟的物理系统的一些概率密度函数（ PDF)
?描述整个系统在空间、能量、时间或多维相空间中的
发展和演化；
? Monte Carlo模拟的主要任务,
? 通过对这些概率密度函数的随机抽样来模拟物理系统的状
态 ;
? 为描述系统的演化所必需的一些附加运算,
? 物理过程的描述 ?从描述物理系统的 pdf出发，随机抽取系统的
可能状态。
第八章 从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
? 本章介绍从一个任意的 pdf获取样本的抽样方法。
1,直接抽样法（反函数法）
2,变换抽样法 ?直接抽样法的一般形式
3,舍选抽样法
4,复合分布的抽样方法
5,混合抽样法
6,近似抽样法（列表法 )
7,多维分布的抽样
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.1 等价的连续概率密度函数
8.1 等价的连续概率密度函数
? 随机变量：连续型、分离型
? 概率密度函数：连续分布、分离分布
? 利用 ?函数，可将分离型的 pdf用连续型的 pdf 描述
? 用同样的方式来讨论分离型和连续型随机变量的抽样方法
8.1 等价的连续概率密度函数
?
?
??
N
i
ii xxpxf
1
)()( ?
?? ?
? ??
??
?
??
?
??
?
?
??
???
?
?
?
?
?
?
???
N
i
ii
N
i
ii
N
i
ii
pxdxxxxp
dxxxpxdxxxfXE
11
1
)(
)()()(
?
?
?? ?
? ??
??
?
??
?
??
?
?
??
?????
?
?
?
?
?
?
?????
N
i
ii
N
i
ii
N
i
ii
pXExdxxxpXEx
dxxxpXExdxxfXExXV
1
2
1
2
1
22
)]([)())((
)()]([)())(()(
?
?
? 已知分离型 pdf,{Pi} ? 分离型随机变量 X的取值为 xi的概率
? 定义一个等价的连续型 pdf,
?
?
?
??
?
??
??
??
)()()(
1)(
ii
i
xfdxxxxf
dxxx
?
?
? 利用与连续型随机变量相同的方式计算分离型随机变量的期望值
和方差,
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.2 pdf的变换
8.2 pdf的变换
x,连续型的随机变量,PDF,f(x)
y = y(x),x的函数,也是随机变量, 求 y(x)的概率密度函数 g(x)
1、若随机变量 x和 y是一一对应的,
2,若随机变量 x和 y不是一一对应的,
[x,x+dx]?[y,y+dy]
X的取值在 [x,x+dx]的概率 ==Y的取值
在 [y,y+dy]的概率,
dy
dxxfyg )()( ?
取绝对值是为了保证 g(y)是非负的
f(x)dx=g(y)dy ?
即有 n个区间 [x,x+dx]?[y,y+dy]
?? dydxxfyg )()(
需要对这 n个区间求和
8.2 pdf的变换
3、推广到 n个随机变量的情况,
Jacobian行列式 ?
nixyy
yyyyxxxx
ii
nn
,,2,1),(
),,(),,( 2121
??
????
??
???
n
nnn
n
y
x
y
x
y
x
y
x
y
x
y
x
n
n
yyy
xxx
J
Jxfyg
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
21
1
2
1
1
1
,
,
),,(
),,(
)()(
21
21
4、特例：如果 y(x)是 x的累积分布函数 (cdf)
)(
1
)(
11
)()()(
xfxF
dx
dydy
dx
xdxfxFxyy
x
?
?
??
????? ?
??
10,1)( ???? yyg
即,y在 [0,1]区间上均匀分布 ? 不管 f(x)取何种形式,累积分布
函数总是在 [0,1]区间上均匀分布
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
)(1 yFx ??
原理 （ 注意,pdf f(x)必须是归一化的 ）,
? 设 y=F(x)为随机变量 x的累积分布函数 ? x和 y是一一对应的
? 先随机抽取 y，然后通过求 F(x)的反函数 F-1(y)得到随机变量 x
的值
? 随机变量 y在区间 [0,1]上均匀分布 ? 利用 [0,1]区间上均匀分布
随机数产生器抽取
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
方法,U[0,1],[0,1]区间上均匀分布的随机数
1,从 U[0,1]抽取随机数 ?；
2,令 F(x)= ?;
3,解方程得 x,)(1 ??? Fx
注：需要知道累积分布函数的解析表达式，且累积分布函数的
反函数存在
分离型随机变量的抽样
直接抽样法适应于分离型的随机变量
?
? ?
?
?
??
?
??
?
??
?
?
?
?
?
???
????
K
i
i
x
N
i
ii
x
p
xdxxp
xdxfxFy
1
1
)(
)()(
?
?
?
?? N
i
ii xxpxf
1
)()( ?
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
方法,
1,计算 yk = yk-1 + pk,k = 2,3,…,N,y 1 = p1
2,从 U[0,1]抽取随机数 ?；
3,求满足 yk-1 < ? < yk 的 K值 ;
4,随机变量的第 k个取值即为欲抽取的值。
0 ? 1
p a pb pc pd
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
p3=0.2 b3+c3
p2=0.3 b2+c2
p1=0.5 b1+c1 a?
13.05.0
3.05.05.0
5.00
???
???
??
?
?
?
33
22
11
cb
cb
cb
?
?
?
例 1、粒子衰变末态的随机抽样
设粒子 a有三种衰变方式，其分支比如下
随机选取每次衰变的衰变方式（衰变道） ?直接抽样法
? ?U[0,1]
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
??
?
?
??
?
?
?
?
?
???
?
?
??
?
?
????
?
?
??
?
?
? ?
rn
n
rnr
n
r
n
nrpp
r
n
pnrB rnr
)!(!
!
,2,1,0,)1(),;( ?例 2、二项式分布的抽样
方法 1,利用上面介绍的直接抽样法，需计算累积分布函数，
当 n很大时，求和计算困难 ;
方法 2,利用二项式分布的定义
1,产生 n个 ? i?U[0,1];
2,统计满足条件 ? i <p（表示成功）的 ? i的数目 r，则 r
表示在 n次实验中成功的次数 ?r即为二项式分布的
抽样值
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
?0,1,2,=r,e r!1 =μ) P ( r ; -r ??
????? nppn,0,
例 3、泊松分布的抽样
方法 1,利用直接抽样法，但计算累积分布函数时非常复杂
方法 2,利用泊松分布的定义：二项式分布的极限形式
1,选取足够大的 n，使 p=?/n相当小，例如,p=0.1
2,产生 n个 ? i?U[0,1];
3,统计满足条件 ? i <p（表示成功）的 ? i的数目 r，则 r
表示在 n次实验中成功的次数 ?r即为泊松分布的抽
样值的近似值,n越大，近似程度越好
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
??
???
??
????
bxax
bxaabxf
,,0
,1)(
)0()( ?? ? xexf x??
例 4、连续型随机变量的直接抽样
1,求区间 [a,b]上均匀分布的随机数 x,
? ?????? xa ab axxdxf )(?
aabx ??? )(?
?产生 ? ?U[0,1];
?
?
2,指数分布
xx x exde ???? ??? ???? ? 1
0
???? ln1)1l n (1 ?????x
?产生 ? ?U[0,1];
?
? ?和（ 1-?）都是 U[0,1]
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.4 变换抽样法 ?直接抽样法的一般形式
8.4 变换抽样法 ?直接抽样法的一般形式
思路,
dy
dxxfyg )()( ?
? 随机变量 y,pdf,g(y) ? 不易进行抽样
? 随机变量 x,pdf,f(x) ? 容易进行抽样
? 如果能够找到 x和 y之间的一个一一对应的变换关系，
y=y(x),使得 g(y)和 f(x)满足关系
则可先由 f(x)分布抽取 x的值 ?，再由变换关系得到 y的值
???满足分布 g(y)
?变换抽样法
8.4 变换抽样法 ?直接抽样法的一般形式
变换抽样法,
?找出 y与 x间的变换关系,y=y(x),f(x)与 g(y)满足关系,
?由 f(x)分布抽取 x的值 ?；
?随机变量 y的取值,?=y(?)
dy
dxxfyg )()( ?
直接抽样法是变换抽样法的一个特殊形式
?X满足 U[0,1],f(x)=1;
?y与 x间的变换关系,y=G-1(x) ?y的累积分布函数的反函
数
? ?? ???? y ydygyGx )()(
8.4 变换抽样法 ?直接抽样法的一般形式
nixyy
yyyyxxxx
ii
nn
,,2,1),(
),,(),,( 2121
??
????
??
???
n
nnn
n
y
x
y
x
y
x
y
x
y
x
y
x
n
n
yyy
xxx
J
Jxfyg
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
21
1
2
1
1
1
,
,
),,(
),,(
)()(
21
21
推广到 n维随机向量的情况,
)(xf ? x? },,,{ 21 n???? ?? ?
)(?? ?ii y?
–由联合概率密度函数 抽取随机向量 的值
–yi的值,
8.4 变换抽样法 ?直接抽样法的一般形式
2
2
2
)(
2
1),( ??
????
??
?
x
eN
???? xy
?? ?? yx
例：高斯分布的抽样方法
进行变量变换,
)1,0(),( NN ???? ?标准正态分布
由 N(0,1)分布抽样得到 y,
如何抽取服从 N(0,1)分布的随机变量？
8.4 变换抽样法 ?直接抽样法的一般形式
???
?
???
??
?
??
?
? ?? ???
???
n
i
i
n
i
i
n
i
ixy
1
2
11
??
1,利用中心极限定理
设 x1,x2,…,x n是一组 n个独立的随机变量,xi的平均值和方差
分别为 μi和 бi，则当 n→∞ 时，变量
服从标准正态分布 N(0,1)
设 x 是在 [0,1]之间均匀分布的随机数,对 n个 x的取值 xi
121221 )()( ???? xVxE ii ?? ?
????????????? ?? ?
? 1221
nnxy n
i i
在 n→∞ 时,y服从 N(0,1)，在实际应用时，可取 n=12
抽样方法,
1、产生 12个 U[0,1]的随机数 ?i
612
1
?? ?
?i
iy ?
2,
8.4 变换抽样法 ?直接抽样法的一般形式
22
2121
2221
2
1
2
1)()(),( yy eeygygyyg ?? ????
??
1),(:],1,0(,
)2s in ()ln2(
)2c o s ()ln2(
2121
212
211
??
??
??
xxfpdfUxx
xxy
xxy
?
?
1
21
2
2
1
2
1
2
2
2
1
y
ytgx
ex
yy
?
??
?
?
?
2,利用变换抽样法
y1和 y2是两个相互独立的、服从标准正态分布的随机数
变换,反变换,
抽样方法,
)2s in ()ln2(
)2c o s ()ln2(
2122
2111
????
????
???
???
y
y
1）产生一对 [0,1]区间均匀分布的随机数 ?1和 ?2；
2）
8.4 变换抽样法 ?直接抽样法的一般形式
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
2
2
2
2
2
1
2
2
2
2
1
2
2
2
1
2
2
2
1
2
2
1
2
2
1
1
1
2
1
2
1
2
1
1
1
2
1
1
2
1
,
,
),(
),(
2
1
2
1
2
2
2
1
2
2
1
2
2
1
21
21
yyyy
eee
y
y
y
ey
y
y
y
y
ey
yy
xx
J
yyyy
y
x
y
x
y
x
y
x
???
??
证明用此方法抽取的 y1,y2满足上面的联合概率分布
雅可比行列式,
JJxxfyyg ??? ),(),( 2121
8.4 变换抽样法 ?直接抽样法的一般形式
222
2222
2
1
2
1
2
1)()(),( yxyx eeeyfxfyxf ???? ?????
???
x
y
yxr
?
???
?
?
t a n
22
222
??
??
s i n2
c o s2
?
?
y
x
?
?
???
?
???
?
??
?? ???
?
?
?
?
? eyxfyxfyxf
yx
g
2
1
),(),(
c o s2,s i n
2
1
s i n2,c o s
2
1
),(
),(
),(
),(
3,采用极坐标变换
x和 y是两个相互独立的、服从标准正态分布的随机数
变换,反变换,
则变量 ?和 ?的联合概率密度函数为
即 ?在 [0,2?] 区间上均匀分布,?服从 ? = 1的 指数分布
8.4 变换抽样法 ?直接抽样法的一般形式
抽样方法,
–产生 ?, ?= 2? r,r?U[0,1]
–产生 ?, ? = -ln r,r?U[0,1]
– x= (2?)1/2 cos?, y=(2?)1/2 sin?
8.4 变换抽样法 ?直接抽样法的一般形式
void test() {
c1 = new TCanvas("c1","Histogram Drawing Options",200,10,700,900);
c1->Divide(1,2);
TH1F * h1 = new TH1F("h1","h1",100,-5.0,5.0);
TH1F * h2 = new TH1F("h2","h2",100,-5.0,5.0);
SetSeed(9,11)
Const Pi = 3.1415926;
for(int i=0; i < 50000; i++) {
double x = randac();
double y = randac();
float y1 = sqrt(-2.0*log(x))*cos(2.0*Pi*y);
float y2 = sqrt(-2.0*log(x))*sin(2.0*Pi*y);
h1->Fill(y1); h2->Fill(y2);
}
c1->cd(1);h1->Draw(); c1->cd(2);h2->Draw();
}
8.4 变换抽样法 ?直接抽样法的一般形式
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.5 舍选抽样法
（ acceptance-rejection sampling)
8.5 舍选抽样法
（ acceptance-rejection sampling)
直接抽样法的困难,
? 许多随机变量的累积分布函数无法用解析函数给出 ;
? 有些随机变量的累积分布函数的反函数不存在或难以求出；
? 即使反函数存在，但计算困难
舍选抽样法,
抽取随机变量 x的一个随机序列 xi,i=1,2,…,按一定的舍选规则
从中选出一个子序列，使其满足给定的概率分布,
假定,
? 随机变量 x的值域为 [a,b];
? x的概率密度函数,f(x)=P*(x)/Z,（其中 Z为归一化因子） ?
难以直接抽样 ;
? Q(x)=Q*(x)/ZQ 是另外一个较为简单的函数 (ZQ为归一化因子）
?可用简单的方法进行抽样 ;
? 在 x的整个取值范围内,cQ*(x) > P*(x),其中 c为一常数
8.5 舍选抽样法
（ acceptance-rejection sampling)
抽样方法,
1,产生两个随机数
? 从 Q(x)分布抽取 ?, ??[a,b];
? 由 [0,cQ*(?)] 区间上的均匀分布抽取 u,u= cQ*(?) ?,
??U[0,1],
2,接收或舍弃取样值 x,
? 如果 u > P*(x)，舍弃，返回到 1,重复上述过程 ;
? 否则，接受 ;
8.5 舍选抽样法
（ acceptance-rejection sampling)
几何解释,
? 在二维图上，随机选取位于曲线 cQ*(x)下的点 [x,u]；
? 选取位于曲线 P*(x)下的那些点，则这些点将服从概率密度为 f(x)的分布
8.5 舍选抽样法
（ acceptance-rejection sampling)
常数 c的选取
? 常数 c应尽可能地小,因为抽样效率与 c成反比 ;
? c=max{P*(x)/Q*(x)},x ?[a,b]
特例,
如果取 Q(x)=1,x?[a,b],即均匀分布,则
? X的抽样,x=(b-a)r+a,r ?U[0,1]
? c可取为 f(x)在 [a,b]区间上的极大值
a b x
f(x)
c
抽取 r1,r2 ?U[0,1]
? = a + (b-a)r1
r2 ? f(?)/c
x = ?
? >
8.5 舍选抽样法
（ acceptance-rejection sampling)
2/2
2
1)(* xexp ??
?
21
11)(*
xxQ ?? ?
? ? 52.12)1(2m a x)(*/)(*m a x 2/2 2 ??
??
?
??
? ??? ?
eexxQxPc
x ??
? ?)a r c t a n (a r c t a n
a r c t a n2
1
)()(
1
1
a r c t a n2
1
)(*/)(*)(
2
ax
a
xdxQxF
xa
dxxQxQxQ
x
aQ
a
a
??????
?
??
?
?
?
?
例 1:标准正态分布的抽样,x?[-a,a]
无法用直接抽样法,累积分布函数无解析
表达式
Breit-wigner or Cauchy分布
8.5 舍选抽样法
（ acceptance-rejection sampling)
? ?)a r c t a n (a r c t a n2t a n
]1,0[
aax
U
???
?
?
?
Qu
U
x
xcQQ
?
?
?
?
?
?
??
]1,0[
1
152.1
)(*
2
由 Q(x)抽取 x ? 直接抽样法
抽取 u
计算 P*(x),如果 u<= P*(x),接受 x
8.5 舍选抽样法
（ acceptance-rejection sampling)
float gaussian_reject(double a)
{
const float c = 1.52;
while(true) {
float eta = randac();
float x = tan(eta * 2.0 * atan(a)+atan(-a));
float q = c * 1/3.1415926*1.0/(1+x*x);
float ksi = randac();
float u = ksi*q;
float p = 1/sqrt(2*3.1415926)*exp(-x*x/2.0);
if(u <= p) break;
}
return x;
}
8.5 舍选抽样法
（ acceptance-rejection sampling)
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,-5.0,5.0);
for(int i=0; i < 5000; i++) {
double x = gaussian_reject(5.0);
h1->Fill(x);
}
c1->cd(2);h1->Draw();
}
8.5 舍选抽样法
（ acceptance-rejection sampling)
2c o s2s i n2s i n2s i n2c o sc o s 22
?????? ????? SC
A
B ?/2
2222
22
2222
2
s i n,c o s
2
s i n,
2
c o s
BA
AB
S
BA
BA
C
BA
B
BA
A
?
??
?
?
???
?
?
?
?
??
??
例 2：利用舍选法产生随机数 C=cos?,S=sin?,其中 ?为 [0,2?]区间内
均匀分布的随机数
方法 1：先产生 [0,2?]间均匀分布的随机数,? = 2? r,
r?U[0,1],然后直接计算 C和 S ? 因需要计算三角函数，
故此方法运算速度慢
方法 2：利用舍选法可避免三角函数运算
令 A和 B为单位圆内直角三角形的两个边，则有
8.5 舍选抽样法
（ acceptance-rejection sampling)
因此，只要产生单位圆内的随机坐标 A和 B,就可用代数运算
求出 C和 S,算法为
1,产生两个 [0,1]区间上均匀分布的随机数 u1和 u2;
2,令 v1=2u1-1,v2=u2,则 v1?U[-1,1],v2?U[0,1];
3,计算 r2=v12+v22,如果 r2 > 1,转到 1，重新产生；
4,令 A=v1,B=v2,计算 C和 S
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.6 复合分布的抽样方法
（ composition method)
8.6 复合分布的抽样方法
（ composition method)
??
?
???
?
1)(;1;0
)()(
dxxfpp
xfpxf
k
k
kk
k
kk
??
?
?
?
??? k
i
i
k
i
i prpUr
1
1
1
];1,0[
1961年由 Marsaglia提出的方法
设随机变量 X的概率密度函数 f(x)可写成一些 PDF的线性叠加,
抽样方法,
1,利用离散型的随机变量的抽样方法抽取序号 k;
2,由 fk(x)抽取 x
8.6 复合分布的抽样方法
（ composition method)
)0()()0()(5.0
0
01
)(
0
01
)(
)(5.0)(5.0)(
5.0)(:
2121
)0[)0,(
),0[)0,(
???????
?
?
? ???
?
?
?
? ????
?
??
?
?
???
?
?
??
?
xexfxexfpp
x
xI
x
xI
xIexIexf
exfP D F
xx
xx
x
其他其他
，
例：用复合法产生双指数分布随机数
– 产生两个 [0,1]区间均匀分布的随机数 r1和 r2;
– 如果 r1 ? 0.5,按 f1(x)抽样；
– 如果 r1 > 0.5,按 f2(x)抽样；
– f1(x)和 f2(x)都可用直接抽样法
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.7 混合抽样法
（ mixed method)
8.7 混合抽样法
(mixed method)
?
混合法,
直接抽样法
舍选法 ?亦称乘抽样法
适用的抽样场合,
直接抽样法不可用、舍选抽样法效率低
–概率密度函数难以积分 ?无累积分布函数的解析表达式；
–累积分布函数的反函数的解析表达式不存在；
–概率密度函数存在尖峰（ spiky);
8.7 混合抽样法
(mixed method)
)()()( xgxfxp ?
假定：概率密度函数可写成下面的因子化形式
其中,
– f(x)包含了 p(x)的峰值部分且可用直接抽样法进行抽样
– g(x)是一个相对变化平缓的函数，包含了 p(x)函数的大
部分的数学复杂性；
8.7 混合抽样法
(mixed method)
)(1)(;)( xfCxfCdxxf
f
b
a f ???
],[,1)();(1)( baxxgxgMxg
g
???
抽样方法,
1,将 f(x)归一化,
2,令
Mg为 g(x)在区间 [a,b]上的极大值
3,利用直接抽样法由 抽取 x; )(xf
4, 抽取 [0,1]区间均匀分布的随机数 r,如果,则接
受 x,否则，返回到 3重新抽样
rxg ?)(
8.7 混合抽样法
(mixed method)
)()()( xgxfxp i
i i?
?
?? ??? i iii
i
iii xpCC xgxfCxp )()()()(
)(xpi?
)(xpi?
设概率密度函数可写成如下的形式,
推广的形式,
?复合抽样法 +混合抽样法 ==〉 乘加抽样法
抽样方法,
1,采用复合抽样法，先确定 p(x)的随机数应由哪一个 抽取；
2,在按 抽样时，采用上面的混合方法
1)(
)()()(;)()(
??
???
??
?
i
i
b
a
i
ii
i
b
a iii
Cdxxp
C
xgxfxpdxxgxfC
8.7 混合抽样法
(mixed method)
? ?1,
2;
)c o s1(;
1
s i n
1
1
),(
0
0
2
2
0m i n11
0
2
2
01
1
0
2
2
2
0
22
00
10
??
?
?
?
?
??
?
?
?
?
?
???
??
??
?
?
?
?
?
?
?
??
?
?
?
?
?
???
Ecm
cm
EE
Ecm
cm
EE
E
E
E
cmrnX
EE
e
e
e
e
e
例,compton散射
微分截面,
如何抽取散射光子的能量？ ==〉 乘加抽样法
8.7 混合抽样法
(mixed method)
? ?
1)(0:]1,[s in
1
1)(
)(
1
)()()(;
2
)1(
ln
1
)()(
1
)(;ln
1
ln
)()()()()(
1
s in
1
1
)(
0
2
2
2
0
2
2
22
2
2
2
0
2
01
11
1
10
0
1
22112
2
0
0
????
?
?
?
?
?
?
?
??
??????
?
?
???????
?
?
?
?
?
?
?
?
?
???
?
?
?
?
?
?
?
?
?
?
?
?
?
?
???
?
?
????
?
?
?
??
?
???
?
?
?
?
?
?
?
?
???
??
??
?
?
???????
?
??
?
?
?
?
?
?
?
gg
dfFf
dfFf
gffgf
8.7 混合抽样法
(mixed method)
抽样方法,r1,r2,r3是三个在 [0,1]区间上均匀分布的随机数
1,确定由哪一个 f函数来抽取 ?：如果 r1<?1/(?1+?2)，选择 f1(?),
否则选 f2(?)；
2,根据 f1或 f2抽取 ?：直接抽样法
22
2
0212
0211
)(:
)(,21
rrFf
erFf r
????
??? ?
????
???
3,计算 sin2?,
?
??? ?????? 1)c o s1()2(s i n
0
22
E
cmttt e
4.计算 g(?),如果 接受 ?,否则返回到第一步 3)( rg ??
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.8 近似抽样法（列表法 )
8.8 近似抽样法（列表法 )
近似抽样法,
用近似的分布函数取代欲抽取的概率密度函数，一般是采用
列表的形式将连续型的概率分布变成分离型的概率分布。
使用的场合,
– 概率密度函数的形式非常复杂，在模拟过程中进行计算
需花费相当多的 CPU时间 ;
– 概率密度函数无解析形式，只能用数值或曲线的形式表
示。
8.8 近似抽样法（列表法 )
bxxxxa n ?????? ?210
nixff ii ?,2,1),( ??
基本方法,
设概率密度函数,f(x),x?[a,b]
1,将区间 [a,b]分成 n个子区间，分点为
2,分点对应的函数值为,
3,对于每一个小区间，利用一简单的函数 fa(x)来近似地
表示原概率分布函数，并使 fa(x)在该区间内的积分与
f(x)在该区间内的积分相等，即俩者的概率相等。
4,利用 f(x)计算每一个子区间的概率值 并计
算 f(x)的累积分布函数在 xi处的值,并
将 F(xi)和 xi存入数据表中 ;
? ?? iixxi dxxfp 1 )(
?? ?? ij ji nipxF 1,,2,1,)( ?
8.8 近似抽样法（列表法 )
)()( 1 ii xFrxF ???
5,抽样方法,
– 随机选择子区间,选取 [0,1]区间内均匀分布的随机数
r,找出满足 的子区间 [xi-1,xi];
– 根据 fa(x)的特点,确定欲抽取的随机数 ?；
几点说明,
1,分点的选取,
f0,f1,...,fn应能充分反映 f(x)的变化状况,即,在 f(x) 变化迅速
的区域分点密一点,变化缓慢的区域分点稀一点,
2,fa(x)的选取,
? 阶梯近似,
? 线性近似
? 二次曲线近似
? 样条函数近似
8.8 近似抽样法（列表法 )
nipxF
dxxfdxxfdxxfp
pdx
xx
p
dxxf
ni
xx
p
xx
dxxf
xf
i
j
ii
b
a
n
i
n
i
x
x
x
x
i
x
x
i
ii
i
x
x
a
ii
i
ii
x
x
a
ni
i
i
i
i
i
i
i
,,2,1)(
1)()()(
)(
,,2,1,
)(
)(
1
1 1
1
11
01
11
1
?
?
??
????
?
?
?
?
?
?
?
?
?
?? ? ??
??
?
?
? ?
?
??
?
??
?
阶梯近似
将 fa(x) 取为阶梯函数，在每一个子区间中 fa(x)都是均匀分布
8.8 近似抽样法（列表法 )
)()( 1 jij xFrxF ???
)()(
)()(
1
1
11
?
?
?? ?
????
jj
ji
jjji xFxF
xFrxxx?
抽样方法,
1,选取 [0,1]区间均匀分布的随机数 ri;
2,找出满足下式的分点 xj-1和 xj,
3,欲抽取的随机数为,
8.8 近似抽样法（列表法 )
nixxxffxx xxfCxf iiii
ii
i
ia,,2,1],,[)()( 11
1
1
1 ?????
?
??
? ?
?
???
??
?
?
?
nidxxfdxxfi
i
i
i
x
x
x
xa,,2,1,)()(1 1 ???? ?? ?
线性近似
利用一系列的折线来近似原分布 f(x),即将 fa(x)取为
其中 C为归一化因子，使得每一子区间内原分布和近似分布的积
分概率相等
8.8 近似抽样法（列表法 )
抽取 r1,r2
找出满足下式要求的 xj-1,xj
? ? 2)()( 112111 ???? ????? jjjjjja xxfrffxFr
211 )( rxxx jjji ?? ???? )1)((
211 rxxx jjji ???? ???
? >
抽样方法,
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.9 多维分布的抽样
8.9 多维分布的抽样
? 将一维分布的抽样方法推广到多维分布
设 n维概率密度为, 的取值域为 n维长方体,)(xf? ],[ ba ??x?
},,,{},,,,{ 2121 nn bbbbaaaa ???? ??
)(xf? 的极大值为 M,则 的随机数 的舍选法 )(xf? },,,{ 21 n???? ?? ?
抽取 r′,ri ?U[0,1],i=1,2,…,n
?i = ai + (bi-ai)ri,i=1,2,…,n
? >
Mfr /)(????
??? ?x
舍选法用于多维抽样
8.9 多维分布的抽样
混合抽样法用于多维分布的抽样
设多维随机变量概率密度函数 可表示为 )(xf?
)()()( xgxHxf ??? ?
式中 为任意多维概率密度函数，其随机数为 )(xg ?
},,,{ 21 n???? ?? ?
H( )?0的极大值为 M，则 f( )的随机数 可
由下列方式抽取,
x? x? },,,{ 21 n???? ?? ?
抽取 r和 ??
MHr /)(???
?? ???
?
>
8.9 多维分布的抽样
?条件密度法,
利用条件概率密度将多维模拟转化为一维模拟问题,
任意 n维概率密度函数 fn(x1,x2,…,xn)可表示为
fn(x1,x2,…,xn) = f(xn| x1,x2,…,xn-1)· fn-1(x1,x2,…,xn-1)
条件概率 x
1,x2,…,xn-1联合概
率密度
???? nnnnn dxxxxfxxxf ),,,(),,,( 211211 ??
fn(x1,x2,…,xn) = f(xn| x1,x2,…,xn-1)·
f(xn-1| x1,x2,…,xn-2)·
…
f(x2| x1)·
f(x1) ?
n个一维概率密度的
乘积
),.,,,,(
),.,,,,(),.,,,,|(
1211
21
121
??
? ?
nn
nn
nn xxxf
xxxfxxxxf
8.9 多维分布的抽样
?条件密度法,
fn(x1,x2,…,xn) = f(xn| x1,x2,…,xn-1)·
f(xn-1| x1,x2,…,xn-2)·
…
f(x2| x1)·
f(x1)
抽样方法,
由 f(x1)抽样得到 x1的随机数 ?1，将 ?1代入 f(x2| x1)中得到
x2的一维分布 f(x2| ?1)；
由 f(x2| ?1)抽样得到 x2的随机数 ?2,…
第二部分,Monte Carlo模拟
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
第八章 从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
? Monte Carlo算法的一个重要组成部分,
描述所要模拟的物理系统的一些概率密度函数（ PDF)
?描述整个系统在空间、能量、时间或多维相空间中的
发展和演化；
? Monte Carlo模拟的主要任务,
? 通过对这些概率密度函数的随机抽样来模拟物理系统的状
态 ;
? 为描述系统的演化所必需的一些附加运算,
? 物理过程的描述 ?从描述物理系统的 pdf出发，随机抽取系统的
可能状态。
第八章 从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
? 本章介绍从一个任意的 pdf获取样本的抽样方法。
1,直接抽样法（反函数法）
2,变换抽样法 ?直接抽样法的一般形式
3,舍选抽样法
4,复合分布的抽样方法
5,混合抽样法
6,近似抽样法（列表法 )
7,多维分布的抽样
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.1 等价的连续概率密度函数
8.1 等价的连续概率密度函数
? 随机变量：连续型、分离型
? 概率密度函数：连续分布、分离分布
? 利用 ?函数，可将分离型的 pdf用连续型的 pdf 描述
? 用同样的方式来讨论分离型和连续型随机变量的抽样方法
8.1 等价的连续概率密度函数
?
?
??
N
i
ii xxpxf
1
)()( ?
?? ?
? ??
??
?
??
?
??
?
?
??
???
?
?
?
?
?
?
???
N
i
ii
N
i
ii
N
i
ii
pxdxxxxp
dxxxpxdxxxfXE
11
1
)(
)()()(
?
?
?? ?
? ??
??
?
??
?
??
?
?
??
?????
?
?
?
?
?
?
?????
N
i
ii
N
i
ii
N
i
ii
pXExdxxxpXEx
dxxxpXExdxxfXExXV
1
2
1
2
1
22
)]([)())((
)()]([)())(()(
?
?
? 已知分离型 pdf,{Pi} ? 分离型随机变量 X的取值为 xi的概率
? 定义一个等价的连续型 pdf,
?
?
?
??
?
??
??
??
)()()(
1)(
ii
i
xfdxxxxf
dxxx
?
?
? 利用与连续型随机变量相同的方式计算分离型随机变量的期望值
和方差,
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.2 pdf的变换
8.2 pdf的变换
x,连续型的随机变量,PDF,f(x)
y = y(x),x的函数,也是随机变量, 求 y(x)的概率密度函数 g(x)
1、若随机变量 x和 y是一一对应的,
2,若随机变量 x和 y不是一一对应的,
[x,x+dx]?[y,y+dy]
X的取值在 [x,x+dx]的概率 ==Y的取值
在 [y,y+dy]的概率,
dy
dxxfyg )()( ?
取绝对值是为了保证 g(y)是非负的
f(x)dx=g(y)dy ?
即有 n个区间 [x,x+dx]?[y,y+dy]
?? dydxxfyg )()(
需要对这 n个区间求和
8.2 pdf的变换
3、推广到 n个随机变量的情况,
Jacobian行列式 ?
nixyy
yyyyxxxx
ii
nn
,,2,1),(
),,(),,( 2121
??
????
??
???
n
nnn
n
y
x
y
x
y
x
y
x
y
x
y
x
n
n
yyy
xxx
J
Jxfyg
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
21
1
2
1
1
1
,
,
),,(
),,(
)()(
21
21
4、特例：如果 y(x)是 x的累积分布函数 (cdf)
)(
1
)(
11
)()()(
xfxF
dx
dydy
dx
xdxfxFxyy
x
?
?
??
????? ?
??
10,1)( ???? yyg
即,y在 [0,1]区间上均匀分布 ? 不管 f(x)取何种形式,累积分布
函数总是在 [0,1]区间上均匀分布
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
)(1 yFx ??
原理 （ 注意,pdf f(x)必须是归一化的 ）,
? 设 y=F(x)为随机变量 x的累积分布函数 ? x和 y是一一对应的
? 先随机抽取 y，然后通过求 F(x)的反函数 F-1(y)得到随机变量 x
的值
? 随机变量 y在区间 [0,1]上均匀分布 ? 利用 [0,1]区间上均匀分布
随机数产生器抽取
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
方法,U[0,1],[0,1]区间上均匀分布的随机数
1,从 U[0,1]抽取随机数 ?；
2,令 F(x)= ?;
3,解方程得 x,)(1 ??? Fx
注：需要知道累积分布函数的解析表达式，且累积分布函数的
反函数存在
分离型随机变量的抽样
直接抽样法适应于分离型的随机变量
?
? ?
?
?
??
?
??
?
??
?
?
?
?
?
???
????
K
i
i
x
N
i
ii
x
p
xdxxp
xdxfxFy
1
1
)(
)()(
?
?
?
?? N
i
ii xxpxf
1
)()( ?
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
方法,
1,计算 yk = yk-1 + pk,k = 2,3,…,N,y 1 = p1
2,从 U[0,1]抽取随机数 ?；
3,求满足 yk-1 < ? < yk 的 K值 ;
4,随机变量的第 k个取值即为欲抽取的值。
0 ? 1
p a pb pc pd
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
p3=0.2 b3+c3
p2=0.3 b2+c2
p1=0.5 b1+c1 a?
13.05.0
3.05.05.0
5.00
???
???
??
?
?
?
33
22
11
cb
cb
cb
?
?
?
例 1、粒子衰变末态的随机抽样
设粒子 a有三种衰变方式，其分支比如下
随机选取每次衰变的衰变方式（衰变道） ?直接抽样法
? ?U[0,1]
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
??
?
?
??
?
?
?
?
?
???
?
?
??
?
?
????
?
?
??
?
?
? ?
rn
n
rnr
n
r
n
nrpp
r
n
pnrB rnr
)!(!
!
,2,1,0,)1(),;( ?例 2、二项式分布的抽样
方法 1,利用上面介绍的直接抽样法，需计算累积分布函数，
当 n很大时，求和计算困难 ;
方法 2,利用二项式分布的定义
1,产生 n个 ? i?U[0,1];
2,统计满足条件 ? i <p（表示成功）的 ? i的数目 r，则 r
表示在 n次实验中成功的次数 ?r即为二项式分布的
抽样值
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
?0,1,2,=r,e r!1 =μ) P ( r ; -r ??
????? nppn,0,
例 3、泊松分布的抽样
方法 1,利用直接抽样法，但计算累积分布函数时非常复杂
方法 2,利用泊松分布的定义：二项式分布的极限形式
1,选取足够大的 n，使 p=?/n相当小，例如,p=0.1
2,产生 n个 ? i?U[0,1];
3,统计满足条件 ? i <p（表示成功）的 ? i的数目 r，则 r
表示在 n次实验中成功的次数 ?r即为泊松分布的抽
样值的近似值,n越大，近似程度越好
8.3 直接抽样法（反函数法）
（ Sampling via Inversion of the cdf）
??
???
??
????
bxax
bxaabxf
,,0
,1)(
)0()( ?? ? xexf x??
例 4、连续型随机变量的直接抽样
1,求区间 [a,b]上均匀分布的随机数 x,
? ?????? xa ab axxdxf )(?
aabx ??? )(?
?产生 ? ?U[0,1];
?
?
2,指数分布
xx x exde ???? ??? ???? ? 1
0
???? ln1)1l n (1 ?????x
?产生 ? ?U[0,1];
?
? ?和（ 1-?）都是 U[0,1]
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.4 变换抽样法 ?直接抽样法的一般形式
8.4 变换抽样法 ?直接抽样法的一般形式
思路,
dy
dxxfyg )()( ?
? 随机变量 y,pdf,g(y) ? 不易进行抽样
? 随机变量 x,pdf,f(x) ? 容易进行抽样
? 如果能够找到 x和 y之间的一个一一对应的变换关系，
y=y(x),使得 g(y)和 f(x)满足关系
则可先由 f(x)分布抽取 x的值 ?，再由变换关系得到 y的值
???满足分布 g(y)
?变换抽样法
8.4 变换抽样法 ?直接抽样法的一般形式
变换抽样法,
?找出 y与 x间的变换关系,y=y(x),f(x)与 g(y)满足关系,
?由 f(x)分布抽取 x的值 ?；
?随机变量 y的取值,?=y(?)
dy
dxxfyg )()( ?
直接抽样法是变换抽样法的一个特殊形式
?X满足 U[0,1],f(x)=1;
?y与 x间的变换关系,y=G-1(x) ?y的累积分布函数的反函
数
? ?? ???? y ydygyGx )()(
8.4 变换抽样法 ?直接抽样法的一般形式
nixyy
yyyyxxxx
ii
nn
,,2,1),(
),,(),,( 2121
??
????
??
???
n
nnn
n
y
x
y
x
y
x
y
x
y
x
y
x
n
n
yyy
xxx
J
Jxfyg
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
21
1
2
1
1
1
,
,
),,(
),,(
)()(
21
21
推广到 n维随机向量的情况,
)(xf ? x? },,,{ 21 n???? ?? ?
)(?? ?ii y?
–由联合概率密度函数 抽取随机向量 的值
–yi的值,
8.4 变换抽样法 ?直接抽样法的一般形式
2
2
2
)(
2
1),( ??
????
??
?
x
eN
???? xy
?? ?? yx
例：高斯分布的抽样方法
进行变量变换,
)1,0(),( NN ???? ?标准正态分布
由 N(0,1)分布抽样得到 y,
如何抽取服从 N(0,1)分布的随机变量？
8.4 变换抽样法 ?直接抽样法的一般形式
???
?
???
??
?
??
?
? ?? ???
???
n
i
i
n
i
i
n
i
ixy
1
2
11
??
1,利用中心极限定理
设 x1,x2,…,x n是一组 n个独立的随机变量,xi的平均值和方差
分别为 μi和 бi，则当 n→∞ 时，变量
服从标准正态分布 N(0,1)
设 x 是在 [0,1]之间均匀分布的随机数,对 n个 x的取值 xi
121221 )()( ???? xVxE ii ?? ?
????????????? ?? ?
? 1221
nnxy n
i i
在 n→∞ 时,y服从 N(0,1)，在实际应用时，可取 n=12
抽样方法,
1、产生 12个 U[0,1]的随机数 ?i
612
1
?? ?
?i
iy ?
2,
8.4 变换抽样法 ?直接抽样法的一般形式
22
2121
2221
2
1
2
1)()(),( yy eeygygyyg ?? ????
??
1),(:],1,0(,
)2s in ()ln2(
)2c o s ()ln2(
2121
212
211
??
??
??
xxfpdfUxx
xxy
xxy
?
?
1
21
2
2
1
2
1
2
2
2
1
y
ytgx
ex
yy
?
??
?
?
?
2,利用变换抽样法
y1和 y2是两个相互独立的、服从标准正态分布的随机数
变换,反变换,
抽样方法,
)2s in ()ln2(
)2c o s ()ln2(
2122
2111
????
????
???
???
y
y
1）产生一对 [0,1]区间均匀分布的随机数 ?1和 ?2；
2）
8.4 变换抽样法 ?直接抽样法的一般形式
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
2
2
2
2
2
1
2
2
2
2
1
2
2
2
1
2
2
2
1
2
2
1
2
2
1
1
1
2
1
2
1
2
1
1
1
2
1
1
2
1
,
,
),(
),(
2
1
2
1
2
2
2
1
2
2
1
2
2
1
21
21
yyyy
eee
y
y
y
ey
y
y
y
y
ey
yy
xx
J
yyyy
y
x
y
x
y
x
y
x
???
??
证明用此方法抽取的 y1,y2满足上面的联合概率分布
雅可比行列式,
JJxxfyyg ??? ),(),( 2121
8.4 变换抽样法 ?直接抽样法的一般形式
222
2222
2
1
2
1
2
1)()(),( yxyx eeeyfxfyxf ???? ?????
???
x
y
yxr
?
???
?
?
t a n
22
222
??
??
s i n2
c o s2
?
?
y
x
?
?
???
?
???
?
??
?? ???
?
?
?
?
? eyxfyxfyxf
yx
g
2
1
),(),(
c o s2,s i n
2
1
s i n2,c o s
2
1
),(
),(
),(
),(
3,采用极坐标变换
x和 y是两个相互独立的、服从标准正态分布的随机数
变换,反变换,
则变量 ?和 ?的联合概率密度函数为
即 ?在 [0,2?] 区间上均匀分布,?服从 ? = 1的 指数分布
8.4 变换抽样法 ?直接抽样法的一般形式
抽样方法,
–产生 ?, ?= 2? r,r?U[0,1]
–产生 ?, ? = -ln r,r?U[0,1]
– x= (2?)1/2 cos?, y=(2?)1/2 sin?
8.4 变换抽样法 ?直接抽样法的一般形式
void test() {
c1 = new TCanvas("c1","Histogram Drawing Options",200,10,700,900);
c1->Divide(1,2);
TH1F * h1 = new TH1F("h1","h1",100,-5.0,5.0);
TH1F * h2 = new TH1F("h2","h2",100,-5.0,5.0);
SetSeed(9,11)
Const Pi = 3.1415926;
for(int i=0; i < 50000; i++) {
double x = randac();
double y = randac();
float y1 = sqrt(-2.0*log(x))*cos(2.0*Pi*y);
float y2 = sqrt(-2.0*log(x))*sin(2.0*Pi*y);
h1->Fill(y1); h2->Fill(y2);
}
c1->cd(1);h1->Draw(); c1->cd(2);h2->Draw();
}
8.4 变换抽样法 ?直接抽样法的一般形式
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.5 舍选抽样法
（ acceptance-rejection sampling)
8.5 舍选抽样法
（ acceptance-rejection sampling)
直接抽样法的困难,
? 许多随机变量的累积分布函数无法用解析函数给出 ;
? 有些随机变量的累积分布函数的反函数不存在或难以求出；
? 即使反函数存在，但计算困难
舍选抽样法,
抽取随机变量 x的一个随机序列 xi,i=1,2,…,按一定的舍选规则
从中选出一个子序列，使其满足给定的概率分布,
假定,
? 随机变量 x的值域为 [a,b];
? x的概率密度函数,f(x)=P*(x)/Z,（其中 Z为归一化因子） ?
难以直接抽样 ;
? Q(x)=Q*(x)/ZQ 是另外一个较为简单的函数 (ZQ为归一化因子）
?可用简单的方法进行抽样 ;
? 在 x的整个取值范围内,cQ*(x) > P*(x),其中 c为一常数
8.5 舍选抽样法
（ acceptance-rejection sampling)
抽样方法,
1,产生两个随机数
? 从 Q(x)分布抽取 ?, ??[a,b];
? 由 [0,cQ*(?)] 区间上的均匀分布抽取 u,u= cQ*(?) ?,
??U[0,1],
2,接收或舍弃取样值 x,
? 如果 u > P*(x)，舍弃，返回到 1,重复上述过程 ;
? 否则，接受 ;
8.5 舍选抽样法
（ acceptance-rejection sampling)
几何解释,
? 在二维图上，随机选取位于曲线 cQ*(x)下的点 [x,u]；
? 选取位于曲线 P*(x)下的那些点，则这些点将服从概率密度为 f(x)的分布
8.5 舍选抽样法
（ acceptance-rejection sampling)
常数 c的选取
? 常数 c应尽可能地小,因为抽样效率与 c成反比 ;
? c=max{P*(x)/Q*(x)},x ?[a,b]
特例,
如果取 Q(x)=1,x?[a,b],即均匀分布,则
? X的抽样,x=(b-a)r+a,r ?U[0,1]
? c可取为 f(x)在 [a,b]区间上的极大值
a b x
f(x)
c
抽取 r1,r2 ?U[0,1]
? = a + (b-a)r1
r2 ? f(?)/c
x = ?
? >
8.5 舍选抽样法
（ acceptance-rejection sampling)
2/2
2
1)(* xexp ??
?
21
11)(*
xxQ ?? ?
? ? 52.12)1(2m a x)(*/)(*m a x 2/2 2 ??
??
?
??
? ??? ?
eexxQxPc
x ??
? ?)a r c t a n (a r c t a n
a r c t a n2
1
)()(
1
1
a r c t a n2
1
)(*/)(*)(
2
ax
a
xdxQxF
xa
dxxQxQxQ
x
aQ
a
a
??????
?
??
?
?
?
?
例 1:标准正态分布的抽样,x?[-a,a]
无法用直接抽样法,累积分布函数无解析
表达式
Breit-wigner or Cauchy分布
8.5 舍选抽样法
（ acceptance-rejection sampling)
? ?)a r c t a n (a r c t a n2t a n
]1,0[
aax
U
???
?
?
?
Qu
U
x
xcQQ
?
?
?
?
?
?
??
]1,0[
1
152.1
)(*
2
由 Q(x)抽取 x ? 直接抽样法
抽取 u
计算 P*(x),如果 u<= P*(x),接受 x
8.5 舍选抽样法
（ acceptance-rejection sampling)
float gaussian_reject(double a)
{
const float c = 1.52;
while(true) {
float eta = randac();
float x = tan(eta * 2.0 * atan(a)+atan(-a));
float q = c * 1/3.1415926*1.0/(1+x*x);
float ksi = randac();
float u = ksi*q;
float p = 1/sqrt(2*3.1415926)*exp(-x*x/2.0);
if(u <= p) break;
}
return x;
}
8.5 舍选抽样法
（ acceptance-rejection sampling)
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,-5.0,5.0);
for(int i=0; i < 5000; i++) {
double x = gaussian_reject(5.0);
h1->Fill(x);
}
c1->cd(2);h1->Draw();
}
8.5 舍选抽样法
（ acceptance-rejection sampling)
2c o s2s i n2s i n2s i n2c o sc o s 22
?????? ????? SC
A
B ?/2
2222
22
2222
2
s i n,c o s
2
s i n,
2
c o s
BA
AB
S
BA
BA
C
BA
B
BA
A
?
??
?
?
???
?
?
?
?
??
??
例 2：利用舍选法产生随机数 C=cos?,S=sin?,其中 ?为 [0,2?]区间内
均匀分布的随机数
方法 1：先产生 [0,2?]间均匀分布的随机数,? = 2? r,
r?U[0,1],然后直接计算 C和 S ? 因需要计算三角函数，
故此方法运算速度慢
方法 2：利用舍选法可避免三角函数运算
令 A和 B为单位圆内直角三角形的两个边，则有
8.5 舍选抽样法
（ acceptance-rejection sampling)
因此，只要产生单位圆内的随机坐标 A和 B,就可用代数运算
求出 C和 S,算法为
1,产生两个 [0,1]区间上均匀分布的随机数 u1和 u2;
2,令 v1=2u1-1,v2=u2,则 v1?U[-1,1],v2?U[0,1];
3,计算 r2=v12+v22,如果 r2 > 1,转到 1，重新产生；
4,令 A=v1,B=v2,计算 C和 S
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.6 复合分布的抽样方法
（ composition method)
8.6 复合分布的抽样方法
（ composition method)
??
?
???
?
1)(;1;0
)()(
dxxfpp
xfpxf
k
k
kk
k
kk
??
?
?
?
??? k
i
i
k
i
i prpUr
1
1
1
];1,0[
1961年由 Marsaglia提出的方法
设随机变量 X的概率密度函数 f(x)可写成一些 PDF的线性叠加,
抽样方法,
1,利用离散型的随机变量的抽样方法抽取序号 k;
2,由 fk(x)抽取 x
8.6 复合分布的抽样方法
（ composition method)
)0()()0()(5.0
0
01
)(
0
01
)(
)(5.0)(5.0)(
5.0)(:
2121
)0[)0,(
),0[)0,(
???????
?
?
? ???
?
?
?
? ????
?
??
?
?
???
?
?
??
?
xexfxexfpp
x
xI
x
xI
xIexIexf
exfP D F
xx
xx
x
其他其他
，
例：用复合法产生双指数分布随机数
– 产生两个 [0,1]区间均匀分布的随机数 r1和 r2;
– 如果 r1 ? 0.5,按 f1(x)抽样；
– 如果 r1 > 0.5,按 f2(x)抽样；
– f1(x)和 f2(x)都可用直接抽样法
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.7 混合抽样法
（ mixed method)
8.7 混合抽样法
(mixed method)
?
混合法,
直接抽样法
舍选法 ?亦称乘抽样法
适用的抽样场合,
直接抽样法不可用、舍选抽样法效率低
–概率密度函数难以积分 ?无累积分布函数的解析表达式；
–累积分布函数的反函数的解析表达式不存在；
–概率密度函数存在尖峰（ spiky);
8.7 混合抽样法
(mixed method)
)()()( xgxfxp ?
假定：概率密度函数可写成下面的因子化形式
其中,
– f(x)包含了 p(x)的峰值部分且可用直接抽样法进行抽样
– g(x)是一个相对变化平缓的函数，包含了 p(x)函数的大
部分的数学复杂性；
8.7 混合抽样法
(mixed method)
)(1)(;)( xfCxfCdxxf
f
b
a f ???
],[,1)();(1)( baxxgxgMxg
g
???
抽样方法,
1,将 f(x)归一化,
2,令
Mg为 g(x)在区间 [a,b]上的极大值
3,利用直接抽样法由 抽取 x; )(xf
4, 抽取 [0,1]区间均匀分布的随机数 r,如果,则接
受 x,否则，返回到 3重新抽样
rxg ?)(
8.7 混合抽样法
(mixed method)
)()()( xgxfxp i
i i?
?
?? ??? i iii
i
iii xpCC xgxfCxp )()()()(
)(xpi?
)(xpi?
设概率密度函数可写成如下的形式,
推广的形式,
?复合抽样法 +混合抽样法 ==〉 乘加抽样法
抽样方法,
1,采用复合抽样法，先确定 p(x)的随机数应由哪一个 抽取；
2,在按 抽样时，采用上面的混合方法
1)(
)()()(;)()(
??
???
??
?
i
i
b
a
i
ii
i
b
a iii
Cdxxp
C
xgxfxpdxxgxfC
8.7 混合抽样法
(mixed method)
? ?1,
2;
)c o s1(;
1
s i n
1
1
),(
0
0
2
2
0m i n11
0
2
2
01
1
0
2
2
2
0
22
00
10
??
?
?
?
?
??
?
?
?
?
?
???
??
??
?
?
?
?
?
?
?
??
?
?
?
?
?
???
Ecm
cm
EE
Ecm
cm
EE
E
E
E
cmrnX
EE
e
e
e
e
e
例,compton散射
微分截面,
如何抽取散射光子的能量？ ==〉 乘加抽样法
8.7 混合抽样法
(mixed method)
? ?
1)(0:]1,[s in
1
1)(
)(
1
)()()(;
2
)1(
ln
1
)()(
1
)(;ln
1
ln
)()()()()(
1
s in
1
1
)(
0
2
2
2
0
2
2
22
2
2
2
0
2
01
11
1
10
0
1
22112
2
0
0
????
?
?
?
?
?
?
?
??
??????
?
?
???????
?
?
?
?
?
?
?
?
?
???
?
?
?
?
?
?
?
?
?
?
?
?
?
?
???
?
?
????
?
?
?
??
?
???
?
?
?
?
?
?
?
?
???
??
??
?
?
???????
?
??
?
?
?
?
?
?
?
gg
dfFf
dfFf
gffgf
8.7 混合抽样法
(mixed method)
抽样方法,r1,r2,r3是三个在 [0,1]区间上均匀分布的随机数
1,确定由哪一个 f函数来抽取 ?：如果 r1<?1/(?1+?2)，选择 f1(?),
否则选 f2(?)；
2,根据 f1或 f2抽取 ?：直接抽样法
22
2
0212
0211
)(:
)(,21
rrFf
erFf r
????
??? ?
????
???
3,计算 sin2?,
?
??? ?????? 1)c o s1()2(s i n
0
22
E
cmttt e
4.计算 g(?),如果 接受 ?,否则返回到第一步 3)( rg ??
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.8 近似抽样法（列表法 )
8.8 近似抽样法（列表法 )
近似抽样法,
用近似的分布函数取代欲抽取的概率密度函数，一般是采用
列表的形式将连续型的概率分布变成分离型的概率分布。
使用的场合,
– 概率密度函数的形式非常复杂，在模拟过程中进行计算
需花费相当多的 CPU时间 ;
– 概率密度函数无解析形式，只能用数值或曲线的形式表
示。
8.8 近似抽样法（列表法 )
bxxxxa n ?????? ?210
nixff ii ?,2,1),( ??
基本方法,
设概率密度函数,f(x),x?[a,b]
1,将区间 [a,b]分成 n个子区间，分点为
2,分点对应的函数值为,
3,对于每一个小区间，利用一简单的函数 fa(x)来近似地
表示原概率分布函数，并使 fa(x)在该区间内的积分与
f(x)在该区间内的积分相等，即俩者的概率相等。
4,利用 f(x)计算每一个子区间的概率值 并计
算 f(x)的累积分布函数在 xi处的值,并
将 F(xi)和 xi存入数据表中 ;
? ?? iixxi dxxfp 1 )(
?? ?? ij ji nipxF 1,,2,1,)( ?
8.8 近似抽样法（列表法 )
)()( 1 ii xFrxF ???
5,抽样方法,
– 随机选择子区间,选取 [0,1]区间内均匀分布的随机数
r,找出满足 的子区间 [xi-1,xi];
– 根据 fa(x)的特点,确定欲抽取的随机数 ?；
几点说明,
1,分点的选取,
f0,f1,...,fn应能充分反映 f(x)的变化状况,即,在 f(x) 变化迅速
的区域分点密一点,变化缓慢的区域分点稀一点,
2,fa(x)的选取,
? 阶梯近似,
? 线性近似
? 二次曲线近似
? 样条函数近似
8.8 近似抽样法（列表法 )
nipxF
dxxfdxxfdxxfp
pdx
xx
p
dxxf
ni
xx
p
xx
dxxf
xf
i
j
ii
b
a
n
i
n
i
x
x
x
x
i
x
x
i
ii
i
x
x
a
ii
i
ii
x
x
a
ni
i
i
i
i
i
i
i
,,2,1)(
1)()()(
)(
,,2,1,
)(
)(
1
1 1
1
11
01
11
1
?
?
??
????
?
?
?
?
?
?
?
?
?
?? ? ??
??
?
?
? ?
?
??
?
??
?
阶梯近似
将 fa(x) 取为阶梯函数，在每一个子区间中 fa(x)都是均匀分布
8.8 近似抽样法（列表法 )
)()( 1 jij xFrxF ???
)()(
)()(
1
1
11
?
?
?? ?
????
jj
ji
jjji xFxF
xFrxxx?
抽样方法,
1,选取 [0,1]区间均匀分布的随机数 ri;
2,找出满足下式的分点 xj-1和 xj,
3,欲抽取的随机数为,
8.8 近似抽样法（列表法 )
nixxxffxx xxfCxf iiii
ii
i
ia,,2,1],,[)()( 11
1
1
1 ?????
?
??
? ?
?
???
??
?
?
?
nidxxfdxxfi
i
i
i
x
x
x
xa,,2,1,)()(1 1 ???? ?? ?
线性近似
利用一系列的折线来近似原分布 f(x),即将 fa(x)取为
其中 C为归一化因子，使得每一子区间内原分布和近似分布的积
分概率相等
8.8 近似抽样法（列表法 )
抽取 r1,r2
找出满足下式要求的 xj-1,xj
? ? 2)()( 112111 ???? ????? jjjjjja xxfrffxFr
211 )( rxxx jjji ?? ???? )1)((
211 rxxx jjji ???? ???
? >
抽样方法,
第八章
从概率分布函数的抽样
(Sampling from Probability Distribution Functions)
8.9 多维分布的抽样
8.9 多维分布的抽样
? 将一维分布的抽样方法推广到多维分布
设 n维概率密度为, 的取值域为 n维长方体,)(xf? ],[ ba ??x?
},,,{},,,,{ 2121 nn bbbbaaaa ???? ??
)(xf? 的极大值为 M,则 的随机数 的舍选法 )(xf? },,,{ 21 n???? ?? ?
抽取 r′,ri ?U[0,1],i=1,2,…,n
?i = ai + (bi-ai)ri,i=1,2,…,n
? >
Mfr /)(????
??? ?x
舍选法用于多维抽样
8.9 多维分布的抽样
混合抽样法用于多维分布的抽样
设多维随机变量概率密度函数 可表示为 )(xf?
)()()( xgxHxf ??? ?
式中 为任意多维概率密度函数，其随机数为 )(xg ?
},,,{ 21 n???? ?? ?
H( )?0的极大值为 M，则 f( )的随机数 可
由下列方式抽取,
x? x? },,,{ 21 n???? ?? ?
抽取 r和 ??
MHr /)(???
?? ???
?
>
8.9 多维分布的抽样
?条件密度法,
利用条件概率密度将多维模拟转化为一维模拟问题,
任意 n维概率密度函数 fn(x1,x2,…,xn)可表示为
fn(x1,x2,…,xn) = f(xn| x1,x2,…,xn-1)· fn-1(x1,x2,…,xn-1)
条件概率 x
1,x2,…,xn-1联合概
率密度
???? nnnnn dxxxxfxxxf ),,,(),,,( 211211 ??
fn(x1,x2,…,xn) = f(xn| x1,x2,…,xn-1)·
f(xn-1| x1,x2,…,xn-2)·
…
f(x2| x1)·
f(x1) ?
n个一维概率密度的
乘积
),.,,,,(
),.,,,,(),.,,,,|(
1211
21
121
??
? ?
nn
nn
nn xxxf
xxxfxxxxf
8.9 多维分布的抽样
?条件密度法,
fn(x1,x2,…,xn) = f(xn| x1,x2,…,xn-1)·
f(xn-1| x1,x2,…,xn-2)·
…
f(x2| x1)·
f(x1)
抽样方法,
由 f(x1)抽样得到 x1的随机数 ?1，将 ?1代入 f(x2| x1)中得到
x2的一维分布 f(x2| ?1)；
由 f(x2| ?1)抽样得到 x2的随机数 ?2,…
