Understanding Interpolation with Cubic Splines

Question

My ultimate question: given data, what function is Interpolation with Method->"Spline" creating?

I recall from my graduate Numerical Analysis class years ago that piecewise spline interpolating functions had certain nice properties. The classic construction of natural splines can be found on Wolfram MathWorld here.

A colleague gave me some chemistry-related data and wanted to fit it with a smooth interpolating function. I was frustrated with how seemingly difficult this was to do with Mathematica. Googling "Mathematica spline interpolation" likely takes one to the MathWorld page above, or to obsolete functions like SplineFit, described here.

I found Interpolation, but for some reason the default settings don't produce a smooth curve. (Isn't that a natural desire?) The MathWorld page references BSplineCurve, which doesn't naturally compute an interpolating function, though the dedicated reader will find an example of interpolation if they read to "Applications/Interpolation" which requires the use of additional, unexplained code to compute knots.

After I computed my own natural cubic spline of the data, I found that Interpolation has a Method option where "Possible settings include "Spline" for spline interpolation...".

All three methods of computing a smooth interpolating "spline" produce different results. With the sample data set pts={{1,2},{2,4},{3,1},{4,3},{5,5},{6,2}}, the blue curve is produced with Interpolation[pts, Method->"Spline", InterpolationOrder->3], the black curve is produced using BSplineCurve and the documentation's Interpolation example, and the red, dashed curve is the classic natural spline per the MathWorld algorithm.

Again, my primary question is: what function is Interpolation with Method->Spline computing? A secondary question would be "Why doesn't Mathematica have a built-in natural spline creating function?

Edit: The BSplineCurve code, taken from the documentation:

dist = Accumulate[Table[EuclideanDistance[pts[[i]], pts[[i+1]]], {i, Length[pts]-1}]];    
param = N[Prepend[dist/Last[dist], 0]];
knots = {0, 0, 0, 0,1/3, 2/3, 1, 1, 1, 1};
m = Table[BSplineBasis[{3, knots}, j - 1, param[[i]]], {i, 6}, {j, 6}];
ctrlpts = LinearSolve[m, pts];
Graphics[BSplineCurve[ctrlpts]]

Also, my cubic spline code isn't concise and isn't worth sharing, but I can confirm it produces the same function as ResourceFunction["CubicSplineInterpolation"][pts,"Natural"], as suggested in the comment.

Answer 1

We can confirm that the Interpolation produces a "not-a-knot" spline when using the "Spline" method by comparing the results of different spline end condition implementations from the CubicSplineInterpolation resource function. That function generates cubic splines in natural, clamped, not-a-knot, and other flavors.

pts = {{1, 2}, {2, 4}, {3, 1}, {4, 3}, {5, 5}, {6, 2}};
builtin = Interpolation[pts, Method -> "Spline", InterpolationOrder -> 3];
resfun = ResourceFunction["CubicSplineInterpolation"][pts, #] &;
notaknot = resfun["NotAKnot"];
Plot[
  Evaluate@Through[{builtin, notaknot}[x]],
  {x, 1, 6},
  PlotLegends -> {"built in", "not-a-knot"}
]

The plot above shows that the two results are identical, confirming that the built-in produces a not-a-knot cubic spline.

---

To attempt to answer your second point, "not-a-knot" splines seem to be generally preferred over natural splines nowadays, as mentioned in comments by @JM. This is reflected in the default choices in e.g. Matlab, Mathematica (as shown above), and SciPy's scipy.interpolate.CubicSpline which all default to not-a-knot end conditions for their cubic splines.

For comparison, here are other kinds of cubic splines obtained from the same data:

natural = resfun["Natural"];
clamped = resfun["Clamped"];
Plot[
  Evaluate@Through[{builtin, notaknot, natural, clamped}[x]],
  {x, 1, 6},
  PlotLegends -> {"built in", "not-a-knot", "natural", "clamped"}
]

---

The question How to make a natural cubic spline was answered by @JM with code for naturalSpline. This was his second, most recent (2020) implementation of natural cubic splines that I could find on this site. He had also contributed an older version (see this answer from 2015).