snpy.spline2 package

Submodules

snpy.dm15temp.dm15temp module

Spline2.py: wrapper for B. Thijsse et al.’s hyper-spline routines.

Yet another spline interpolation routine. The problem: given a set of experimental data with noise, find the spline with the optimal number of knots.

Solution:

They use the usual kind of routines to determine least-squares splines from a given set of knot points. The problem REALLY boils down to: how many knots do you use? There are two extremes: put a knot point on each data point to get an interpolating spline (which sucks for experimental data with noise). The other extreme is to have the minimal set of knots to define a polynomial of order k (e.g., a cubic). This also sucks. Somewhere between the two extremes is a number of knots that optimally recovers the information in the data and smooths out the noise.

spline2 starts with a large number of knots (interpolating spline) and iteratively removes knots until a figure of merit reaches some prescribed value. In this case, this figure of merit is the Durbin-Watson statistic, which measures the auto- correlation between the residuals of the spline fit.

For more details, see:
  • Barend J. Thijsse et al., “A Practical Algorithm for Least-Squares spline Approximation of Data Containing Noise”, Computers in Physics, vol 12 no. 4 July 1998
  • http://structureandchange.3me.tudelft.nl/
snpy.spline2.spline2.eval_extrema(tck)

Attempts to find the extrema of the spline (where S’(x)==0).

Parameters:tck (3-tuple) – tuple of (knots, coefficients, k) that are returned as the first output of spline2.
Returns:(xextr, yextr, signs) where
  • xextr are the x-positions of the extrema,
  • yextr are S(extr), and
  • signs provice the sign of the 2nd derivative S’‘(extr):
    • signs[i] < 0 –> maximum,
    • signs[i] > 0 –> minimum,
    • signs[i] close to 0 –> inflection point
Return type:3-tuple
snpy.spline2.spline2.eval_inflect(tck)

Attempts to find the inflection points of the spline (where S’‘(x)==0).

Parameters:tck (3-tuple) – tuple of (knots, coefficients, k) that are returned as the first output of spline2.
Returns:(xinflect, yinflect, dyinflect) where
  • xinflect are the x-positions of the inflection points,
  • yminflect are S(xinflect), and
  • dyinflect are S’(xinflect).
Return type:3-tuple
snpy.spline2.spline2.eval_integ(x0, x1, tck)

Evaluates the integral from x0 to x1 of the spline defined by tck.

Parameters:
  • x0 (float) – lower limit of integration
  • x1 (float) – upper limit of integration
  • tck (3-tuple) – tuple of (knots, coefficients, k) that are returned as the first output of spline2.
Returns:

the integration.
Return type:

float
snpy.spline2.spline2.eval_x(value, tck)

Attempts to find the roots of the equation S(x) = value for the spline defined by tck.

Parameters:
  • value (float) – value to solve the root for
  • tck (3-tuple) – tuple of (knots, coefficients, k) that are returned as the first output of spline2.
Returns:

roots of the equation.
Return type:

float array
snpy.spline2.spline2.evalsp(x, tck, deriv=0)

Evaluate a spline computed by spline2.

Parameters:
  • x (float array or scalar) – The spline is evaluated at the points x.
  • tck (3-tuple) – a tuple of (knots, coefficents, k) that are returned as the first output of spline2.
  • deriv (int) – if > 0, compute the deriv-th derivative of the spline
Returns:

evaluated interpolant or derivative thereof.
Return type:

float array
snpy.spline2.spline2.spline2(x, y, w=None, sigma=None, rsigma=None, xrange=None, degree=3, acfsearch=0, acffunc='exp', ksi=None, n=None, allownonopt=1, lopt=None, rejlev=0.05, xlog=0, interactive=0, verbose=0, full_output=0)

Solve for the optimal spline given a set of (x,y) data.

Parameters:
  • w (float array) – specify weights (1/sigma) for each data point
  • sigma (float) – specify an absolute sigma for all data points: w[i] = 1/sigma this will then force a regular chi-square fit instead of DW
  • rsigma (float) – specify a relative sigma for all data ponts: w[i] = 1/(y*rsigma)
  • xrange (2-tuple) – Only use data in interval (xrange[0], xrange[1])
  • degree (int) – degree of the spline (default: 3 for cubic spline)
  • acfsearch (bool) – perform an automated search for autocorrelation.
  • acffunc (str) – Use acffunc as the autocorrelation function. can be one of ‘exp’,’gauss’,’linear’, or ‘sinc’. Default: ‘exp’
  • ksi (float) – Use a specific autocorrelation length equal to ksi
  • n (int) – Only search for autocorrelation on index interval n
  • allownonopt (bool) – Allow splines with non-optimized breakpoints (default True)
  • lopt (int) – Force knot optimization to start with lpot knots.
  • rejlev (float) – Use rejection level on statistical tests of rejlev (default 0.05)
  • xlog (bool) – Take the log10(x) before spline fitting. Default: False
  • verbose (bool) – Lots of output to stdout. Default: False
  • interactive (bool) – Allows the user to choose the optimal spline manually.
  • full_output (bool) – along with tck, return the following statistics: rms, dws (Durbin-Watson statistic), lfin (final number of knot points), ksi (computed auto-correlation length), acffit (indicator of how well the assumed auto- correlation function represents the data),
Returns:

(t, c, k) if full_output=0 ((t,c,k), rms, dws, lfin, ksi, acffit) if full_output=1

  • t: array of lfin+1 knot points

  • c: array of lfin+k-1 spline coefficients

  • k: order of the spline (note: order = degree+1, so this is 4

    for a cubic spline!)

The tuple (t,c,k) can be input to routines such as evalsp().

  • rms: dispersion of the fitted spline
  • dws: Durbin-Watson statistic
  • lfin: final number of knot points
  • ksi: computed auto-correlation length
  • acffit: how well the correlations agree with assumed funcion.
Return type:

(tuple)