I’ve seen several people ask lately how to compute the distribution (CDF) function for a standard normal random variable, often denoted Φ(x). They want to know how to compute it in Java, or Python, or C++, etc. Every language has its own standard libraries, and in general I recommend using standard libraries. However, sometimes you want to minimize dependencies. Or maybe you want more transparency than your library allows. The code given here is in Python, but it is so compact that it could easily be ported to any other language.

I just posted Python code for computing the error function, erf(x). The normal density Φ(x) is a simple transformation of erf(x). Given code for erf(x), here’s code for Φ(x).

def phi(x): return 0.5*( 1.0 + erf(x/math.sqrt(2)) )

After deriving the transformations between erf(x) and Φ(x) several times, including their complements and inverses, I wrote them down to save. See the PDF file Relating Φ and erf.

See also stand alone code for computing the inverse of the standard normal CDF.

If high accuracy isn’t essential, there’s also Polya’s approximation to the standard Gaussian c.d.f. which Aludaat and Alotdat have recently improved, at least in a maximum error sense (Applied Mathematical Sciences, Vol. 2, 2008, no. 9, 425 – 429). Polya’s is still of interest, at least to me, because its peak deviations are farther out in number of standard errors than Aludaat & Alodat’s. Both have error which approaches zero as the number of standard errors increases in magnitude beyond two. Polya’s maximum deviation from Gaussian c.d.f. is a tad over 0.003 at just over 1.5 standard errors.

The Polya approximation is: 0.5 +- 0.5 sqrt(1-exp(-2z^2/pi))

The inverse is available by algebra.