Результаты SciPy interp1d отличаются от результатов MatLab interp1 - PullRequest
9 голосов
/ 21 ноября 2011

Я конвертирую программу 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:

-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, которую я конвертирую, должна была использовать линейную интерполяцию, поэтому вопрос, вероятно, спорный.

Ответы [ 2 ]

11 голосов
/ 21 ноября 2011

Основной метод интерполяции, который scipy.interpolate.interp1d и interp1 различен.Scipy использует процедуры netlib fitpack, которые выдают стандартные непрерывные кубические сплайны C2.«Кубический» аргумент в interp1 использует кусочно-кубические интерполяционные полиномы Эрмита, которые не являются непрерывными в C2.См. здесь для объяснения того, что делает Matlab.

Я подозреваю, что это источник разницы, которую вы видите.

0 голосов
/ 22 января 2014

При текущем использовании scipy http://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.PchipInterpolator.html Это создаст монотонную кубическую интерполяцию пройденного y = f (x) и использует алгоритм pchip для определения уклонов в точках.

Таким образом, для каждого раздела (x, y), который вы передаете, алгоритм pchip вычислит (x, dy / dx), и в этих точках будет только кубическое прохождение 2 точек с известной производной. Для каждой конструкции он будет непрерывным с непрерывной первой производной.

...