Unit 10, lesson 124

Finite Differences and Derivative Approximation

How to estimate a derivative from function values or measured data with forward, backward and central differences, and why making the step size smaller helps only up to a point.

  • Theory part 247
  • Practice part 248
  • Numerical Differentiation

A sensor records the height of a drone every half second. Nobody hands you a formula, yet you want its velocity. The derivative is a limit of difference quotients, so the natural move is to stop before the limit and use a difference quotient with a small but fixed step. That simple idea works well, and it fails in two interesting ways, one from stopping early and one from the computer’s rounding.

Key terms

Forward, backward and central differences

Definition Let ff be defined on an open interval containing xx and let h>0h > 0 be a step size small enough that x−hx - h and x+hx + h lie in that interval. The forward, backward and central difference approximations to f′(x)f'(x) are

D+f(x)=f(x+h)−f(x)h,D−f(x)=f(x)−f(x−h)h,D0f(x)=f(x+h)−f(x−h)2h.D_+ f(x) = \frac{f(x + h) - f(x)}{h}, \qquad D_- f(x) = \frac{f(x) - f(x - h)}{h}, \qquad D_0 f(x) = \frac{f(x + h) - f(x - h)}{2h}.

If ff is differentiable at xx, each of them tends to f′(x)f'(x) as h→0+h \to 0^+. The central difference is the average of the other two.

In plain words These are slopes of secant lines. The forward one looks ahead, the backward one looks behind, and the central one looks both ways from the point.

Example Take f(x)=exf(x) = e^x at x=0x = 0, where f′(0)=1f'(0) = 1, and use h=0.1h = 0.1.

ApproximationValueError
D+f(0)=e0.1−10.1D_+ f(0) = \frac{e^{0.1} - 1}{0.1}1.0517091.0517090.05170.0517
D−f(0)=1−e−0.10.1D_- f(0) = \frac{1 - e^{-0.1}}{0.1}0.9516260.951626−0.0484-0.0484
D0f(0)=e0.1−e−0.10.2D_0 f(0) = \frac{e^{0.1} - e^{-0.1}}{0.2}1.0016681.0016680.00170.0017

The central difference is about thirty times more accurate here at the same cost of two function values.

Truncation error of a difference formula

Definition The truncation error of a difference formula is the difference between the formula, computed in exact arithmetic, and f′(x)f'(x). Taylor’s theorem gives it. If f′′f'' is continuous on [x,x+h][x, x + h], then for some ξ\xi in (x,x+h)(x, x + h)

D+f(x)−f′(x)=h2f′′(ξ).D_+ f(x) - f'(x) = \frac{h}{2} f''(\xi).

In the same way, if f′′f'' is continuous on [x−h,x][x - h, x], then for some η\eta in (x−h,x)(x - h, x)

D−f(x)−f′(x)=−h2f′′(η).D_- f(x) - f'(x) = -\frac{h}{2} f''(\eta).

If f′′′f''' is continuous on [x−h,x+h][x - h, x + h], then for some ζ\zeta in (x−h,x+h)(x - h, x + h)

D0f(x)−f′(x)=h26f′′′(ζ).D_0 f(x) - f'(x) = \frac{h^2}{6} f'''(\zeta).

So the one-sided formulas are first order, with error proportional to hh, and the central formula is second order, with error proportional to h2h^2.

In plain words Halving hh roughly halves the error of a one-sided difference but cuts the error of a central difference by about four.

Example For exe^x at 00 every derivative is close to 11 near 00, so the predicted errors are about h2\frac{h}{2} and h26\frac{h^2}{6}. With h=0.1h = 0.1 these are 0.050.05 and 0.001670.00167, matching the table above. With h=0.01h = 0.01 the actual errors are 0.005020.00502 and 0.00001670.0000167, again as predicted.

Cancellation

