Ordinary least squares fits a single global function, one straight line or one polynomial of fixed degree, to an entire dataset. When the underlying relationship changes shape across the range of \(x\), no single global polynomial degree describes it well everywhere at once: a degree low enough to avoid oscillating between distant regions of the data is too rigid locally, and a degree high enough to bend to a local feature over-fits everywhere else. Locally weighted regression sidesteps the choice of a global polynomial entirely by fitting a new, low-degree model separately near every point of interest, using mostly nearby observations and all but ignoring everything far away.
This is the idea behind LOWESS (locally weighted scatterplot smoothing) and its generalization LOESS (locally estimated scatterplot smoothing), introduced by William Cleveland in 1979 and extended together with Susan Devlin in 1988. Both methods are built directly on top of weighted least squares, a small generalization of the ordinary normal equation already familiar from multiple linear regression.
From Ordinary to Weighted Least Squares
Recall the multiple linear regression setup: a design matrix \(X\in\mathbb{R}^{n\times p}\), a target vector \(y\in\mathbb{R}^n\), and a parameter vector \(\beta\in\mathbb{R}^p\), fit by minimizing the sum of squared residuals \(\lVert y-X\beta\rVert^2\), which leads to the normal equation \(X^TX\beta=X^Ty\). That treats every observation as equally important. Weighted least squares relaxes this by attaching a nonnegative weight \(w_i\) to each observation, penalizing residual \(i\) by \(w_i\) instead of by \(1\):
\[ \mathcal{L}(\beta) = \sum_{i=1}^n w_i\bigl(y_i-x_i^T\beta\bigr)^2 = (y-X\beta)^TW(y-X\beta), \qquad W=\operatorname{diag}(w_1,\dots,w_n). \]
Differentiating with respect to \(\beta\) and setting the gradient to zero,
\[ \nabla_\beta\mathcal{L}(\beta) = -2X^TW(y-X\beta) = 0 \quad\Longrightarrow\quad X^TWX\,\beta = X^TWy, \]
gives the weighted normal equation \(\beta=(X^TWX)^{-1}X^TWy\), identical in form to ordinary least squares with \(X^TX\) and \(X^Ty\) replaced by their weighted counterparts. Setting every \(w_i=1\) recovers the ordinary case exactly; a large weight forces the fitted line to pass close to that observation, and a weight of zero removes it from the fit entirely.
Building a Local Neighborhood
Weighted least squares becomes local regression once the weights are made to depend on distance from a chosen query point \(x_0\): observations near \(x_0\) get a large weight, observations far from it get a small or zero weight. A separate weighted regression is then solved for every \(x_0\) at which a fitted value is needed, each time with a freshly computed set of weights and, in general, a different fitted line.
Two choices control this: how many neighbors to include, and how quickly their influence should fall off with distance.
Span and Neighborhood Size
Rather than fixing an absolute distance, Cleveland's original formulation fixes a span \(0<f\le1\), the fraction of the \(n\) data points to use, giving \(q=\lceil fn\rceil\) neighbors. For a query point \(x_0\), \(d_{\max}(x_0)\) denotes the distance to the \(q\)-th closest of the \(n\) points. Because \(f\) is a fraction rather than an absolute width, the neighborhood automatically shrinks in dense regions of the data and widens in sparse ones. A small span produces a fitted curve that follows the data closely, low bias but high variance; a large span produces a smoother curve, high bias but low variance, the same bias-variance trade-off already familiar from choosing a polynomial degree or a ridge-regression penalty.
The Tricube Weight Function
Within the neighborhood, weights should decrease smoothly from \(1\) at \(x_0\) to \(0\) at the boundary rather than dropping abruptly. Cleveland's choice is the tricube function,
\[ w_i(x_0) = \begin{cases} \bigl(1-u_i^3\bigr)^3, & u_i < 1 \\[4pt] 0, & u_i \ge 1 \end{cases} \qquad\text{where}\qquad u_i = \frac{|x_i-x_0|}{d_{\max}(x_0)}. \]
It is designed to have zero slope at both \(u=0\) and \(u=1\), so the weight profile is smooth rather than kinked at the point of maximum weight or at the neighborhood boundary. Any observation at exactly the \(q\)-th smallest distance receives weight exactly zero, so a neighborhood built from the \(q\) nearest points typically contributes only \(q-1\) points with strictly positive weight.
LOWESS: Locally Weighted Linear Fits
The original 1979 method, LOWESS, fits a weighted linear (degree \(1\)) local model at every data point \(x_i\): design matrix rows \((1,\,x_j)\), weights \(w_j(x_i)\) from the tricube function centered at \(x_i\), solved with the weighted normal equation above. The fitted value at \(x_i\) is the local line evaluated at \(x_i\) itself, \(\hat y_i=\hat\beta_0(x_i)+\hat\beta_1(x_i)\,x_i\); connecting the \(n\) fitted values \((x_i,\hat y_i)\) produces the smoothed curve. A fitted value at any other point is produced the same way, by centering a fresh neighborhood there and solving again.
LOESS: Local Polynomial Regression
Cleveland and Devlin's 1988 generalization, LOESS, allows the local model to be linear or quadratic, and extends the same weighted-least-squares construction to multiple predictors by replacing the scalar distance \(|x_i-x_0|\) with a suitable multivariate distance. In practice the two names are now used almost interchangeably for the univariate case; the description above under LOWESS is exactly the special case of LOESS with local polynomial degree \(1\). Raising the local degree to \(2\) adds a quadratic term \(x^2\) to the design matrix and reduces the bias a straight local line shows near a peak or a trough, at the cost of estimating one additional parameter in every neighborhood, and therefore higher variance.
Robustness: Downweighting Outliers
A single outlier can dominate a local fit exactly as it would dominate an ordinary least-squares fit globally, and a local neighborhood, containing far fewer points than the full dataset, is even more exposed to this. To guard against it, Cleveland's procedure fits the weighted regression once, computes the residuals \(e_i=y_i-\hat y_i\), and derives a second, robustness weight from each residual using Tukey's bisquare function,
\[ \delta_i = B\!\left(\frac{e_i}{6\,\operatorname{median}|e|}\right), \qquad B(u) = \begin{cases}(1-u^2)^2, & |u|<1 \\ 0, & |u|\ge1.\end{cases} \]
Scaling by six times the median absolute residual makes the cutoff adapt to the overall noise level instead of using a fixed threshold. Points fit poorly by the first pass get a small or zero \(\delta_i\); the regression is then solved again using the combined weight \(w_i\delta_i\) in place of \(w_i\), and the whole process is repeated a small, fixed number of times, Cleveland found two robustifying iterations enough in practice. This is an instance of iteratively reweighted least squares, the same general strategy behind fitting robust regression models and, in a different guise, logistic regression.
Worked Example
Consider twenty points \(x_i=0,1,\dots,19\) sampled from a noisy oscillation, together with a query point \(x_0=10\) and a span of \(f=0.4\), so \(q=\lceil0.4\cdot20\rceil=8\) neighbors are used. Sorting the twenty distances \(|x_i-10|\) gives \(d_{\max}(10)=4\); of the eight points within that radius, the two at exactly distance \(4\) receive weight zero, leaving seven with positive tricube weight:
| \(x_i\) | \(y_i\) | \(|x_i-x_0|\) | \(u_i\) | \(w_i\) |
|---|---|---|---|---|
| 7 | -0.23 | 3 | 0.75 | 0.1932 |
| 8 | -1.69 | 2 | 0.50 | 0.6699 |
| 9 | -2.17 | 1 | 0.25 | 0.9539 |
| 10 | -1.75 | 0 | 0.00 | 1.0000 |
| 11 | -1.29 | 1 | 0.25 | 0.9539 |
| 12 | -0.52 | 2 | 0.50 | 0.6699 |
| 13 | 0.10 | 3 | 0.75 | 0.1932 |
Plugging these seven weighted observations into the weighted normal equation gives the sums \(\sum w_i\approx4.634\), \(\sum w_ix_i\approx46.34\), \(\sum w_iy_i\approx-6.556\), \(\sum w_ix_i^2\approx474.15\), \(\sum w_ix_iy_i\approx-62.96\), and solving the resulting \(2\times2\) system for the local slope and intercept yields \(\hat\beta_1\approx0.242\), \(\hat\beta_0\approx-3.833\), so
\[ \hat y(10) = \hat\beta_0+\hat\beta_1\cdot10 \approx -1.41. \]
A single global least-squares line through all twenty points, by contrast, is nearly flat (\(\hat y=-0.018x+0.44\)) because the oscillation averages out over the full range, predicting \(\hat y(10)\approx0.27\), far from the actual local trend of the data near \(x=10\). This is exactly the situation local regression is designed for: the global model is not wrong on average, but it is a poor description of the neighborhood around any particular point when the true relationship bends back and forth across the domain.
When Local Regression Helps, and When It Doesn't
- Local regression trades a single closed-form model for a family of small ones, so there is no fixed set of coefficients to report or interpret; it is primarily a tool for visualizing and smoothing a trend, not for stating a compact predictive equation.
- Every fitted value requires solving its own weighted least-squares problem, so evaluating the curve at \(m\) query points costs roughly \(m\) times as much as the single solve an ordinary or ridge regression needs.
- Near the edges of the data, a neighborhood can only extend in one direction, so the local fit there relies on fewer, one-sided points and is noticeably less reliable than in the interior, the same edge effect seen when extrapolating any regression model beyond its training range.
- A closely related but distinct nonparametric approach is kernel regression, such as the Nadaraya-Watson estimator, which uses a distance-based weight function much like the tricube weight above but averages the \(y_i\) directly instead of fitting a local weighted line or polynomial through them.