= Solution
Assume initially that the design points are distinct and ordered on $[a,b]$. The <cubic smoothing spline> minimizes
$$
\boxed{\sum_{i=1}^n\{Y_i-g(X_i)\}^2+\lambda\int_a^b\{g''(t)\}^2\,dt,\qquad\lambda>0,}
$$
over functions with an absolutely continuous first derivative and square-integrable second derivative. Dividing the residual term by $n$ simply changes the scaling assigned to $\lambda$. The <minimum roughness property of the natural cubic spline interpolant> reduces the problem to <natural cubic splines> with knots at the design points: among functions with any prescribed fitted values, that interpolant has the smallest integrated squared curvature. The minimizer is therefore a <natural cubic smoothing spline>, with natural second-derivative boundary conditions and linear tails when extended.
Choose a full basis $B_1,\ldots,B_n$ for this natural-spline space. Put $B_{ij}=B_j(X_i)$ and $\Omega_{jk}=\int_a^b B_j''(t)B_k''(t)\,dt$. Writing $g(t)=\sum_jc_jB_j(t)$ changes the objective to
$$
\|Y-Bc\|^2+\lambda c^T\Omega c.
$$
Differentiation gives $(B^TB+\lambda\Omega)\widehat c=B^TY$. Although the curvature penalty has an affine nullspace, the full basis at distinct design points makes the combined matrix <positive-definite matrix>[positive definite]. Hence
$$
\boxed{\widehat Y=S_\lambda Y,\qquad
S_\lambda=B(B^TB+\lambda\Omega)^{-1}B^T.}
$$
Equivalently, use the cardinal basis that interpolates the observations, making $B=I$ and $S_\lambda=(I+\lambda\Omega)^{-1}$. These are <smoothing matrices>, symmetric but generally not idempotent. Their <effective degrees of freedom> are $\operatorname{tr}(S_\lambda)$. As $\lambda\downarrow0$ the fit interpolates; as $\lambda\uparrow\infty$ it approaches <ordinary least squares> on the affine functions. With replicated predictors, use a design matrix for distinct knots and observation weights; the basis form remains valid when identifiable.
Useful choices include <Leave-one-out cross-validation>, <generalized cross-validation>, and estimation of smoothness by <restricted maximum likelihood> through a mixed-model representation. For fixed $\lambda$, the <leave-one-out residual identity for a linear smoother> gives
$$
\mathrm{CV}(\lambda)=\frac1n\sum_i\left[\frac{Y_i-\widehat Y_i}{1-(S_\lambda)_{ii}}\right]^2,
\qquad
\mathrm{GCV}(\lambda)=\frac{\|Y-\widehat Y\|^2/n}{[1-\operatorname{tr}(S_\lambda)/n]^2}.
$$
Minimize these over a suitable range. Choose smoothing by predictive performance or a variance-component likelihood, not by arbitrarily forcing a visually appealing curve.
Back to article page
