Nonparametric regression

Nonparametric regression estimates a regression function without restricting it to a fixed finite-dimensional family.

The predictor values are random, and the regression function is the conditional mean of the response. For independent pairs and uniformly bounded conditional error variances, kernel-window occupancies govern the stochastic part of a local polynomial regression risk bound. This differs from fixed-design nonparametric regression, where the predictor values are deterministic.

A smoothing parameter controls the amount of averaging or regularization in a fitted function. Increasing the penalty multiplier $\lambda$ in a cubic smoothing spline discourages curvature and generally reduces effective degrees of freedom. Increasing a smoothing bandwidth enlarges local averaging neighbourhoods. Software may use a monotone transformed parameter such as spar, so its numerical value should not be confused with the raw penalty multiplier.

A penalized function estimate balancing a residual loss against an integrated derivative penalty. For squared-error loss and the second derivative roughness penalty, the estimate is a cubic smoothing spline. Higher derivative orders give other spline degrees.

Let $L_{i}g$ be the average of $g$ on successive knot intervals. Minimizing the displayed loss over $C^1[a,b]$ gives a quadratic spline with constant tails outside the first and last knots. To prove the reduction, choose the unique constant-tail quadratic spline $q$ with the same interval averages as a candidate $g$, and put $v=g-q$. Since $q^{\prime \prime }$ is constant on every knot interval, $q^{\prime }$ is continuous and vanishes on both tails, integration by parts gives $\int q^{\prime }{v}^{\prime }=-{\sum }_i{q}^{\prime \prime }{|}_{I_i}{\int }_{I_i}v=0$. Consequently the data loss is unchanged and $\int (g^{\prime }{)}^2=\int (q^{\prime }{)}^2+\int (v^{\prime }{)}^2$. The penalty decreases strictly unless $v$ is constant; its zero interval averages then force $v=0$. The finite-dimensional quadratic optimization on this spline space also gives a unique minimizer because the interval-average map is injective.

For $\lambda >0$ and at least two distinct design points, the unique $C^2[a,b]$ minimizer of squared residuals plus $\lambda J(g)$ is $N(I+\lambda \Gamma {)}^{-1}\mathbf{Y}$. The minimum roughness property of the natural cubic spline interpolant reduces the optimization to the value vector. The matrix $I+\lambda \Gamma$ is a positive-definite matrix, so the normal equations have a unique solution; the interpolation equality case gives uniqueness as a function. Linear tails are part of the natural boundary convention.

For distinct ordered observations $x_i$, let $s_u$ be their natural cubic spline interpolant with values $u_i$. Its roughness is a quadratic form, $\int (s_u^{\prime \prime }{)}^2=u^{T}Ku$, where $K$ is positive semidefinite and its nullspace consists of vectors sampled from affine functions. Minimizing penalized least squares $\Vert y-u{\Vert }^2+\lambda u^{T}Ku$ gives $\hat{u}=(I+\lambda K{)}^{-1}y$. Thus the smoothing matrix is $A_{\lambda }=(I+\lambda K{)}^{-1}$ and its effective degrees of freedom are $\mathrm{tr}(A_{\lambda })$. As $\lambda \downarrow 0$ the fit approaches interpolation and the trace approaches $n$; as $\lambda \to \infty$ it approaches the least-squares affine line and the trace approaches $2$.

The fixed-design nonparametric regression model observes $Y_i=m(x_i)+\sigma {ϵ}_i$ at deterministic design points $x_i$. For independent and identically distributed random variables ${ϵ}_i$ with the standard normal distribution, comparing two mean functions has Kullback-Leibler divergence ${\sum }_i(m_1(x_i)-m_0(x_i){)}^2/(2{\sigma }^2)$.

At design sites $i/n$, let ${\tilde{Y}}_n$ be the step function taking response $Y_i$ on $((i-1)/n,i/n]$. Project it onto an interval-adapted wavelet basis through resolution $J$. Each coefficient is ${\hat{c}}_{\lambda }={\sum }_i{Y}_i{\int }_{I_i}{b}_{\lambda }(t)dt$. Its expectation is the corresponding coefficient of the step function of sampled regression means. For a Haar wavelet basis, dyadic sample size and resolution at most the sampling scale, this is exactly $n^{-1}{\sum }_i{Y}_i{b}_{\lambda }(i/n)$ with right-closed pointwise versions. Independent errors of variance ${\sigma }^2$ give coefficient variance at most ${\sigma }^2/n$: Cauchy-Schwarz inequality on each sampling cell yields ${\sum }_i({\int }_{I_i}{b}_{\lambda }{)}^2\le \Vert b_{\lambda }{\Vert }_2^2/n$.

