Интегральная функция интерполяционного полинома в питоне

Question

У меня есть набор точек, и в этих точках я использовал scipy для вычисления интерполяционного полинома. Я хочу иметь примитив этой функции

self.p=interpolate.interp1d(self.__discreteDistribution['x'], self.__discreteDistribution['y'], kind='cubic')

Я могу легко использовать scipy для вычисления значения интеграла за интервал, используя

integrate.quad(self.p, 0, max)

Вместо этого я хочу иметь примитив self.p(). Я нашел симпатию, но у меня нет аналитической версии моего интерполирующего полинома.

Что бы вы сделали по этому поводу?

Answer 1

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

Способ 1: UnivariateSpline.

In [1]: import numpy as np

In [2]: x = np.arange(8)

In [3]: y = x

In [4]: from scipy.interpolate import interp1d

In [5]: from scipy.interpolate import interp1d, UnivariateSpline

In [6]: spl = UnivariateSpline(x, y, s=0)

In [7]: spl.<TAB>
spl.antiderivative        spl.get_coeffs            spl.roots
spl.derivative            spl.get_knots             spl.set_smoothing_factor
spl.derivatives           spl.get_residual          
spl.ext                   spl.integral              

In [8]: spl.integral(0, 1)
Out[8]: 0.5000000000000001

Две причуды UnivariateSpline: во-первых, использовать s=0 для интерполяции (в отличие от подгонки по методу наименьших квадратов). Во-вторых, остерегайтесь экстраполяции за пределы игры. По умолчанию, UnivariateSpline экстраполирует значения за пределами (это можно контролировать в конструкторе), но .integral предполагает, что сплайн равен нулю.

In [9]: spl.integral(-1, 1)
Out[9]: 0.5000000000000001

Способ 2: сплев, сплит и шина.

In [13]: from scipy.interpolate import splev, splint, splrep

In [14]: tck = splrep(x, y, s=0)

In [15]: splint(0, 1, tck)
Out[15]: 0.5000000000000001

Это эквивалентно использованию UnivariateSpline, только интерфейс немного отличается. Смотрите документы для деталей.

Способ 3: interp1d.

Под капотом, interp1d также использует b-сплайны (если вы не запрашиваете kind='linear' или 'near'), но процедуры оценки отличаются.interp1d создает вызываемый объект, который затем может быть передан интегратору общего назначения.

In [18]: from scipy.interpolate import interp1d

In [19]: interp = interp1d(x, y, kind='cubic')

In [20]: from scipy.integrate import quad

In [21]: quad(interp, 0, 1)
Out[21]: (0.5000000000000024, 5.5511151231258095e-15)

Опять же, следите за выходом за пределы: поведение результата, созданного interp1d, не очень полезно (даже если оно в некоторой степени управляемо).

Answer 2

Тот interp1d method кажется, просто дает вам функциональный объект для функции интерполяции, который позволяет вам оценивать произвольные значения функции.

Глядя на документацию, я не вижу интерфейса с внутренним представлением. Я предполагаю, что это сумма кусочно определенных кубических полиномов. В этом случае ваш примитив будет суммой кусочно определенных квадратичных полиномов. Это действительно было бы полезно для вас?

Помимо непосредственного вычисления сплайнов (как предложено в комментариях Askewchan), вы можете попытаться использовать значения функции с approximate_taylor_polynomial ( doc), чтобы пересмотреть и получить poly1d (doc) объекты на каждом подинтервале, а затем использовать poly1d.integ ( doc) для каждого подинтервала, чтобы получить коэффициенты примитива.