Solved – How to make Marsaglia polar method return values [0, 1)

normal distributionrandom-generation

I've implemented the Marsaglia polar method to generate random numbers that are normally distributed. Unfortunately the code as shown on this website returns numbers that are not within the range [0, 1).

Code:

     double x1, y1, w;

     do {
             x1 = 2.0 * ranf() - 1.0;
             y1 = 2.0 * ranf() - 1.0;
             w = x1 * x1 + y1 * y1;
     } while ( w >= 1.0 || w == 0 );

     w = sqrt( (-2.0 * ln( w ) ) / w );
     return x1 * w;

where ranf() returns values [0, 1).

For example, if (x1, y1) was (-0.43458512585358, -0.07521858582050478), this returns a value of -1.7830255550765148. Obviously not a value within my expected range.

I've seen some implementations that multiplies this by the standard deviation and adds the mean. But if I want to get numbers that range from [0, 1), what should I use as my standard deviation and mean? currently, I'm using a standard deviation of 1/8.125 ~= 0.121212… and 0.5 as the mean, but I've only stumbled on this by experimentation. What is the official way for me to get this appropriately ranged?

To answer the question in the comments below: I want to generate normal random variables within the range 0 and 1. Obviously I'd like the mean to be at 0.5 and values no be distributed around this mean.

More for my enlightenment, will the polar method really return numbers over the entire real line? If so, then it should be easy for me to scale the result down to [0, 1), assuming I don't hit floating point underflow due to the division to scale things down.

New edit:
What I'm trying to do is take Aniko's advise in his answer to my other Power Factor Problem question. So my first thought was that simulate my distribution of velocities, was to implement my on GaussianRandom class by deriving from Random and overriding Random.Sample() and the 3 other methods as recommended by MSDN. Since Random.Sample() is suppose to return a double in the range [0, 1), then I've hit this stumbling block.

I've come to the same conclusion as mbq's answer after staring at the ln(x) graph where x is in (0, 1] and sqrt(x) where x is in (0, 2].

From a practical/programming standpoint, I can't override Random.Sample(), because later when I actually call GaussianRandom.NextDouble(), there is no way for me to pass in the actual standard deviations for my velocities. It looks like I'll have to use the Marsaglia polar method in the raw instead of wrapped within the Random class. My "cartridge" generator will have to return x1 * w * stdDev + mean.

Best Answer

It is just impossible, since normal distribution has nonzero probability for the whole $R$. Of course because of the numerical issues it would be not the case for Marsaglia method (this obviously means that it is not creating normally distributed numbers), yet this numerical range is so huge that it is nonsense to scale it to [0,1].

EDIT: So, as long as I understand, you want this number to be in [0,1) to overload .net Random; you don't need to! Write your own class, say NormalRandom which will hold one Random object and will use its output to generate (using Marsaglia) normal random numbers from $-\infty$ to $\infty$ and just use them. Or rather use something already done; Google found http://www.codeproject.com/KB/recipes/Random.aspx .
On margin, remember that when it is said that something which cannot be negative is "normally distributed" (for instance mass) it does mean that it is not from normal but from "almost-normal" distribution which looks pretty normal, but still has one tail cut out (nevertheless simple rejection of negative outputs from generator is a fair way to deal with it).

Related Question