SciPy interp1d results are different than MatLab interp1

Question

I'm converting a MatLab program to Python, and I'm having problems understanding why scipy.interpolate.interp1d is giving different results than MatLab interp1.

In MatLab the usage is slightly different:

yi = interp1(x,Y,xi,'cubic')

SciPy:

f = interp1d(x,Y,kind='cubic')
yi = f(xi)

For a trivial example the results are the same: MatLab:

interp1([0 1 2 3 4], [0 1 2 3 4],[1.5 2.5 3.5],'cubic')
  1.5000 2.5000 3.5000

Python:

interp1d([1,2,3,4],[1,2,3,4],kind='cubic')([1.5,2.5,3.5])
  array([ 1.5,  2.5,  3.5])

But for a real-world example they are not the same:

x =   0.0000e+000  2.1333e+001  3.2000e+001  1.6000e+004  2.1333e+004  2.3994e+004
Y =   -6   -6   20   20   -6   -6
xi =  0.00000 11.72161 23.44322 35.16484...  (2048 data points)

Matlab:

-6.0000e+000
-1.2330e+001
-3.7384e+000
  ...
 7.0235e+000
 7.0028e+000
 6.9821e+000

SciPy:

array([[ -6.00000000e+00],
       [ -1.56304101e+01],
       [ -2.04908267e+00],
       ..., 
       [  1.64475576e+05],
       [  8.28360759e+04],
       [ -5.99999999e+00]])

Any thoughts as to how I can get results that are consistent with MatLab?

Edit: I understand that there is some latitude in implementation for cubic interpolation algorithms which probably accounts for the differences I'm seeing. It also seems that the original MatLab program that I am converting should have used linear interpolation, so the question is probably moot.

Answer 1

The underlying interpolation method that scipy.interpolate.interp1d and interp1 are different. Scipy uses the netlib fitpack routines, which yields standard, C2 continuous cubic splines. The "cubic" argument in interp1 uses piecewise cubic hermite interpolating polynomials, which are not C2 continuous. See here for an explanation of what Matlab does.

I suspect that is the source of the difference you are seeing.