SSJ
V. 2.6.

umontreal.iro.lecuyer.functionfit
Class SmoothingCubicSpline

java.lang.Object
  extended by umontreal.iro.lecuyer.functionfit.SmoothingCubicSpline
All Implemented Interfaces:
MathFunction, MathFunctionWithDerivative, MathFunctionWithFirstDerivative, MathFunctionWithIntegral

public class SmoothingCubicSpline
extends Object
implements MathFunction, MathFunctionWithFirstDerivative, MathFunctionWithDerivative, MathFunctionWithIntegral

Represents a cubic spline with nodes at (xi, yi) computed with the smoothing cubic spline algorithm of Schoenberg. A smoothing cubic spline is made of n + 1 cubic polynomials. The ith polynomial of such a spline, for i = 1,…, n - 1, is defined as Si(x) while the complete spline is defined as

S(x) = Si(x),         for x∈[xi-1, xi].

For x < x0 and x > xn-1, the spline is not precisely defined, but this class performs extrapolation by using S0 and Sn linear polynomials. The algorithm which calculates the smoothing spline is a generalization of the algorithm for an interpolating spline. Si is linked to Si+1 at xi+1 and keeps continuity properties for first and second derivatives at this point, therefore Si(xi+1) = Si+1(xi+1), S'i(xi+1) = S'i+1(xi+1) and S''i(xi+1) = S''i+1(xi+1).

The spline is computed with a smoothing parameter ρ∈[0, 1] which represents its accuracy with respect to the initial (xi, yi) nodes. The smoothing spline minimizes

L = ρi=0n-1wi(yi-Si(xi))2 + (1 - ρ)∫x0xn-1(S''(x))2dx

In fact, by setting ρ = 1, we obtain the interpolating spline; and we obtain a linear function by setting ρ = 0. The weights wi > 0, which default to 1, can be used to change the contribution of each point in the error term. A large value wi will give a large weight to the ith point, so the spline will pass closer to it. Here is a small example that uses smoothing splines:


   int n;
   double[] X = new double[n];
   double[] Y = new double[n];
   // here, fill arrays X and Y with n data points (x_i, y_i)
   // The points must be sorted with respect to x_i.

   double rho = 0.1;
   SmoothingCubicSpline fit = new SmoothingCubicSpline(X, Y, rho);

   int m = 40;
   double[] Xp = new double[m+1];       // Xp, Yp are spline points
   double[] Yp = new double[m+1];
   double h = (X[n-1] - X[0]) / m;      // step

   for (int i = 0; i <= m; i++) {
      double z = X[0] + i * h;
      Xp[i] = z;
      Yp[i] = fit.evaluate(z);          // evaluate spline at z
   }


Constructor Summary
SmoothingCubicSpline(double[] x, double[] y, double rho)
          Constructs a spline with nodes at (xi, yi), with weights = 1 and smoothing factor ρ = rho.
SmoothingCubicSpline(double[] x, double[] y, double[] w, double rho)
          Constructs a spline with nodes at (xi, yi), with weights wi and smoothing factor ρ = rho.
 
Method Summary
 double derivative(double z)
          Evaluates and returns the value of the first derivative of the spline at z.
 double derivative(double z, int n)
          Evaluates and returns the value of the n-th derivative of the spline at z.
 double evaluate(double z)
          Evaluates and returns the value of the spline at z.
 int getFitPolynomialIndex(double x)
          Returns the index of P, the Polynomial instance used to evaluate x, in an ArrayList table instance returned by getSplinePolynomials().
 double getRho()
          Returns the smoothing factor used to construct the spline.
 Polynomial[] getSplinePolynomials()
          Returns a table containing all fitting polynomials.
 double[] getWeights()
          Returns the weights of the points.
 double[] getX()
          Returns the xi coordinates for this spline.
 double[] getY()
          Returns the yi coordinates for this spline.
 double integral(double a, double b)
          Evaluates and returns the value of the integral of the spline from a to b.
 
Methods inherited from class java.lang.Object
equals, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
 

Constructor Detail

SmoothingCubicSpline

public SmoothingCubicSpline(double[] x,
                            double[] y,
                            double[] w,
                            double rho)
Constructs a spline with nodes at (xi, yi), with weights wi and smoothing factor ρ = rho. The xi must be sorted in increasing order.

Parameters:
x - the xi coordinates.
y - the yi coordinates.
w - the weight for each point, must be > 0.
rho - the smoothing parameter
Throws:
IllegalArgumentException - if x, y and z do not have the same length, if rho has wrong value, or if the spline cannot be calculated.

SmoothingCubicSpline

public SmoothingCubicSpline(double[] x,
                            double[] y,
                            double rho)
Constructs a spline with nodes at (xi, yi), with weights = 1 and smoothing factor ρ = rho. The xi must be sorted in increasing order.

Parameters:
x - the xi coordinates.
y - the yi coordinates.
rho - the smoothing parameter
Throws:
IllegalArgumentException - if x and y do not have the same length, if rho has wrong value, or if the spline cannot be calculated.
Method Detail

evaluate

public double evaluate(double z)
Evaluates and returns the value of the spline at z.

Specified by:
evaluate in interface MathFunction
Parameters:
z - argument of the spline.
Returns:
value of spline.

integral

public double integral(double a,
                       double b)
Evaluates and returns the value of the integral of the spline from a to b.

Specified by:
integral in interface MathFunctionWithIntegral
Parameters:
a - lower limit of integral.
b - upper limit of integral.
Returns:
value of integral.

derivative

public double derivative(double z)
Evaluates and returns the value of the first derivative of the spline at z.

Specified by:
derivative in interface MathFunctionWithFirstDerivative
Parameters:
z - argument of the spline.
Returns:
value of first derivative.

derivative

public double derivative(double z,
                         int n)
Evaluates and returns the value of the n-th derivative of the spline at z.

Specified by:
derivative in interface MathFunctionWithDerivative
Parameters:
z - argument of the spline.
n - order of the derivative.
Returns:
value of n-th derivative.

getX

public double[] getX()
Returns the xi coordinates for this spline.

Returns:
the xi coordinates.

getY

public double[] getY()
Returns the yi coordinates for this spline.

Returns:
the yi coordinates.

getWeights

public double[] getWeights()
Returns the weights of the points.

Returns:
the weights of the points.

getRho

public double getRho()
Returns the smoothing factor used to construct the spline.

Returns:
the smoothing factor.

getSplinePolynomials

public Polynomial[] getSplinePolynomials()
Returns a table containing all fitting polynomials.

Returns:
Table containing the fitting polynomials.

getFitPolynomialIndex

public int getFitPolynomialIndex(double x)
Returns the index of P, the Polynomial instance used to evaluate x, in an ArrayList table instance returned by getSplinePolynomials(). This index k gives also the interval in table X which contains the value x (i.e. such that xk < x <= xk+1).

Returns:
Index of the polynomial check with x in the Polynomial list returned by methodgetSplinePolynomials

SSJ
V. 2.6.

To submit a bug or ask questions, send an e-mail to Pierre L'Ecuyer.