# How to calculate the integral in normal distribution?

The factory is making products with this normal distribution: $\mathcal{N}(0, 25)$. What should be the maximum error accepted with the probability of 0.90? [Result is 8.225 millimetre]

How will I calculate it? How to integrate: $\exp\left(- \frac{x^2}{2} \right)$ ?

-
What have you tried? Do you know how to set up the problem? –  Sasha May 14 '12 at 17:31
Are you supposed to use a table? Software? A numerical quadrature method like Simpson's rule? –  Ｊ. Ｍ. May 14 '12 at 17:33
Can a table help me in that? I would be very happy if I knew how to calculate this integral... Then I would be able to calculate Phi(x/5) - Phi(-x/5) = 0.90 and then find x that solve the equation... I am not allowed to use any software like Matlab, Mathematica or any other application, I can only use a calculator. –  user1111261 May 14 '12 at 17:42
There should be tables for the CDF of the (standard) normal distribution in the usual statistics textbooks... –  Ｊ. Ｍ. May 14 '12 at 18:07
Definite integrals of that function are found by numerical methods rather than by finding a closed-form antiderivative. In exercises of this kind usually one gets the value of the integral either from software or from a table in the back of the book. –  Michael Hardy May 14 '12 at 18:20

There are several methods to approach this, but I am going to use one that meets your requirement (clarified in a comment) that one must forego the use of computational engines like Mathematica, instead opting for a calculator. Additionally, I feel that using a table of normal distribution values is cheating, so I will be foregoing their use as well.

First, we need the equation for $\mathcal{N}(0,25)$, which, by definition, is: \begin{align*} f(x) &= \mathcal{N}(\mu,\sigma ^2)\\ &= \mathcal{N}(0,25)\\ &= \frac{1}{\sigma\sqrt{2\pi}}\,e^{ -\frac{(x-\mu)^2}{2\sigma^2} }\\ &= \frac{1}{5 \sqrt{2 \pi }}\,e^{-\frac{x^2}{50}} \end{align*} Now, we simply need to integrate this from $-x$ to $x$, set it equal to $.90$, and solve for $x$ (our answer): $$F(x) = \frac{1}{5 \sqrt{2 \pi }}\int_{-x}^x e^{-\frac{x^2}{50}} \,\mathrm{d}x=0.9$$ However, we run into problems when we realize that this isn't a simple integral to take! Since we already know the answer, we can actually exploit this to our advantage by finding an elementary function that estimates $F(x)$ with an adequately small margin of error after the integral. One simple function you can use to estimate $F(x)$ is a Taylor series.

We first need to find a Taylor series for $f(x)$ using the formula for a Taylor series: $$\sum\limits_{n=0}^{\infty } \frac {f^{(n)}(a)}{n!} \, (x-a)^{n}$$ One can easily recognize the pattern for our function when $a=0$ (the center of our bell curve) to generate this series: $$f(x) = \frac{1}{5 \sqrt{2 \pi}}\sum\limits_{k=0}^{\infty} \frac{(-1)^k x^{2k}}{50^k k!}$$ We can now trivially take the integral of this series where we otherwise wouldn't have been able to: $$F(x) = \frac{1}{5 \sqrt{2 \pi}}\sum\limits_{k=0}^{\infty} \frac{(-1)^k x^{2 k+1}}{50^k k! (2 k+1)}$$ Since we can't use our pencil and paper to evaluate the series all the way to infinity, we need to figure out how many terms we need to go out to get an acceptably accurate answer. This can be accomplished by figuring out the error caused by not going out to infinity before we actually solve the equation. By plugging in our already-known value of $8.225$ for $x$ in the $k$th term of the series and evaluating just the series part from $k$ to $\infty$, we get the error that the series will have from $0$ to $k-1$. Since we, again, can't go all the way out to infinity, we can get a slight, but adequate, underestimate of the error by simply evaluating the series from $k$ to $k$ since each subsequent term in the series is substantially smaller.

I started by plugging in $8.225$ for $x$ when $k=7$ and got this (do this on your calculator): $$\sum\limits_{k=7}^{7 } \frac{(-1)^k (8.225)^{2 k+1}}{50^k k! (2 k+1)} = -\frac{(8.225)^{15}}{59062500000000000} \approx -0.000903081$$ Since our integrated series is the equivalent of $F$ in $F(x) - F(-x) = \int_{-x}^x f(x)\,\mathrm{d}x$ where $f(x)$ is the equation for $\mathcal{N}(0,25)$, we need to multiply our answer by $2$ to get a rough estimate of the error our final answer will have: $-0.000903081 \times 2 = -0.001806162$.

