Nonparametric Regression and Local Polynomial Estimators
Kernel regression fits a constant near each evaluation point. A local polynomial fit also estimates slope and higher-order terms, allowing its bias to reflect additional smoothness. We derive the weights and their exact polynomial reproduction, then prove pointwise, integrated, and uniform risk bounds. The uniform bound requires Gaussian noise and control between evaluation points.
For a random pair \((X,Y)\) with \(\mathbb E|Y|<\infty\), let \(f(x)=\mathbb E[Y\mid X=x]\) be a version of the regression function. Write \(Y_i=f(X_i)+\xi_i\). Independent identically distributed pairs give \(\mathbb E[\xi_i\mid X_i]=0\); the residual need not be independent of its design value. Below, risk calculations use a deterministic design \(X_i\in[0,1]\) and independent mean-zero errors with the stated moment bounds.
The Nadaraya–Watson estimate is \(\widehat f_n^{\rm NW}(x)=\sum_iY_iK((X_i-x)/h)/\sum_iK((X_i-x)/h)\) when the denominator is nonzero; define it as zero otherwise. Thus it is linear in the responses, with weights depending only on the design and \(x\). A zero denominator is handled by this piecewise definition, not by multiplying an undefined quotient by zero.
There is also a density plug-in interpretation. If the pair has joint density \(p(x,y)\) and marginal \(p_X(x)>0\), then \(f(x)=\int yp(x,y)dy/p_X(x)\). For a kernel with \(\int K=1\), \(\int uK(u)du=0\), and finite absolute first moment, set
The change of variable \(u=(Y_i-y)/h\) gives \(\int K((Y_i-y)/h)dy=h\) and
Therefore \(\int\widehat p_n(x,y)dy=(nh)^{-1}\sum_iK((X_i-x)/h)\) and \(\int y\widehat p_n(x,y)dy=(nh)^{-1}\sum_iY_iK((X_i-x)/h)\). Their ratio is precisely the Nadaraya–Watson estimator.
1. Local polynomial estimators
For a nonnegative kernel with positive local weight sum, the Nadaraya–Watson estimate uniquely minimizes \(\sum_i(Y_i-t)^2K((X_i-x)/h)\). Replace the constant \(t\) by a degree-\(\ell\) polynomial in local coordinates. Put
Define \(\widehat\theta_n(x)\) as the minimizer over \(t\in\mathbb R^{\ell+1}\) of
For fixed \(x\), the fitted function of the nearby coordinate is a polynomial. The overall map \(x\mapsto\widehat f_n(x)\) generally is not a single polynomial, because the weights and fit change with \(x\). The local estimator remains linear in the response vector whenever the local matrix is invertible.
1.1. Pointwise Risk
Let \(u_i=(X_i-x)/h\) and define
After discarding the term independent of \(t\), the least-squares criterion divided by \(nh\) is \(-2t^\top a_{nx}+t^\top B_{nx}t\). Its normal equation is \(B_{nx}\widehat\theta_n=a_{nx}\). Positive definiteness makes the criterion strictly convex and gives a unique minimizer. The following assumptions refer to the bandwidth sequence under consideration.
Assumption 1. For some \(\lambda_0>0\), every sufficiently large \(n\) satisfies \(\inf_{x\in[0,1]}\lambda_{\min}(B_{nx})\ge\lambda_0\).
Assumption 2. For some \(a_0>0\), every interval \(A\subset[0,1]\) satisfies \(n^{-1}\sum_i\mathbf1\{X_i\in A\}\le a_0\max(|A|,1/n)\).
Assumption 3. The kernel is nonnegative, supported on \([-1,1]\), and bounded by \(K_{\max}\). The bounds below also hold for signed kernels if the matrix in Assumption 1 remains positive definite, but nonnegativity suffices for the stated least-squares construction.
Solving the normal equations gives the weights
The inverse is essential. Since \(\|U(0)\|=1\), \(\|B_{nx}^{-1}\|\le\lambda_0^{-1}\), and \(\sup_{|u|\le1}\|U(u)\|\le(\sum_{j\ge0}1/(j!)^2)^{1/2}<2\), we have \(|W_{ni}^*|\le2K_{\max}/(nh\lambda_0)\). Only observations in \([x-h,x+h]\cap[0,1]\) contribute. Their number is at most \(a_0\max(2nh,1)=2a_0nh\) when \(h\ge1/(2n)\). Therefore, with \(C_* =\max(2K_{\max}/\lambda_0,4a_0K_{\max}/\lambda_0)\),
The other essential fact is exact polynomial reproduction. For the coordinate vector \(e_k\), indexed from \(k=0\),
Thus the weights sum to one and annihilate the centered monomials of degrees \(1,\ldots,\ell\), even near the boundary.
Suppose \(f\in\Sigma(\beta,L)\) on \([0,1]\), with \(\ell=\lceil\beta\rceil-1\). The Taylor remainder after degree \(\ell\) is bounded by \(L|z-x|^\beta/\ell!\). For \(\ell=0\) this is the Hölder condition; for \(\ell\ge1\), the integral remainder from the preceding kernel-density note gives this bound since its integral constant is at most \(1/\ell!\). Polynomial reproduction cancels all Taylor terms except \(f(x)\), so
If the errors are independent, centered, and satisfy \(\mathbb E\xi_i^2\le\sigma_{\max}^2\), then
Hence \(\operatorname{MSE}(x)\le q_1^2h^{2\beta}+q_2/(nh)\), uniformly in \(x\) and \(f\) in the class. The minimizing bandwidth for this upper bound is \(h_*=[q_2/(2\beta q_1^2)]^{1/(2\beta+1)}n^{-1/(2\beta+1)}\). Any fixed positive multiple of this rate gives \(\sup_f\sup_x\mathbb E_f|\widehat f_n(x)-f(x)|^2=O(n^{-2\beta/(2\beta+1)})\). Integrating this uniform bound over \([0,1]\), by Tonelli, gives the same rate for expected squared \(L^2[0,1]\) error. Random designs require conditional versions of the assumptions and control of exceptional design events before these become unconditional risk bounds.
Remark 4. If \(K\ge0\) is positive on a set of positive measure, then \(B=\int U(u)U(u)^\top K(u)du\) is positive definite. Indeed, \(v^\top Bv=\int(v^\top U(u))^2K(u)du\). A nonzero vector \(v\) defines a nonzero polynomial, which has only finitely many roots and cannot vanish on a positive-measure set. This proves strict positivity. It does not alone prove Assumption 1: the discrete design must approximate the appropriate integral, and at the boundary that integral is one-sided.
Remark 5. Here is a complete verification for the equally spaced design \(X_i=i/n\). Suppose \(K(u)\ge K_{\min}\mathbf1\{|u|\le\Delta\}\), with \(K_{\min},\Delta>0\), \(\Delta\le1\), and \(h\to0\), \(nh\to\infty\). Then the rescaled design \(u_i=(i/n-x)/h\) has mesh \(d_n=1/(nh)\). For every \(x\le1-h\Delta\), it covers \([0,\Delta]\) up to an endpoint discrepancy of at most \(d_n\). If \(x>1-h\Delta\), it covers \([-\Delta,0]\) up to the same discrepancy for all large \(n\), since \(1-2h\Delta>1/n\).
On either interval, \(q_v(u)=(v^\top U(u))^2\) and \(q_v'\) are bounded uniformly over \(\|v\|=1\). On each mesh cell the Riemann-sum error is at most \(d_n\) times the oscillation, hence \(O(d_n^2)\); summing cells and endpoint errors gives \(O(d_n)\), uniformly over both \(x\) and unit \(v\). Consequently
Each integral is positive by the polynomial argument. Its minimum over the compact unit sphere is strictly positive by continuity, and there are just two intervals. Taking half this minimum proves Assumption 1. Also the grid has at most \(n|A|+1\le2\max(n|A|,1)\) points in an interval, proving Assumption 2 with \(a_0=2\). The requirement is \(nh\to\infty\), not \(nh\to0\).
1.2. Uniform Risk
Now assume independent Gaussian errors \(\xi_i\sim \mathcal{N}(0,\sigma_\xi^2)\), a deterministic design satisfying Assumptions 1–2, and a kernel satisfying Assumption 3 that is globally Lipschitz. All constants may depend on the smoothness, design, noise and kernel bounds, but not on \(n\) or on \(f\) within its class. Lipschitz regularity of the kernel does not restrict the regression smoothness \(\beta\) to one.
We first prove the required maximum bounds. If \(P(|Z_i|>t)\le K_i e^{-D_it^p}\), \(p\ge1\), set \(C_i=((1+K_i)/D_i)^{1/p}\). Tail integration gives
With \(C=\max_iC_i\), the maximum of these nonnegative exponential transforms is at most their sum. Jensen, using convexity of \(z\mapsto e^{z^p}-1\), gives \(\mathbb E\max_{i\le M}|Z_i|\le C[\log(1+M)]^{1/p}\). Independence is unnecessary for this bound.
For centered Gaussian vectors \(\eta_j\in\mathbb R^d\), suppose each coordinate variance is at most \(v^2>0\). Gaussian integration gives \(\mathbb E\exp(\eta_{jk}^2/(4v^2))\le\sqrt2\). Hence Jensen and \(\exp(\max a_{jk})\le\sum\exp(a_{jk})\) give
Zero variance is trivial. The vectors may be dependent.
Let \(S_i(x)=U((X_i-x)/h)K((X_i-x)/h)\). The compactly supported map \(u\mapsto U(u)K(u)\) is globally Lipschitz, with some constant \(L_S\): on \([-1,1]\) use the product rule bounds, and across the support endpoints use Lipschitz continuity and the fact that \(K\) vanishes there. Thus \(\|S_i(x)-S_i(z)\|\le L_S|x-z|/h\). Define
Use the grid \(x_j=j/n^2\), \(j=0,\ldots,n^2\), including both endpoints. Every \(x\) is within \(n^{-2}\) of a grid point. Put \(\eta_j=(nh)^{-1/2}\sum_i\xi_iS_i(x_j)\), of dimension \(d=\ell+1\). The grid and increment bounds give
The normalization of \(\eta_j\) produces \(1/(nh)\), not \(1/(nh)^2\), in the first term. Each coordinate variance is bounded by
Furthermore, Cauchy–Schwarz gives \(\mathbb E(\sum_i|\xi_i|)^2\le n\sum_i\mathbb E\xi_i^2=n^2\sigma_\xi^2\). Substitution in the Gaussian maximum inequality proves
If \(nh\to\infty\), the remainder \((nh)^{-4}=o((nh)^{-1})\). Continuity of the weights under these assumptions also makes the supremum measurable, since it equals a supremum over rational evaluation points. Combining with \(\|b\|_\infty\le q_1h^\beta\), we obtain
Balancing these terms gives \(h=\alpha(\log n/n)^{1/(2\beta+1)}\) and
Thus the root-risk normalization is \(\psi_n=(\log n/n)^{\beta/(2\beta+1)}\); the bandwidth has a different exponent.
Why the logarithm is necessary in general. For the equally spaced Gaussian design, a matching lower bound follows from a finite testing construction. Choose a fixed smooth nonnegative bump \(g\) supported inside \((0,1)\), with maximum one. On \(M=\lfloor1/h\rfloor\) disjoint intervals of length \(h\), place functions \(f_j(x)=a h^\beta g((x-t_j)/h)\), and let \(f_0=0\). A sufficiently small fixed \(a>0\) makes every \(f_j\) belong to \(\Sigma(\beta,L)\): differentiating \(\ell\) times scales its Hölder seminorm by exactly \(h^{\beta-\ell}h^{-(\beta-\ell)}=1\). Two distinct bumps have sup-norm distance \(a h^\beta\).
The Gaussian sample laws satisfy
Take \(h=c(\log n/n)^{1/(2\beta+1)}\), with \(c>0\) small. Then the right side is at most \(\frac14\log M\) for all large \(n\), since \(\log M\sim\log n/(2\beta+1)\). To see the testing consequence explicitly, let \(J\) be uniform on \(\{1,\ldots,M\}\), and draw data from \(P_J\). Its mutual information is at most \(M^{-1}\sum_jD(P_j\|P_0)\): subtracting the mutual information from that average gives the nonnegative divergence of the mixture from \(P_0\). For any classifier with average error probability \(P_e\), conditioning also on its error indicator gives \(H(J\mid\text{data})\le\log2+P_e\log M\). Since \(H(J)=\log M\),
for large \(n\). Convert any function estimator into a classifier by selecting the nearest bump in sup norm, breaking ties deterministically. If its error is below \(a h^\beta/2\), the triangle inequality forces correct classification. Hence its maximal expected squared sup-norm error is at least \((a h^\beta/2)^2P_e\ge a^2h^{2\beta}/8\). This proves the matching order \((\log n/n)^{2\beta/(2\beta+1)}\) for this design, and justifies the minimax claim under a concrete model. ◻