10.3.4 The Error Function

In Chapter 9 we introduced the error function:

erf ( t ) = 2 π ∫ 0 t     e - s 2 d s .

According to the Fundamental Theorem of Calculus, `text[erf (] t text[)]` is a constant multiple of a particular antiderivative of the function `e^(-t^2)`. We have seen that `text[erf (] t text[)]` can be used to calculate values of the standard normal distribution function. If you ask a computer algebra system to integrate `e^(-t^2text[/]2)`, it will probably respond with an answer that involves the error function.

Also in Chapter 9 we cited computer evidence that the integral definition is not an efficient way to evaluate `text[erf (] t text[)]` — nor is it the way the computer program itself evaluates this function. This suggests the possibility of a formula for `text[erf (] t text[)]` that requires much less arithmetic than would be required for numerical evaluation of an integral. We raised a similar issue with respect to the elementary functions `e^x`, `sin x`, `cos x`, `ln x`, and `tan^(-1) x` (and all the other functions in your computer algebra system or on your calculator): How can these functions be calculated quickly by a device that does nothing but arithmetic? For the five functions just mentioned, we have an answer of sorts in their Taylor series formulas, at least for values of `x` that are not too far from zero. Our answer is, "Add up enough terms of the series to get a partial sum that is acceptably close to the total sum."

However, finding a Taylor series for `text[erf (] t text[)]` is not a trivial task. We know its first derivative is `2e^(-t^2)`/`sqrt(pi)`, but subsequent derivatives get increasingly difficult to calculate. The next three (computed with help from a computer algebra system) are

- 4 ⁢ t ⁢ e - t 2 π , - 4 ⁢   e - t 2 + 8 ⁢ t 2 ⁢ ⁢ e - t 2 π , and 24 ⁢   t ⁢   e - t 2 - 16 ⁢ ⁢ t 3 ⁢ ⁢ e - t 2 π .

All we really need to know is the value of each of these functions at `t=0` and enough of a pattern to figure out as many Taylor coefficients as we need for any given calculation.

Checkpoint 2

It appears from Checkpoint 2 that direct calculation of Taylor coefficients is a lot of work for relatively little payoff. From the function and four derivatives, we have only two nonzero terms of the Taylor series and not enough information to guess the pattern of subsequent coefficients. We need to work smarter, not harder!

We will see in Section 10.6 that the exponential series,

e x ⁢ = 1 + x + x 2 2 + x 3 6 + ⋅⋅⋅ + x k k ! + ⋅⋅⋅ ,

converges for all real numbers `x`. This means that any real-number expression can be substituted for `x` and the resulting series will still converge. In particular, if we substitute `x=-t^2`, we find

e - t 2	= 1 + ( - t 2 ) + ( - t 2 ) 2 2 + ( - t 2 ) 3 6 + ⋅⋅⋅ + ( - t 2 ) k k ! + ⋅⋅⋅
	= 1 - t 2 + t 4 2 - t 6 6 + ⋅⋅⋅ + ( - 1 ) k   t 2 k k ! + ⋅⋅⋅ .

Now we integrate both sides from, say, `0` to `x`:

∫ 0 x e - t 2 d t	= ∫ 0 x 1 d t - ∫ 0 x t 2 d t + ∫ 0 x t 4 2 d t - ∫ 0 x t 6 6 d t + ⋅⋅⋅ + ∫ 0 x ( - 1 ) k   t 2 k k ! d t + ⋅⋅⋅
	= t | 0 x - t 3 3 | 0 x + t 5 2 ! ⋅ 5 | 0 x - t 7 3 ! ⋅ 7 | 0 x + ⋅⋅⋅ + ( - 1 ) k   t 2 k + 1 k ! ⋅ ( 2 k + 1 ) | 0 x + ⋅⋅⋅
	= x - x 3 3 + x 5 2 ! ⋅ 5 - x 7 3 ! ⋅ 7 + ⋅⋅⋅ + ( - 1 ) k   x 2 k + 1 k ! ⋅ ( 2 k + 1 ) + ⋅⋅⋅ .

Thus `text[erf (] t text[)]` can be evaluated by substituting a value for `x` in this series, calculating enough terms to achieve the desired accuracy, and multiplying the result by `2`/`sqrt(pi)`.

Is this correct? A computer algebra system tells us that `text[erf (]1text[)]=0.84270079`, to the eight decimal places shown. When we divide this by `2text[/]sqrt(pi)`, we get `0.74682413`, which should be the result of substituting `x=1` in the series we just calculated. Let's make the substitution and start calculating partial sums:

`1`	`=1`
`1-1/3`	`=0.66666667`
`1-1/3+1/10`	`=0.76666667`
`1-1/3+1/10-1/42`	`=0.74285714`
`1-1/3+1/10-1/42+1/216`	`=0.74748667`
`1-1/3+1/10-1/42+1/216-1/1320`	`=0.74672920`
`1-1/3+1/10-1/42+1/216-1/1320+1/9360`	`=0.74683603`
`1-1/3+1/10-1/42+1/216-1/1320+1/9360-1/75600`	`=0.74682281`
`1-1/3+1/10-1/42+1/216-1/1320+1/9360-1/75600+1/685440`	`=0.74682427`
`1-1/3+1/10-1/42+1/216-1/1320+1/9360-1/75600+1/685440-1/6894720`	`=0.74682412`

Thus with only ten terms of the series we have seven-place accuracy in the number we expected to calculate. It should not be hard to see now how values of `text[erf (] t text[)]` — or of the normal distribution function — can be calculated by simple arithmetic. This is much less arithmetic than would be required for a numerical integration involving many values of an exponential function.

Checkpoint 3