SciPy的interp1d结果与MatLab的interp1不同。

Question

我正在将一个MatLab程序转换成Python，但是我不理解为什么scipy.interpolate.interp1d给出的结果与MatLab interp1不同。

在MatLab中使用略有不同:

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

SciPy：

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

举个简单的例子，结果是相同的： 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])

但对于一个现实世界的例子，它们并不相同：

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:

->

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

你认为我如何才能得到与MatLab一致的结果？

编辑：我明白三次样条插值算法的实现可能存在一些差异，这可能解释了我看到的差异。另外，原始的MatLab程序应该使用线性插值，所以这个问题可能没有意义。

- LMO

Answer 1

scipy.interpolate.interp1d 和 interp1 使用的基本插值方法是不同的。Scipy 使用 netlib 的 fitpack 程序，其产生标准的 C2 连续的三次样条。 interp1 中的 "cubic" 参数使用分段三次 Hermite 插值多项式，这些多项式不是 C2 连续的。请查看此处以了解 Matlab 所做的内容。我怀疑这就是你看到的差异的来源。

- talonmies

Answer 2

在当前的scipy中使用http://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.PchipInterpolator.html。这将创建通过传递的y=f(x)的单调立方插值，并使用pchip算法确定点上的斜率。

因此，对于每个部分，您传递了（x，y），pchip算法将计算（x，dy/dx），并且只有一个立方体通过这些点，这些点具有已知导数。根据构造，它将具有连续的一阶导数。

- Benn Malengier