Hard thresholding keeps detail coefficients whose magnitude exceeds a prescribed threshold and sets the rest to zero, while retaining the coarse approximation coefficients. Soft thresholding instead uses $\mathrm{sgn}(\hat{b})(|\hat{b}|-\tau {)}_{+}$. The simultaneous wavelet coefficient noise bound motivates thresholds of order $\sigma \sqrt{\log n/n}$ in a normalized Gaussian regression model with at most order $n$ details. This is a nonlinear selection rule: coefficients can be removed individually rather than discarding an entire finest resolution.

For coefficient weights $w_{\lambda i}$ and independent mean-zero sub-Gaussian random variables of parameter $\sigma$, put $q_{\lambda }^2={\sum }_i{w}_{\lambda i}^2$. Multiplying moment-generating-function bounds shows that the weighted noise has parameter $\sigma q_{\lambda }$. The exponential Markov inequality, optimized in its parameter, gives tail bound $2\exp (-t^2/(2{\sigma }^2{q}_{\lambda }^2))$. At the displayed thresholds, a union bound over $N$ coefficients gives simultaneous control with probability at least $1-\delta$. No coefficient independence is required. Integrated wavelet weights satisfy $q_{\lambda }\le n^{-1/2}$. Finite variance alone only supports the weaker Chebyshev choice ${\tau }_{\lambda }=\sigma q_{\lambda }\sqrt{N/\delta }$.

Let $x_i=i/n$, $1\le i\le n$, and let $N_x$ count sites within distance $h/2$ of $x\in [0,1]$. If $nh\ge 2$, a nearest site lies in the closed window, so $N_x\ge 1$. A clipped window has length at least $\min (h,1)/2$, and a grid count gives $N_x\ge n\min (h,1)/2-1$. Combining these two bounds yields $N_x\ge n\min (h,1)/4$. For $0<h\le 2$, the stronger clipped-length bound $h/2$ directly gives $N_x\ge \lfloor nh/2\rfloor \ge nh/4$. These estimates include both endpoints of the domain.

At an interior target point in regular fixed-design nonparametric regression, mean functions with Lipschitz bound $L$ have pointwise minimax risk of order $(L{\sigma }^2/n{)}^{2/3}$ for fixed positive $L,\sigma$ as $n\to \infty$. A smoothing bandwidth of order $({\sigma }^2/(n{L}^2){)}^{1/3}$ balances squared bias of an estimator and variance; the two-point lower bound with a triangular bump gives the matching lower rate.

A linear estimator at $x_0$ has the form $\hat{m}(x_0)={\sum }_i{W}_i(x_0)Y_i$, where the weights depend on the design and tuning parameters but not on the observed responses.

A smoothing matrix maps an observation vector linearly to fitted values, $\hat{Y}=SY$. Its diagonal entries quantify each observation's direct contribution to its own fitted value.

Kernel regression estimates a regression function using localized kernel weights rather than a fixed finite-dimensional functional form. The Nadaraya–Watson estimator is the local-constant version; local polynomial regression improves its local approximation by fitting a polynomial.

Local polynomial regression fits a low-degree polynomial by weighted least squares near each target point.

For equally spaced fixed-design nonparametric regression, independent mean-zero uniformly finite-variance errors, a bounded compactly supported kernel, and a uniformly invertible local Gram matrix, degree $\ell \ge s-1$ gives the displayed bound for $s$ bounded derivatives. Kernel weights have absolute sum $O(1)$ and squared sum $O((nh{)}^{-1})$. Polynomial reproduction removes a degree-$s-1$ Taylor polynomial. Without the degree requirement the general bias of an estimator is only $O(h^{\min (s,\ell +1)})$.

Local linear regression fits $a+b(X_i-x)$ to responses near each target $x$ by weighted least squares, then returns the fitted intercept. With weights $w_i=K((X_i-x)/h)$, put $S_k={\sum }_i{w}_i(X_i-x{)}^k$ and $T_k={\sum }_i{w}_i(X_i-x{)}^k{Y}_i$. If $S_0{S}_2-S_1^2>0$, then $\hat{f}(x)=(S_2{T}_0-S_1{T}_1)/(S_0{S}_2-S_1^2)$. Its polynomial reproduction property of local polynomial regression reproduces affine functions, reducing boundary bias relative to a local constant fit.