Since we are looking for a three decimal accuracy, we should proceed and try the series when $k=8$: $$\sum\limits_{k=8}^{8 } \frac{(-1)^k (8.225)^{2 k+1}}{50^k k! (2 k+1)} = \frac{(8.225)^{17}}{26775000000000000000} \approx 0.000134766$$ Finally, $0.000134766 \times 2 = 0.000269532$. This is of adequate accuracy for our final calculation; since the error calculated when $k=8$ is acceptable, we will evaluate the terms of the series from $k=0$ to $7$ ($k$ used in error calculation minus $1$): $$F(x) \approx \frac{1}{5 \sqrt{2 \pi}}\sum\limits_{k=0}^{7} \frac{(-1)^k x^{2 k+1}}{50^k k! (2 k+1)}$$ This comes out to be: $$\frac{x-\frac{x^3}{150}+\frac{x^5}{25000}-\frac{x^7}{5250000}+\frac{x^9}{1350000000}-\frac{x^{11}}{412500000000}+\frac{x^{13}}{146250000000000}-\frac{x^{15}}{59062500000000000}}{5 \sqrt{2 \pi }}$$ Since this is our estimate of $F(x)$, and this is an odd function (a function is odd if $-f(x) = f(-x)$), we simply need to multiply this function by $2$ to get $F(x) - F(-x)$ (what we need to set to $0.90$ and solve): $$\frac{1}{5} \sqrt{\frac{2}{\pi}} (x-\frac{x^3}{150}+\frac{x^5}{25000}-\frac{x^7}{5250000}+\frac{x^9}{1350000000}-\frac{x^{11}}{412500000000}+\frac{x^{13}}{146250000000000}-\frac{x^{15}}{59062500000000000})$$

