在C / C ++中以正态分布生成随机数
有谁知道我可以轻松生成随机数在C / C + +的正常分布?
http://www.mathworks.com/access/helpdesk/help/toolbox/stats/normrnd.html
我不想使用任何Boost。
我知道Knuth详细地谈了这个,但是我现在没有他的书。
Box-Muller变换是常用的。 这正确地产生具有正态分布的值。
http://en.wikipedia.org/wiki/Normal_distribution#Generating_values_from_normal_distribution
http://en.wikipedia.org/wiki/Box_Muller_transform
math很容易。 你生成两个统一的数字,从那里你得到两个正态分布的数字。 返回一个,另存为另一个随机数的下一个请求。
编辑:自2011年8月12日以来,我们有C + + 11直接提供了std::normal_distribution
,这是我今天的方式。
这是原来的答案:
这里有一些解决scheme是按上升的复杂度sorting的
-
从0到1 加上12个统一的数字并减去6.这将匹配正常variables的均值和标准偏差。 一个明显的缺点是,范围限于+/- 6 – 不像真正的正态分布。
-
上面列出了Box-Muller变换 ,实现相对简单。 如果您需要非常精确的样本,请注意Box-Muller变换与一些统一的发生器相结合,会产生称为Neave效应的exception。
HR Neave,“使用乘法同余伪随机数发生器的Box-Muller变换”,Applied Statistics,22,92-97,1973
-
为了获得最佳精度,我build议绘制制服,并应用逆累积正态分布来得到正态分布的variables。 你可以find一个非常好的逆累积正态分布algorithm
https://web.archive.org/web/20151030215612/http://home.online.no/~pjacklam/notes/invnorm/
希望有所帮助
彼得
一个简单而快速的方法就是求和一些均匀分布的随机数并取其平均值。 请参阅中心极限定理以获得充分的解释。
这里有一个基于一些引用的C ++示例。 这是快速和肮脏的,你最好不要重新发明和使用boost库。
#include "math.h" // for RAND, and rand double sampleNormal() { double u = ((double) rand() / (RAND_MAX)) * 2 - 1; double v = ((double) rand() / (RAND_MAX)) * 2 - 1; double r = u * u + v * v; if (r == 0 || r > 1) return sampleNormal(); double c = sqrt(-2 * log(r) / r); return u * c; }
你可以使用一个QQ图来检查结果,看看它接近一个真正的正态分布(排列你的样本1..x,把排名变成总数的比例,即多less个样本,得到z值并绘制它们,向上的直线就是期望的结果)。
我创build了一个用于正态分布随机数生成基准的C ++开源项目 。
它比较了几种algorithm,其中包括
- 中心极限定理方法
- Box-Muller变换
- 马尔加利亚极地法
- Zigguratalgorithm
- 逆变换抽样方法。
-
cpp11random
使用C ++ 11std::normal_distribution
和std::minstd_rand
(实际上它是在clang中的Box-Muller变换)。
iMac Corei5-3330S@2.70GHz上的单精度( float
)版本的结果,铛6.1,64位:
为了正确,程序validation了样本的均值,标准差,偏度和峰度。 与其他方法相比,通过将4,8或16个统一数相加的CLT方法没有良好的峰度。
Zigguratalgorithm比其他algorithm具有更好的性能。 但是,它不适用于SIMD并行,因为它需要查表和分支。 具有SSE2 / AVX指令集的Box-Muller比zigguratalgorithm的非SIMD版本要快得多(x1.79,x2.99)。
因此,我build议使用Box-Muller作为SIMD指令集的体系结构,否则可能是Ziggurat。
PS基准使用最简单的LCG PRNG来生成均匀分布的随机数。 所以对某些应用来说可能是不够的。 但是性能比较应该是公平的,因为所有的实现都使用相同的PRNG,所以基准testing主要是testing转换的性能。
使用std::tr1::normal_distribution
。
std :: tr1命名空间不是boost的一部分。 这是包含来自C ++技术报告1的库添加的名称空间,并且可以使用最新的Microsoft编译器和gcc,与boost无关。
这是如何生成现代C ++编译器的样本。
#include <random> ... std::mt19937 generator; double mean = 0.0; double stddev = 1.0; std::normal_distribution<double> normal(mean, stddev); cerr << "Normal: " << normal(generator) << endl;
您可以使用GSL 。 给出了一些完整的例子来演示如何使用它。
看看: http : //www.cplusplus.com/reference/random/normal_distribution/ 。 这是产生正态分布的最简单的方法。
如果你使用C ++ 11,你可以使用std::normal_distribution
:
#include <random> std::default_random_engine generator; std::normal_distribution<double> distribution(/*mean=*/0.0, /*stddev=*/1.0); double randomNumber = distribution(generator);
还有很多其他的发行版可以用来转换随机数引擎的输出。
我遵循http://www.mathworks.com/help/stats/normal-distribution.html中给出的PDF的定义,并提出了这个问题:;
const double DBL_EPS_COMP = 1 - DBL_EPSILON; // DBL_EPSILON is defined in <limits.h>. inline double RandU() { return DBL_EPSILON + ((double) rand()/RAND_MAX); } inline double RandN2(double mu, double sigma) { return mu + (rand()%2 ? -1.0 : 1.0)*sigma*pow(-log(DBL_EPS_COMP*RandU()), 0.5); } inline double RandN() { return RandN2(0, 1.0); }
这可能不是最好的方法,但它很简单。
看看我发现的。
这个库使用Zigguratalgorithm。
comp.lang.c常见问题清单共享三种不同的方式,可以轻松生成高斯分布的随机数字。
你可以看看它: http : //c-faq.com/lib/gaussian.html
存在用于逆累积正态分布的各种algorithm。 量化金融最stream行的是在http://chasethedevil.github.io/post/monte-carlo-inverse-cumulative-normal-distribution/上进行testing
此外,它显示了Ziggurat类似的方法的缺点。
计算机是确定性的设备。 计算中没有随机性。 另外,CPU中的运算器可以通过一些有限的整数(在有限域中进行评估)和有限的实有理数集来评估。 并且还执行按位操作。 math与更多伟大的集合,如[0.0,1.0]与无限的点数交易。
你可以用一些控制器听电脑里面的电线,但是它的分布是否均匀? 我不知道。 但如果假设它的信号是累加值大量的独立随机variables的结果,那么你将接收到近似正态分布的随机variables(这已经在概率论中certificate了)
存在称为 – 伪随机生成器的algorithm。 正如我所说的,伪随机生成器的目的是模拟随机性。 而美好的标准是: – 经验分布(从某种意义上说 – 逐点,一致,L2)收敛到理论值 – 你从随机发生器收到似乎是相互依赖的。 当然,从“真实的angular度来看”并不是真的,但我们假设它是真实的。
一种stream行的方法 – 你可以用均匀分布来求解12个IRV。但是,在推导中值极限定理和傅里叶变换的帮助下,泰勒级数,要求有n – > + inf几次的假设。 所以举例来说理论上 – 我个人并不认为人们是如何以均匀的分布来执行12个IRV的。
我在大学有过可变性理论。 尤其对我来说,这只是一个math问题。 在大学我看到了以下模型:
double generateUniform(double a, double b) { return uniformGen.generateReal(a, b); } double generateRelei(double sigma) { return sigma * sqrt(-2 * log(1.0 - uniformGen.generateReal(0.0, 1.0 -kEps))); } double generateNorm(double m, double sigma) { double y2 = generateUniform(0.0, 2 * kPi); double y1 = generateRelei(1.0); double x1 = y1 * cos(y2); return sigma*x1 + m; }
这样的方式如何仅仅是一个例子,我想它存在另一种实现它的方法。
在这本书“莫斯科,BMSTU,2004:XVI概率论,例6.12,第246-247页”中可以find正确的证据,Krishchenko Alexander Petrovich ISBN 5-7038-2485-0
不幸的是我不知道这本书翻译成英文的存在。