A nonnegative regression kernel supplies weights $K((x_i-x)/h)$ in a local weighted least squares fit. Unlike a probability density function, its integral need not equal one, because rescaling every weight does not change the fitted coefficients.

A smoothing bandwidth $h>0$ determines the scale over which a regression kernel averages data, or a kernel for density estimation spreads each observation. It has the same units as the observations. In kernel density estimation, changing $h$ changes the visible detail while preserving the total integral one. For a Gaussian density kernel, each bump has standard deviation $h$ and its individual contribution has peak height $1/(nh\sqrt{2\pi })$.

If the smoothing bandwidth is too small relative to the spacings in a statistical sample, the kernel density estimator shows narrow spikes tied to individual observations: apparent peaks may be sampling accidents. If it is too large relative to features of the unknown probability density function, it blurs distinct groups together and spreads fitted probability far from the observations. An intermediate smoothing bandwidth can reveal broad structure, but visual smoothness alone cannot establish accuracy.

The first three panels use exactly $(-1.8,-1.4,-1.1,-0.7,0.7,1.1,1.4,1.8)$ with smoothing bandwidths $0.08$, $0.35$ and $1.2$, respectively. Every curve has area one; the two groups in the observed statistical sample become less distinct as $h$ increases. No true probability density function is specified in this illustration, so the middle panel is illustrative, not a certified optimal fit. The fourth panel uses $(0.1,0.4,0.9)$ and $h=0.5$ to illustrate boundary bias of a kernel density estimator.

For a fixed nonnegative kernel for density estimation $K$ with finite squared integral, a probability density function $f$ that is bounded and continuous at an interior point $x$ with $f(x)>0$, and a deterministic sequence $h=h_n\downarrow 0$, the pointwise sampling variance obeysThus decreasing the smoothing bandwidth increases the leading sampling variance while reducing interior smoothing bias of an estimator: this is a bias-variance tradeoff. The conditions $h_n\to 0$ and $n{h}_n\to \infty$ express two goals, increasingly fine resolution and enough observations contributing locally. Practical methods such as plug-in bandwidth selection use the data to choose $h$; pointwise optimal kernel bandwidth depends on local features of the unknown probability density function. The smoothing bandwidth remains a tuning choice in nonparametric statistics, and knowing the support of a probability distribution can require more than a width adjustment near a boundary.

Estimate the unknown curvature in the AMISE minimizer and substitute it into that formula. A consistent positive curvature estimate makes the selected bandwidth ratio-consistent. Pilot smoothing is needed because differentiating an undersmoothed density estimator magnifies variance. A normal-reference pilot can initialize an estimate of a higher derivative functional; the final density-curvature estimate is computed from the data rather than assuming the density is normal.

Twice integrating by parts gives $R(f^{\prime \prime })=\int f{f}^{(4)}$. Estimate this by averaging a fourth-derivative kernel density estimator at the observations. Its expectation contains diagonal term $K^{(4)}(0)/(n{g}^5)$ and leading smoothing bias $-g^2{\mu }_2(K)R(f^{\prime \prime \prime })/2$. Balancing them yields pilot scale $g\asymp R(f^{\prime \prime \prime }{)}^{-1/7}{n}^{-1/7}$. With suitable derivatives and tails, the estimation error is bounded by $O_p(g^2+(n{g}^5{)}^{-1}+n^{-1/2}+n^{-1}{g}^{-9/2})$, giving an $O_p(n^{-2/7})$ relative-bandwidth bound after substitution. More carefully calibrated bias cancellation can improve that bound.

// Target: statistical-inference.bigb

The uniform smoothing kernel assigns equal weight inside the window $[x-h,x+h]$ and zero weight outside. Its Nadaraya–Watson estimator is the average of observations in that window.

The unit-width centred box kernel is $K(u)=1_{[-1/2,1/2]}(u)$. It has integral one and squared $L^2$ norm one. Its rescaling $K_h(u)=h^{-1}K(u/h)$ averages a function on intervals of length $h$. It is obtained from the uniform smoothing kernel with support $[-1,1]$ by rescaling, and is useful in both kernel density estimation and local polynomial regression.

The local polynomial Gram matrix is the weighted matrix of inner products of the local monomial basis. Its invertibility makes the local fit unique.

