Divided differences

Notation

Note that the divided difference ${\displaystyle [y_{k},\ldots ,y_{k+j}]}$ depends on the values ${\displaystyle x_{k},\ldots ,x_{k+j}}$ and ${\displaystyle y_{k},\ldots ,y_{k+j}}$, but the notation hides the dependency on the x-values. If the data points are given by a function f, $${\displaystyle (x_{0},y_{0}),\ldots ,(x_{n},y_{n})=(x_{0},f(x_{0})),\ldots ,(x_{n},f(x_{n}))}$$ one sometimes writes the divided difference in the notation $${\displaystyle f[x_{k},\ldots ,x_{k+j}]\ {\stackrel {\text{def}}{=}}\ [f(x_{k}),\ldots ,f(x_{k+j})]=[y_{k},\ldots ,y_{k+j}].}$$Other notations for the divided difference of the function ƒ on the nodes x₀, ..., x_n are: $${\displaystyle f[x_{k},\ldots ,x_{k+j}]={\mathopen {[}}x_{0},\ldots ,x_{n}]f={\mathopen {[}}x_{0},\ldots ,x_{n};f]=D[x_{0},\ldots ,x_{n}]f.}$$

Polynomials and power series

The matrix $${\displaystyle J={\begin{pmatrix}x_{0}&1&0&0&\cdots &0\\0&x_{1}&1&0&\cdots &0\\0&0&x_{2}&1&&0\\\vdots &\vdots &&\ddots &\ddots &\\0&0&0&0&\;\ddots &1\\0&0&0&0&&x_{n}\end{pmatrix}}}$$ contains the divided difference scheme for the identity function with respect to the nodes ${\displaystyle x_{0},\dots ,x_{n}}$, thus ${\displaystyle J^{m}}$ contains the divided differences for the power function with exponent ${\displaystyle m}$. Consequently, you can obtain the divided differences for a polynomial function ${\displaystyle p}$ by applying ${\displaystyle p}$ to the matrix ${\displaystyle J}$: If $${\displaystyle p(\xi )=a_{0}+a_{1}\cdot \xi +\dots +a_{m}\cdot \xi ^{m}}$$ and $${\displaystyle p(J)=a_{0}+a_{1}\cdot J+\dots +a_{m}\cdot J^{m}}$$ then $${\displaystyle T_{p}(x)=p(J).}$$ This is known as Opitz' formula.

Now consider increasing the degree of ${\displaystyle p}$ to infinity, i.e. turn the Taylor polynomial into a Taylor series. Let ${\displaystyle f}$ be a function which corresponds to a power series. You can compute the divided difference scheme for ${\displaystyle f}$ by applying the corresponding matrix series to ${\displaystyle J}$: If $${\displaystyle f(\xi )=\sum _{k=0}^{\infty }a_{k}\xi ^{k}}$$ and $${\displaystyle f(J)=\sum _{k=0}^{\infty }a_{k}J^{k}}$$ then $${\displaystyle T_{f}(x)=f(J).}$$

Expanded form

$${\displaystyle {\begin{aligned}f[x_{0}]&=f(x_{0})\\f[x_{0},x_{1}]&={\frac {f(x_{0})}{(x_{0}-x_{1})}}+{\frac {f(x_{1})}{(x_{1}-x_{0})}}\\f[x_{0},x_{1},x_{2}]&={\frac {f(x_{0})}{(x_{0}-x_{1})\cdot (x_{0}-x_{2})}}+{\frac {f(x_{1})}{(x_{1}-x_{0})\cdot (x_{1}-x_{2})}}+{\frac {f(x_{2})}{(x_{2}-x_{0})\cdot (x_{2}-x_{1})}}\\f[x_{0},x_{1},x_{2},x_{3}]&={\frac {f(x_{0})}{(x_{0}-x_{1})\cdot (x_{0}-x_{2})\cdot (x_{0}-x_{3})}}+{\frac {f(x_{1})}{(x_{1}-x_{0})\cdot (x_{1}-x_{2})\cdot (x_{1}-x_{3})}}+\\&\quad \quad {\frac {f(x_{2})}{(x_{2}-x_{0})\cdot (x_{2}-x_{1})\cdot (x_{2}-x_{3})}}+{\frac {f(x_{3})}{(x_{3}-x_{0})\cdot (x_{3}-x_{1})\cdot (x_{3}-x_{2})}}\\f[x_{0},\dots ,x_{n}]&=\sum _{j=0}^{n}{\frac {f(x_{j})}{\prod _{k\in \{0,\dots ,n\}\setminus \{j\}}(x_{j}-x_{k})}}\end{aligned}}}$$

