interpolate._cubic

Interpolation algorithms using piecewise cubic polynomials.

Module Contents

Classes

PchipInterpolator(self,x,y,axis=0,extrapolate=None) rPCHIP 1-d monotonic cubic interpolation.
Akima1DInterpolator(self,x,y,axis=0) Akima interpolator
CubicSpline(self,x,y,axis=0,bc_type=”not-a-knot”,extrapolate=None) Cubic spline data interpolator.

Functions

pchip_interpolate(xi,yi,x,der=0,axis=0) Convenience function for pchip interpolation.
class PchipInterpolator(x, y, axis=0, extrapolate=None)

rPCHIP 1-d monotonic cubic interpolation.

x and y are arrays of values used to approximate some function f, with y = f(x). The interpolant uses monotonic cubic splines to find the value of new points. (PCHIP stands for Piecewise Cubic Hermite Interpolating Polynomial).

x : ndarray
A 1-D array of monotonically increasing real values. x cannot include duplicate values (otherwise f is overspecified)
y : ndarray
A 1-D array of real values. y’s length along the interpolation axis must be equal to the length of x. If N-D array, use axis parameter to select correct axis.
axis : int, optional
Axis in the y array corresponding to the x-coordinate values.
extrapolate : bool, optional
Whether to extrapolate to out-of-bounds points based on first and last intervals, or to return NaNs.

__call__ derivative antiderivative roots

Akima1DInterpolator CubicSpline BPoly

The interpolator preserves monotonicity in the interpolation data and does not overshoot if the data is not smooth.

The first derivatives are guaranteed to be continuous, but the second derivatives may jump at .

Determines the derivatives at the points , , by using PCHIP algorithm [1]_.

Let , and are the slopes at internal points . If the signs of and are different or either of them equals zero, then . Otherwise, it is given by the weighted harmonic mean

where and .

The end slopes are set using a one-sided scheme [2]_.

[1]F. N. Fritsch and R. E. Carlson, Monotone Piecewise Cubic Interpolation, SIAM J. Numer. Anal., 17(2), 238 (1980). :doi:`10.1137/0717021`.
[2]see, e.g., C. Moler, Numerical Computing with Matlab, 2004. :doi:`10.1137/1.9780898717952`
__init__(x, y, axis=0, extrapolate=None)
roots()

Return the roots of the interpolated function.

_edge_case(h1, m0, m1)
_find_derivatives(y)
pchip_interpolate(xi, yi, x, der=0, axis=0)

Convenience function for pchip interpolation. xi and yi are arrays of values used to approximate some function f, with yi = f(xi). The interpolant uses monotonic cubic splines to find the value of new points x and the derivatives there.

See PchipInterpolator for details.

xi : array_like
A sorted list of x-coordinates, of length N.
yi : array_like
A 1-D array of real values. yi’s length along the interpolation axis must be equal to the length of xi. If N-D array, use axis parameter to select correct axis.
x : scalar or array_like
Of length M.
der : int or list, optional
Derivatives to extract. The 0-th derivative can be included to return the function value.
axis : int, optional
Axis in the yi array corresponding to the x-coordinate values.

PchipInterpolator

y : scalar or array_like
The result, of length R or length M or M by R,
class Akima1DInterpolator(x, y, axis=0)

Akima interpolator

Fit piecewise cubic polynomials, given vectors x and y. The interpolation method by Akima uses a continuously differentiable sub-spline built from piecewise cubic polynomials. The resultant curve passes through the given data points and will appear smooth and natural.

x : ndarray, shape (m, )
1-D array of monotonically increasing real values.
y : ndarray, shape (m, …)
N-D array of real values. The length of y along the first axis must be equal to the length of x.
axis : int, optional
Specifies the axis of y along which to interpolate. Interpolation defaults to the first axis of y.

__call__ derivative antiderivative roots

PchipInterpolator CubicSpline PPoly

New in version 0.14.

Use only for precise data, as the fitted curve passes through the given points exactly. This routine is useful for plotting a pleasingly smooth curve through a few given points for purposes of plotting.

[1] A new method of interpolation and smooth curve fitting based
on local procedures. Hiroshi Akima, J. ACM, October 1970, 17(4), 589-602.
__init__(x, y, axis=0)
extend(c, x, right=True)
from_spline(tck, extrapolate=None)
from_bernstein_basis(bp, extrapolate=None)
class CubicSpline(x, y, axis=0, bc_type="not-a-knot", extrapolate=None)

Cubic spline data interpolator.

Interpolate data with a piecewise cubic polynomial which is twice continuously differentiable [1]_. The result is represented as a PPoly instance with breakpoints matching the given data.

x : array_like, shape (n,)
1-d array containing values of the independent variable. Values must be real, finite and in strictly increasing order.
y : array_like
Array containing values of the dependent variable. It can have arbitrary number of dimensions, but the length along axis (see below) must match the length of x. Values must be finite.
axis : int, optional
Axis along which y is assumed to be varying. Meaning that for x[i] the corresponding values are np.take(y, i, axis=axis). Default is 0.
bc_type : string or 2-tuple, optional

Boundary condition type. Two additional equations, given by the boundary conditions, are required to determine all coefficients of polynomials on each segment [2]_.

