Main Content

Least-squares spline approximation

returns the B-form of the spline `spline`

= spap2(`knots`

,`k`

,`x`

,`y`

) *f* of order
`k`

with the given knot sequence `knots`

for which

(*) y(:,j) = f(x(j)), all j

in the weighted mean-square sense, meaning that the sum

$$\sum _{j}w(j)|y(:,j)-f\left(x(j)\right){|}^{2}$$

is minimized, with default weights equal to 1. The data
values `y(:,j)`

can be scalars, vectors, matrices, or
ND-arrays, and |*z*|^{2} is the sum of
the squares of all the entries of *z*. Data points with the
same site are replaced by their average.

If the sites `x`

satisfy the Schoenberg-Whitney
conditions

$$\begin{array}{l}\text{knots}(j)x(j)\text{knots}(j+k)\\ (**)\text{}j=1,\mathrm{...},\text{length}(x)=\text{length(knots)}-k\end{array}$$

then there is a unique spline of the given order and knot sequence satisfying
(*) exactly. No spline is returned unless (**) is satisfied for some subsequence
of `x`

.

`spap2(`

,
with `l`

,`k`

,`x`

,`y`

) `l`

a positive integer, returns the B-form of a
least-squares spline approximant, but with the knot sequence chosen for you. The
knot sequence is obtained by applying `aptknt`

to an appropriate
subsequence of `x`

. The resulting piecewise-polynomial consists
of `l`

polynomial pieces and has `k-2`

continuous derivatives. If you feel that a different distribution of the
interior knots might do a better job, follow this up
with

sp1 = spap2(newknt(spline),k,x,y));

`spap2({knorl1,...,knorlm},k,{x1,...,xm},y) `

provides a least-squares spline approximation to *gridded*
data. Here, each `knorli`

is either a knot sequence or a
positive integer. Further, `k`

must be an
`m`

-vector, and `y`

must be an
(`r+m`

)-dimensional array, with
`y(:,i1,...,im)`

the datum to be fitted at the
`site`

`[x{1}(i1),...,x{m}(im)]`

, all `i1`

, ...,
`im`

. However, if the spline is to be scalar-valued, then,
in contrast to the univariate case, `y`

is permitted to be an
`m`

-dimensional array, in which case
`y(i1,...,im)`

is the datum to be fitted at the
`site`

`[x{1}(i1),...,x{m}(im)]`

, all `i1`

, ...,
`im`

.

`spap2({knorl1,...,knorlm},k,{x1,...,xm},y,w) `

also lets you specify the weights. In this `m`

-variate case,
`w`

must be a cell array with `m`

entries,
with `w{i}`

a nonnegative vector of the same size as
`xi`

, or else `w{i}`

must be empty, in
which case the default weights are used in the `i`

th
variable.

`spcol`

is called on to provide the
almost block-diagonal collocation matrix
(*B _{j}*,

`slvblk`

solves the linear system
(*) in the (weighted) least-squares sense, using a block QR factorization. Gridded data are fitted, in tensor-product fashion, one variable at a time, taking advantage of the fact that a univariate weighted least-squares fit depends linearly on the values being fitted.