Newton polynomial

Given a set of k + 1 data points

${\displaystyle (x_{0},y_{0}),\ldots ,(x_{j},y_{j}),\ldots ,(x_{k},y_{k})}$

where no two x_j are the same, the Newton interpolation polynomial is a linear combination of Newton basis polynomials

${\displaystyle N(x):=\sum _{j=0}^{k}a_{j}n_{j}(x)}$

with the Newton basis polynomials defined as

${\displaystyle n_{j}(x):=\prod _{i=0}^{j-1}(x-x_{i})}$

for j > 0 and ${\displaystyle n_{0}(x)\equiv 1}$.

The coefficients are defined as

${\displaystyle a_{j}:=[y_{0},\ldots ,y_{j}]}$

where

${\displaystyle [y_{0},\ldots ,y_{j}]}$

is the notation for divided differences.

Thus the Newton polynomial can be written as

${\displaystyle N(x)=[y_{0}]+[y_{0},y_{1}](x-x_{0})+\cdots +[y_{0},\ldots ,y_{k}](x-x_{0})(x-x_{1})\cdots (x-x_{k-1}).}$

Newton's formula is of interest because it is the straightforward and natural differences-version of Taylor's polynomial. Taylor's polynomial tells where a function will go, based on its y value, and its derivatives (its rate of change, and the rate of change of its rate of change, etc.) at one particular x value. Newton's formula is Taylor's polynomial based on finite differences instead of instantaneous rates of change.

As with other difference formulas, the degree of a Newton interpolating polynomial can be increased by adding more terms and points without discarding existing ones. Newton's form has the simplicity that the new points are always added at one end: Newton's forward formula can add new points to the right, and Newton's backward formula can add new points to the left.

The accuracy of polynomial interpolation depends on how close the interpolated point is to the middle of the x values of the set of points used. Obviously, as new points are added at one end, that middle becomes farther and farther from the first data point. Therefore, if it isn't known how many points will be needed for the desired accuracy, the middle of the x-values might be far from where the interpolation is done.

Gauss, Stirling, and Bessel all developed formulae to remedy that problem.[2]

Gauss's formula alternately adds new points at the left and right ends, thereby keeping the set of points centered near the same place (near the evaluated point). When so doing, it uses terms from Newton's formula, with data points and x values renamed in keeping with one's choice of what data point is designated as the x₀ data point.

Stirling's formula remains centered about a particular data point, for use when the evaluated point is nearer to a data point than to a middle of two data points.

Bessel's formula remains centered about a particular middle between two data points, for use when the evaluated point is nearer to a middle than to a data point.

Bessel and Stirling achieve that by sometimes using the average of two differences, and sometimes using the average of two products of binomials in x, where Newton's or Gauss's would use just one difference or product. Stirling's uses an average difference in odd-degree terms (whose difference uses an even number of data points); Bessel's uses an average difference in even-degree terms (whose difference uses an odd number of data points).

For any given finite set of data points, there is only one polynomial of least possible degree that passes through all of them. Thus, it is appropriate to speak of the "Newton form", or Lagrange form, etc., of the interpolation polynomial. However, the way the polynomial is obtained matters. There are several similar methods, such as those of Gauss, Bessel and Stirling. They can be derived from Newton's by renaming the x-values of the data points, but in practice they are important.

For the special case of x_i = i, there is a closely related set of polynomials, also called the Newton polynomials, that are simply the binomial coefficients for general argument. That is, one also has the Newton polynomials ${\displaystyle p_{n}(z)}$ given by

${\displaystyle p_{n}(z)={z \choose n}={\frac {z(z-1)\cdots (z-n+1)}{n!}}}$

In this form, the Newton polynomials generate the Newton series. These are in turn a special case of the general difference polynomials which allow the representation of analytic functions through generalized difference equations.

Solving an interpolation problem leads to a problem in linear algebra where we have to solve a system of linear equations. Using a standard monomial basis for our interpolation polynomial we get the very complicated Vandermonde matrix. By choosing another basis, the Newton basis, we get a system of linear equations with a much simpler lower triangular matrix which can be solved faster.

