### Numerical differentiation: Difference formulas for the first derivative

### General difference formulas

The difference formulas presented earlier have been introduced via Taylor's theorem. In some applications it is necessary to use difference formulas with a higher accuracy (think of rounding errors) or with grid points that are not equidistantly distributed. Here we will present a method to derive new difference formulas. The 3-point and 5-point central difference formulas are special cases of derivatives computed through a **Savitsky-Golay filter**.

The general approach for \(n+1\) different grid points \(x_0, x_1, \ldots x_n\), which differ from each other by a multiple of the step size \(h\), is to find coefficients \(c_0, c_1, \ldots, c_n\) such that the expression \[Q(h)\stackrel{\mathrm{def}}{=}\sum_{i=0}^n c_if(x_i)\] approximates the derivative \(f'(x_0)\) with the largest possible order of the truncation error. We find the coefficients \(c_0, c_1, \ldots, c_n\) by developing each \(f(x_i)\) into a Taylor polynomial about \(x_0\) of sufficiently high degree and then by constructing and solving equations for the unknown \(c_i\) 's.

The description above is rather cryptic, but examples make this systematic approach clear.

The 3-point central difference formula We consider 3 grid points \(x_0-h\), \(x_0\) and \(x_0+h\). We define the formula \[Q(h)=A\cdot f(x_0-h)+B\cdot f(x_0)+C\cdot f(x_0+h)\] and determine the Taylor approximation of \(Q(h)\). We know that \[\begin{aligned} f(x_0-h) &= f(x_0)-f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2+O(h^3)\\ f(x_0+h) &= f(x_0)+f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2+O(h^3)\end{aligned}\] and therefore \[\begin{aligned}Q(h) &= A\cdot ( f(x_0)-f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2)+B\cdot f(x_0)\\ &\phantom{=}{}+C\cdot (f(x_0)+f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2) + O(h^3)\\ \\ &= (A+B+C)f(x_0)+(C-A)f'(x_0)h+\tfrac{1}{2}(A+C)f''(x_0)h^2+O(h^3)\end{aligned}\] The equations that \(A, B, C\) must satisfy in order for \(Q(h)\) to maximise the order of the error are \[A+B+C = 0,\quad (C-A)h = 1,\quad A+C = 0\] These equations are easy to solve: \[A=-\frac{1}{2h}, B=0, C=\frac{1}{2h}\] So we have found \[Q(h)=\frac{f(x_0+h)-f(x_0-h)}{2h}\] and \[f'(x_0)=\frac{f(x_0+h)-f(x_0-h)}{2h}+O(h^2)\] This is equal to the central difference formula.

One-sided 3-point difference formulas The left-sided 3-point difference formula: \[f'(x_0)\approx\frac{-3f(x_0)+4f(x_0+h)-f(x_0+2h)}{2h}\] The right-sided 3-point difference formula: \[f'(x_0)\approx\frac{3f(x_0)-4f(x_0-h)+f(x_0-2h)}{2h}\] One-sided 3-point difference formulas are useful for approximating derivatives at the edges of a finite discrete signal.

We take 3 grid points \(x_0\), \(x_0+h\) and \(x_0+2h\). We define the formula \[Q(h)=A\cdot f(x_0)+B\cdot f(x_0+h)+C\cdot f(x_0+2h)\] and determine the Taylor approximation of \(Q(h)\). We know that \[\begin{aligned} f(x_0+h) &= f(x_0)+f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2+O(h^3)\\ f(x_0+2h) &= f(x_0)+2f'(x_0)h+\tfrac{1}{2}f''(x_0)(2h)^2+O(h^3)\end{aligned}\] and therefore \[\begin{aligned}Q(h) &=A\cdot f(x_0) + B\cdot ( f(x_0)+f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2)\\ &\phantom{=}{}+C\cdot (f(x_0)+2f'(x_0)h+\tfrac{4}{2}f''(x_0)h^2) + O(h^3)\\ \\ &= (A+B+C)f(x_0)+(B+2C)f'(x_0)h+\tfrac{1}{2}(B+4C)f''(x_0)h^2+O(h^3)\end{aligned}\] The equations that \(A, B, C\) must satisfy in order for \(Q(h)\) to maximise the order of the error are \[A+B+C = 0,\quad (B+2C)h = 1,\quad B+4C = 0\] These equations are easy to solve: \[A=-\frac{3}{2h}, B=\frac{2}{h}, C=-\frac{1}{2h}\] So we have found \[Q(h)=\frac{-3f(x_0)+4f(x_0+h)-f(x_0+2h)}{2h}\] and \[f'(x_0)=\frac{-3f(x_0)+4f(x_0+h)-f(x_0+2h)}{2h}+O(h^2)\] This is a left-sided 3-point difference formula.

