DTU immoptibox

 

3.1. Cubic Splines
 
The toolbox contains a number of functions that can be used to determine a cubic spline, which either interpolates given points or is found as a weighted least squares fit to a set of given data points. In both cases it is possible to apply quite general boundary conditions. The theory is described in [L5].
Physical spline

The function splinefit is used to determine the spline s so that it fits to (or interpolates) given data points, possibly with constraints of the form,
  d11s(a) + d12s'(a) + d13s"(a)   =   d14
d21s(b) + d22s'(b) + d23s"(a)   =   d24
where [a,b] is the range of the spline knots.

splineval and splinedif can be used to evaluate the spline and derivatives of it, respectively.

The spline s is represented by a struct  S  with fields
fail Performance indicator. The other fields are only significant if  S.fail = 0 .
  fail = 0 :   No problems.
  1 :   is not a real valued vector of length at least 2.
  2 :   Knots are not in increasing order.
  3 :   Some data abscissa is not real or is outside  [x(1), x(end)].
  4 :   No data ordinates are given.
  5 :   Too few active data points (ie points with strictly positive weights).
  6 :   The Schoenberg-Whitney conditions are not satisfied.
  7 :   When given, the boundary condition matrix D must be a 2*4 matrix.
x Knots.
c Coefficients in the B-spline representation of s.
pp Piecewise polynomial representation of s.
sdy Estimated standard deviation of data.
sdc Estimated standard deviation of the coefficients  c .
 
Go to the top Return to immoptibox main page

 

3.1.1. User's guide to   splinefit
This function determines a cubic spline as a weighted least squares fit to given data points, possibly with given boundary conditions. The algorithm is described in Section 5.1 of [L5]. If the number of data points is equal to the number of degrees of freedom, one gets the interpolating spline.

Typical calls are
    S = splfit(tyw,x)
    S = splfit(tyw,x,D)

Input parameters
tyw Data points and weights. Array with 2 or 3 columns,
  tyw(:,1) Data abscissas.
  tyw(:,2) Data ordinates.
  tyw(:,3) Weights; only points with positive weights are included.
If tyw holds less than 3 columns, then all weights are set to 1.
x Knots. Must satisfy
    x(1) < x(2) <= ... <= x(end-1) < x(end) .
D 2*4 array with the dij presented in the introduction. a = x(1), b = x(end).
 
Output parameters
S Struct representing the spline, as described in the introduction.
 
Go to the top Return to immoptibox main page

 

3.1.2. User's guide to   splineval
This function can be used to evaluate a cubic spline, computed by splinefit. The evaluation is done by means of the piecewise 3rd order polynomial provided in S.pp.

For arguments outside the knot range we use the 2nd order polynomial obtained by continuation of the spline beyond its end knots.

Typical calls are
  f = splval(S,t)
f = splval(t,S)
The alternative formulation is introduced in order that splineval can be used in conjunction with fminbnd, fzero, quad, and other standard MATLAB function functions.
 
Input parameters
S Struct representing the spline, as described in the introduction.
t Vector with arguments for the spline.
 
Output parameters
f Vector of the same type as t, with fi = s(ti) .
 
Go to the top Return to immoptibox main page

 

3.1.3. User's guide to   splinedif
This function can be used to evaluate derivatives of a cubic spline, computed by splinefit. The evaluation is done by means of a differentiated version of the piecewise 3rd order polynomial provided in S.pp.

Arguments outside the knot range are not allowed.

Typical calls are
  f = splinedif(S,t)
f = splinedif(S,t,d)
 
Input parameters
S Struct representing the spline, as described in the introduction.
t Vector with arguments for the spline.
If any t(i) is outside the knot range, you get an error return.
d Differentiation order. Default d = 1.
 
Output parameters
f Vector of the same type as t, with fi = s(d)(t i) .
 
Go to the top Return to immoptibox main page

 

3.1.4. Example
The data set wild.dat provided in the toolbox was generated by adding "noise" to a function defined on [-1, 1] with slope approximately equal to 0.2 at both end points. This is used at the right hand end in the following MATLAB program, while we use natural spline condition at the other end.

The 6 interior knots are distributed so that they are closest where the function varies most, cf Section 5.2 in [L5].

load wild.dat
D = [0 0 1 0; 0 1 0 .2];
x = [-1 -.25 0 .2 .3 .4 .6 1];
S = splinefit(wild,x,D);
t = linspace(-1,1,201);
subplot(211)
plot(wild(:,1),wild(:,2),'.b',...
  t,splineval(S,t),'-r')
title('Data and fit. n = 7')
subplot(223)
plot(t,splinedif(S,t),'-r'), grid on
xlabel('First derivative')
subplot(224)
plot(t,splinedif(S,t,2),'-r'), grid on
xlabel('Second derivative')
Fit of 'wild.dat'
Note how the position of the knots is seen clearly in the piecewise linear s''.
 
Go to the top Return to immoptibox main page