For k + 1 data points we construct the Newton basis as

${\displaystyle n_{j}(x):=\prod _{i=0}^{j-1}(x-x_{i})\qquad j=1,\ldots ,k.}$

Using these polynomials as a basis for ${\displaystyle \Pi _{k}}$ we have to solve

${\displaystyle {\begin{bmatrix}1&&\ldots &&0\\1&x_{1}-x_{0}&&&\\1&x_{2}-x_{0}&(x_{2}-x_{0})(x_{2}-x_{1})&&\vdots \\\vdots &\vdots &&\ddots &\\1&x_{k}-x_{0}&\ldots &\ldots &\prod _{j=0}^{k-1}(x_{k}-x_{j})\end{bmatrix}}{\begin{bmatrix}a_{0}\\\\\vdots \\\\a_{k}\end{bmatrix}}={\begin{bmatrix}y_{0}\\\\\vdots \\\\y_{k}\end{bmatrix}}}$

to solve the polynomial interpolation problem.

This system of equations can be solved iteratively by solving

${\displaystyle \sum _{i=0}^{j}a_{i}n_{i}(x_{j})=y_{j}\qquad j=0,\dots ,k.}$

While the interpolation formula can be found by solving a linear system of equations, there is a loss of intuition in what the formula is showing and why Newton's interpolation formula works is not readily apparent. To begin, we will need the following result:

${\displaystyle [y_{1},\ldots ,y_{n}]=[y_{n},\ldots ,y_{1}]}$ . This equality means that reversing the terms of the divided difference has no effect on the end result. We will prove this result with induction.

Basis: ${\textstyle [y_{1},y_{2}]={\frac {[y_{2}]-[y_{1}]}{x_{2}-x_{1}}}={\frac {[y_{1}]-[y_{2}]}{x_{1}-x_{2}}}=[y_{2},y_{1}]}$

Induction: Suppose the result holds for a divided difference involving fewer than ${\displaystyle n+1}$ terms. Using the induction hypothesis, we see that ${\displaystyle [y_{1},\ldots ,y_{n+1}]={\frac {[y_{2},\ldots ,y_{n+1}]-[y_{1},\ldots ,y_{n}]}{x_{n+1}-x_{1}}}={\frac {[y_{n},\ldots ,y_{1}]-[y_{n+1},\ldots ,y_{2}]}{x_{1}-x_{n+1}}}=[y_{n+1},\ldots ,y_{1}]}$ where the induction hypothesis was used at the 2nd equality.

In order to derive the interpolation formula, we will now use the following result which will also be proven with induction:

${\displaystyle [y_{1},\ldots ,y_{n}](x_{n}-x_{1})\cdot \ldots \cdot (x_{n}-x_{n-1})=y_{n}-P(x_{n})}$where ${\displaystyle P(x)}$ is the unique polynomial of degree (at most) ${\displaystyle n-2}$ passing through the points ${\displaystyle (x_{1},y_{1}),\ldots ,(x_{n-1},y_{n-1})}$. With this result, we can now exactly quantify the error between the interpolation polynomial ${\displaystyle P(x)}$ at ${\displaystyle (x_{n},y_{n})}$ and the true value ${\displaystyle y_{n}}$.

Basis: ${\displaystyle [y_{1},y_{2}](x_{2}-x_{1})={\frac {y_{2}-y_{1}}{x_{2}-x_{1}}}(x_{2}-x_{1})=y_{2}-y_{1}=y_{2}-P(x_{2})}$ where ${\displaystyle P(x)=y_{1}}$ is the unique polynomial of degree 0 passing through ${\displaystyle (x_{1},y_{1})}$.

Induction: Suppose the result holds for when there are fewer than ${\displaystyle n+1}$ points. Let ${\displaystyle P(x)}$ be the polynomial of degree (at most) ${\displaystyle n-1}$ passing through ${\displaystyle (x_{1},y_{1}),\ldots ,(x_{n},y_{n}).}$

