[SciPy-User] interpolate.interp1d - constructing a cubic interpolator is Slllooooooow.

[Evgeni Burovski]

Hi,

First, the "why" question: interp1d(..., kind='cubic') constructs a
continuously differentiable cubic spline. This requires solving an
N-by-N linear system of equations, where N=1000 in your example.

As an implementation detail, cubic spline interpolation can be made to
use banded matrices, but
`interp1d(..., kind='cubic')` uses full matrices under the hood, with
corresponding consequences for the time/memory footprint.

On the other hand, UnivariateSpline (and its equivalent splrep/splev
combo) uses a different, more efficient, formulation of the linear
algebraic problem. This explains what Eraldo noted, that the
UnivariateSpline is faster to construct (and likely faster to evaluate
as well). Notice that to get an interpolating spline from FITPACK
wrappers, you can use `s=0` like so:
>>> UnivariateSpline(x, y, k=3, s=0)

Now, depending on your use case you might not need a continuous
differentiability. If this is the case, you do not really need a
global spline fitting process. For local interpolating schemes, you
can use Pchip or Akima1DInterpolator, or you can build your own cubic
Hermite interpolator using PPoly or BPoly.from_derivatives.
(do *not* use PiecewisePolynomial, if you care about speed at all).

Finally, If you feel adventurous, you can try
https://github.com/scipy/scipy/pull/3174
which replaces the implementation of `interp1d(..., kind='cubic')`,
and adds several new spline constructors for interpolating and LSQ
splines. That PR is not yet merged, so comments welcome.


Cheers,

Evgeni



On Thu, Jan 15, 2015 at 9:49 AM, Eraldo Pomponi
<eraldo.pomponi at gmail.com> wrote:
> Dear Andrew,
>
> On Thu, Jan 15, 2015 at 9:46 AM, Andrew Nelson <andyfaff at gmail.com> wrote:
>>
>> I'm intending to use interpolation in a curvefitting function.  So have
>> been investigating the use of interpolate.interp1d.
>> I'd prefer to use cubic interpolation but it seems to take ages:
>>
>> import numpy as np
>>
>> from scipy.interpolate import interp1d
>>
>> a = np.linspace(-2 * np.pi, 2 * np.pi, 1000)
>>
>> b = np.cos(a)
>>
>> %timeit interp1d(a, b)
>>
>> 10000 loops, best of 3: 71.6 µs per loop
>>
>> %timeit interp1d(a, b, kind='cubic')
>>
>> 1 loops, best of 3: 5.15 s per loop
>>
>>
>> I'm wondering why it takes 5 orders of magnitude (x72000) longer to
>> calculate a cubic interpolator than a linear interpolator?
>
>
> I can reproduce your results but I cannot comment on the reason why there
> exist a so big difference between the two cases. On the other hand,
> following the documentation, I would go for the use of the more recent
> UnivariateSpline
> (http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.interpolate.UnivariateSpline.html#scipy.interpolate.UnivariateSpline)
> class that doesn't have this drawback. On my system, following your example,
> I get:
>
> %timeit spl = UnivariateSpline(a,b,k=3)
>
> 1000 loops, best of 3: 271 µs per loop
>
>
> Cheers,
>
> Eraldo
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>