Normal distribution random numbers

Normal distribution random numbers

[Zelphir Kaltstahl]

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

Re: Normal distribution random numbers

[Zelphir Kaltstahl]

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
>
>

Re: Normal distribution random numbers

[Arne Babenhauserheide]

Hi Zelphir,

Zelphir Kaltstahl <zelphirkaltstahl@posteo.de> writes:
> 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 …
…
> 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.

I don’t know how random:normal does it, just that I used it.

See https://www.gnu.org/software/guile/manual/html_node/Random.html

If you want to add defined covariance, you can use cholesky
decomposition: https://hg.sr.ht/~arnebab/wisp/browse/examples/cholesky.scm?rev=91ec8dc32652

Best wishes,
Arne
-- 
Unpolitisch sein
heißt politisch sein
ohne es zu merken

Re: Normal distribution random numbers

[Zelphir Kaltstahl]

Hi Arne!

Thanks for the pointers!

On 5/30/20 11:30 PM, Arne Babenhauserheide wrote:
> Hi Zelphir,
>
> Zelphir Kaltstahl <zelphirkaltstahl@posteo.de> writes:
>> 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 …
> …
>> 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.
> I don’t know how random:normal does it, just that I used it.
>
> See https://www.gnu.org/software/guile/manual/html_node/Random.html
>
> If you want to add defined covariance, you can use cholesky
> decomposition: https://hg.sr.ht/~arnebab/wisp/browse/examples/cholesky.scm?rev=91ec8dc32652
>
> Best wishes,
> Arne

Re: Normal distribution random numbers

[tantalum]

surely not the ideal way to generate numbers with a normal distribution, 
but there is a way to use custom probabilities from a list, which i 
think is nice to know.
it works like this:

have a list of probabilities. can be as long as you want. i think that 
is called probability density.
->
(1 0 3 1)

create the cumulative sums for this list. that is, sums like this
   (a b c ...) -> (a (+ a b) (+ a b c) ...)
i think that is called cumulative distribution.
->
(1 1 4 5)

create a random number up to the largest sum
->
(random 5)

return the first index of the list of cumulative sums that is greater 
than the random number.
given the distribution above, you would see index 2 a lot, never index 
1, and index 0 and 3 rarely.

~~~
(use-modules ((srfi srfi-1) #:select (last)))

(define (cusum a . b)
   "calculate cumulative sums from the given numbers.
    (a b c ...) -> (a (+ a b) (+ a b c) ...)"
   (cons a (if (null? b) (list) (apply cusum (+ a (car b)) (cdr b)))))

(define* (random-discrete-f probabilities #:optional (state 
*random-state*))
   "(real ...) [random-state] -> procedure:{-> integer}"
   (let* ((cuprob (apply cusum probabilities)) (sum (last cuprob)))
     (lambda ()
       (let ((deviate (random sum state)))
         (let loop ((a 0) (b cuprob))
           (if (null? b) a (if (< deviate (car b)) a (loop (+ 1 a) (cdr 
b)))))))))

(define random* (random-discrete-f (list 1 0 3 1) 
(random-state-from-platform)))
(display (random*))
~~~

Re: Normal distribution random numbers

[Zelphir Kaltstahl]

Hi!

Interesting technique, that when thinking about it, intuitively makes
sense to me. Thanks for sharing!

Regards,
Zelphir   

On 31.05.20 17:12, tantalum wrote:
> surely not the ideal way to generate numbers with a normal
> distribution, but there is a way to use custom probabilities from a
> list, which i think is nice to know.
> it works like this:
>
> have a list of probabilities. can be as long as you want. i think that
> is called probability density.
> ->
> (1 0 3 1)
>
> create the cumulative sums for this list. that is, sums like this
>   (a b c ...) -> (a (+ a b) (+ a b c) ...)
> i think that is called cumulative distribution.
> ->
> (1 1 4 5)
>
> create a random number up to the largest sum
> ->
> (random 5)
>
> return the first index of the list of cumulative sums that is greater
> than the random number.
> given the distribution above, you would see index 2 a lot, never index
> 1, and index 0 and 3 rarely.
>
> ~~~
> (use-modules ((srfi srfi-1) #:select (last)))
>
> (define (cusum a . b)
>   "calculate cumulative sums from the given numbers.
>    (a b c ...) -> (a (+ a b) (+ a b c) ...)"
>   (cons a (if (null? b) (list) (apply cusum (+ a (car b)) (cdr b)))))
>
> (define* (random-discrete-f probabilities #:optional (state
> *random-state*))
>   "(real ...) [random-state] -> procedure:{-> integer}"
>   (let* ((cuprob (apply cusum probabilities)) (sum (last cuprob)))
>     (lambda ()
>       (let ((deviate (random sum state)))
>         (let loop ((a 0) (b cuprob))
>           (if (null? b) a (if (< deviate (car b)) a (loop (+ 1 a) (cdr
> b)))))))))
>
> (define random* (random-discrete-f (list 1 0 3 1)
> (random-state-from-platform)))
> (display (random*))
> ~~~
>
>

Re: Normal distribution random numbers

[Mikael Djurfeldt]

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
> >
> >
>
>

Re: Normal distribution random numbers

[Zelphir Kaltstahl]

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
> <mailto: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
>     >
>     >
>

Re: Normal distribution random numbers

[Mikael Djurfeldt]

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
>> >
>> >
>>
>>

Re: Normal distribution random numbers

[Mikael Djurfeldt]

(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
>>> >
>>> >
>>>
>>>