[SciPy-Dev] Splines

[Charles R Harris]

On Thu, Feb 21, 2013 at 2:41 AM, Pauli Virtanen <pav at iki.fi> wrote:

> Hi,
>
> Charles R Harris <charlesr.harris <at> gmail.com> writes:
> > There have been several threads on the list about splines
> > and consolidation of splines. For instance, there are several
> > uniform spline implementations for images and signal processing,
> > various low level functions in fitlib that are unexposed, and
> > perhaps useful altenatives to b-splines for some applications
> > like straight interpolation. For myself, I've started implementing
> > several functions in pure Python with an eye to converting them
> > to Cython once the interface and documentation is in place,
> > mostly for doing things that fitpack doesn't do because it is
> > very integrated, as opposed to supplying a basic toolset.
> > As part of this project, I'd like to get some feedback
> > on which functions people use most and what features they
> > would like to see. I'm not interested in the high level classes
> > at this point, either the current classes or combo functions
> > like fpcurf, but rather a collection of good lower level function
> > that could be used to implement the higher level functions
> > in a more flexible way. Thoughts?
>
> Great! I was going to start on this for 0.13.0, but this should
> speed things up considerably :)
>
>
Heh. The time has come :)


> Overall, I think what we need is are (i) a well-specified spline
> format, and (ii) solid basic functions for evaluating and
> manipulating them, (iii) ditto for tensor product splines.
>

We probably need more than one, I'd guess at least two, b-splines and
piecewise polynomials, with conversion from b-splines to piecewise.


> How to construct the splines (interpolation, smoothing, etc.) should
> be considered as a separate problem. We can turn to FITPACK for
> smoothing splines, but it should not be used for interpolating
> splines.
>
>
Agree. In particular, the construction of splines with optimal knot
positions is an ongoing research area with no generally useful solution as
far as I know, although there are some interesting papers out there.


>     ***
>
> Some misc thoughts on this:
>
> * The spline data format should be documented and set in
> stone as a first step. Users (and future developers) will
> want to toy around with it.
>

We will probably want a least two formats, maybe three: b-splines, uniform
b-splines, and piecewise polynomials. The format should contain the minimal
data needed to specify the spline. For b-splines I don't see a problem with
the (t,c,k) format for 1-D, but for gridded data in higher dimensions we
might want t and k to be lists and c to be a multidimensional array. That
format can also be adapted for uniform splines, but as the latter have
enough simplifying properties that it might be better to have a separate
format for them.


>
> Also, the data format for tensor product N-dim splines needs
> to be set. They are what we are missing, and what people are
> constantly asking for. We don't want them turn to
> `ndimage.map_coordinates`, which is clunky to use.
>
> The Fitpack tck format looks like this:
>
> https://github.com/pv/scipy-work/commits/spline-unify
>
>
That class is already pretty high level ;)


> Currently, there's also a second B-spline data format used in
> scipy.interpolate with different padding, but we may want to
> stick with the FITPACK one, it's probably as good as any.
>
> * Functions for splines with varying x-coordinates are needed.
> Uniform grid splines would be a nice bonus as a speed-gain
> optimization.
>
> * The 1-D spline routines should be able to work over an
> arbitrary axis of multidimensional data. Even better if this
> can be done without reshaping and copying the input data
> (e.g. with Numpy iterators).
>
> This sounds like providing strided 1-D loops for heavy lifting,
> and bolting array iterators on top.
>
> * Functions for integration + differentiation of splines
> as as abstract objects would be useful. For efficiency,
> evaluation of derivatives & integrals probably might need
> to be provided separately.
>

Agree, but I'd make all these functions separate and also have a conversion
function from b-splines to piecewise polynomial representation.


>
> * For tensor product splines, evaluation on a scattered point
> set + on a grid would be useful.
>
> * It will probably be easiest to start from a clean slate,
> rather than trying to reuse scipy.interpolate.
>
>
Yes.


> * FITPACK should not be used as a basis for this work, no
> need to use ancient FORTRAN 77 code for simple stuff. We can
> use its routines for generating smoothing splines, though.
>
>
Not to mention it is hard to read these days with the goto's, minimal
indentation, and portmanteau functions.


> * Routines for constructing interpolating splines --- I think most
> of the time people will use these for simple gridded data interpolation
> rather than smoothing. FITPACK's knot selection is nice when it works,
> but often it doesn't (or requires careful fiddling), so we should have
> something simple and robust as a default algorithm.
>
> * I'm not sure what to do with the various boundary conditions
> and out-of-bounds value handling. It's probably best to leave room
> for various ways to do this...
>
> * Scipy also has two implementations of piecewise polynomials
> --- these should be consolidated into one, too.
>
>
Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20130221/5752e88b/attachment.html>