Taylor Series Approximation for the Normal Distribution

An incredibly accurate and efficient method of approximating the integral of the standard normal distribution.

For a while now, I’ve been thinking about finding better ways of estimating the integral of the probability density function (for the standard normal distribution), as shown below.

f(z)=12πez22

I started with analyzing different types of Riemann approximations, but soon realized that it could be quite inaccuate, since by using about a ten million subintervals, the approximate would only be accurate to around 6-7 decimal places. In addition, it would also be resource-intensive for an ordinary computer to execute.

More recently, I thought about perhaps instead using a power series in order to approximate the integral. I realized that it would make sense to begin with using a Taylor series for ex centered about 0, since the probability function is e-based and is symmetric about the y-axis. Then, I proceeded level by level in order to create a Taylor series for the probability function, as shown below.

ez=1+z+z22!+z33!+z44!+...ez22=1+(z22)+(z22)22!+(z22)33!+(z22)44!+...12πez22=12π(1+(z22)+(z22)22!+(z22)33!+(z22)44!+...)=12π(1z221!+z442!z683!+z8164!...)

Now that we have the Taylor series for the probability function, we can now easily integrate this (what is essentially a) polynomial.

12πez22dz=12π(z110!z3321!+z5542!z7783!+z99164!...)=12πn=0(1)nz2n+1(2n+1)2nn!

Also, I found the interval of convergence for this infinite series in order to check that the approximation of the integral would apply for zR, which is the domain of the probability function. The following limit confirmed this, since as n, the limit goes to 0, whose absolute value is always less than 1. Therefore, this statement is always true, and thus the interval of convergence is z(,).

limn|(1)n+1z2n+3(2n+3)2n+1(n+1)!(2n+1)2nn!(1)nz2n+1|

Hypothetically, if we were to find an exact value for this integral (assuming we take a definite integral), then we would need to evaluate the infinite series. We could find the exact sum for an infinite series if it was geometric or telescoping, but it doesn’t appear to be the case. So, I thought of further decomposing the Taylor series for the integral of the probability function into roots by taking each term of the series and setting it equal to dndzznn!, where n is the indexing value of the series. Obviously, the even values of n would correspond to a dndz value of 0, since that is the only possible numeric component of the term that could “vary.”

Next, I wanted to analyze the progression of the absolute value of dndz, since I already knew that we were dealing with an alternating series. Here below is a table that displays the values of the absolute value of dndz for the first 6 odd values of n.

n
1
3
5
7
9
11
|dndz||z=0
1
1
3
15
105
945

Ignoring the first term, the absolute value of dndz seems to first be multiplied by 3, then by 5, then by 7, and so on. I tried coming up with a sequence generator for the absolute value of dndz, but it turned out to be quite challenging. After collaborating with a few students in my math class, we determined the sequence generator to be (2n)!n!2n, assuming that n starts at 1.

After substituting this into the sigma notation for the Taylor series for the integral of the probability function, here is the result:

12πez22dz=12π(z+n=1(1)nz2n+1n!2n(2n+1))

The best part about using this series is that when we set the upper bound of the sum to just 5, we get 12 decimal places of accuracy. This makes for a much more efficient computation for a computer to handle, instead of using a Riemann sum of ten million subintervals. My next step is to find a way to find the exact value of the infinite series - no approximations. I’ll save that for another post :)