Code covered by the BSD License  

Download apps, toolboxes, and other File Exchange content using Add-On Explorer in MATLAB.

» Watch video

Highlights from
Free-knot spline approximation

4.8 | 5 ratings Rate this file 46 Downloads (last 30 days) File Size: 102 KB File ID: #25872 Version: 1.19
image thumbnail

Free-knot spline approximation


Bruno Luong (view profile)


17 Nov 2009 (Updated )

Least squares approximation of 1D data using free-knots spline

| Watch this File

File Information

The purpose of this function is to provide a flexible and robust fit to one-dimensional data using free-knot splines. The knots are free and able to cope with rapid change in the underlying model. Knot removal strategy is used to fit with only a small number of knots.

Optional L2-regularization on the derivative of the spline function can be used to enforce the smoothness.

Shape preserving approximation can be enforced by specifying the lower and upper bounds of the derivative(s) of the spline function on sub-intervals. Furthermore specific values of the spline function and its derivative can be specified on a set of discrete data points.

I did not test QUADPROG engine, but I have implemented it. Any feedback is welcome.


Pseudo Inverse, Multiple Same Size Linear Solver, and Min/Max Filter inspired this file.

MATLAB release MATLAB 7.9 (R2009b)
Other requirements optimization toolbox or QP solvers available at: (QPC) QPC solver is strongly recommended.
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (18)
04 Jul 2016 Elhassan Abdou -ETRO-

Dear Bruno,

I want to force the spline to be straight line between the first fixed knot and the second one.
Is this possible?

Comment only
28 Jan 2015 Gaute Hope

How do you prefer this work to be cited?

Comment only
28 Jan 2015 Gaute Hope  
04 Nov 2014 fasfa

fasfa (view profile)

Thank you, any possible modification to make it work for k=1? if not, do you know any other method for constant piece wise approximation?

Comment only
30 Oct 2014 Bruno Luong

Bruno Luong (view profile)

To Fasta:

Yes for k=1; the gradient wrt knots are dirac like, and the gradient method used by here cannot handle it correctly.

Yes, it's a major work to extend to 2D.

Comment only
30 Oct 2014 fasfa

fasfa (view profile)

I think there is an error when you use piecewise constant (k=1). The results are not the expected. Is it ok?

Also I would like to know if it would be very difficult to modify it and use it for surfaces instead of lines (1D -> 2D).

17 Jun 2014 Chi-Fu

Chi-Fu (view profile)

05 Dec 2013 Markus

Markus (view profile)


I´m always getting the following warning message:

Warning: Options LargeScale = 'off' and Algorithm = 'trust-region-reflective'
conflict. Ignoring Algorithm and running active-set algorithm. To run
trust-region-reflective, set LargeScale = 'on'. To run active-set without this
warning, set Algorithm = 'active-set'.

What am I doing wrong?

Comment only
05 Jul 2012 Bruno Luong

Bruno Luong (view profile)


not exactly like you want but you can enforce the y-value between two knots to be zero:

y = sin(x);
y = y + 0.1*randn(size(y));

nknots = 5;
lo = -inf(1,nknots);
up = +inf(1,nknots);
lo(3) = 0;
up(3) = 0;
shape = struct('p',0,'lo',lo,'up',up);
options = struct('shape', shape,'animation', 1, 'knotremoval','none');

Comment only
03 Jul 2012 Richard

Brilliant package Bruno. Quick question.

I know you can fix any given knot, but is it possible to fix a given not to a y-value but let the least squares find the best x-value for it?

I'm using 4 knots/3 lines to represent my data but I would always like the 3rd knot to have y=0. While this condition is met some of the time (by chance) ideally I would like to enforce it.

Comment only
26 Jan 2012 Bruno Luong

Bruno Luong (view profile)

Nima, it is not rotational invariant. Because the fit is carried out using the least-squares to the ordinate data (y) only.

Comment only
23 Jan 2012 Nima Zarif

