我想在Python中指定一个分布的概率密度函数,然后从该分布中随机选取N个数字。我该如何做呢?
我想在Python中指定一个分布的概率密度函数,然后从该分布中随机选取N个数字。我该如何做呢?
通常,您需要使用累积概率密度函数的倒数。一旦掌握了这一点,那么在分布上生成随机数就很简单:
import random
def sample(n):
return [ icdf(random.random()) for _ in range(n) ]
或者,如果你使用NumPy:
import numpy as np
def sample(n):
return icdf(np.random.random(n))
在这两种情况下,icdf
是一个反累积分布函数,它接受介于0和1之间的值,并输出相应的分布值。
为了说明icdf
的性质,我们将以一个简单的均匀分布为例,该分布在10和12之间:
概率分布函数在10到12之间为0.5,在其他地方为零
累积分布函数在10以下为0(没有样本低于10),在12以上为1(没有样本高于12),并在线性增加值之间(PDF的积分)
反向累积分布函数仅在0到1之间定义。 在0处,它是10,在12处,它是1,并在线性变化值之间变化
当然,困难的部分是获得反向累积密度函数。 这实际上取决于您的分布,有时您可能拥有解析函数,有时您可能希望求助于插值。 数值方法可能很有用,因为可以使用数值积分来创建CDF,并可以使用插值来反转它。
"""
Draws a random number from given probability density function.
Parameters
----------
pdf -- the function pointer to a probability density function of form P = pdf(x)
interval -- the resulting random number is restricted to this interval
pdfmax -- the maximum of the probability density function
integers -- boolean, indicating if the result is desired as integer
max_iterations -- maximum number of 'tries' to find a combination of random numbers (rand_x, rand_y) located below the function value calc_y = pdf(rand_x).
returns a single random number according the pdf distribution.
"""
def draw_random_number_from_pdf(pdf, interval, pdfmax = 1, integers = False, max_iterations = 10000):
for i in range(max_iterations):
if integers == True:
rand_x = np.random.randint(interval[0], interval[1])
else:
rand_x = (interval[1] - interval[0]) * np.random.random(1) + interval[0] #(b - a) * random_sample() + a
rand_y = pdfmax * np.random.random(1)
calc_y = pdf(rand_x)
if(rand_y <= calc_y ):
return rand_x
raise Exception("Could not find a matching random number within pdf in " + max_iterations + " iterations.")
在我看来,如果你不需要检索大量的随机变量,这种解决方案比其他解决方案表现更好。另一个好处是你只需要PDF,避免计算CDF、反向CDF或权重。