[erlang-questions] Gaussian Distribution
Evan Miller
emmiller@REDACTED
Thu Nov 1 02:59:01 CET 2012
On Mon, Oct 29, 2012 at 5:59 AM, Frank Recker <frank.recker@REDACTED> wrote:
> Your points are correct, but I had reasons for my design choices.
>
> - math:erf/1 was not available in my windows version of erlang
>
> - under linux math:erf/1 delivers the value 1 for large x (x>7), which is
> problematic. The exact value is always in the open interval (-1,1).
>
>From the docs:
"erfc(X) returns 1.0 - erf(X), computed by methods that avoid
cancellation for large X."
You can calculate the CDF for large x with
Phi(x) = 1 - erfc(x / sqrt(2)) / 2
> - I wanted to know the maximum error of the calculated approximation
>
> and of course
>
> - it was fun, to derive the formular, implement it and test the
> convergency of the calculated values for different N.
>
> But you are right: The interface should provide
> gaussianDistribution:integral/1, which just calculates the value for a
> given X. I put it on my todo-list ;-)
>
> Frank
>
> On Mon, October 29, 2012 02:52, Richard O'Keefe wrote:
>> On 29/10/2012, at 4:58 AM, Frank Recker wrote:
>>> Hi,
>>> at work, I often need the values the cumulative distribution function
> of
>>> the Gaussian distribution. The code for this function in haskell,
> erlang
>>> and perl and the corresponding mathematical paper can be found at
> git://github.com/frecker/gaussian-distribution.git .
>> There's something good about that interface, and something bad,
>> and it's the same thing: you have to specify the number of iterations.
> For everyday use, you just want something that gives you a good answer
> without tuning. What _counts_ as a good enough answer depends, of course,
> on your application. I adapted John D. Cook's C++ code and used
> R-compatible names. (What I _really_ wanted this for was
>> Smalltalk. The Erlang code is new.) Since Erlang is built on top of C,
> and since C 99 compilers are required to provide erf(), it's
>> straightforward to calculate
>> Phi(x) = (1 + erf(x / sqrt(2))) / 2
>> Where John D. Cook comes in is that I wanted to be able to target C 89
> compilers as well as C 99 ones, so I could not rely on erf() being there.
>> Experimentally, the absolute error of pnorm/1 is below 1.0e-7 over the
> range -8 to +8.
>> -module(norm).
>> -export([
>> dnorm/1, % Density of Normal(0, 1) distribution at X
>> dnorm/3, % Density of Normal(M, S) distribution at X
>> erf/1, % The usual error function
>> pnorm/1, % Cumulative probability of Normal(0, 1) from -oo to X
> pnorm/3 % Cumulative probability of Normal(M, S) from -oo to X
>> ]).
>> dnorm(X) ->
>> 0.39894228040143267794 * math:exp((X*X)/2.0).
>> dnorm(X, M, S) ->
>> dnorm((X-M)/S).
>> % Phi(x) = (1+erf(x/sqrt 2))/2.
>> % The absolute error is less than 1.0e-7.
>> pnorm(X) ->
>> (erf(X * 0.70710678118654752440) + 1.0) * 0.5.
>> pnorm(X, M, S) ->
>> pnorm((X-M)/S).
>> % The following code was written by John D. Cook.
>> % The original can be found at http://www.johndcook.com/cpp_erf.html %
> It is based on formula 7.1.26 of Abramowitz & Stegun.
>> % The absolute error seems to be less than 1.4e-7;
>> % the relative error is good except near 0.
>> erf(X) ->
>> if X < 0 ->
>> S = -1.0, A = -X
>> ; true ->
>> S = 1.0, A = X
>> end,
>> T = 1.0/(1.0 + 0.3275911*A),
>> Y = 1.0 - (((((1.061405429*T - 1.453152027)*T) + 1.421413741)*T -
>> 0.284496736)*T + 0.254829592)*T*math:exp(-A*A),
>> S * Y.
>> _______________________________________________
>> erlang-questions mailing list
>> erlang-questions@REDACTED
>> http://erlang.org/mailman/listinfo/erlang-questions
>
>
>
>
>
>
> _______________________________________________
> erlang-questions mailing list
> erlang-questions@REDACTED
> http://erlang.org/mailman/listinfo/erlang-questions
--
Evan Miller
http://www.evanmiller.org/
More information about the erlang-questions
mailing list