Effective kernel weights are the coefficients $W_i(x_0)$ for which a local polynomial estimate can be written $\hat{m}(x_0)={\sum }_i{W}_i(x_0)Y_i$.

Local polynomial weights of degree $p$ reproduce every multivariate polynomial of total degree at most $p$: applying the estimator to its design values returns its value at the target point.

If a degree-$p$ local polynomial is fitted in the scaled coordinate $(x_i-x)/h$, its coefficient of that coordinate divided by $h$ estimates the first derivative. Its effective weights differentiate every polynomial of degree at most $p$ exactly at the target.

This local-constant estimator is one particular kernel regression method.

The Nadaraya–Watson estimator is the degree-zero local polynomial estimateAt a support boundary it generally has first-order bias because it reproduces constants but not linear functions.

For the unit-width symmetric box kernel for density estimation, suppose the design density is positive and differentiable at $x$, the regression function is twice differentiable there, and the conditional variance and design density are bounded near $x$. Write $G_n=(nh{)}^{-1}{\sum }_i(Y_i-m(x))1_{\{|X_i-x|\le h/2\}}$. Symmetry cancels its linear Taylor term, giving $\mathbb{E}{G}_n=O(h^2)$; independence gives $\mathrm{Var}{G}_n=O((nh{)}^{-1})$. The absolute error bound for a ratio estimator with a controlled denominator gives the displayed rate if the small-denominator contribution is negligible. The balance $h\asymp n^{-1/5}$ gives $O(n^{-2/5})$. Interior design positivity matters: support boundaries can prevent the linear cancellation.

A local constant average reproduces constants but need not reproduce a linear regression function near an endpoint. For $x_i=i/n$, $K=1_{[-1,1]}$, $m(x)=x$, and $nh\to \infty$, its deterministic estimate at zero is $(\lfloor nh\rfloor +1)/(2n)\sim h/2$. A local linear fit removes this first-order boundary bias of an estimator.

For fixed $\alpha \in [0,1)$, design $x_i=i/n$, $m\in C^2[0,1]$, a nonzero nonnegative symmetric piecewise $C^1$ regression kernel supported on $[-1,1]$, and $h\to 0$, $n{h}^2\to \infty$, write ${\mu }_{r,\alpha }={\int }_{-\alpha }^1{u}^{r}K(u)du$ and ${\beta }_r={\mu }_{r,\alpha }/{\mu }_{0,\alpha }$. On the rescaled grid of spacing $1/(nh)$, the Riemann sum error for $u^{r}K(u)$ is $O((nh{)}^{-1})$, by bounded variation. The normalized first and second local moments are therefore $h{\beta }_1+O(n^{-1})$ and $h^2{\beta }_2+O(h/n)$. A uniform second-order Taylor expansion gives the displayed bias of an estimator, since $n^{-1}=o(h^2)$. Equivalently it is $h{\beta }_1{m}^{\prime }(0)+h^2(\alpha {\beta }_1+{\beta }_2/2)m^{\prime \prime }(0)+o(h^2)$. The first-order term explains the generic boundary bias of a local constant estimator.

For an equally spaced design with independent errors of variance at most ${\sigma }^2$ and a regression function with Lipschitz constant $L$, the unit-width box kernel estimator averages the $N_x$ observations in the window. Its bias is at most $Lh/2$, and its variance is at most ${\sigma }^2/N_x$. The window occupancy for an equally spaced regression design gives, whenever $nh\ge 2$,In the usual range $h\le 1$, this is $O(\sigma /\sqrt{nh}+Lh)$ uniformly up to the boundary. For arbitrarily large bandwidth, the stochastic error cannot continue to decrease as $1/\sqrt{nh}$: once all observations are included the estimator is simply their mean.

At an interior point where the regression function is differentiable, a symmetric compactly supported regression kernel makes the first-order bias of an estimator vanish asymptotically. In regular fixed-design nonparametric regression with $h\to 0$ and $nh\to \infty$, the local constant estimator then has squared bias of an estimator $o(h^2)$ and variance $R(K)/(nh)+o((nh{)}^{-1})$ for a unit-integral regression kernel. The differentiability remainder is pointwise and need not hold uniformly over a function class.

Isotonic regression estimates a sequence or function subject to an order constraint such as ${\theta }_1\le \cdots \le {\theta }_n$.

For bounded monotone sequences observed with independent regular location noise, the normalized squared-error minimax risk is of order $n^{-2/3}$.