Hey Bruno thanks for the awesome package, it works great!
I just have one problem, is this approximation rotation invariant? Since when I intentionally rotate my data points the knots are completely different with increased fit error compared to the baseline profile! please let me know if I am doing sth wrong

Comment only
23 Jun 2011 Doug Weathers

Hi Bruno,

Nice package, many great features!

I am having a problem with the 'startingknots' parameter.

When I type this:


I get this:

BSFK starts
??? Error using ==> sparse
Sparse matrix sizes must be non-negative integers less than MAXSIZE as defined by
COMPUTER. Use HELP COMPUTER for more details.

Error in ==> BSFK>BuildDineqMat at 1623
D = sparse(row,col,val,m,n);

Error in ==> BSFK>UpdateConstraints at 2039
[D LU X] = BuildDineqMat(t, knotidx, k, shape);

Error in ==> BSFK>InitPenalization at 2077
smoothing = UpdateConstraints(smoothing, t, shape, pntcon, periodic);

Error in ==> BSFK at 393
smoothing = InitPenalization(y, t, k, d, lambda, p, regmethod, ...

However, if I use 'chebyschev' instead of a vector of starting knots, it works. But it misses some important knots.

Any ideas?



Comment only
03 Oct 2010 Bruno Luong

Bruno Luong (view profile)

Michael, I'll reply in the appropriate place (minmaxfilt)

Comment only
30 Sep 2010 Michael Teo

Does it handle NaN data?

ePeriod = 3;
aData = [ 5;1;3;NaN;8;2;3;NaN;1;9 ];
minmaxfilt(eData, ePeriod, 'max', 'valid')];

Actual output:

If we take NaN as a empty data, the expected output is:

08 Jun 2010 Bruno Luong

Bruno Luong (view profile)

Just discover an issue with continuous regularization. In the mean time, please use the discrete regularization

Comment only
01 Jun 2010 Bruno Luong

Bruno Luong (view profile)

Periodic spline is now available

Comment only
09 Jan 2010 x yuij

x yuij (view profile)

18 Nov 2009 1.1

Update description, more options added to control the fit, discrete regularization

19 Nov 2009 1.2

Remove NaN data before fitting, change TRY/CATCH ME syntax for better compatibility (tested under 2006B), estimate automatic of the noise standard deviation

28 Nov 2009 1.3

Change title and description

04 Dec 2009 1.4

A major enhancement with shape preserving splines

08 Dec 2009 1.5

Point-wise constraints. Discover an error of the Jacobian formula in [Schutze/Schwetlick 97] paper, modify the calculation accordingly. This concern only the constrained fitting.

09 Dec 2009 1.6

Correct another bug in the Jacobian calculation (constrained case)

09 Dec 2009 1.7

Change the description.

10 Dec 2009 1.8

Singular constraints will issue a warning (instead of an error). Refine the Gauss-Newton direction. Fix few minor bugs.

14 Dec 2009 1.9

Correct a bug in UpdateConstraints that did not update the knot positions. Precasting data to double. Update more frequently the scaling matrix. Reduce the Lagrange's tolerance to detect active set of QPC solver

18 Mar 2010 1.10

fixed small bug when calling QP engine minqdef

10 Apr 2010 1.11

fix a bug with parsing k and nknots
Spline order can be as low as k=1 (piecewise constant fit)

31 May 2010 1.12

New feature: Periodic spline

01 Jun 2010 1.13

Remove some redundant code, modify test program

07 Jun 2010 1.14

A more robust conversion in pp form is implemented

07 Jun 2010 1.15

Fix a small bug (eigs with 'sa' option requires true symmetric matrix, which is now always the case by symmetrizing)

10 Jun 2010 1.16

Fix the bug for continuous regularization

05 Jul 2011 1.18

Fix the bug when starting knots are provided

06 Jul 2012 1.19

Fix a bug when checking for knot collision

Contact us