If bc_type is a string, then the specified condition will be applied at both ends of a spline. Available conditions are:

  • ‘not-a-knot’ (default): The first and second segment at a curve end are the same polynomial. It is a good default when there is no information on boundary conditions.
  • ‘periodic’: The interpolated functions is assumed to be periodic of period x[-1] - x[0]. The first and last value of y must be identical: y[0] == y[-1]. This boundary condition will result in y'[0] == y'[-1] and y''[0] == y''[-1].
  • ‘clamped’: The first derivative at curves ends are zero. Assuming a 1D y, bc_type=((1, 0.0), (1, 0.0)) is the same condition.
  • ‘natural’: The second derivative at curve ends are zero. Assuming a 1D y, bc_type=((2, 0.0), (2, 0.0)) is the same condition.

If bc_type is a 2-tuple, the first and the second value will be applied at the curve start and end respectively. The tuple values can be one of the previously mentioned strings (except ‘periodic’) or a tuple (order, deriv_values) allowing to specify arbitrary derivatives at curve ends:

  • order: the derivative order, 1 or 2.
  • deriv_value: array_like containing derivative values, shape must be the same as y, excluding axis dimension. For example, if y is 1D, then deriv_value must be a scalar. If y is 3D with the shape (n0, n1, n2) and axis=2, then deriv_value must be 2D and have the shape (n0, n1).
extrapolate : {bool, ‘periodic’, None}, optional
If bool, determines whether to extrapolate to out-of-bounds points based on first and last intervals, or to return NaNs. If ‘periodic’, periodic extrapolation is used. If None (default), extrapolate is set to ‘periodic’ for bc_type='periodic' and to True otherwise.
x : ndarray, shape (n,)
Breakpoints. The same x which was passed to the constructor.
c : ndarray, shape (4, n-1, …)
Coefficients of the polynomials on each segment. The trailing dimensions match the dimensions of y, excluding axis. For example, if y is 1-d, then c[k, i] is a coefficient for (x-x[i])**(3-k) on the segment between x[i] and x[i+1].
axis : int
Interpolation axis. The same axis which was passed to the constructor.

__call__ derivative antiderivative integrate roots

Akima1DInterpolator PchipInterpolator PPoly

Parameters bc_type and interpolate work independently, i.e. the former controls only construction of a spline, and the latter only evaluation.

When a boundary condition is ‘not-a-knot’ and n = 2, it is replaced by a condition that the first derivative is equal to the linear interpolant slope. When both boundary conditions are ‘not-a-knot’ and n = 3, the solution is sought as a parabola passing through given points.

When ‘not-a-knot’ boundary conditions is applied to both ends, the resulting spline will be the same as returned by splrep (with s=0) and InterpolatedUnivariateSpline, but these two methods use a representation in B-spline basis.

New in version 0.18.0.

In this example the cubic spline is used to interpolate a sampled sinusoid. You can see that the spline continuity property holds for the first and second derivatives and violates only for the third derivative.

>>> from scipy.interpolate import CubicSpline
>>> import matplotlib.pyplot as plt
>>> x = np.arange(10)
>>> y = np.sin(x)
>>> cs = CubicSpline(x, y)
>>> xs = np.arange(-0.5, 9.6, 0.1)
>>> plt.figure(figsize=(6.5, 4))
>>> plt.plot(x, y, 'o', label='data')
>>> plt.plot(xs, np.sin(xs), label='true')
>>> plt.plot(xs, cs(xs), label="S")
>>> plt.plot(xs, cs(xs, 1), label="S'")
>>> plt.plot(xs, cs(xs, 2), label="S''")
>>> plt.plot(xs, cs(xs, 3), label="S")
>>> plt.xlim(-0.5, 9.5)
>>> plt.legend(loc='lower left', ncol=2)
>>> plt.show()

In the second example, the unit circle is interpolated with a spline. A periodic boundary condition is used. You can see that the first derivative values, ds/dx=0, ds/dy=1 at the periodic point (1, 0) are correctly computed. Note that a circle cannot be exactly represented by a cubic spline. To increase precision, more breakpoints would be required.

>>> theta = 2 * np.pi * np.linspace(0, 1, 5)
>>> y = np.c_[np.cos(theta), np.sin(theta)]
>>> cs = CubicSpline(theta, y, bc_type='periodic')
>>> print("ds/dx={:.1f} ds/dy={:.1f}".format(cs(0, 1)[0], cs(0, 1)[1]))
ds/dx=0.0 ds/dy=1.0
>>> xs = 2 * np.pi * np.linspace(0, 1, 100)
>>> plt.figure(figsize=(6.5, 4))
>>> plt.plot(y[:, 0], y[:, 1], 'o', label='data')
>>> plt.plot(np.cos(xs), np.sin(xs), label='true')
>>> plt.plot(cs(xs)[:, 0], cs(xs)[:, 1], label='spline')
>>> plt.axes().set_aspect('equal')
>>> plt.legend(loc='center')
>>> plt.show()

The third example is the interpolation of a polynomial y = x**3 on the interval 0 <= x<= 1. A cubic spline can represent this function exactly. To achieve that we need to specify values and first derivatives at endpoints of the interval. Note that y’ = 3 * x**2 and thus y’(0) = 0 and y’(1) = 3.

>>> cs = CubicSpline([0, 1], [0, 1], bc_type=((1, 0), (1, 3)))
>>> x = np.linspace(0, 1)
>>> np.allclose(x**3, cs(x))
True
[1]Cubic Spline Interpolation on Wikiversity.
[2]Carl de Boor, “A Practical Guide to Splines”, Springer-Verlag, 1978.
__init__(x, y, axis=0, bc_type="not-a-knot", extrapolate=None)
_validate_bc(y, expected_deriv_shape, axis)

Validate and prepare boundary conditions.

validated_bc : 2-tuple
Boundary conditions for a curve start and end.
y : ndarray
y casted to complex dtype if one of the boundary conditions has complex dtype.