Re: Normal distribution random numbers

[Mikael Djurfeldt <mikael@djurfeldt.com>]

(Maybe I was sloppy now. An anomaly might actually be detected sooner due
to increasing distances in the discretization of values in the tail, but I
still don't think this is something to worry about.)

Den tors 4 juni 2020 17:11Mikael Djurfeldt <mikael@djurfeldt.com> skrev:

> Yes
>
> Den tors 4 juni 2020 17:08Zelphir Kaltstahl <zelphirkaltstahl@posteo.de>
> skrev:
>
>> Hi Mikael!
>>
>> Thanks for putting that into perspective and giving some numbers!
>>
>> When I looked at the code of Guile for random:normal, I also guessed,
>> that it makes use of that Box-Muller-transform, but wasn't sure, so thanks
>> for confirming that as well.
>>
>> So basically the tails are wrong, but to draw a number in the area where
>> the tails are wrong is so unlikely, that it would take that much time, as
>> stated in your number example, if I understand this correctly(?/.)
>>
>> Regards,
>>
>> Zelphir
>> On 04.06.20 17:03, Mikael Djurfeldt wrote:
>>
>> Hi Zelphir,
>>
>> random:normal actually uses the Box-Muller-transform. But since it uses
>> 64 bits, we only loose values that would be generated once in 2*10^20. That
>> is, if we could draw one billion numbers per second, such values would be
>> drawn once in 7000 years. So, we would start noticing an anomaly after
>> maybe 100000 years or so.
>>
>> But maybe we should replace this with some more correct and efficient
>> algorithm at some point.
>>
>> Best regards,
>> Mikael
>>
>> Den lör 30 maj 2020 22:43Zelphir Kaltstahl <zelphirkaltstahl@posteo.de>
>> skrev:
>>
>>> I just realized, that I did not check what Guile implements as
>>> non-SRFIs. I found:
>>> https://www.gnu.org/software/guile/manual/html_node/Random.html which
>>> has `random:normal`! I should have checked that first. Still good to
>>> know, what a can of worms normal distribution implementation can be.
>>>
>>> On 30.05.20 22:21, Zelphir Kaltstahl wrote:
>>> > Hi Guile Users!
>>> >
>>> > I recently wrote a little program involving lots of uniformly
>>> > distributed random integers. For that I used SRFI-27 and it works fine.
>>> >
>>> > Then I thought: How would I get normal distributed random numbers? I
>>> > don't have a project or program in mind for this, but it struck me,
>>> that
>>> > I do not know, how to get a normal distribution from a uniform
>>> > distribution. So I dug into the matter …
>>> >
>>> > Turns out the math is not really my friend:
>>> >
>>> > * https://stackoverflow.com/a/3265174 – OK, if that's true, then don't
>>> > use Box-Muller-Transform
>>> > * https://stackoverflow.com/a/86885 – The what? I need to somehow
>>> > inverse the Gaussian distribution to get a function to calculate normal
>>> > distributed values from uniformly distributed values? Something like
>>> > that. Safe to say it is above my current math skills.
>>> > * The wiki page also does not help me much:
>>> > https://en.wikipedia.org/wiki/Inverse_transform_sampling Seems too
>>> > complicated.
>>> >
>>> > So I thought: "OK, maybe I can simply copy, how other languages
>>> > implement it!" The wiki page mentions, that R actually makes use of the
>>> > inverse thingy. So I set out to look at R source code:
>>> >
>>> > * https://github.com/wch/r-source/blob/master/src/nmath/rnorm.c – OK,
>>> > looks simple enough … Lets see what `norm_rand` is …
>>> > * https://github.com/wch/r-source/blob/master/src/nmath/snorm.c#L62 –
>>> > yeah … well … I'm not gonna implement _that_ pile of … Just look at the
>>> > lines
>>> >
>>> https://github.com/wch/r-source/blob/master/src/nmath/snorm.c#L135-L196
>>> > what a mess! Not a single comment to help understanding in it. Such a
>>> > disappointment.
>>> > * Python also seems to only use an approximation with magic constants:
>>> > https://github.com/python/cpython/blob/3.8/Lib/random.py#L443
>>> >
>>> > So it seems, that there is no easy way to implement it properly with
>>> > correct tails to the left and right side of the distribution, something
>>> > clean and not made with mathematical traps built-in. Or is there?
>>> >
>>> > I found a post about using 2 normal distributions to do
>>> > Box-Muller-transform:
>>> >
>>> https://www.alanzucconi.com/2015/09/16/how-to-sample-from-a-gaussian-distribution/
>>> >
>>> > However, it seems to require a uniform float not integer and it is the
>>> > Box-Muller-transform, which is said to clamp between -6 and 6 according
>>> > to the people writing the answers on stackoverflow.
>>> >
>>> > So my question is: Is there a good implementation in the Guile universe
>>> > already? (Or a simple way to implement it?) I don't really need it
>>> right
>>> > now, but I think this thing could be an obstacle for many people
>>> without
>>> > serious math knowledge and it would be good to know, where to find it,
>>> > should one have need for normal distributed random numbers.
>>> >
>>> > Regards,
>>> > Zelphir
>>> >
>>> >
>>>
>>>