Коэффициенты получения кубической сплайна get

Question

Скажем, у меня есть два массива на питоне, и я хочу получить (и фактически использовать) кубическую сплайн-интерполяцию между этими точками. (IE: Я хочу интегрировать функцию). Я бы предпочел способ использовать numpy scipy.

Я знаю о scipy.interpolate.interp1d. Однако это только позволяет мне вычисляться точки, для сказать очень просто функция:

Теперь я мог бы сделать что-то просто нравится:

import numpy as np 
import scipy.interpolate 
import matplotlib.pyplot as plt 

y = np.array([0,2,3,4,8,10,12,12,12,10,9,8,7,1,0,0,1,2]) 
x = np.array(range(len(y))) 
xvals = np.linspace(0, len(y)-1, len(y)*100, endpoint = False) 
func = scipy.interpolate.interp1d(x, y, kind = "cubic") 
yvals = func(xvals) 
plt.plot(xvals,yvals) 
plt.plot(x,y, "o") 

Однако я хотел бы дополнительно обработать этот кубический сплайн (т.е. мне нужно чтобы получить интеграцию) .. для ручного делать вещи, мне нужно, чтобы получить факторы, так:

a_i * x^3 + b_i * x^2 + c_i * x + d_i where i goes from 0 to n/3 

(п = число elemetns - это просто определение г-го кубического)

Поэтому я ожидаю список кортежей (или 2d-массива), описывающих все сплайны. - Или способ получить i-й кубик, и действительно, действительно хотел бы получить удобство «x-to-i», чтобы найти, в какой сплайне я сейчас.

(Хотя, конечно, эта последняя проблема - простой поиск первого значения больше, чем ссылка в отсортированном списке - я мог бы сделать это вручную, если понадобится).

Answer 1

Для интерполяции вы можете использовать scipy.interpolate.UnivariateSpline(..., s=0).

Он имеет, помимо прочего, способ integrate.

EDIT: s=0 параметр для конструктора UnivariateSpline заставляет сплайн проходить через все точки данных. Результат в базисе B-сплайнов, вы можете получить узлы и коэффициенты с методами get_coefs() и get_knots(). Формат такой же, как в FITPACK, dierckx, на netlib. Обратите внимание, что формат tck, используемый внутри, на interp1d (который опирается на splmake ATM) и UnivariateSpline (в качестве альтернативы, splrep/splev) несовместим.

EDIT2: вы можете получить кусочно-полиномиальное представление сплайна с PPoly.from_spline --- но это не в соответствии с iterp1d. Используйте splrep (..., s = 0), чтобы получить интерполирующий сплайн, а затем преобразуйте результат.

Answer 2

Просто мысль, которая может помочь. From the docs, вы можете получить кубический сплайн-интерполяции по-другому, которые могут помочь вам:

import numpy as np 
import scipy.interpolate 
import matplotlib.pyplot as plt 

y = np.array([0,2,3,4,8,10,12,12,12,10,9,8,7,1,0,0,1,2]) 
x = np.array(range(len(y))) 
xvals = np.linspace(0, len(y)-1, len(y)*100, endpoint = False) 
func = scipy.interpolate.splrep(x, y, s=0) 
yvals = scipy.interpolate.splev(xvals, func, der=0) 

# display original vs cubic spline representation for security... 
plt.figure() 
plt.plot(x, y, 'x', xvals, yvals, x, y, 'b') 
plt.legend(['Linear', 'Cubic Spline']) 
plt.axis([-0.05, 20, -2, 20]) 
plt.title('Cubic-spline interpolation') 
plt.show() 

Это дает вам доступ к коэффициентам, если вы хотите с помощью

pp = scipy.interpolate.spltopp(func[0][1:-1],func[1],func[2]) 

#Print the coefficient arrays, one for cubed terms, one for squared etc 
print(pp.coeffs) 

Это также дает пример на странице о том, как интегрировать с помощью этого кубического сплайна представления (изменено константы, чтобы удовлетворить вашу ситуацию, я надеюсь - ваш пробег может варьироваться):

def integ(x, tck, constant=0): 
    x = np.atleast_1d(x) 
    out = np.zeros(x.shape, dtype=x.dtype) 
    for n in xrange(len(out)): 
     out[n] = scipy.interpolate.splint(0, x[n], tck) 
    out += constant 
    return out 

const_of_integration = 0 
yint = integ(xvals, func, const_of_integration)