生成泊松变量++
我实现这个功能，以生成泊松随机变量生成泊松变量++
typedef long unsigned int luint;
luint poisson(luint lambda) {
double L = exp(-double(lambda));
luint k = 0;
double p = 1;
do {
k++;
p *= mrand.rand();
} while(p > L);
return (k-1);
}
其中mrand是梅森旋转算法的随机数发生器。我发现，随着我增加lambda，预期的分布将会出现错误，平均值达到750左右。这是由于数值逼近还是我犯了错误？
exp（-750）是一个非常小的数字，非常接近最小可能的双倍，所以你的问题是数值。无论如何，你的复杂度在lambda中是线性的，所以对于高lambda来说算法不是很有效率。除非你有很好的理由来自己编写代码，否则使用现有的库实现可能是有意义的，因为这些数值算法对于你遇到的精确问题而言往往会很敏感。
我想我会使用正常的近似值，因为在我的情况下，lambda总是一个很大的数值。 – Bob 2011-04-14 03:28:29
由于在表达式(p>L)中只使用L，因此本质上是测试(log(p) > -lambda)。这不是一个非常有用的转变。当然，你不需要exp（-750），但是你只需要溢出p。
现在，p只是Π（mrand.rand（）），而log（p）是log（Π（mrand.rand（）））是Σ（log（mrand.rand（））。必要的改造：
double logp = 0;
do {
k++;
logp += log(mrand.rand());
} while(logp > -lambda);
double只有11指数位，但52位尾数因此这是数值稳定大规模增加支付的价格是您在每次迭代需要log，而不是。一个exp前面
在这样的情况下，你不需要多次调用随机数发生器。 u需要是累积概率的一个表：
double c[k] = // the probability that X <= k (k = 0,...)
然后生成随机数0 <= r < 1，并采取第一整数X使得c[X] > r。您可以通过二进制搜索找到这个X。
要生成此表，我们需要的个体概率
p[k] = lambda^k/(k! e^lambda) // // the probability that X = k
如果lambda很大，这成为非常不准确，因为你已经找到。但是我们可以在这里使用一个技巧：从最大值开始（或接近），用k = floor[lambda]，假装p[k]等于1。然后用递推关系
p[i+1] = (p[i]*lambda)/(i+1)
和i < k使用
p[i-1] = (p[i]*i)/lambda
这确保了最大的概率有最大可能的精确计算p[i]为i > k。
现在只需使用c[i+1] = c[i] + p[i+1]计算c[i]，直到c[i+1]与c[i]相同。然后你可以通过除以这个极限值c[i]来标准化数组;或者您可以保持原样，并使用随机数0 <= r < c[i]。
如果你走“现有的库”路线，你的编译器可能已经支持C++ 11 std :: random包。这里是你如何使用它：
#include <random>
#include <ctime>
#include <iostream>
std::mt19937 mrand(std::time(0)); // seed however you want
typedef long unsigned int luint;
luint poisson(luint lambda)
{
std::poisson_distribution<luint> d(lambda);
return d(mrand);
}
int main()
{
std::cout << poisson(750) << '\n';
std::poisson_distribution<luint> d(750);
std::cout << d(mrand) << '\n';
std::cout << d(mrand) << '\n';
}
我已经使用了两种方法上面：
我试图模仿现有的接口。
如果用平均值创建std :: poisson_distribution，则反复使用该分配的效果相同（如main（）中所做的那样）。
下面是示例输出对我来说：
751
730
779
IIRC，泊松变量具有指数分布。因此，这是http://*.com/questions/2106503/pseudorandom-number-generator-exponential-distribution的精确副本。但即使我错了，那里给出的方法也应该起作用。 – MSalters 2011-04-14 09:09:18
@ MSalters：泊松分布是离散的 - 它只需要整数值。指数分布是连续的。所以他们不一样（虽然他们是相关的）。 – TonyK 2011-04-14 11:16:20
*的权利：“如果给定时间间隔[0，t]内的到达数量遵循泊松分布，平均值=λt，那么到达间隔时间的长度将遵循指数分布，平均值为1/λ“。这是两者之间的有效转换，在结构上类似于我在下面提出的算法。 – MSalters 2011-04-14 11:57:12