Page 1 of 2

Smoothing splines (clamped spline)

Posted: July 17th, 2019, 7:15 pm
by Alan
There are a zillion resources on splines, but I can't seem to find this. 
I need an algorithm for a smoothing (cubic) spline where the endpoint derivatives are specified (the so-called 'clamped spline').

With no smoothing (so just ordinary interpolation), there are many treatments.
With so-called 'natural' endpoints (second derivative is zero) + smoothing, there are many treatments.

Trying to not reinvent the wheel.
Does anyone know a link or a cite?

Re: Smoothing splines (clamped spline)

Posted: July 17th, 2019, 9:32 pm
by Cuchulainn
Alan. 
My 2018 C++ book has a chapter on this + code. I remember a number of boundary conditons.
I don't have the book here. 1200 pages is a bit too much across the sea :-)

BTW what do you mean by 'smoothing', precisely?

enum CubicSplineBC {SecondDeriv, FirstDeriv};
I have just sent you the code

The best references are Bulirch and Stoer and the magical Ahlberg, Nielsen and Walsh (in my bibliography). I have implemented them in my C++ code.

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 12:47 am
by Alan
Hi Daniel,

Thank you for the code. I haven't opened it yet. Let me know if I will find the solution to the following problem, which will answer your question.

I want to find the piecewise cubic polynomial [$]g(t)[$] that solves the following problem:

(*) minimize [$] \left( \sum_{i=1}^n \{ y_i - g(t_i) \}^2 + \lambda \int_a^b \{ g''(u) \}^2 \, du \right)[$], 

subject to 
(A) given values for [$]g'(a)[$] and [$]g'(b)[$], where [$] a < t_1 < t_2 < \cdots < t_n < b[$], 
(B) given data [$]\{t_i,y_i\}[$] and
(C) given [$]\lambda > 0[$], the 'smoothing' parameter.

Without (A), this is the standard smoothing spline problem, solved by Reinsch (1967) (fulltext here) and with many subsequent treatments. 

BTW, I've got your 2018 book. Correct me if I'm wrong, but it looks like the cubic splines discussed there are just *interpolating* splines; i.e., the solution to (*) when [$]\lambda=0[$]?
.

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 10:23 am
by Cuchulainn
That's clear, Alan. Thanks

So, you want SSSP + (A) together? In that case I reckon that the algorithm can be modified as an optimisation problem. But instead of solving a tridiagonal matrix system we must minimise a nonlinear function for the coefficients of the cubic spline. Looks benign enough especially since the integral of the 2nd derivative can be computed exactly.
The 'easiest' would be a derivative-free method, e.g. DE.

Does this make sense?

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 1:19 pm
by Cuchulainn
Here is Python code (MM can also do it for sure). The bigger the smooth the better the fit.
import numpy as np
import matplotlib.pyplot as plt

import csaps

np.random.seed(1234)

x = np.linspace(-5., 5., 25)
y = np.exp(-(x/2.5)**2) + (np.random.rand(25) - 0.2) * 0.3

sp = csaps.UnivariateCubicSmoothingSpline(x, y, smooth=0.85)

xs = np.linspace(x[0], x[-1], 150)
ys = sp(xs)

plt.plot(x, y, 'o', xs, ys, '-')
plt.show()

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 1:36 pm
by FaridMoussaoui
Alan, check "Numerical Analysis" by Burden, Faires & Reynods (1981). There is something about a clamped cubic spline algorithm.

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 3:21 pm
by rmodafferi
Hi Alan, you can find the theory in Green & Silverman (1994), and Turlach (2005) the theory and a very readable implementation in Matlab in Martinez & Martinez (2016)

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 5:30 pm
by Alan
Thanks, guys but it looks like you are all citing either (i) cubic spline interpolation (the interpolant goes through the data points) or (ii) standard smoothing splines (Green & Silverman, for example, which have a so-called "natural" boundary condition [$]g''(a) = g''(b) = 0[$]).

I am asking about what is apparently a non-standard case: smoothing splines with a first derivative condition at endpoints. It's problem (*) above with conditions (A) + (B) + (C).

