Smoothing noisy data with spline functions

Smoothing splines are well known to provide nice curves which smooth discrete, noisy data. We obtain a practical, effective method for estimating the optimum amount of smoothing from the data. Derivatives can be estimated from the data by differentiating the resulting (nearly) optimally smoothed spline.

We consider the model y _i ( t _i )+ _i , i =1, 2, ..., n , t _i [0, 1], where g W ₂ ^{( m )} ={ f : f , f , ..., f ^{( m 1)} abs. cont., f ^{( m )} ₂[0,1]}, and the { _i } are random errors with E _i =0, E _{i j} = ² _ij . The error variance ² may be unknown. As an estimate of g we take the solution g _n, to the problem: Find f W ₂ ^(m) to minimize $$\frac{1}{n}\sum\limits_{j = 1}^n {(f(t_j ) - y_j )^2 + \lambda \int\limits_0^1 {(f^{(m)} (u))^2 du} }$$ . The function g _n, is a smoothing polynomial spline of degree 2 m 1. The parameter controls the tradeoff between the "roughness" of the solution, as measured by $$\int\limits_0^1 {[f^{(m)} (u)]^2 du}$$ , and the infidelity to the data as measured by $$\frac{1}{n}\sum\limits_{j = 1}^n {(f(t_j ) - y_j )^2 }$$ , and so governs the average square error R( ; g)=R( ) defined by $$R(\lambda ) = \frac{1}{n}\sum\limits_{j = 1}^n {(g_{n,\lambda } (t_j ) - g(t_j ))^2 }$$ . We provide an estimate $$\hat \lambda$$ , called the generalized cross-validation estimate, for the minimizer of R( ) . The estimate $$\hat \lambda$$ is the minimizer of V ( ) defined by $$V(\lambda ) = \frac{1}{n}\parallel (I - A(\lambda ))y\parallel ^2 /\left[ {\frac{1}{n}{\text{Trace(}}I - A(\lambda ))} \right]^2$$ , where y=(y ₁, ..., y _n)^t and A ( ) is the n n matrix satisfying (g _n, ( t ₁), ..., g _n, ( t _n))^t= A ( ) y . We prove that there exist a sequence of minimizers $$\tilde \lambda = \tilde \lambda (n)$$ of EV( ) , such that as the (regular) mesh {t _i} _i=1 ⁿ becomes finer, $$\mathop {\lim }\limits_{n \to \infty } ER(\tilde \lambda )/\mathop {\min }\limits_\lambda ER(\lambda ) \downarrow 1$$ . A Monte Carlo experiment with several smooth g 's was tried with m =2, n =50 and several values of ², and typical values of $$R(\hat \lambda )/\mathop {\min }\limits_\lambda R(\lambda )$$ were found to be in the range 1.01---1.4. The derivative g of g can be estimated by $$g'_{n,\hat \lambda } (t)$$ . In the Monte Carlo examples tried, the minimizer of $$R_D (\lambda ) = \frac{1}{n}\sum\limits_{j = 1}^n {(g'_{n,\lambda } (t_j ) - } g'(t_j ))$$ tended to be close to the minimizer of R(¿) , so that $$\hat \lambda$$ was also a good value of the smoothing parameter for estimating the derivative.