Definition Subtractive cancellation is the loss of significant digits when two nearly equal numbers, each known only to a fixed number of significant digits, are subtracted. If aa and bb agree in their first kk significant digits, the computed a−ba - b has about kk fewer correct significant digits than aa and bb had.

In plain words The digits the two numbers share cancel out, and only the less reliable trailing digits are left.

Example Work with 8 significant digits. Then e0.001e^{0.001} is stored as 1.00100051.0010005 and e0e^0 as 1.00000001.0000000. Their difference 0.00100050.0010005 has only 5 significant digits. Dividing by h=0.001h = 0.001 gives D+f(0)≈1.0005D_+ f(0) \approx 1.0005, and the three lost digits can never be recovered. A smaller hh makes f(x+h)f(x + h) and f(x)f(x) agree in even more digits, so the loss gets worse.

Roundoff in difference quotients

Definition Suppose each computed value of ff has a relative error of at most ε\varepsilon, about 1.1×10−161.1 \times 10^{-16} in double precision, and ∣f∣≤F\lvert f \rvert \le F near xx. Then the computed numerator of D+f(x)D_+ f(x) can be off by about 2εF2\varepsilon F, so the computed quotient carries a roundoff error of up to about

2εFh.\frac{2\varepsilon F}{h}.

For the central difference the numerator is divided by 2h2h, so the corresponding term is about εFh\frac{\varepsilon F}{h}.

In plain words Rounding errors in ff are tiny, but the formula divides them by hh. A small hh makes them large.

Example For exe^x at 00 with h=10−12h = 10^{-12}, the roundoff term is about

2⋅1.1×10−1610−12≈2×10−4.\frac{2 \cdot 1.1 \times 10^{-16}}{10^{-12}} \approx 2 \times 10^{-4}.

The computed forward difference is in fact off by about 9×10−59 \times 10^{-5}. That is worse than its error of 5×10−55 \times 10^{-5} at the much larger step h=10−4h = 10^{-4}, and more than ten thousand times worse than its error of 6×10−96 \times 10^{-9} at h=10−8h = 10^{-8}.

Choosing the step size

Definition For the forward difference, with ∣f′′∣≤M\lvert f'' \rvert \le M near xx, the total error is bounded by about

E(h)=Mh2+2εFh.E(h) = \frac{Mh}{2} + \frac{2\varepsilon F}{h}.

This is smallest at

h∗=2εFM,h^* = 2\sqrt{\frac{\varepsilon F}{M}},

