Last modified: July 26, 2024

This article is written in: 🇺🇸

Cubic Spline

Cubic spline interpolation is a refined mathematical tool frequently used within numerical analysis. It's an approximation technique that employs piecewise cubic polynomials, collectively forming a cubic spline. These cubic polynomials are specifically engineered to pass through a defined set of data points, hence striking a balance between overly simple (like linear) and overly intricate (like high-degree polynomial) interpolations.

Conceptual Illustration:

Imagine you have a set of points on a 2D plane. A cubic spline will "thread" through these points, forming a smooth curve that does not exhibit sudden bends or wiggles:

The spline is constructed from piecewise cubic segments that meet at the data points with continuous first and second derivatives, ensuring a visually natural and mathematically smooth interpolation.

Mathematical Formulation

We start with a set of n+1$n+1$ data points:

(x0,y0),(x1,y1),…,(xn,yn)$$(x_0,y_0),(x_1,y_1),\dots ,(x_n,y_n)$$

with x0<x1<⋯<xn$x_0<x_1<\cdots <x_n$.

A cubic spline is a function S(x)$S\left(x\right)$ defined piecewise on each interval [xi,xi+1]$[x_i,x_{i+1}]$, i=0,1,…,n−1$i=0,1,\dots ,n-1$:

Si(x)=ai+bi(x−xi)+ci(x−xi)2+di(x−xi)3$$S_i\left(x\right)=a_i+b_i(x-x_i)+c_i(x-x_i{)}^2+d_i(x-x_i{)}^3$$

where ai,bi,ci,di$a_i,b_i,c_i,d_i$ are the unknown coefficients for the cubic polynomial on the interval [xi,xi+1]$[x_i,x_{i+1}]$.

Key Requirements:

I. Interpolation Condition:

Each segment must pass through the given data points:

Si(xi)=yi,Si(xi+1)=yi+1.$$S_i\left(x_i\right)=y_i,S_i\left(x_{i+1}\right)=y_{i+1}.$$

II. Continuity of the Spline:

The function S(x)$S\left(x\right)$ must be continuous at the interior points:

Si(xi+1)=Si+1(xi+1).$$S_i\left(x_{i+1}\right)=S_{i+1}\left(x_{i+1}\right).$$

III. Continuity of First Derivative:

The first derivatives of adjacent segments must match:

S′i(xi+1)=S′i+1(xi+1).

IV. Continuity of Second Derivative:

Similarly, the second derivatives must also be continuous:

S″i(xi+1)=S″i+1(xi+1).

V. Boundary Conditions:

For a natural cubic spline, the second derivatives at the boundaries are set to zero:

S″0(x0)=0,S″n−1(xn)=0.

These conditions lead to a system of linear equations that determine the coefficients ai,bi,ci, and di.

Derivation of the Coefficients

I. Divided Differences and Step Sizes:

Let the spacing between consecutive points be:

hi=xi+1−xi,i=0,1,…,n−1

II. Formulation in Terms of Second Derivatives:

Let Mi=S″(xi). If we can determine Mi for all i, we can write each piece Si(x) as:

$$S_i(x) = M_{i}\frac{(x_{i+1}-x)^3}{6h_i} + M_{i+1}\frac{(x - x_i)^3}{6h_i}

• \left(y_i - \frac{M_i h_i^2}{6}\right)\frac{x_{i+1}-x}{h_i}

• \left(y_{i+1} - \frac{M_{i+1}h_i^2}{6}\right)\frac{x - x_i}{h_i}.$$

This form ensures the correct boundary conditions and continuity of derivatives once Mi are found.

III. System of Equations for Mi:

Using the conditions for continuity of first and second derivatives, one obtains a tridiagonal system of equations in terms of the unknown second derivatives Mi. For natural splines, we have:

M0=0,Mn=0.

The remaining Mi's (for i=1,…,n−1) are found by solving this system:

μiMi−1+2Mi+λiMi+1=di, where μi,λi,di depend on the data points and hi.

IV. Once Mi are determined, the coefficients ai,bi,ci,di of the cubic spline in the standard form:

Si(x)=ai+bi(x−xi)+ci(x−xi)2+di(x−xi)3

can be computed directly:

ai=yi

bi=yi+1−yihi−hi3(2Mi+Mi+1)

ci=Mi

di=Mi+1−Mi3hi.

Algorithm Steps

I. Data Preparation:

• Sort the data points (xi,yi) in ascending order by xi.
• Compute hi=xi+1−xi for i=0,…,n−1.

II. Form the Equations:

Set up the linear system to solve for the second derivatives Mi. For a natural spline:

M0=0,Mn=0

For i=1,…,n−1:

hi−1Mi−1+2(hi−1+hi)Mi+hiMi+1=3(yi+1−yihi−yi−yi−1hi−1)

III. Solve the Tridiagonal System:

Use an efficient method (like the Thomas algorithm) to solve for M1,M2,…,Mn−1.

IV. Calculate the Coefficients:

With Mi known, compute:

ai=yi,bi,ci=Mi,di

as described above.

V. Interpolation:

For a given x, find the interval [xi,xi+1] such that xi≤x≤xi+1 and evaluate:

Si(x)=ai+bi(x−xi)+ci(x−xi)2+di(x−xi)3

Example

Consider three points: (0,0),(1,0.5),(2,0).

• h0=1,h1=1.
• Boundary conditions: M0=M2=0.

We form the equation for i=1:

1⋅M0+2(1+1)M1+1⋅M2=3((0−0.5)/1−(0.5−0)/1)=3(−0.5−0.5)=−3.

Since M0=M2=0:

4M1=−3⟹M1=−0.75

Then:

a0=y0=0,c0=M0=0

b0=(y1−y0)/h0−h0(2M0+M1)/3=0.5/1−(2⋅0+(−0.75))/3=0.5+0.25=0.75

d0=(M1−M0)/(3h0)=(−0.75−0)/3=−0.25

And similarly for i=1.

This gives piecewise polynomials smoothly interpolating the given points.

Advantages

• Splines ensure smoothness by producing curves with continuous first and second derivatives, leading to a natural and visually appealing fit.
• The method provides controlled behavior, avoiding the large oscillations often observed in high-degree polynomial interpolation.
• Local control means that changing a single data point only impacts the neighboring spline segments, rather than the entire curve.

Limitations

• The method involves complexity, as it requires setting up and solving a linear system of equations, making it more complicated than simpler interpolation techniques.
• Computational cost is higher than methods like linear interpolation, particularly for large datasets with many control points.
• Boundary conditions must be specified (e.g., natural, clamped) to define the behavior at the endpoints, which can influence the overall fit of the spline.