Relay-Version: version B 2.10 5/3/83; site utzoo.UUCP
Posting-Version: version B 2.10.3 4.3bsd-beta 6/6/85; site ucbvax.BERKELEY.EDU
Path: utzoo!watmath!clyde!burl!ulysses!ucbvax!brahms!weemba
From: weemba@brahms.BERKELEY.EDU (Matthew P. Wiener)
Newsgroups: net.math
Subject: Re: Series for normal distribution function
Message-ID: <12180@ucbvax.BERKELEY.EDU>
Date: Wed, 5-Mar-86 00:41:42 EST
Article-I.D.: ucbvax.12180
Posted: Wed Mar  5 00:41:42 1986
Date-Received: Fri, 7-Mar-86 04:10:41 EST
References: <443@cubsvax.UUCP>
Sender: usenet@ucbvax.BERKELEY.EDU
Reply-To: weemba@brahms.UUCP (Matthew P. Wiener)
Organization: University of California, Berkeley
Lines: 76

In article <443@cubsvax.UUCP> winston@cubsvax.UUCP (Ken Winston) writes:
>The standard normal distribution function is
>
>	N(t) = 1/sqrt(2*pi) * integral(-infinity to t) [exp(-z^2/2)]*dz.
>
>I have seen the following approximation a few different times now as a way
>of finding N(t) on a computer. It seems to be accurate in the usual range, for
>t<=0:
>
>    N(t) ~= x*exp(-t^2/2)*(.436184 - .120167*x + .937298*x^2)/sqrt(2*pi),
>	where x = 1/(1+.3327*|t|)
>
>For t>0 use N(t) = 1-N(-t).
>
>Does anybody
>know where this comes from? What is the general method that gives rise to this
>series?

This particular representation is due to C Hastings _Approximations for
Digital Computers_.

First, get the function to look as geometrically like a polynomial (or other
easy to compute function as possible).  Asymptotes should be factored out
first, for example.  Second, economize the approximation.  This consists of
fudging the coefficients so that the error is spread around more uniformly,
usually by replacing the highest terms by known lower degree approximants.
There is a whole theory on how to do this, starting with the Chebyshev
polynomials.

In the particular case,	one starts with the asymptotic expansion

                       2
               2     -t    1             1       1.3       1.3.5
    N(t) = -------- e     ---  [ 1  -  -----  + ------  -  -----  + ... ]
           sqrt(pi)       2 t              2        2 2       2 3
                                       (2 t )   (2 t )    (2 t )

easily derived by repeatedly integrating by parts.  It diverges, but the
error is less than the first omitted term.  It can be used as is for x>5
or so.  (Depends on how much accuracy you need.)  A Taylor series, either
for N(t) or N(t)/exp(-t^2) converges too slowly for t>.7 or so.

Where Hastings came up with the change of variables x=1/(1 + c t) is beyond
me.  In terms of the geometry of the error function, the resulting change
gives one a much more polynomic geometry, and therefore a much better chance
at getting good approximation.

One then has to economize the coefficients.  Consider T_4(x)=8x^4-8x^2+1.
Pretty innocuos looking, isn't it?  Actually, it satisfies the important
property that |T_4(x)|<=1 for |x|<=1.  Thus x^4 is approximately x^2-1/8,
with error at most 1/8, on the interval [-1,1].  In general, T_n(x) is a
polynomial with leading term 2^(n-1)*x^n, and can be used to simulate x^n
with a lower degree polynomial by committing an error of 2^(1-n).  For
other situations and intervals, similar methods have been found.

Suppose for example, you want a good degree 7 approximation to sine(x)
on some interval.  One takes the Taylor series out to degree 11, commit
a Taylor series error by stopping, and then commit two further errors by
economizing out the x^9 and x^11 terms.  These three errors together are
much smaller than the degree 7 Taylor series error, and are more uniformly
distributed.

A very readable account is in F Acton _Numerical Methods That Work_,
chapters 1 and 12.  Not only is it a good book (although be warned that
some of his *particular* comments about certain methods are no longer
valid--old methods get revived a lot in numerical analyis), it is the
only numerical analysis book that makes it clear that THINKING is the
most important technique for solving numerical problems on computers.
It is also the FUNNIEST math book ever written.

"A favorite form of lunacy among engineers ..."
"If your author were writing another typical book on Numerical Analysis ..."
"reducing ... costly n-dimensional wanderings through foggy Numberland."
"... you end up by getting angry and throwing the guy out of your office."

ucbvax!brahms!weemba	Matthew P Wiener/UCB Math Dept/Berkeley CA 94720