which is of order ε≈10−8\sqrt{\varepsilon} \approx 10^{-8} when FF and MM are of moderate size. For the central difference, with ∣f′′′∣≤M3\lvert f''' \rvert \le M_3 near xx, the bound M3h26+εFh\frac{M_3 h^2}{6} + \frac{\varepsilon F}{h} is smallest at

h∗=(3εFM3)1/3,h^* = \left( \frac{3\varepsilon F}{M_3} \right)^{1/3},

which is of order ε1/3≈5×10−6\varepsilon^{1/3} \approx 5 \times 10^{-6}.

In plain words Truncation wants a small step, roundoff wants a big one. The best step sits in between, and below it the answers get worse, not better.

Example Here are the sizes of the actual errors of both formulas for exe^x at 00, computed in double precision.

hhForward errorCentral error
10−110^{-1}5.2×10−25.2 \times 10^{-2}1.7×10−31.7 \times 10^{-3}
10−210^{-2}5.0×10−35.0 \times 10^{-3}1.7×10−51.7 \times 10^{-5}
10−410^{-4}5.0×10−55.0 \times 10^{-5}1.7×10−91.7 \times 10^{-9}
10−510^{-5}5.0×10−65.0 \times 10^{-6}1.2×10−111.2 \times 10^{-11}
10−810^{-8}6.1×10−96.1 \times 10^{-9}6.1×10−96.1 \times 10^{-9}
10−1210^{-12}8.9×10−58.9 \times 10^{-5}3.3×10−53.3 \times 10^{-5}

The forward difference is best near h=10−8h = 10^{-8} and the central one near h=10−5h = 10^{-5}, just as the formulas for h∗h^* predict.

Derivatives from sampled data

Definition Suppose we know only the values yk=f(tk)y_k = f(t_k) at equally spaced points tk=t0+kht_k = t_0 + kh for k=0,1,…,Nk = 0, 1, \dots, N. At an interior point the derivative is estimated by the central difference, and at the first and last points by a forward and a backward difference,

f′(tk)≈yk+1−yk−12h,f′(t0)≈y1−y0h,f′(tN)≈yN−yN−1h.f'(t_k) \approx \frac{y_{k+1} - y_{k-1}}{2h}, \qquad f'(t_0) \approx \frac{y_1 - y_0}{h}, \qquad f'(t_N) \approx \frac{y_N - y_{N-1}}{h}.

The second derivative at an interior point is estimated by the second central difference

f′′(tk)≈yk+1−2yk+yk−1h2.f''(t_k) \approx \frac{y_{k+1} - 2y_k + y_{k-1}}{h^2}.

If f(4)f^{(4)} is continuous on [tk−1,tk+1][t_{k-1}, t_{k+1}], this estimate minus f′′(tk)f''(t_k) equals h212f(4)(ξ)\frac{h^2}{12} f^{(4)}(\xi) for some ξ\xi in (tk−1,tk+1)(t_{k-1}, t_{k+1}). Measurement errors in the data are divided by hh or h2h^2, so they act like roundoff.

In plain words With a table of numbers, use the neighbors on both sides when you have them and a one-sided difference at the ends.

Example A drone’s height in meters is recorded every half second.

tt000.50.51.01.01.51.52.02.0
height002.32.34.14.15.45.46.26.2

At t=1.0t = 1.0 the velocity is about 5.4−2.31.0=3.1\frac{5.4 - 2.3}{1.0} = 3.1 meters per second. At t=0t = 0 the forward difference gives 2.3−00.5=4.6\frac{2.3 - 0}{0.5} = 4.6, and at t=2.0t = 2.0 the backward difference gives 6.2−5.40.5=1.6\frac{6.2 - 5.4}{0.5} = 1.6. The acceleration at t=1.0t = 1.0 is about

5.4−2⋅4.1+2.30.25=−2.0\frac{5.4 - 2 \cdot 4.1 + 2.3}{0.25} = -2.0

meters per second squared, so the drone is slowing its climb.

Check yourself

  1. Compute the central difference for f(x)=x3f(x) = x^3 at x=1x = 1 with h=0.1h = 0.1, and compare it with the error formula.
Show the answer

The central difference is

1.13−0.930.2=1.331−0.7290.2=3.01.\frac{1.1^3 - 0.9^3}{0.2} = \frac{1.331 - 0.729}{0.2} = 3.01.

The exact derivative is 33, so the error is 0.010.01. The formula gives h26f′′′=0.016⋅6=0.01\frac{h^2}{6} f''' = \frac{0.01}{6} \cdot 6 = 0.01, exact here because f′′′=6f''' = 6 is constant.

  1. A forward difference gives an error of 0.0040.004 with h=0.02h = 0.02. Roughly what error do you expect with h=0.01h = 0.01, assuming roundoff is still negligible?
Show the answer

The forward difference is first order, so halving hh roughly halves the error, giving about 0.0020.002.

  1. Why is h=10−15h = 10^{-15} a bad step for a forward difference in double precision?
Show the answer

The roundoff term is about

2εFh≈2.2×10−16F10−15≈0.2F,\frac{2\varepsilon F}{h} \approx \frac{2.2 \times 10^{-16} F}{10^{-15}} \approx 0.2F,

so the result can be wrong in the first decimal place. The best step is near 10−810^{-8}.

How this lesson connects

It also uses finite precision.