准确性:C++11的binomial_distribution<int>与R返回的值不一致

Question

我需要在C++中生成遵循超几何分布的样本。但是，对于我的情况，我可以用二项分布来近似它，没有任何问题。

因此，我想使用C++11中的std实现，如果我在计算概率时生成许多样本，我得到的值与R告诉我的值不同。更重要的是，随着样本数量的增加，差异并不会变小。R和C++的参数相同。

因此，问题是:为什么我得不到相同的结果，我可以做什么/我应该相信哪一个？

参见下面的R和C++代码。C++程序计算R值的差值。即使我让程序运行很长一段时间，这个数字也不会变小，而是在E-5，E-6，E-7数量级上下摇摆。

R：

dbinom(0:2, 2, 0.48645948945615974379)
#0.26372385596962805154 0.49963330914842424280 0.23664283488194759464

C++：

#include <iostream>
#include <iomanip>
#include <random>

using namespace std;

class Generator {
public:
    Generator();
    virtual ~Generator();
    int binom();
private:
    std::random_device randev;
    std::mt19937_64 gen;
    std::binomial_distribution<int> dist;
};
Generator::Generator() : randev(), gen(randev()), dist(2,0.48645948945615974379) { }
Generator::~Generator() {}
int Generator::binom() { return dist(gen); }

int main() {
    Generator rd;
    const double nrolls = 10000000; // number of experiments
    double p[3]={};
    for (int k=1; k<100; ++k) {
        for (int i=0; i<nrolls; ++i) {
            int number = rd.binom();
            ++p[number];
        }

        cout << "Samples=" << setw(8) << nrolls*k <<
            "   dP(0)="<<setw(13)<<p[0]/(nrolls*k)-0.26372385596962805154<<
            "   dP(1)="<<setw(13)<<p[1]/(nrolls*k)-0.49963330914842424280<<
            "   dP(2)="<<setw(13)<<p[2]/(nrolls*k)-0.23664283488194759464<<endl;
    }
    cout<<"end";
    return 0;
}

选择性输出：

Samples=   1e+07   dP(0)=  -2.0056e-05   dP(1)=  9.49909e-05   dP(2)= -7.49349e-05
Samples=   1e+08   dP(0)=   1.5064e-05   dP(1)=  3.43609e-05   dP(2)= -4.94249e-05
Samples= 9.9e+08   dP(0)= -2.06449e-05   dP(1)=  5.93429e-06   dP(2)=  1.47106e-05

Answer 1

这真的应该是一个评论。

我看不出你的数据有什么问题。您正在进行10**9次重复。因此，根据中心极限定理，您应该看到精度约为10**(-4.5)。这确实是你所看到的。dP(0)和dP(2)的符号波动是另一个好兆头。如果您多次运行您的程序，那么最后一行上的符号是否始终显示相同的模式。如果不是，这是另一个好兆头。

顺便说一句，在我看来，R给了你太多的数字。对于doubles，你只有大约15位的准确度。