Nonparametric Regression

Throughout our journey in regression, we've fit global models—functions defined by a single set of parameters that apply uniformly across the entire input domain. Linear regression uses one slope and intercept everywhere. Even polynomial regression, despite its flexibility, commits to a single polynomial that governs all predictions.\n\nBut what if the relationship between variables changes fundamentally across different regions of the input space? What if the slope is steep for small values of $x$ but nearly flat for large values? What if there are local patterns that no single polynomial can capture?\n\nLocal regression addresses these challenges by fitting separate models in different neighborhoods of the input space. Rather than asking 'What single function best describes all the data?', it asks 'What function best describes the data near this specific point?'\n\nThis deceptively simple shift in perspective unlocks remarkable flexibility—and introduces fascinating new challenges in bias-variance tradeoffs, computational complexity, and the curse of dimensionality.

From Global to Local:\n\nConsider regression as answering the question: Given a new point $x_0$, what should we predict for $y$?\n\nGlobal approach: Fit a model $f(x; \boldsymbol{\beta})$ using all training data, then evaluate $\hat{y} = f(x_0; \hat{\boldsymbol{\beta}})$.\n\nLocal approach: Fit a model using only data points near $x_0$, giving more weight to closer points.\n\nThe key insight is that locally, even complex global relationships often appear simple. A sinusoidal curve looks linear when you zoom in enough. A complicated economic trend might be well-approximated by a line within any small time window.\n\nThe Classical Motivation:\n\nLocal regression emerged from exploratory data analysis in the 1970s-1980s. Researchers like Cleveland (1979) and Cleveland & Devlin (1988) developed LOESS (LOcally Estimated Scatterplot Smoothing) as a way to visualize trends in noisy data without committing to a parametric form.

Intuition Through Visualization:\n\nImagine you want to estimate the trend at point $x_0$. The local regression approach:\n\n1. Define a neighborhood around $x_0$ containing some fraction of the data\n2. Weight the data so nearby points count more than distant ones\n3. Fit a simple model (line or low-degree polynomial) using weighted least squares\n4. Extract the prediction at $x_0$ from this local model\n5. Repeat for every point where you want a prediction\n\nThis sounds computationally expensive (fit a model for every prediction point!)—and it is! But the resulting flexibility is remarkable, and modern computing makes it practical for many applications.

The Weighted Least Squares Framework:\n\nLet $(x_1, y_1), (x_2, y_2), \ldots, (x_n, y_n)$ be our training data. To predict at a target point $x_0$, we solve a weighted least squares problem:\n\n$$\min_{\beta_0, \beta_1, \ldots, \beta_p} \sum_{i=1}^{n} w_i(x_0) \left[ y_i - \beta_0 - \beta_1 (x_i - x_0) - \ldots - \beta_p (x_i - x_0)^p \right]^2$$\n\nwhere $w_i(x_0)$ is the weight assigned to observation $i$ based on its distance from $x_0$.\n\nKey components:\n- Centering at $x_0$: We use $(x_i - x_0)$ rather than $x_i$ so that $\beta_0$ directly gives us the prediction at $x_0$\n- Polynomial degree $p$: Controls local model complexity (typically $p = 0, 1, \text{ or } 2$)\n- Weight function $w_i(x_0)$: Determines how locality is defined\n\nThe prediction at $x_0$ is simply $\hat{y}(x_0) = \hat{\beta}_0$.

Matrix Formulation:\n\nFor local polynomial regression of degree $p$, define:\n\n$$\mathbf{X}_0 = \begin{bmatrix} 1 & (x_1 - x_0) & (x_1 - x_0)^2 & \cdots & (x_1 - x_0)^p \\ 1 & (x_2 - x_0) & (x_2 - x_0)^2 & \cdots & (x_2 - x_0)^p \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & (x_n - x_0) & (x_n - x_0)^2 & \cdots & (x_n - x_0)^p \end{bmatrix}$$\n\n$$\mathbf{W}_0 = \text{diag}(w_1(x_0), w_2(x_0), \ldots, w_n(x_0))$$\n\nThe weighted least squares solution is:\n\n$$\hat{\boldsymbol{\beta}}_0 = (\mathbf{X}_0^T \mathbf{W}_0 \mathbf{X}_0)^{-1} \mathbf{X}_0^T \mathbf{W}_0 \mathbf{y}$$\n\nAnd the prediction at $x_0$ is:\n\n$$\hat{y}(x_0) = \mathbf{e}_1^T \hat{\boldsymbol{\beta}}_0 = \mathbf{e}_1^T (\mathbf{X}_0^T \mathbf{W}_0 \mathbf{X}_0)^{-1} \mathbf{X}_0^T \mathbf{W}_0 \mathbf{y}$$\n\nwhere $\mathbf{e}_1 = [1, 0, 0, \ldots, 0]^T$ extracts the intercept.

