如何将均匀分布(如大多数随机数发生器产生的,例如0.0到1.0之间)转换为正态分布? 如果我想要选择均值和标准差怎么办?
Ziggurat算法对此非常有效,尽管Box-Muller变换更容易从头开始实现(而不是疯狂慢)。
有很多方法:
-
请勿使用Box Muller。特别是当您绘制许多高斯数时。 Box Muller得出的结果被限制在-6和6之间(假设双精度。浮点会使情况更糟)。而且它确实比其他可用方法效率低。
-
Ziggurat很好,但需要进行表查找(由于缓存大小问题,需要进行一些特定于平台的调整)
-
均匀率是我的最爱,只有几次加法/乘法和对数的1/50(例如,看那里)。
-
反转CDF是高效的(并且被忽略了,为什么?),如果您搜索google,就可以快速实现它。准随机数是必需的。
将任何函数的分布更改为另一个函数都涉及使用所需函数的逆函数。
换句话说,如果您针对特定的概率函数p(x),则可以通过对其积分-> d(x)=积分(p(x))并使用其反函数来获得分布:Inv(d(x)) 。现在使用随机概率函数(具有均匀分布)并通过函数Inv(d(x))转换结果值。您应该根据选择的功能获得带有分布的随机值。
这是通用的数学方法-通过使用它,您现在可以选择具有逆或良好逆近似的任何概率或分布函数。
希望这对您有所帮助,并感谢您对使用分布的简短评论,而不是概率本身。
这是使用Box-Muller变换的极坐标形式的javascript实现。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
| /*
* Returns member of set with a given mean and standard deviation
* mean: mean
* standard deviation: std_dev
*/
function createMemberInNormalDistribution(mean,std_dev){
return mean + (gaussRandom()*std_dev);
}
/*
* Returns random number in normal distribution centering on 0.
* ~95% of numbers returned should fall between -2 and 2
* ie within two standard deviations
*/
function gaussRandom() {
var u = 2*Math.random()-1;
var v = 2*Math.random()-1;
var r = u*u + v*v;
/*if outside interval [0,1] start over*/
if(r == 0 || r >= 1) return gaussRandom();
var c = Math.sqrt(-2*Math.log(r)/r);
return u*c;
/* todo: optimize this algorithm by caching (v*c)
* and returning next time gaussRandom() is called.
* left out for simplicity */
} |
使用中央极限定理Wikipedia条目mathworld条目可发挥您的优势。
生成n个均匀分布的数字,将它们相加,减去n * 0.5,您将得到近似正态分布的输出,均值等于0,方差等于(1/12) * (1/sqrt(N))(关于最后一个均匀分布,请参见Wikipedia)
n = 10可以使您快得快一半。如果您想要的东西超过一半,请寻求轮胎解决方案(正态分布在Wikipedia条目中已指出)
我会使用Box-Muller。关于以下两点:
您每次迭代都会得到两个值
通常,您缓存一个值,然后返回另一个值。在下一次调用样本时,您返回缓存的值。
Box-Muller给出Z得分
然后,您必须按标准偏差缩放Z分数,并添加均值以获得正态分布中的完整值。
其中R1,R2是随机统一数字:
正态分布,SD为1:sqrt(-2 * log(R1))* cos(2 * pi * R2)
这是正确的……不需要执行所有这些慢循环!
八年后我可以添加一些东西似乎令人难以置信,但是对于Java,我想向读者介绍Random.nextGaussian()方法,该方法为您生成均值0.0和标准差1.0的高斯分布。
简单的加法和/或乘法将改变均值和标准差,以满足您的需求。
标准的Python库模块random具有您想要的:
normalvariate(mu, sigma)
Normal distribution. mu is the mean, and sigma is the standard deviation.
对于算法本身,请查看Python库中random.py中的函数。
手动输入在这里
问:如何将均匀分布(如大多数随机数生成器产生的,例如0.0到1.0之间)转换为正态分布?
对于软件实现,我知道几个随机生成器名称,这些名称在[0,1]中为您提供了伪统一的随机序列(Mersenne Twister,线性一致生成器)。我们称它为U(x)
存在于数学领域中称为概率论。
第一件事:如果您想对r.v建模积分分布为F的情况下,您可以尝试仅求F ^ -1(U(x))的值。在理论上证明了这样的r.v.将具有积分分布F。
当F ^ -1可以无问题地解析得出时,第2步可适用于生成r.v.?F,而无需使用任何计数方法。 (例如exp.distribution)
为了建模正态分布,您可以计算y1 * cos(y2),其中y1?[2pi]是均匀的。 y2是相对分布。
问:如果我要选择平均值和标准差怎么办?
您可以计算sigma * N(0,1)+ m。
可以证明,这种移位和缩放导致N(m,sigma)
我想您应该在EXCEL:=norminv(rand();0;1)中尝试此操作。这将产生应该以零均值和单位方差正态分布的随机数。可以提供任何值" 0",以便数字具有期望的均值,并且通过更改" 1",您将获得与输入平方相等的方差。
例如:=norminv(rand();50;3)将产生为MEAN = 50 VARIANCE = 9的正态分布数字。
使用实现的函数rnorm()也比创建正态分布的随机数生成器要快,因为它比写随机数生成器要快。参见以下代码作为证明
1 2 3 4 5
| n <- length(z)
t0 <- Sys.time()
z <- rnorm(n)
t1 <- Sys.time()
t1-t0 |
这是使用Box-Muller变换的极坐标形式的Matlab实现:
功能randn_box_muller.m:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
| function [values] = randn_box_muller(n, mean, std_dev)
if nargin == 1
mean = 0;
std_dev = 1;
end
r = gaussRandomN(n);
values = r.*std_dev - mean;
end
function [values] = gaussRandomN(n)
[u, v, r] = gaussRandomNValid(n);
c = sqrt(-2*log(r)./r);
values = u.*c;
end
function [u, v, r] = gaussRandomNValid(n)
r = zeros(n, 1);
u = zeros(n, 1);
v = zeros(n, 1);
filter = r==0 | r>=1;
% if outside interval [0,1] start over
while n ~= 0
u(filter) = 2*rand(n, 1)-1;
v(filter) = 2*rand(n, 1)-1;
r(filter) = u(filter).*u(filter) + v(filter).*v(filter);
filter = r==0 | r>=1;
n = size(r(filter),1);
end
end |
并调用histfit(randn_box_muller(10000000),100);这是结果:
显然,与Matlab内置randn相比,它确实效率很低。
我有以下代码可能会有所帮助:
1 2 3 4 5 6 7 8
| set.seed(123)
n <- 1000
u <- runif(n) #creates U
x <- -log(u)
y <- runif(n, max=u*sqrt((2*exp(1))/pi)) #create Y
z <- ifelse (y < dnorm(x)/2, -x, NA)
z <- ifelse ((y > dnorm(x)/2) & (y < dnorm(x)), x, z)
z <- z[!is.na(z)] |
1 2 3 4 5 6
| function distRandom(){
do{
x=random(DISTRIBUTION_DOMAIN);
}while(random(DISTRIBUTION_RANGE)>=distributionFunction(x));
return x;
} |
近似:
1 2 3
| function rnd_snd() {
return (Math.random()*2-1)+(Math.random()*2-1)+(Math.random()*2-1);
} |
参见http://www.protonfish.com/random.shtml