逆（CDF）变换采样的错误分布

Question

我试图模拟一个几何分布，使用逆CDF方法，但我得到的结果略有错误，我不知道为什么。你知道吗

更具体地说，形状系数p=0.8的几何分布应具有以下特征：

mean: 1.25 
variance: 0.31

但是，运行下面的代码，我得到：

mean: 0.6224363901913519
var: 0.391813011265263
[Finished in 0.3s]

如你所见，我得到的平均值与预期值相差很大。你知道吗

np.log(uniform[i])/np.log(1-p) is the result of solving the equation: F(X) = R for X in terms of R, F(X) = CDF of geometric distribution = 1 - (1 - p)^k.

R是区间（0,1）上的均匀分布。你知道吗

因此，解决它会产生以下结果：

X = ln(1-R)/ln(1-p)

然而，由于1-R和R都均匀分布在（0,1）上，我们可以做以下简化：

X = ln(R)/ln(1-p)

上面的公式是正确的，应该得到一个几何分布样本。你知道吗

import numpy as np

n = 10000
p = 0.8
geo_dist = np.zeros(n,dtype = np.float64)
uniform = np.random.uniform(0, 1, n)
for i in range(n):
    geo_dist[i] = np.log(uniform[i])/np.log(1-p)
print("mean: " +str(geo_dist.mean()))
print("var: " +str(geo_dist.var())) 

我试着通过使用第64页在一个绝望的尝试，以修复什么应该是一个琐碎的脚本，没有用。你知道吗

我还尝试使用scipy生成均匀分布制服.rvs（）而不是np.统一问题依然存在。你知道吗

如果p=0.5：

expected mean: 2
expected variance : 2

但是，我编写的代码有以下结果：

mean: 1.4440009653569306
var: 2.0421079966161093
[Finished in 0.3s]

有人知道为什么这样不行吗？ 非常感谢。你知道吗

Answer 1

实际上，您正在对lambda等于-1/ln（1-p）的连续exponential distribution进行采样

好的，这是正确采样的代码，上限应用于指数输出

import numpy as np

N = 100000
p = 0.8

q = np.random.random(N)
g = np.ceil(np.log(1.0 - q)/np.log(1.0-p))

print(np.mean(g))
print(np.var(g))

哪个指纹

1.25055
0.3146946975

请注意：

1. 最好使用无显式循环的NumPy矢量化功能

2. 用(1-R) -> R替换从U（0,1）采样的R对NumPy RNG不起作用-它返回半封闭范围[0…1]中的值，这意味着您可以不时获得log（0）和FP exception。