What algorithm used in interp1d function in scipy.interpolate

Question

So i was writing a python program for my numerical course, and I had to code a cubic spline program. So i implement the formula for cubic spline given in books like Numerical methods by Chapra and canale and Numerical mathematics by chenny and kincaid.

so my data is

x=[1.0,3.0,4.0,7.0]
y=[1.5,4.5,9.0,25.5]

Using this data and applying cubic spline I get for x=1.5 , y=1.79122340426

While using this same data but using the scipy function gives:

  >>> scipy.interpolate.interp1d(x, y, kind='cubic')(1.5)
array(1.265624999999932)

So, why is that difference in result? It is obvious that they are not using the same formula. What is the cubic spline formula used in that scipy function? Is it a natural cubic spline formula or something improved? Note: The value 1.2656 is more acurate.

Trusting

Answer 1

EDIT: @ev-br in the comments for this answer provided an important correction to my answer. In fact interp1D spline is not FITPACK based. Check the comments for the link provided by @ev-br.

The Scipy functions for curve fitting are based on FITPACK. Try seeing the documentation on the functions you are using and you'll be able to see a "References" chapter where something like this will appear:

Notes
-----
See splev for evaluation of the spline and its derivatives. Uses the
FORTRAN routine curfit from FITPACK.
If provided, knots `t` must satisfy the Schoenberg-Whitney conditions,
i.e., there must be a subset of data points ``x[j]`` such that
``t[j] < x[j] < t[j+k+1]``, for ``j=0, 1,...,n-k-2``.
References
----------
Based on algorithms described in [1]_, [2]_, [3]_, and [4]_:
.. [1] P. Dierckx, "An algorithm for smoothing, differentiation and
   integration of experimental data using spline functions",
   J.Comp.Appl.Maths 1 (1975) 165-184.
.. [2] P. Dierckx, "A fast algorithm for smoothing data on a rectangular
   grid while using spline functions", SIAM J.Numer.Anal. 19 (1982)
   1286-1304.
.. [3] P. Dierckx, "An improved algorithm for curve fitting with spline
   functions", report tw54, Dept. Computer Science,K.U. Leuven, 1981.
.. [4] P. Dierckx, "Curve and surface fitting with splines", Monographs on
   Numerical Analysis, Oxford University Press, 1993.

These references in particular where taken from the source of fitpack.py on the function "splrep". If you need to do a very thorough comparison between your algorithm and the spline from interp1D just go to the docs:

scipy.interpolate.interp1d

And you'll see a link called [source] right after the definition of the name of the function (so: scipy.interpolate.interp1D [source]). Remember that there's a lot of routine handlers on these functions so be patient while navigating the source.

Ogren

Generally correct, splrep/splev and UnivariateSpline uses FITPACK, but interp1d is not. – Hubris

If it does not than they chose some strange names for their sub-modules because inside the interp1D source you can clearly see the call for "splmake" (if interpolation is spline) which in turn has B = _fitpack._bsplmat(order, xk). I would interpret this as using FITPACK. – Ogren

No, _bspldismat is not using FITPACK the Dierckx fortran library :-). – Hubris

Here's the source: github.com/scipy/scipy/blob/master/scipy/interpolate/src/… – Hubris

Well in that case I'm at a loss. I truly thought scipy splines where, generically speaking, re-implementations or calls to the FITPACK routines. I'll defer to more informed opinion. I'll edit my answer to call into attention your commentary. Sorry for the misinterpretation. – Ogren

yes, it's a bit of a mess, sadly. Most of them are, as you rightly say (BTW, +1 from me). interp1d is an odd one out. – Hubris

Same here, and better safe than sorry. I'm not sure it's relevant to the poster but I would hate to think someone took some wrong assumptions based on my answer. Thanks for the insight. – Ogren