${\displaystyle {\begin{aligned}&[y_{1},\ldots ,y_{n+1}](x_{n+1}-x_{1})\cdot \ldots \cdot (x_{n+1}-x_{n})\\&={\frac {[y_{2},\ldots ,y_{n+1}]-[y_{1},\ldots ,y_{n}]}{x_{n+1}-x_{1}}}(x_{n+1}-x_{1})\cdot \ldots \cdot (x_{n+1}-x_{n})\\&=\left([y_{2},\ldots ,y_{n+1}]-[y_{1},\ldots ,y_{n}]\right)(x_{n+1}-x_{2})\cdot \ldots \cdot (x_{n+1}-x_{n})\\&=[y_{2},\ldots ,y_{n+1}](x_{n+1}-x_{2})\cdot \ldots \cdot (x_{n+1}-x_{n})-[y_{1},\ldots ,y_{n}](x_{n+1}-x_{2})\cdot \ldots \cdot (x_{n+1}-x_{n})\\&=(y_{n+1}-Q(x_{n+1}))-[y_{1},\ldots ,y_{n}](x_{n+1}-x_{2})\cdot \ldots \cdot (x_{n+1}-x_{n})\\&=y_{n+1}-(Q(x_{n+1})+[y_{1},\ldots ,y_{n}](x_{n+1}-x_{2})\cdot \ldots \cdot (x_{n+1}-x_{n}))\end{aligned}}}$

where ${\displaystyle Q(x)}$ is the unique polynomial of degree (at most) ${\displaystyle n-2}$ passing through the points ${\displaystyle (x_{2},y_{2}),\ldots ,(x_{n},y_{n})}$. The second to last equality comes from the induction hypothesis as ${\displaystyle [y_{2},\ldots ,y_{n+1}]}$ involves ${\displaystyle n}$ points and thus ${\displaystyle [y_{2},\ldots ,y_{n+1}](x_{n+1}-x_{2})\cdot \ldots \cdot (x_{n+1}-x_{n})=y_{n+1}-Q(x_{n+1})}$. Getting close to the desired result, we now claim that ${\displaystyle Q(x)+[y_{1},\ldots ,y_{n}](x-x_{2})\cdot \ldots \cdot (x-x_{n})=P(x)}$as both polynomials pass through ${\displaystyle (x_{1},y_{1}),\ldots ,(x_{n},y_{n})}$ and are of degree (at most) ${\displaystyle n-1}$. Both of these criteria uniquely define a polynomial. The fact that the left hand side passes through ${\displaystyle (x_{2},y_{2}),\ldots ,(x_{n},y_{n})}$ is readily apparent by how ${\displaystyle Q(x)}$ is defined to be. To demonstrate the left hand side passes through ${\displaystyle (x_{1},y_{1})}$, we will use the first result proved above along with the induction hypothesis:

${\displaystyle {\begin{aligned}&Q(x_{1})+[y_{1},\ldots ,y_{n}](x_{1}-x_{2})\cdot \ldots \cdot (x_{1}-x_{n})\\&=Q(x_{1})+[y_{n},\ldots ,y_{1}](x_{1}-x_{n})\cdot \ldots \cdot (x_{1}-x_{2})\\&=Q(x_{1})+y_{1}-Q(x_{1})\\&=y_{1}\\&=P(x_{1})\\\end{aligned}}}$

where the 2nd equality follows from the fact that ${\displaystyle Q(x)}$ is the polynomial of degree (at most) ${\displaystyle n-2}$ passing through ${\displaystyle (x_{n},y_{n}),\ldots ,(x_{2},y_{2})}$ satisfying the induction hypothesis. Continuing the induction step above, we now see that ${\displaystyle [y_{1},\ldots ,y_{n+1}](x_{n+1}-x_{1})\cdot \ldots \cdot (x_{n+1}-x_{n})=y_{n+1}-P(x_{n+1})}$where ${\displaystyle P(x)}$ is the polynomial of degree ${\displaystyle n-1}$ passing through ${\displaystyle (x_{1},y_{1}),\ldots ,(x_{n},y_{n}).}$ Thus the proof is complete.