The Weight Function:\n\nThe weights $w_i(x_0)$ determine how 'local' the fit is. They should satisfy:\n1. Locality: Large for points near $x_0$, small or zero for distant points\n2. Smoothness: Continuous decay prevents abrupt changes in predictions\n3. Symmetry: $w(u) = w(-u)$ treats left and right neighbors equally\n\nThe Tri-cube Weight Function (Cleveland's standard choice):\n\n$$W(u) = \begin{cases} (1 - |u|^3)^3 & \text{if } |u| < 1 \\ 0 & \text{otherwise} \end{cases}$$\n\nThis function has appealing properties:\n- Compact support: Zero outside $[-1, 1]$, so distant points are completely ignored\n- Smooth at boundaries: Both the function and its first derivative are continuous at $|u| = 1$\n- Peaked at center: Maximum weight at $u = 0$, smooth decay outward

Bandwidth and the Span Parameter:\n\nThe key hyperparameter in LOESS is the bandwidth $h$ (also called the span or $\alpha$), which controls the size of the local neighborhood.\n\nTwo ways to specify bandwidth:\n\n1. Fixed bandwidth $h$: The neighborhood has a fixed width. Points within distance $h$ of $x_0$ receive non-zero weight:\n$$w_i(x_0) = W\left( \frac{x_i - x_0}{h} \right)$$\n\n2. Nearest-neighbor fraction $\alpha$ (LOESS standard): The neighborhood includes a fraction $\alpha$ of the data. The bandwidth $h(x_0)$ is the distance to the $k$-th nearest neighbor, where $k = \lfloor \alpha n \rfloor$:\n$$h(x_0) = |x_{(k)} - x_0|$$\nwhere $x_{(k)}$ is the $k$-th nearest neighbor of $x_0$.\n\nThe nearest-neighbor approach is adaptive: neighborhoods are smaller in dense regions and larger in sparse regions, promoting consistent local sample sizes.

Choosing the Local Polynomial Degree:\n\nThe degree $p$ of the local polynomial affects both fitting and boundary behavior.\n\nLocal constant ($p = 0$): Nadaraya-Watson Estimator\n\nThe simplest case fits a horizontal line (constant) at each point:\n$$\hat{y}(x_0) = \frac{\sum_{i=1}^{n} w_i(x_0) y_i}{\sum_{i=1}^{n} w_i(x_0)}$$\n\nThis is simply a weighted average of nearby $y$ values. Fast to compute but suffers from boundary bias—the estimate is pulled toward the center of the local neighborhood.\n\nLocal linear ($p = 1$): The Standard Choice\n\nFits a line at each point. The key advantage: automatic boundary correction. At boundaries, the local line automatically adjusts its intercept to minimize bias.\n\nLocal quadratic ($p = 2$):\n\nFits a parabola at each point. Better captures curvature in the true function, but requires more data and is more variable.

Boundary Bias Illustrated:\n\nConsider estimating $f(x)$ at $x_0$ near the left boundary. With local constant fitting:\n- Most data points are to the right of $x_0$\n- The weighted average is pulled toward the right\n- If the function is increasing, we underestimate $f(x_0)$\n\nWith local linear fitting:\n- We fit a line $\hat{y} = \beta_0 + \beta_1 (x - x_0)$\n- The line extrapolates the local trend\n- The intercept $\beta_0$ is a bias-corrected estimate at $x_0$\n\nThis boundary correction is one of the most important reasons to use local linear regression over kernel smoothing with local constant fits.

Algorithm: LOESS with Local Linear Regression\n\nGiven training data $(x_1, y_1), \ldots, (x_n, y_n)$, span parameter $\alpha$, and evaluation points $x_0^{(1)}, \ldots, x_0^{(m)}$:\n\nFor each evaluation point $x_0$:\n\n1. Find $k$-nearest neighbors: $k = \lfloor \alpha n \rfloor$\n2. Compute bandwidth: $h(x_0) = $ distance to $k$-th nearest neighbor\n3. Compute weights: $w_i = W\left( \frac{x_i - x_0}{h(x_0)} \right)$ using tri-cube kernel\n4. Solve weighted least squares:\n $$\min_{\beta_0, \beta_1} \sum_{i=1}^{n} w_i [y_i - \beta_0 - \beta_1(x_i - x_0)]^2$$\n5. Extract prediction: $\hat{y}(x_0) = \hat{\beta}_0$

The Outlier Problem:\n\nStandard LOESS uses least squares, which is notoriously sensitive to outliers. A single aberrant point can substantially distort the local fit, affecting predictions in its neighborhood.\n\nCleveland's Robustifying Procedure:\n\nCleveland (1979) proposed an iterative reweighting scheme that down-weights observations with large residuals:\n\n1. Initial fit: Compute $\hat{y}_i$ using standard LOESS\n2. Compute residuals: $r_i = y_i - \hat{y}_i$\n3. Scale residuals: $u_i = r_i / (6 \cdot \text{median}|r_i|)$\n4. Compute robustness weights: $\delta_i = B(u_i)$ where $B$ is the bisquare function\n5. Refit: Multiply original weights by $\delta_i$ and recompute LOESS\n6. Iterate steps 2-5 (typically 2-3 times)

The Bisquare Weight Function:\n\n$$B(u) = \begin{cases} (1 - u^2)^2 & \text{if } |u| < 1 \\ 0 & \text{otherwise} \end{cases}$$\n\nThis function:\n- Gives weight 1 to residuals much smaller than $6 \cdot \text{median}|r|$\n- Gradually down-weights larger residuals\n- Completely ignores residuals larger than $6 \cdot \text{median}|r|$\n\nThe factor of 6 is chosen so that under Gaussian errors, roughly 99.7% of observations get non-zero robustness weights in the initial iteration.

Bias-Variance Tradeoff in Local Regression:\n\nFor local linear regression with bandwidth $h$, the asymptotic bias and variance at a point $x_0$ are:\n\nBias:\n$$\text{Bias}[\hat{f}(x_0)] \approx \frac{h^2}{2} f''(x_0) \mu_2$$\n\nwhere $\mu_2 = \int u^2 K(u) du$ for kernel $K$.\n\nVariance:\n$$\text{Var}[\hat{f}(x_0)] \approx \frac{\sigma^2}{n h p(x_0)} \nu_0$$\n\nwhere $p(x_0)$ is the density of $x$ at $x_0$, and $\nu_0 = \int K(u)^2 du$.\n\nKey insights:\n- Bias increases with $h^2$: larger bandwidth → more smoothing → more bias\n- Bias depends on $f''(x_0)$: flatter regions can tolerate larger bandwidth\n- Variance decreases with $nh$: larger bandwidth → more data → less variance\n- Variance depends inversely on local data density $p(x_0)$

Optimal Bandwidth:\n\nMinimizing the asymptotic mean squared error (MSE) yields the optimal bandwidth:\n\n$$h_{\text{opt}}(x_0) \propto \left( \frac{\sigma^2}{n p(x_0) [f''(x_0)]^2} \right)^{1/5}$$\n\nConvergence rate:\n$$\text{MSE}[\hat{f}(x_0)] = O(n^{-4/5})$$\n\nThis is slower than the parametric rate of $O(n^{-1})$ but faster than the local constant rate of $O(n^{-2/3})$ (for local linear vs. local constant).\n\nPractical implication: The $n^{-4/5}$ rate means you need more data to achieve the same accuracy as a correctly specified parametric model—the price of flexibility.

We've established the complete framework for local regression, from intuition to rigorous implementation. Let's consolidate the key concepts:

What's Next:\n\nLOESS is one approach to local fitting. The next page explores kernel smoothing—a closely related technique that uses fixed bandwidth and has a cleaner theoretical foundation. We'll see how kernel smoothers relate to LOESS and when each is preferred.