Similarly, you can also find the right-sided 3-point difference formula: \[f'(x_0)=\frac{3f(x_0)-4f(x_0-h)+f(x_0-2h)}{2h}+O(h^2)\]

The 5-point central difference formula \[f'(x_0)\approx\frac{f(x_0-2h)-8f(x_0-h)+8f(x_0+h)-f(x_0+2h)}{12h}\] The truncation error is generally smaller than with the 3-point central difference formula.

We consider 5 grid points \(x_0-2h\), \(x_0-h\), \(x_0\), \(x_0+h\) and \(x_0+2h\). We define the formula \[\begin{aligned}Q(h)&=A\cdot f(x_0-2h)+B\cdot f(x_0-h)+C\cdot f(x_0)\\ &\phantom{=}{}+ D\cdot f(x_0+h)+E\cdot f(x_0+2h)\end{aligned} \] and determine the Taylor approximation of \(Q(h)\). We know that \[\begin{aligned} f(x_0-2h) &= f(x_0)-2f'(x_0)h+2f''(x_0)h^2-\tfrac{4}{3}f'''(x_0)h^3+\tfrac{2}{3}f''''(x_0)h^4+O(h^5)\\ f(x_0-h) &= f(x_0)-f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2-\tfrac{1}{6}f'''(x_0)h^3+\tfrac{1}{24}f''''(x_0)h^4+O(h^5)\\ f(x_0+h) &= f(x_0)+f'(x_0)h+\tfrac{1}{2}f''(x_0)h^2+\tfrac{1}{6}f'''(x_0)h^3+\tfrac{1}{24}f''''(x_0)h^4+O(h^5)\\ f(x_0+2h) &= f(x_0)+2f'(x_0)h+2f''(x_0)h^2+\tfrac{4}{3}f'''(x_0)h^3+\tfrac{2}{3}f''''(x_0)h^4+O(h^5)\end{aligned}\] and therefore \[\begin{aligned}Q(h) &= (A+B+C+D+E)f(x_0)+(-2A-B+D+2E)f'(x_0)h\\&\phantom{=}{}+\tfrac{1}{2}(4A+B+D+4E)f''(x_0)h^2+\tfrac{1}{6}(-8A-B+D+8E)f'''(x_0)h^3\\&\phantom{=}{}+\tfrac{1}{24}(16A+B+D+16E)f''''(x_0)h^4+O(h^5)\end{aligned}\] The equations that \(A, B, C, D, E\) must satisfy in order for \(Q(h)\) to maximise the order of the error are \[\begin{aligned}A+B+C +D+\;\;E&= 0\\ -2A-B+\phantom{+C}D+\;2E &= \frac{1}{h}\\ 4A+B+\phantom{+C}D+\;4E &= 0\\ -8A-B+\phantom{+C}D+\;8E &= 0\\ 16A+B+\phantom{+C}D+16E &= 0\end{aligned}\] These equations are relatively easy to solve: \[A=\frac{1}{12h}, B=-\frac{2}{3h},C=0, D=\frac{2}{3h}, E=-\frac{1}{12h}\] So we have found \[Q(h)=\frac{f(x_0-2h)-8f(x_0-h)+8f(x_0+h)-f(x_0+2h)}{12h}\] and \[f'(x_0)\approx\frac{f(x_0-2h)-8f(x_0-h)+8f(x_0+h)-f(x_0+2h)}{12h}\] This is the 5-point central difference formula. It can be verified that the truncation error has order \(h^4\).