I will guess this non-standard case has been solved, but if I can't find a reference in a week or so, I'll try it from scratch.  

 

Re: Smoothing splines (clamped spline)

Posted: July 18th, 2019, 8:40 pm
by FaridMoussaoui
Well, as your are looking for [$]N[$] locally cubic spline, each spline is defined by [$]4[$] parameters. (e.g [$]4 N[$] unknowns).
The cubics need to match at the knots points ([$] 2 N[$] equations). Make also the first derivative & second derivatives continuous at the interior knots points ([$] N -1[$] equations each). That  [$] 4 N - 2[$] equations.
The "natural" boundary conditions or your set first derivatives provide two more equations. That's [$]4 N [$] equations.

Re: Smoothing splines (clamped spline)

Posted: July 19th, 2019, 3:00 pm
by Cuchulainn
My feeling is that I understand the requirement * + A + B + C.
The solver compute the coefficients of the cubic spline. Natural BCs are not very difficult.

So put Farid's 4N unknowns in a (some) NL solver, yes? I don't think it is all that difficult but you probably need to do preprocessing of the input.

Re: Smoothing splines (clamped spline)

Posted: July 19th, 2019, 4:30 pm
by FaridMoussaoui
For Alan specific problem with the smoothing parameter, the problem was already solved in
Carl de Boor (1978), A Practical Guide to Splines (Chapter XIV). While it is solved with "natural" BC
(allowing to have a tridiagonal matrix and adding 2 equations to the linear system), solving for your specified equations adds also 2 equations and now a band matrix.

Just let me know if you want to receive a copy of the book at your optioncity email.

Re: Smoothing splines (clamped spline)

Posted: July 19th, 2019, 4:39 pm
by Cuchulainn
Well, as your are looking for [$]N[$] locally cubic spline, each spline is defined by [$]4[$] parameters. (e.g [$]4 N[$] unknowns).
The cubics need to match at the knots points ([$] 2 N[$] equations). Make also the first derivative & second derivatives continuous at the interior knots points ([$] N -1[$] equations each). That  [$] 4 N - 2[$] equations.
The "natural" boundary conditions or your set first derivatives provide two more equations. That's [$]4 N [$] equations.
I get (N = number of intervals)

4N unknowns

N+1 interpolation points
3(N-1) cnty conditions 
2 BCs

So, the same.

Re: Smoothing splines (clamped spline)

Posted: July 19th, 2019, 4:41 pm
by Cuchulainn
Well, as your are looking for [$]N[$] locally cubic spline, each spline is defined by [$]4[$] parameters. (e.g [$]4 N[$] unknowns).
The cubics need to match at the knots points ([$] 2 N[$] equations). Make also the first derivative & second derivatives continuous at the interior knots points ([$] N -1[$] equations each). That  [$] 4 N - 2[$] equations.
The "natural" boundary conditions or your set first derivatives provide two more equations. That's [$]4 N [$] equations.
I get (N = number of intervals)

4N unknowns

N+1 interpolation points
3(N-1) cnty conditions 
2 BCs

So, the same.

Re: Smoothing splines (clamped spline)

Posted: July 19th, 2019, 4:46 pm
by Alan
For Alan specific problem with the smoothing parameter, the problem was already solved in
Carl de Boor (1978), A Practical Guide to Splines (Chapter XIV). While it is solved with "natural" BC
(allowing to have a tridiagonal matrix and adding 2 equations to the linear system), solving for your specified equations adds also 2 equations and now a band matrix.

Just let me know if you want to receive a copy of the book at your optioncity email.

Thanks, Farid. I've got that book and will (re)-read that chapter with your comment in mind. I bought these books (de Boor and Green & Silverman) many years ago for a project, but now need to revisit them for this variant problem. 

Re: Smoothing splines (clamped spline)

Posted: July 21st, 2019, 10:13 pm
by jherekhealy
Yes, the BFIT routine of de Boor will do what you're after. FYI, this is part of the NSWC Fortran library https://raw.githubusercontent.com/jacob ... src/nswc.f