Now, by setting this equal to $.90$, rearranging the equation as a polynomial, and using a method of our choice to solve polynomials on a calculator (like Newton's method to converge on the answer), we find that the relevant solution is: $$x \approx 8.22473 \approx 8.225$$

-

There are several ways to compute the cumulative normal distribution.

Simple Series Integration

First of all, we can start with $$e^{-x^2/2}=1-\frac{x^2}{2^1\cdot1!}+\frac{x^4}{2^2\cdot2!}-\frac{x^6}{2^3\cdot3!}+\dots$$ and integrate to get \begin{align} \frac1{\sqrt{2\pi}}\int_0^xe^{-t^2/2}\,\mathrm{d}t &=\frac1{\sqrt{2\pi}}\left(x-\frac{x^3}{3\cdot2^1\cdot1!}+\frac{x^5}{5\cdot2^2\cdot2!}-\frac{x^7}{7\cdot2^3\cdot3!}+\dots\right)\\ &=\frac1{\sqrt{2\pi}}\sum_{k=0}^\infty(-1)^k\frac{x^{2k+1}}{(2k+1)2^kk!} \end{align}

Unilateral Power Series

To get a series with no sign changes, write $$\Omega(x)=e^{x^2/2}\int_0^xe^{-t^2/2}\,\mathrm{d}t$$ and note that $$\Omega'(x)=1+x\Omega(x)$$ which leads to the following recursion for the coefficients $$a_{n+2}=\frac{a_n}{n+2}$$ Because $\Omega(0)=0$ and $\Omega'(0)=1$, we get \begin{align} \frac1{\sqrt{2\pi}}\int_0^xe^{-t^2/2}\,\mathrm{d}t &=\frac1{\sqrt{2\pi}}e^{-x^2/2}\Omega(x)\\ &=\frac1{\sqrt{2\pi}}e^{-x^2/2}\sum_{k=0}^\infty \frac{x^{2k+1}}{(2k+1)!!} \end{align}

Asymptotic Expansion

The terms of the convergent series above can get quite big as $x$ gets big, so let's compute an asymptotic expansion for the tail.

Consider \begin{align} \int_x^\infty e^{-t^2/2}\,\mathrm{d}t &=\int_0^\infty e^{-(t+x)^2/2}\,\mathrm{d}t\\ &=e^{-x^2/2}\int_0^\infty e^{-xt-t^2/2}\,\mathrm{d}t\\ &=\frac{e^{-x^2/2}}{x}\int_0^\infty e^{-u-u^2/(2x^2)}\,\mathrm{d}u\\ &=\frac{e^{-x^2/2}}{x}\int_0^\infty e^{-v}u'\,\mathrm{d}v\\ \end{align} where $v=u+\dfrac{u^2}{2x^2}$; that is, $u=x^2(\sqrt{1+2v/x^2}-1)$. Using the binomial theorem, to expand $\sqrt{1+2v/x^2}$, we get $$\frac1{\sqrt{2\pi}}\int_x^\infty e^{-t^2/2}\,\mathrm{d}t \sim\frac1{\sqrt{2\pi}}e^{-x^2/2}\sum_{k=0}^\infty\frac{(-1)^k(2k-1)!!}{x^{2k+1}}$$ where $(-1)!!=1$.

Note that this series is not convergent, but asymptotic.

-
A good summary of some of the simplest approximations. –  Ｊ. Ｍ. Aug 2 '12 at 12:32
Could you please change your $\Phi(x)$ notation to some other notation? $\Phi(x)$ is commonly used to denote the cumulative probability distribution function of the standard normal random variable, and your choice might create unnecessary confusion. Thanks. –  Dilip Sarwate Aug 2 '12 at 12:39
@DilipSarwate: Does $\Omega$ work better? Since I defined what I meant by $\Phi$, and it only appeared as an auxiliary function in a small section, I didn't think it would much matter what I called it. –  robjohn Aug 2 '12 at 13:30

There's tables that usually accompany probability books that give you the solution over a certain interval, but the integral of the normal distribution (the Gaussian function) is known as the error function $$\frac{1}{\sqrt{2\pi}}\int{e^{-\frac{x^2}{2}}}dx=\frac{1}{2}\mathrm{erf}(\frac{x}{\sqrt{2}})+C$$

The erf function is equal to -1 at negative infinity, so the CDF of the standard normal distribution (σ = 1, μ = 0) is: $$\mathrm{\Phi}(a)=\frac{1}{2}\mathrm{erf}(\frac{a}{\sqrt{2}})+\frac{1}2$$

-
"the integral of the normal distribution (the Gaussian function) is known as the error function (sqrt(pi)/2)*erfi(x)" Well, that depends on what you call the normal distribution, what integral you are talking about, and what you mean by erfi(x). If by integral you mean the cumulative distribution function $\Phi(x)$ mentioned in the comments by the OP, then your assertion is incorrect. –  Dilip Sarwate Jul 15 '12 at 20:17
I clarified my question based on your feedback. Thanks –  Ataraxia Jul 30 '12 at 19:07
You still need limits on the integral, e.g. $\displaystyle \int_0^x$ instead of just $\displaystyle \int$, and you really should use a different variable of integration than the one used in the limit. –  Dilip Sarwate Jul 30 '12 at 19:54
Other than the missing constant of integration, why are limits necessary? It goes without saying that if you're trying to find a CDF, you need to add limits and evaluate the definite integral. In the second equation you'll notice that I used "a" as the (upper) limit variable. And the question is talking about the CDF, so the lower limit is negative infinity. –  Ataraxia Jul 30 '12 at 20:04
The point you are missing is that the definition of $\text{erf}(x)$ is $$\text{erf}(x) = \frac{2}{\sqrt{\pi}}\int_0^x \exp(-t^2)\,\mathrm dt.$$ Writing that $$\frac{1}{\sqrt{2\pi}}\int{e^{-\frac{x^2}{2}}}dx=\frac{1}{2}\mathrm{erf}(\frac{‌​x}{\sqrt{2}})+C$$ or equivalently, $$\frac{2}{\sqrt{\pi}}\int \exp(-x^2)\,\mathrm dx = \text{erf}(x)+C$$ as you have done is not useful since $C$ is arbitrary. What, for example, is $\text{erf}(0)$? $0$ according to the definition but $$\frac{2}{\sqrt{\pi}}\int \exp(-x^2)\,\mathrm dx - C$$ according to you. –  Dilip Sarwate Jul 31 '12 at 1:11