Newton forward and backward difference interpolation

Newton’s forward and backward difference formulas interpolate in a table whose x values are equally spaced, such as a printed function table or readings taken at a fixed time step. Both build the same polynomial as Lagrange or divided differences, but the equal spacing lets you skip every division and read the coefficients straight off a difference table. You can check every step below with the forward and backward difference calculator.

The difference table

With points x₀, x₁ = x₀ + h, …, xₙ and values y₀ … yₙ, the forward differences are

Δyᵢ = yᵢ₊₁ − yᵢ,   Δ²yᵢ = Δyᵢ₊₁ − Δyᵢ,   Δᵏyᵢ = Δᵏ⁻¹yᵢ₊₁ − Δᵏ⁻¹yᵢ

Each column is the difference of neighbouring entries in the column before it, always lower minus upper. The examples below use √x at x = 1.0 … 1.4, rounded to four decimals as in a printed table (h = 0.1):

The forward formula uses the top diagonal (bold). The backward formula uses the last entry of each column: 1.1832, 0.0430, −0.0018, 0.0000, −0.0004. These are the backward differences ∇ᵏy₄, because ∇ᵏyₙ = Δᵏyₙ₋ₖ. There is only one table; “forward” and “backward” just name which diagonal you read.

The two formulas

Forward, for x near the start of the table, with p = (x − x₀)/h:

y = y₀ + p·Δy₀ + p(p − 1)/2!·Δ²y₀ + p(p − 1)(p − 2)/3!·Δ³y₀ + …

Backward, for x near the end of the table, with p = (x − xₙ)/h (so p is negative inside the table):

y = yₙ + p·∇yₙ + p(p + 1)/2!·∇²yₙ + p(p + 1)(p + 2)/3!·∇³yₙ + …

If you keep every term, both formulas give exactly the same polynomial and the same answer. The choice only matters when you stop early, which is usually the point: the formula that starts next to x gets most of the accuracy from its first two or three terms.

Worked example 1: forward, √1.05

x = 1.05 is near the top, so use the forward formula. p = (1.05 − 1.0)/0.1 = 0.5.

The true value is √1.05 = 1.024695. Two terms give 1.0244 (linear interpolation), three give 1.024675 and four give 1.024700, within 5 × 10⁻⁶. The fifth term makes the answer worse (1.024716). That is not a mistake: the third and fourth differences here are mostly rounding noise from the four-decimal table, so adding them adds noise.

Worked example 2: backward, √1.37

x = 1.37 is near the bottom, so use the backward formula. p = (1.37 − 1.4)/0.1 = −0.3.

• y₄ = 1.1832
• p·∇y₄ = −0.3 × 0.0430 = −0.0129 → 1.1703
• p(p + 1)/2!·∇²y₄ = (−0.3)(0.7)/2 × (−0.0018) = −0.105 × (−0.0018) = 0.000189 → 1.170489
• ∇³y₄ = 0, so the cubic term adds nothing.

The true value is √1.37 = 1.170470, so the quadratic estimate is off by 2 × 10⁻⁵, which is the size of the rounding in the table. Running the forward formula to the same degree would start from 1.0 and need p = 3.7, giving a much larger truncation error.

How many terms to keep

The first term you leave out is a good estimate of your error. In theory the error after the term in Δⁿ is

p(p − 1)…(p − n) / (n + 1)! · hⁿ⁺¹ · f⁽ⁿ⁺¹⁾(ξ)

for some ξ in the range of the points used, and hⁿ⁺¹f⁽ⁿ⁺¹⁾ ≈ Δⁿ⁺¹y₀. In practice: keep adding terms while they keep shrinking, and stop at the column where the differences become small and irregular (sign flips, values of a few units in the last decimal). That column is noise. A constant column means the data is exactly polynomial and every later term is zero.

Doing it in Python

numpy.diff(y, k) returns the k-th forward differences, so the top diagonal is its first element. This reproduces example 1:

import numpy as np
from math import factorial

x = np.array([1.0, 1.1, 1.2, 1.3, 1.4])
y = np.array([1.0, 1.0488, 1.0954, 1.1402, 1.1832])
h = x[1] - x[0]

# Top diagonal: Δ^k y0 for k = 0..n
top = [np.diff(y, k)[0] for k in range(len(y))]

def newton_forward(xq, degree=len(y) - 1):
    p = (xq - x[0]) / h
    total, factor = 0.0, 1.0
    for k in range(degree + 1):
        total += factor * top[k] / factorial(k)
        factor *= p - k
    return total

print(round(newton_forward(1.05, 3), 6))   # 1.0247

For production code, numpy.interp or SciPy’s splines are usually better choices; see interpolation in Python.

Common mistakes

• Unequal spacing. The formulas assume a constant h. Check every step before you start; if one differs, use divided differences instead.
• Wrong sign in the backward factors. Forward uses p(p − 1)(p − 2)…, backward uses p(p + 1)(p + 2)…. Mixing them up is the most common exam error.
• Positive p in the backward formula. p = (x − xₙ)/h is negative for any x inside the table. A positive p means you used x₀ instead of xₙ.
• Subtracting the wrong way. Every difference is lower minus upper. One reversed entry flips the sign of everything to its right.
• Forgetting the factorials. The third term divides by 2!, the fourth by 3!, and so on.
• Using every column blindly. High differences of rounded data are noise, as example 1 shows. More terms are not automatically better.
• Extrapolating. Both formulas will happily return a value outside the table, but the error grows fast. See choosing an interpolation method.

More guides

• What is interpolation?
• How to interpolate steam tables
• How to interpolate in Excel and Google Sheets
• Interpolation in Python
• Choosing an interpolation method