6
$\begingroup$

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.

interpolating splines of the same data

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.

Improve this question
$\endgroup$
5
  • $\begingroup$ Please include the data and the code you used that generated the non-smooth result you mention. Also please include the code you used to generate the "classic natural spline" from your sample data. $\endgroup$
    – MarcoB
    May 16, 2022 at 0:13
  • $\begingroup$ You might be interested in the resource function CubicSplineInterpolation. It has options for natural, clamped, not-a-knot flavors of cubic splines, and it returns an InterpolatingFunction $\endgroup$
    – MarcoB
    May 16, 2022 at 0:35
  • $\begingroup$ I've updated with the BSplineCurve code. Also, the CubicSplineInterpolation ResourceFunction does what I want. I stumbled on it as I was researching this question, though didn't try it. I don't understand the role of these resource functions; they seem like user-contributed functions, correct? They don't seem like a regular part of the Wolfram Language. $\endgroup$
    – GregH
    May 16, 2022 at 0:52
  • $\begingroup$ OP, the purpose of the (mostly) user contributed resource functions are that they add non-traditional functionality to the WL. They are well-vetted by WRI, and serve to meet (seemingly) niche use cases such as these. You might also be interested in any of the other interpolation resource functions like AkimaSpline or AkimaInterpolation. $\endgroup$
    – CA Trevillian
    May 16, 2022 at 1:20
  • 1
    $\begingroup$ I don't know if it's documented, but it produces a not-a-knot spline, like Matlab's csapi(x, y). $\endgroup$
    – Michael E2
    May 16, 2022 at 1:48

1 Answer 1

Reset to default
7
$\begingroup$

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"}
]

plots of built-in and not-a-knot interpolation results overlap perfectly

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"}
]

plot with the four kinds of splines

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).

Improve this answer
$\endgroup$
4
  • $\begingroup$ A very nice answer. $\endgroup$
    – Michael E2
    May 16, 2022 at 17:13
  • $\begingroup$ @MichaelE2 Thank you! $\endgroup$
    – MarcoB
    May 17, 2022 at 2:36
  • $\begingroup$ Agree with @MichaelE2; well done, @MarcoB. I'll have to read up on not-a-knot splines to see why they are preferable over other splines. $\endgroup$
    – GregH
    May 17, 2022 at 17:44
  • 2
    $\begingroup$ @Greg, mostly, this is because when you approximate a smooth function with cubic splines, the not-a-knot condition is generally fourth-order accurate, while the natural condition is only second-order accurate. (Numerical analysis textbooks that discuss splines will have a more precise statement of this.) $\endgroup$
    – J. M.'s eventual burnout ♦
    May 17, 2022 at 18:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.