问题
已提供一个Rand7()的API可以随机生成1到7的数字,使用Rand7实现Rand10,Rand10可以随机生成1到10的数字。
思路
简单说:
(1)通过(Rand N - 1) % 10 + 1的方法,可以求出Rand10,当N是10的倍数的时候。
(2)用( Rand7 - 1 ) * 7 + Rand7可以随机生成1-49,记作Rand49。
(3)如果可以通过Rand49计算出Rand40,即随机生成1-40,就可以通过Rand40 % 10来取得Rand10。
(4)如何通过Rand49计算出Rand40呢?可以通过拒绝采样的方法来计算:用Rand49生成一个数,如果它位于41-49,则丢弃,继续生成,当生成一个1-40的数时返回。这样子做可以近似看成随机生成1-40之间的数。
(5)到这里就可以把代码写出来了,但是为什么这样子做是有用的呢?下面阐述一下从拒绝采样到本题的一个思路历程。
拒绝采样(Reject Sampling)
首先了解几个知识点:
连续分布中,累积分布函数CDF为\(F(x) = P \{X < x\}\),表示随机变量X落在区间\((-\infty, x)\)的概率。
概率密度函数PDF的意义在积分的时候得到体现,$ \int_{a}^{b} f(x) dx = F(b) - F(a) = P { a \le X \lt b }$。
离散分布中,累积分布函数CDF一样有意义,\(F(x) = P \{X < x\}\)表示小于x的那些\(x_k\)处的\(p_k\)之和。
概率质量函数PMF表示取得单个点的概率,\(f(x) = P \{ X = x\}\)。
拒绝采样基于这样一个前提: 对一个随机变量取样, 等价于从这个变量所服从的密度函数下方的区域均匀取样。
以下说的p(x),q(x)指的是连续分布中的概率密度函数。我们用概率密度函数来表示一个概率分布。
假设我们想对p(x)分布采样,但是由于种种原因,我们很难对p(x)采样。但是另一个分布q(x)是可以采样的。那么我们可以根据q(x)的采样来得到一个在p(x)上的近似采样,具体做法如下。
(1)将q(x)和一个大于1的常数M相乘,使得M * q(x)可以把p(x)罩住,如下图所示(图源见参考文献)。
因为概率密度函数面积为1,所以q(x)肯定是有些地方比p(x)大,有些地方比q(x)小,如果全都比p(x)大,面积就超过1了。为了罩住p(x),也就是处处比p(x)大,需要乘以一个大于1的常数M,这就是M * q(x)。
(2)从q(x)中获得一个采样点x(i),对于x(i)计算一个接受概率 $ \alpha = \frac{p(x_i)}{Mq(x_i)}\(,从Uniform(0,1)随机生成一个值\)\mu\(,如果\)\alpha \ge \mu$,则接受x(i)作为一个来自p(x)的采样点。否则继续采样。
因为从q(x)分布取得这个点的概率是q(x),并不是p(x),(实际上对于连续分布,随机变量为单个数值的概率为0,但是这里我们讨论的是采样,从密度函数下方的区域均匀采样,这个场景下讨论时我们认为这个概率是有意义的值的而不是0),要把它变成p(x)需要多一个接受概率来进行选择,理想情况是选择p(x) / q(x)作为接收概率,因为q(x)乘以这个接受概率刚好就等于p(x)。
但是p(x) / q(x)是可能大于1的,正如前面所说p(x)有些地方比q(x)高,所以对q(x)乘上一个M,得到一个处于(0,1)的p(x) / (M * q(x)),可以作为接受概率,这个接受概率可以从一个(0,1)的均匀分布上取得。从图中可以看出,p(x)越接近M * q(x)的时候,接受概率越大,这很符合我们的采样直觉,因为M * q(x)正比于q(x)而M又是已经给出,所以p(x)越接近M * q(x)也表示p(x)越接近q(x),这个时候它们的概率比较接近,当然要接受概率大一些。反过来,当p(x)和M * q(x)差很多时,p(x)和q(x)的差距也很大,这个时候根据q(x)采样到的样本很大概率要被拒绝。
最终采样的概率为\(\frac{p(x)}{M q(x)} q(x) = \frac {1}{M} p(x)\),其中\(\frac{p(x)}{M q(x)}\)可以从一个均匀分布采样,q(x)是已知的可以采样的分布,通过这两个分布采样可以得到一个近似p(x)的采样。显然,当我们取的M越接近1的时候,这个近似采样越接近p(x)上的采样。如果p(x)和q(x)完全一样,M就是1,在q上采样就等于在p上采样。
回到本题
Rand7以随机生成1-7,那么用( Rand7 - 1 ) * 7 + Rand7就可以随机生成1-49,记作Rand49()。
如果能够计算出Rand40,即随机生成1-40,就可以通过(Rand40 - 1) % 10 + 1的方法来取得Rand10。
那么知道Rand49的情况下如何计算出Rand40呢?这个场景就类似上面拒绝采样的场景:我们可以从q(x)分布上采样,如何模拟出一个p(x)上的近似采样。对离散分布当然也可以用拒绝采样,把上述的概率密度函数写成概率质量函数即可。
把Rand40对应的分布当作上述的p(x),然后把Rand49对应的分布当作上述q(x),然后M取49/40。
我们在q(x)上采样相当于用Rand49生成数字,每个数生成概率为 1/49,乘上M就是 1/40,对于每个随机变量,M * q(x)为1/40。p(x)是Rand40对应的分布,随机变量为1-40时,p(x) = 1/ 40,那么接受概率为p(x) / (M * q(x))就是1。随机变量为41-49时,p(x)为0,所以接受概率为0。
利用被拒绝的数
在Rand49取到41-49时,会被拒绝,那么可不可以在取到41-49时不拒绝,而是把这些数利用起来呢?
现在考虑取到的数是41-49之间的情况,减掉40规整为1-9,假设这个值为val,val处于1-9之间,考虑再用rand7(),构建一个[1, 63]的均匀分布,也就是 (val - 1) * rand7() + 7。再根据拒绝采样在[1, 63]的均匀分布中采1-60之间的值,然后用%10的方法得到Rand10。
接着上面的操作,如果采的数是61-63时,先减掉60规整为1-3,假设这个值为val,val处于1-3之间,考虑再用rand7(),构建出[1, 21]的均匀分布,也就是 (val - 1) * rand7() + 7。再根据拒绝采样采1-20之间的值,用%10的方法得到Rand10。
接着上面的操作,如果采的数是21时,只有一个数,这个时候就直接拒绝,继续按上面的一系列操作采样。
代码
# The rand7() API is already defined for you.
# def rand7():
# @return a random integer in the range 1 to 7
class Solution(object):
def rand10(self):
"""
:rtype: int
"""
while True:
val = (rand7() - 1) * 7 + rand7()
if(val <= 40):
break
return (val-1) % 10 + 1
class Solution(object):
def rand10(self):
"""
:rtype: int
"""
while True:
val = (rand7() - 1) * 7 + rand7()
if val <= 40:
return (val-1) % 10 + 1
val -= 40
val = (val - 1) * 7 + rand7()
if val <= 60:
return (val-1) % 10 + 1
val -= 60
val = (val - 1) * 7 + rand7()
if val <= 20:
return (val-1) % 10 + 1
参考资料
原文地址:https://www.cnblogs.com/liaohuiqiang/p/9857273.html