Spline regularization

Question

I am fitting some B-splines to data, but the data has a "gap" region where the spline is less constrained by the data. I want to devise a regularization scheme to help prevent the spline from deviating too much from where it should be.

One standard approach to spline regularization, is to minimize the curvature of the spline, approximated by the second derivative integrated over the full interval. For example, if we have a spline $f(t)$, defined by $$ f(t) = \sum_i x_i N_i(t) $$ where $N_i(t)$ are the normalized B-spline basis functions defined with some knot vector. We can then define: \begin{align} \textrm{curvature} &\approx \int dt \left| \frac{d^2 f(t)}{dt^2} \right|^2 \\ &= x^T G x \end{align} with $$ G_{ij} = \int dt \frac{d^2 N_i(t)}{dt^2} \frac{ d^2 N_j(t)}{dt^2} $$ Then we can define a regularized least squares problem as, $$ \min_x || b - A x ||^2 + \lambda^2 x^T G x $$ and solve this, for example with a normal equations approach. The matrix $G$ is not positive definite, so I'm not sure if this can be solved with QR or SVD.

In any case, here are two examples of my results in using this regularization.

The period between 2011 and 2014 is the data gap region, and you can see the unregularized spline solutions (green) have a large deviation from their "baseline". The blue curves show my attempt to regularize by imposing the curvature regularization. However, you can see that the result still shows a bump in the data gap region, and also kills the high-frequency oscillations, which I would like to preserve.

So my question is: can anyone think of a better way to regularize these splines to prevent the large deviation from the baseline in the gap period (2011-2014), but also preserve the high-frequency oscillations?

One possibility is to do a linear fit to the spline in the range 2014-2019, and then minimize the distance to this linear fit in the 2011-2014 gap period. But I don't know how to put such a condition into the usual form $x^T \Lambda x$ for some matrix $\Lambda$.

Another possibility is to somehow high-pass filter the spline, though I don't know how to implement that as a regularization condition either.

Answer 1

As your problem is a local regression problem, I would not use a spline fit, but LO(W)ESS. This estimates $f(x)$ at a sample point $x$ by a weighted least squares fit to the points $x_i$ that are the $k$ nearest neighbors of $x$. For details, see

Cleveland: "Robust locally weighted regression and smoothing scatterplots." Journal of the American Statistical Association 74.368, pp. 829-836, 1979

In R, it is out-of-the box available through the function loess, as in the following example using the R builtin data set mtcars:

alpha <- 0.75  # fraction of points used as neighbors
smoothed <- loess(y ~ x, data=data.frame(x=mtcars$disp, y=mtcars$mpg), span=alpha, family="gaussian", degree=2)
x <- seq(min(mtcars$disp), max(mtcars$disp), by=10)
plot(mtcars$disp, mtcars$mpg)
lines(x, predict(smoothed, x), col="red")