With the help of the polynomial function ${\displaystyle \omega (\xi )=(\xi -x_{0})\cdots (\xi -x_{n})}$ this can be written as $${\displaystyle f[x_{0},\dots ,x_{n}]=\sum _{j=0}^{n}{\frac {f(x_{j})}{\omega '(x_{j})}}.}$$

Peano form

If ${\displaystyle x_{0}<x_{1}<\cdots <x_{n}}$ and ${\displaystyle n\geq 1}$, the divided differences can be expressed as $${\displaystyle f[x_{0},\ldots ,x_{n}]={\frac {1}{(n-1)!}}\int _{x_{0}}^{x_{n}}f^{(n)}(t)\;B_{n-1}(t)\,dt}$$ where ${\displaystyle f^{(n)}}$ is the ${\displaystyle n}$-th derivative of the function ${\displaystyle f}$ and ${\displaystyle B_{n-1}}$ is a certain B-spline of degree ${\displaystyle n-1}$ for the data points ${\displaystyle x_{0},\dots ,x_{n}}$, given by the formula $${\displaystyle B_{n-1}(t)=\sum _{k=0}^{n}{\frac {(\max(0,x_{k}-t))^{n-1}}{\omega '(x_{k})}}}$$

This is a consequence of the Peano kernel theorem; it is called the Peano form of the divided differences and ${\displaystyle B_{n-1}}$ is the Peano kernel for the divided differences, all named after Giuseppe Peano.

Hermite–Genocchi formula

When ⁠${\displaystyle f}$⁠ is ⁠${\displaystyle n}$⁠ times continuously differentiable on a convex domain containing the nodes ⁠${\displaystyle x_{0},\dots ,x_{n}}$⁠, including the confluent case of repeated nodes in arbitrary order (see § Generalization), the divided difference can be expressed as an integral over a ⁠${\displaystyle n}$⁠-simplex ⁠${\displaystyle \Delta _{n}}$⁠: $${\displaystyle f[x_{0},\dots ,x_{n}]=\int _{\Delta _{n}}f^{(n)}({\textstyle \sum _{k=0}^{n}t_{k}x_{k}})\,dt_{0}\cdots dt_{n},}$$ where $${\displaystyle \Delta _{n}=\{(t_{0},\dots ,t_{n})\in \mathbb {R} ^{n}\mid t_{i}\geq 0,\,{\textstyle \sum _{k=0}^{n}t_{k}=1}\}}$$ is the ⁠${\displaystyle n}$⁠-simplex.

Cauchy's integral formula

When ⁠${\displaystyle f}$⁠ is analytic on a simply connected domain ⁠${\displaystyle \Omega \subset \mathbb {C} }$⁠ containing the nodes ⁠${\displaystyle x_{0},\dots ,x_{n}}$⁠, including the confluent case of repeated nodes in arbitrary order (see § Generalization), the divided difference admits a contour integral representation via Cauchy's integral formula: $${\displaystyle f[x_{0},\dots ,x_{n}]={\frac {1}{2\pi i}}\oint _{\Gamma }{\frac {f(z)}{\prod _{k=0}^{n}(z-x_{k})}}\,dz,}$$ where ⁠${\displaystyle \Gamma }$⁠ is a positively oriented, simple closed contour in ⁠${\displaystyle \Omega }$⁠ enclosing all the nodes ⁠${\displaystyle x_{0},\dots ,x_{n}}$⁠.

Forward and backward differences

When the data points are equidistantly distributed we get the special case called forward differences. They are easier to calculate than the more general divided differences.

Given n+1 data points $${\displaystyle (x_{0},y_{0}),\ldots ,(x_{n},y_{n})}$$ with $${\displaystyle x_{k}=x_{0}+kh,\ {\text{ for }}\ k=0,\ldots ,n{\text{ and fixed }}h>0}$$ the forward differences are defined as: $${\displaystyle {\begin{aligned}\Delta ^{(0)}y_{k}&:=y_{k},\qquad k=0,\ldots ,n\\\Delta ^{(j)}y_{k}&:=\Delta ^{(j-1)}y_{k+1}-\Delta ^{(j-1)}y_{k},\qquad k=0,\ldots ,n-j,\ j=1,\dots ,n.\end{aligned}}}$$
whereas the backward differences are defined as: $${\displaystyle {\begin{aligned}\nabla ^{(0)}y_{k}&:=y_{k},\qquad k=0,\ldots ,n\\\nabla ^{(j)}y_{k}&:=\nabla ^{(j-1)}y_{k}-\nabla ^{(j-1)}y_{k-1},\qquad k=0,\ldots ,n-j,\ j=1,\dots ,n.\end{aligned}}}$$
Thus the forward difference table is written as: $${\displaystyle {\begin{matrix}y_{0}&&&\\&\Delta y_{0}&&\\y_{1}&&\Delta ^{2}y_{0}&\\&\Delta y_{1}&&\Delta ^{3}y_{0}\\y_{2}&&\Delta ^{2}y_{1}&\\&\Delta y_{2}&&\\y_{3}&&&\\\end{matrix}}}$$
whereas the backwards difference table is written as: $${\displaystyle {\begin{matrix}y_{0}&&&\\&\nabla y_{1}&&\\y_{1}&&\nabla ^{2}y_{2}&\\&\nabla y_{2}&&\nabla ^{3}y_{3}\\y_{2}&&\nabla ^{2}y_{3}&\\&\nabla y_{3}&&\\y_{3}&&&\\\end{matrix}}}$$

The relationship between divided differences and forward differences is $${\displaystyle [y_{j},y_{j+1},\ldots ,y_{j+k}]={\frac {1}{k!h^{k}}}\Delta ^{(k)}y_{j},}$$whereas for backward differences:^{[citation needed]}$${\displaystyle [{y}_{j},y_{j-1},\ldots ,{y}_{j-k}]={\frac {1}{k!h^{k}}}\nabla ^{(k)}y_{j}.}$$

Explicit formula

When the data points are equispaced we can also derive an explicit formula for ${\displaystyle [y_{j},\ldots ,y_{j+k}]}$. For any fixed ${\displaystyle k\leq n}$ and ${\displaystyle j}$ such that ${\displaystyle j+k\leq n}$, $${\displaystyle [y_{j},\ldots ,y_{j+k}]={\frac {(-1)^{j+k}}{k!h^{k}}}\sum _{i=j}^{j+k}(-1)^{i}y_{i}{\binom {k}{i-j}}.}$$

Proof. We prove this by induction on ${\displaystyle k}$:

Base case: Let ${\displaystyle k=0}$ and ${\displaystyle j=0\ldots n}$. Then by definition ${\displaystyle [y_{j}]=y_{j}}$, and $${\displaystyle {\frac {(-1)^{j}}{0!h^{0}}}\sum _{i=j}^{j}(-1)^{i}y_{i}{\binom {0}{i-j}}=(-1)^{2j}y_{j}=y_{j}.}$$

Induction step: Assume the above formula holds until ${\displaystyle k}$ and consider ${\displaystyle j}$ such that ${\displaystyle j+k+1\leq n}$. By the recursive definition we have

$${\displaystyle [y_{j},\ldots ,y_{j+k+1}]={\frac {1}{h(k+1)}}\left([y_{j+1},\ldots ,y_{j+k+1}]-[y_{j},\ldots ,y_{j+k}]\right).}$$

We can use our inductive hypothesis on both members of the left side, and obtain

$${\displaystyle {\frac {1}{h(k+1)}}\left({\frac {(-1)^{j+k+1}}{k!h^{k}}}\sum _{i=j+1}^{j+k+1}(-1)^{i}y_{i}{\binom {k}{i-j}}-{\frac {(-1)^{j+k}}{k!h^{k}}}\sum _{i=j}^{j+k}(-1)^{i}y_{i}{\binom {k}{i-j}}\right).}$$

This can be rearranged as

$${\displaystyle {\frac {(-1)^{j+k+1}}{(k+1)!h^{k+1}}}\sum _{i=j}^{j+k+1}(-1)^{i}y_{i}\left({\binom {k}{i-j}}+{\binom {k}{i-j-1}}\right)}$$

where ${\displaystyle {\binom {k}{s}}=0}$ when ${\displaystyle s<k}$. We obtain our thesis for ${\displaystyle k+1}$ by substituting the identity $${\displaystyle {\binom {k}{i-j}}+{\binom {k}{i-j-1}}={\binom {k+1}{i-j}}.}$$

Generalized case example

To evaluate a confluent divided difference with non-adjacent repeated nodes, such as ⁠${\displaystyle f[z_{1},z_{2},z_{1}]}$⁠, the node sequence is first permuted to group identical nodes contiguously: $${\displaystyle f[z_{1},z_{2},z_{1}]=f[z_{1},z_{1},z_{2}].}$$ The evaluation then proceeds by constructing the divided difference table: $${\displaystyle {\begin{matrix}x_{0}=z_{1}&f(x_{0})&&\\&&{\frac {f'(z_{1})}{1!}}&\\x_{1}=z_{1}&f(x_{1})&&{\frac {{\frac {f(x_{2})-f(x_{1})}{x_{2}-x_{1}}}-{\frac {f'(z_{1})}{1!}}}{x_{2}-x_{0}}}\\&&{\frac {f(x_{2})-f(x_{1})}{x_{2}-x_{1}}}&\\x_{2}=z_{2}&f(x_{2})&&\\\end{matrix}}}$$ Extracting and simplifying the rightmost entry from the resulting divided difference table yields: $${\displaystyle f[z_{1},z_{2},z_{1}]={\frac {f(z_{2})-f(z_{1})}{(z_{2}-z_{1})^{2}}}-{\frac {f'(z_{1})}{z_{2}-z_{1}}}.}$$