An illustrated guide to interpolation methods
Introduction
Interpolation is a way to estimate a function \( f(x) \) from samples. We assume that we are given a set of points \( (x_1, x_2, \dots ,x_n) \) known as knots along with the value of the function at those knots \[ (y_1, y_2, \dots, y_n) = (f(x_1), f(x_2), \dots, f(x_n)). \]
The fundamental building blocks of an interpolation method are the basis functions: a set of functions \(\varphi_k(x)\) that take the value 1 at a particular knot \(x_k\) and the value 0 at every other knot. Once we have basis functions we can interpolate between knots simply by taking the linear combination \[ \hat f(x) = \sum_{k=1}^n y_k \varphi_k(x). \] Different methods of constructing basis functions with this property lead to different interpolation methods. To illustrate some of the choices and trade-offs between interpolation methods, we consider the problem of reconstructing the functions \( f(x) = x^3 - 25x\) and \( g(x) = x - 5 \lfloor \frac{1}{5} x \rfloor \) from uniformly-spaced samples.
Shiftable basis functions
We first look at three methods that use shiftable basis functions. That is, there is a single function \( \phi(x) \) such that all the basis functions can be generated by shifting \( \phi \) to \(x_k\): \[ \varphi_k(x) - \phi(x - x_k). \]
Nearest-neighbor interpolation
This is probably the simplest method of interpolation. To estimate the value of the function at a point \(x\), we just use the value of the function at the knot \(x_k\) that is nearest to \(x\). This leads to rectangular basis functions \(\varphi_k(x)\) that take the value 1 at all point that are nearer to \(x_k\) than any other knot, and take the value 0 elsewhere. These are shiftable basis functions that correspond to the rect filter \[ \phi(x) = \begin{cases} 1 & \textrm{if } |x| \le 0.5 \\ 0 & \textrm{otherwise.} \end{cases} \]
The basis function is shown here:
The result of nearest neighbor interpolation is shown here:
Linear interpolation
Linear interpolation connects each pair of neighboring knots \( (x_k, y_k) \) and \( (x_{k+1}, y_{k+1}) \)with a line segment. This results in interpolation by a piecewise affine function. The shiftable basis function used for linear interpolation is shown here:
The result of linear interpolation is shown here:
Sinc interpolation
Sinc interpolation results from interpolation using shifted sinc basis functions \[ \phi(x) = \operatorname{sinc}(x) = \frac{\sin(\pi x)}{\pi x}. \] This basis function is shown here:
Sinc interpolation is important in the fields of harmonic analysis and signal processing where it is often used to reconstruct band-limited signals from samples. The result of sinc interpolation is shown here:
Non-shiftable basis functions
We now look at some additional interpolation methods that use basis functions whose shape changes as a function of \(k\). Here we will plot both the basis function corresponding to \(x_k = 0\), but also all the basis functions on the same axes.
Polynomial interpolation
Polynomial interpolation, sometimes known as Lagrange interpolation, involves fitting a polynomial to the data points. To fit \(n\) points requires a polynomial of degree \(n-1\). Each of the basis functions are themselves polynomials of degree \(n-1\) and can be computed using the formula \[ \varphi_k(x) = \prod_{\begin{smallmatrix}1 \leq m\leq n\\ m \neq k\end{smallmatrix}}{\frac {x-x_{m}}{x_{k}-x_{m}}}. \]
These basis functions are illustrated here:
The result of nearest polynomial interpolation is shown here:
We see that the polynomial basis functions can easily interpolate the cubic function (since it’s a polynomial itself, the interpolation is exact!) but have a lot of trouble interpolating the sawtooth function. The existence of these artifacts in the sawtooth interpolation is known as Runge’s phenomenon, which is a common problem in polynomial interpolation.
Kriging
Kriging, or Gaussian process regression, is a method that performs an optimal interpolation based on a statistical prior over possible functions. It is often assumed that the values \( y = f(x) \) are produced from a Gaussian process \( Y(x) \) with mean 0 and known covariance \( Cov(Y(x), Y(z)) = c(x, z) \). The optimal basis functions can be computed in terms of the covariance \(c\) using the formula
\[
\begin{bmatrix}
\varphi_1(x) \\ \vdots \\ \varphi_n(x)
\end{bmatrix}
= \begin{bmatrix}
c(x_1,x_1) & \cdots & c(x_1,x_n) \
\vdots & \ddots & \vdots \
c(x_n,x_1) & \cdots & c(x_n,x_n) \end{bmatrix}^{-1}
\begin{bmatrix}
c(x_1,x) \\ \vdots \\ c(x_n,x)
\end{bmatrix}
\]
The basis functions corresponding to the squared exponential covariance function
\[
c(x, z) = e^{-\alpha |x - z|^2}
\]
for \( \alpha = 1/2 \) are illustrated here:
The result of interpolation by Kriging is shown here:
Cubic splines
Cubic spline interpolation is another method that relies on piecewise cubic interpolation between knots. The splines are constrained to be continuous and to have continuous first and second derivatives, which leads to a tridiagonal system of linear equations to solve for the cubic polynomials. The basis functions corresponding to cubic spline interpolation are illustrated here:
The result of cubic spline interpolation is shown here: