假设现在手上只有一个均匀分布的随机数发生器,如何利用这个符合均匀分布的随机数发生器来生成符合任意想要的分布的随机数呢,这里介绍了两种方法,第一种方法方便快捷,但需要目标分布的概率分布累积函数(CDF)可逆;第二种方法更为通用,通过抽样的方法,即使目标分布的cdf不可逆仍然能够work。
根据概率分布的定义,对任意分布的概率密度求积分得1,即 $$ \int_{- \infty}^{+ \infty}f(x)dx=1 \quad\quad\quad\quad(1) $$
定义其CDF为F(x),有\(dF(x) = f(x)dx\) ,令\(y = F(x)\)带入上式:
$$ \int_{0}^{1}dy=1 \quad\quad\quad\quad(2) $$
注意到y在区间[F(a),F(b)]上符合均匀分布。
综上,目标的随机数x服从概率密度f(x),其CDF为\(y=F(x)\),y服从[0, 1]区间上的均匀分布,假设F(x)可逆,有\(x=F^{-1}(y)\)
指数分布的CDF为\(y=F(x)=1-e^{-\lambda x}\)
求逆得\(x=F^{-1}(y)=-\frac{1}{\lambda}ln(1-y)\)
因为y服从[0,1]的均匀分布,所以1-y也服从[0,1]的均匀分布,所以可以简写为 \(x=F^{-1}(y)=-\frac{1}{\lambda}lny\)
上面的方法,可以直接通过公式推到的方式得到目标分布的解析式,通过代入均匀分布上的随机数即可得到目标分布的随机数,非常高效,但问题是有些分布的CDF求逆非常困难,这时就需要一种新的方法可以在不求逆的情况下也可以得到服从概率密度为f(x)的随机数。
这里假设借助第一种方法可以高效的产生服从概率密度为g(x)的随机数,这个g(x)和f(x)有“某种形式的近似”(到底是那种形式的近似,后面可以看到),且\(\sup_x {\frac{f(x)}{g(x)}}<=c\)。
然后执行下面的步骤即可生成服从概率密度为f(x)的随机数:
1 从概率分布G中生成一个随机数Y
2 从[0,1]上的均匀分布生成一个随机数U
3 如果\(U < \frac{f(Y)}{c*g(Y)}\),成功返回X=Y,否则重复第1步
每次通过第3步的判断条件,以一定概率成功,否则会重复执行整个步骤,那么一个很自然的问题就是,平均执行几次才能成功返回?
设\(p = P(U < \frac{f(Y)}{c*g(Y)})\),那么第n次实验成功的概率分布为几何概率分布,\(P(N=n)=(1-p)^{n-1}p\),期望为\(E(n)=\frac{1}{p}\)
$$p = P(U < \frac{f(Y)}{cg(Y)}) = P(U < \frac{f(Y)}{cg(Y)}|Y=y)*P(Y=y) $$
$$ = \int_{- \infty}^{+ \infty}\frac{f(y)}{cg(y)} g(y)dy = \frac{1}{c} \quad\quad\quad\quad(3)$$
因此\(E(n) = c\)
我们当然希望实验次数越少越好(\(c \approx 1\)),当c最小是有\(c=\sup_x {\frac{f(x)}{g(x)}}\),当然如果g(x)和f(x)相等时就有c=1,当然这是不可能的,因为f(x)的概率密度累积函数很难求逆,而g(x)确是通过概率密度累积函数求逆得到的可以用均匀分布高效生成的一种分布,所以在选择g(x)更像是一种艺术,选择的好,可以让我们生成符合概率密度为f(x)的随机数非常容易,否则会是低效的。
我们要证明的目标是\(P(Y=y | U < \frac{f(Y)}{cg(Y)}) = F(y)\)
$$P(Y=y | U < \frac{f(Y)}{cg(Y)}) = P(U < \frac{f(Y)}{cg(Y)}|Y=y) * P(Y=y) / P(U < \frac{f(Y)}{cg(Y)})$$ $$ = \frac{F(y)}{cG(y)} * G(y) / (\frac{1}{c}) = F(y) \quad\quad\quad\quad(4) $$
证毕
选择\(g(x)=e^{-x}\),正态分布\(f(x) = \frac{1}{\sqrt{2\pi}}e^{-\frac{x^{2}}{2} }\),因此\(sup_x\frac{f(x)}{g(x)} =\frac{1}{\sqrt{2\pi}}e^{\frac{x - x^{2}}{2} } \approx 1.32\),即实验成功期望次数是1.32,计算好了\(\frac{f(x)}{g(x)} \),按照上述3步操作即可生成正态分布随机数。
~~~~
Ref:
http://blog.codinglabs.org/articles/methods-for-generating-random-number-distributions.html
http://blog.sciencenet.cn/blog-117333-563857.html
http://www.columbia.edu/~ks20/4703-Sigman/4703-07-Notes-ARM.pdf