All this work now leads to where Newton's interpolation formula comes from. Rearranging the result above, we note that ${\displaystyle P(x)+[y_{1},\ldots ,y_{n+1}](x-x_{1})\cdot \ldots \cdot (x-x_{n})}$is the polynomial of degree (at most) ${\displaystyle n}$ passing through ${\displaystyle (x_{1},y_{1}),\ldots ,(x_{n+1},y_{n+1})}$, and thus we see that "extending" a polynomial ${\displaystyle P(x)}$ to the next point ${\displaystyle (x_{n+1},y_{n+1})}$ requires adding the term ${\displaystyle [y_{1},\ldots ,y_{n+1}](x-x_{1})\cdot \ldots \cdot (x-x_{n})}$ giving us Newton's interpolation formula.

The limit of the Newton polynomial if all nodes coincide is a Taylor polynomial, because the divided differences become derivatives.

${\displaystyle \lim _{(x_{0},\dots ,x_{n})\to (z,\dots ,z)}f[x_{0}]+f[x_{0},x_{1}]\cdot (\xi -x_{0})+\dots +f[x_{0},\dots ,x_{n}]\cdot (\xi -x_{0})\cdot \dots \cdot (\xi -x_{n-1})=}$

${\displaystyle =f(z)+f'(z)\cdot (\xi -z)+\dots +{\frac {f^{(n)}(z)}{n!}}\cdot (\xi -z)^{n}}$

As can be seen from the definition of the divided differences new data points can be added to the data set to create a new interpolation polynomial without recalculating the old coefficients. And when a data point changes we usually do not have to recalculate all coefficients. Furthermore, if the x_i are distributed equidistantly the calculation of the divided differences becomes significantly easier. Therefore, the divided-difference formulas are usually preferred over the Lagrange form for practical purposes.

Examples

The divided differences can be written in the form of a table. For example, for a function f is to be interpolated on points ${\displaystyle x_{0},\ldots ,x_{n}}$. Write

${\displaystyle {\begin{matrix}x_{0}&f(x_{0})&&\\&&{f(x_{1})-f(x_{0}) \over x_{1}-x_{0}}&\\x_{1}&f(x_{1})&&{{f(x_{2})-f(x_{1}) \over x_{2}-x_{1}}-{f(x_{1})-f(x_{0}) \over x_{1}-x_{0}} \over x_{2}-x_{0}}\\&&{f(x_{2})-f(x_{1}) \over x_{2}-x_{1}}&\\x_{2}&f(x_{2})&&\vdots \\&&\vdots &\\\vdots &&&\vdots \\&&\vdots &\\x_{n}&f(x_{n})&&\\\end{matrix}}}$

Then the interpolating polynomial is formed as above using the topmost entries in each column as coefficients.

For example, suppose we are to construct the interpolating polynomial to f(x) = tan(x) using divided differences, at the points

${\displaystyle x_{0}=-{\tfrac {3}{2}}}$	${\displaystyle x_{1}=-{\tfrac {3}{4}}}$	${\displaystyle x_{2}=0}$	${\displaystyle x_{3}={\tfrac {3}{4}}}$	${\displaystyle x_{4}={\tfrac {3}{2}}}$
${\displaystyle f(x_{0})=-14.1014}$	${\displaystyle f(x_{1})=-0.931596}$	${\displaystyle f(x_{2})=0}$	${\displaystyle f(x_{3})=0.931596}$	${\displaystyle f(x_{4})=14.1014}$

Using six digits of accuracy, we construct the table

${\displaystyle {\begin{matrix}-{\tfrac {3}{2}}&-14.1014&&&&\\&&17.5597&&&\\-{\tfrac {3}{4}}&-0.931596&&-10.8784&&\\&&1.24213&&4.83484&\\0&0&&0&&0\\&&1.24213&&4.83484&\\{\tfrac {3}{4}}&0.931596&&10.8784&&\\&&17.5597&&&\\{\tfrac {3}{2}}&14.1014&&&&\\\end{matrix}}}$

Thus, the interpolating polynomial is

${\displaystyle -14.1014+17.5597(x+{\tfrac {3}{2}})-10.8784(x+{\tfrac {3}{2}})(x+{\tfrac {3}{4}})+4.83484(x+{\tfrac {3}{2}})(x+{\tfrac {3}{4}})(x)+0(x+{\tfrac {3}{2}})(x+{\tfrac {3}{4}})(x)(x-{\tfrac {3}{4}})=}$

${\displaystyle =-0.00005-1.4775x-0.00001x^{2}+4.83484x^{3}}$

Given more digits of accuracy in the table, the first and third coefficients will be found to be zero.

Another example:

The sequence ${\displaystyle f_{0}}$ such that ${\displaystyle f_{0}(1)=6,f_{0}(2)=9,f_{0}(3)=2}$ and ${\displaystyle f_{0}(4)=5}$, i.e., they are ${\displaystyle 6,9,2,5}$ from ${\displaystyle x_{0}=1}$ to ${\displaystyle x_{3}=4}$.

You obtain the slope of order ${\displaystyle 1}$ in the following way:

• ${\displaystyle f_{1}(x_{0},x_{1})={\frac {f_{0}(x_{1})-f_{0}(x_{0})}{x_{1}-x_{0}}}={\frac {9-6}{2-1}}=3}$
• ${\displaystyle f_{1}(x_{1},x_{2})={\frac {f_{0}(x_{2})-f_{0}(x_{1})}{x_{2}-x_{1}}}={\frac {2-9}{3-2}}=-7}$
• ${\displaystyle f_{1}(x_{2},x_{3})={\frac {f_{0}(x_{3})-f_{0}(x_{2})}{x_{3}-x_{2}}}={\frac {5-2}{4-3}}=3}$

As we have the slopes of order ${\displaystyle 1}$, it is possible to obtain the next order:

• ${\displaystyle f_{2}(x_{0},x_{1},x_{2})={\frac {f_{1}(x_{1},x_{2})-f_{1}(x_{0},x_{1})}{x_{2}-x_{0}}}={\frac {-7-3}{3-1}}=-5}$
• ${\displaystyle f_{2}(x_{1},x_{2},x_{3})={\frac {f_{1}(x_{2},x_{3})-f_{1}(x_{1},x_{2})}{x_{3}-x_{1}}}={\frac {3-(-7)}{4-2}}=5}$

Finally, we define the slope of order ${\displaystyle 3}$:

• ${\displaystyle f_{3}(x_{0},x_{1},x_{2},x_{3})={\frac {f_{2}(x_{1},x_{2},x_{3})-f_{2}(x_{0},x_{1},x_{2})}{x_{3}-x_{0}}}={\frac {5-(-5)}{4-1}}={\frac {10}{3}}}$

Once we have the slope, we can define the consequent polynomials:

• ${\displaystyle p_{0}(x)=6}$.
• ${\displaystyle p_{1}(x)=6+3(x-1)}$
• ${\displaystyle p_{2}(x)=6+3(x-1)-5(x-1)(x-2)}$.
• ${\displaystyle p_{3}(x)=6+3(x-1)-5(x-1)(x-2)+{\frac {10}{3}}(x-1)(x-2)(x-3)}$

• De numeris triangularibus et inde de progressionibus arithmeticis: Magisteria magna, a work by Thomas Harriot describing similar methods for interpolation, written 50 years earlier than Newton's work but not published until 2009
• Newton series
• Neville's schema
• Polynomial interpolation
• Lagrange form of the interpolation polynomial
• Bernstein form of the interpolation polynomial
• Hermite interpolation
• Carlson's theorem
• Table of Newtonian series

• Module for the Newton Polynomial by John H. Mathews