parv bhadra

An Ode to Splines

roller coaster

Splines don’t try to be everything all at once. They don’t impose some grand, sweeping polynomial narrative onto data. They’re modest, thoughtful, and precise - handling their own little segment of the world with grace, making everything look effortlessly smooth.

They’re the kind of function you bring home to meet your dataset.

You've probably seen them in CAD software, robotics path planning, and numerical solutions to differential equations - but today, I just want to talk about why they they’re the perfect balance of simplicity and power, and honestly, they deserve more appreciation.

The Problem with High-Degree Polynomial Interpolation

Let's say you have a dataset with n+1 points (x0,y0),(x1,y1),,(xn,yn). The naive way to interpolate these points is to fit an nth-degree polynomial that passes through every single one of them. Enter the Vandermonde matrix, an absolute nightmare for numerical stability. This polynomial will:

  1. Oscillate aggressively (Runge's phenomenon) near the edges because high-degree polynomials have a habit of being extra.
  2. Be computationally expensive to evaluate and differentiate.
  3. Not generalize well - if you add an extra point, you basically have to start over.

Example of a Vandermonde Matrix

For a set of points with x-coordinates x0,x1,x2,...,xn, the Vandermonde matrix looks like:

V = | 1    x₀    x₀²   ...   x₀ⁿ |
    | 1    x₁    x₁²   ...   x₁ⁿ |
    | 1    x₂    x₂²   ...   x₂ⁿ |
    | ⋮     ⋮      ⋮    ⋱     ⋮   |
    | 1    xₙ    xₙ²   ...    xₙⁿ |

For instance, with points at x=1,2,...,10, the Vandermonde matrix would be:

V = | 1    1     1     1     1      1       1        1         1         1  |
    | 1    2     4     8    16     32      64      128       256        512 |
    | 1    3     9    27    81    243     729     2187      6561      19683 |
    | 1    4    16    64   256   1024    4096    16384     65536     262144 |
    | ⋮     ⋮     ⋮      ⋮     ⋮       ⋮       ⋮        ⋮          ⋮          ⋮ |
    | 1   10   100  1000 10000 100000 1000000 10000000 100000000 1000000000 |

The condition number of this matrix grows exponentially with its size, making the interpolation problem highly sensitive to small perturbations in the data.

To see this in action, consider fitting a 9th-degree polynomial to 10 data points sampled from Runge's function:

f(x)=11+25x2

You'll notice that the polynomial wildly oscillates near the boundaries, making it unreliable for interpolation. This is where splines shine. They don't try to do too much, they don't overreach—they just handle their own local segment with precision.

What is a Spline?

A spline is a piecewise polynomial function that stitches together multiple lower-degree polynomials while ensuring smoothness at the points where they meet. Formally, a spline of degree k consists of polynomials Si(x) on each subinterval [xi,xi+1] such that:

  1. Si(x) and Si+1(x) match at xi+1.
  2. Their first k1 derivatives also match at xi+1.

Splines are cooperative, working in small sections and blending seamlessly into one another. Unlike high-degree polynomials, which insist on making everything about them, splines are content with solving just their piece of the problem.

Linear Splines: The Basic Form

A linear spline is just piecewise linear interpolation. It's quick and easy, but it's also discontinuous in the first derivative—making it a bad choice for anything that requires smooth motion.

Si(x)=yi+yi+1yixi+1xi(xxi),x[xi,xi+1]

Cubic Splines: The Industry Standard

Cubic splines hit the sweet spot of smoothness (continuous up to the second derivative) and computational efficiency. They take the form:

Si(x)=ai+bi(xxi)+ci(xxi)2+di(xxi)3,x[xi,xi+1]

To solve for the coefficients ai,bi,ci,di, we impose:

  1. Interpolation constraints: Si(xi)=yi, Si(xi+1)=yi+1.
  2. Continuity of first derivative: Si(xi+1)=Si+1(xi+1).
  3. Continuity of second derivative: Si(xi+1)=Si+1(xi+1).
  4. Boundary conditions: Natural splines assume S(x0)=S(xn)=0; clamped splines specify first derivatives at the endpoints.

This leads to a system of linear equations that we solve using tridiagonal matrix methods—much more stable than the Vandermonde matrix from polynomial interpolation.

Solving the Spline System: The Tridiagonal Matrix

When constructing a cubic spline, we generate a system of equations for the second derivatives S(xi) at each knot. The earlier system can be written in matrix form as a tridiagonal system:

| 2h₀      h₀       0      ...       0     |   | M₁ |   | d₁ |
| h₀    2(h₀+h₁)    h₁     ...       0     |   | M₂ |   | d₂ |
| 0        h₁    2(h₁+h₂)  ...       0     | × | M₃ | = | d₃ |
| ⋮        ⋮        ⋮       ⋱         h₈    |   | ⋮   |    | ⋮  |
| 0        0        0      h₈       2h₈    |   |M₉  |   | d₉ |

where hi=xi+1xi are the interval widths, Mi=S(xi) are the second derivatives at the knots, and di are terms computed from the given data points.

This matrix is:

This system is efficiently solved using the Thomas algorithm, a specialized form of Gaussian elimination with O(n) complexity.

Splines make everything feel structured, precise, and under control. Where high-degree polynomials are chaotic, splines are composed. Where brute-force approaches are inefficient, splines find a way to be elegant.

Final Thoughts

Splines provide a structured and computationally stable approach to interpolation. By segmenting the problem into local polynomial fits, they achieve smooth transitions without the numerical pitfalls of high-degree polynomials. That's why I love them. They get the job done with minimal fuss, proving that you don't need to overcomplicate things to achieve elegance.