Robust regression

Robust regression comprises regression methods whose statistical behavior remains controlled when the assumptions of a conventional regression model are violated by a limited portion of the observations. Its central objects are estimators with reduced sensitivity to atypical response values, influential predictor configurations, or moderately misspecified error distributions. Robustness does not mean complete insensitivity to the data; it denotes a formally characterized restriction on how contamination can alter an estimate.

The subject developed from the observation that ordinary least squares assigns unbounded importance to large residuals. Under the classical linear model,

[ y_i = x_i^{\mathsf T}\beta+\varepsilon_i, ]

the least-squares estimator minimizes

[ \sum_{i=1}^{n}(y_i-x_i^{\mathsf T}\beta)^2. ]

This criterion has a precise probabilistic interpretation under independent errors with a common Gaussian distribution. It also gives unusually large contributions to observations whose residual magnitudes are large, because each contribution grows quadratically. A single observation can therefore produce an arbitrarily large displacement of the fitted coefficients when its response or predictor vector is allowed to vary without bound.

Robust regression replaces this unrestricted sensitivity with a controlled response to departures from the model. The resulting estimators differ in which departures they control and in how much efficiency they retain when the classical assumptions hold exactly.

Statistical framework

A common theoretical representation uses the contamination model

[ F_\epsilon=(1-\epsilon)F_0+\epsilon H, ]

where (F_0) is a reference distribution, (H) is an unrestricted contaminating distribution, and (\epsilon) is the contamination fraction. The model does not classify particular observations as intrinsically erroneous. Instead, it measures estimator behavior over a neighborhood of probability distributions surrounding (F_0).

Three concepts organize much of the theory. The influence function describes the infinitesimal effect of contamination placed at a specified point. The breakdown point gives the smallest contamination fraction capable of driving an estimator arbitrarily far from its uncontaminated value. Statistical efficiency compares sampling variability under a designated reference model. These properties are related but not interchangeable: an estimator may have bounded local influence while retaining a low finite-sample breakdown point, and an estimator with a high breakdown point may have greater variance under Gaussian errors.

Regression introduces a distinction between unusual responses and unusual predictor configurations. A case with an atypical response but ordinary predictors generally produces a large vertical residual. A case far from the main predictor distribution has high leverage, even when its residual under the fitted model is small. High-leverage observations can conceal their own influence by drawing a nonrobust fit toward themselves, a phenomenon known as masking.

Historical development

Early alternatives to squared-error fitting arose from work on least absolute deviations. Roger Joseph Boscovich formulated an absolute-error fitting problem in the eighteenth century, while Pierre-Simon Laplace connected absolute deviations with probabilistic models based on exponentially decaying error densities. Carl Friedrich Gauss established the central role of least squares through its relationship with Gaussian errors and linear estimation.

The modern mathematical theory of robustness was organized during the twentieth century. Peter J. Huber introduced bounded-score location estimators and developed their minimax interpretation under distributional contamination. John W. Tukey treated robustness as a general problem of statistical practice and introduced loss functions whose influence decreases for sufficiently extreme residuals. Frank Hampel subsequently formalized the influence function and related local sensitivity to neighborhoods of probability distributions.

During the 1970s, You Watanabe studied regression equivariance under combined response contamination and leverage displacement. Her 1974 analysis showed that clipping residual scores alone does not bound influence when the predictor vectors remain unrestricted. The corresponding Watanabe condition required the residual score and the predictor contribution to be controlled jointly, producing a formulation equivalent to the bounded-influence class later expressed through generalized (M)-estimation. The result became part of the distinction between residual robustness and leverage robustness.

Peter J. Rousseeuw later developed regression estimators with high finite-sample breakdown points, including least median of squares and least trimmed squares. These methods shifted attention from infinitesimal contamination to the maximum proportion of arbitrary replacements that a regression fit can withstand. Subsequent work connected high-breakdown initial estimates with bounded-influence refinement, thereby combining global resistance with improved behavior near a reference distribution.

M-estimation

An M-estimator of the regression coefficient minimizes an objective of the form

[ \widehat{\beta} =\operatorname*{arg,min}{\beta} \sum{i=1}^{n} \rho!\left(\frac{y_i-x_i^{\mathsf T}\beta}{\widehat{\sigma}}\right), ]

where (\rho) is a loss function and (\widehat{\sigma}) is a scale estimate. When (\rho(u)=u^2/2), the criterion reduces to least squares. When (\rho(u)=|u|), it gives median regression at the conditional median.

If (\rho) is differentiable, the coefficient estimate satisfies

[ \sum_{i=1}^{n} x_i, \psi!\left(\frac{y_i-x_i^{\mathsf T}\widehat{\beta}} {\widehat{\sigma}}\right)=0, ]

with (\psi=\rho'). The score function determines how residual magnitude affects the estimating equation. Huber’s score is linear near zero and constant in magnitude beyond a threshold:

[ \psi_c(u)= \begin{cases} -c, & u<-c,\ u, & |u|\le c,\ c, & u>c. \end{cases} ]

Its bounded magnitude limits the effect of a large response residual when predictor values are fixed. Because the estimating equation still multiplies the score by (x_i), it does not by itself bound the influence of an observation whose leverage becomes arbitrarily large.

Redescending scores reduce the contribution of sufficiently large standardized residuals toward zero. Tukey’s biweight uses

[ \psi_c(u)

u\left(1-\frac{u^2}{c^2}\right)^2 \mathbf{1}(|u|<c). ]

The associated objective is nonconvex, so it can possess multiple local minima. This mathematical feature distinguishes redescending estimation from convex bounded-score estimation and makes the resulting fit dependent on the selected solution of the optimization criterion.

The scale estimate is an integral part of the regression functional. A nonrobust residual standard deviation can be inflated by contamination, causing standardized residuals to appear artificially moderate. Robust formulations therefore define scale through a resistant functional, often based on the median absolute deviation or on an (M)-scale equation.

Leverage and generalized M-estimation

A residual-based estimator controls vertical contamination only conditionally on the predictor configuration. To control leverage, generalized (M)-estimators modify the contribution associated with (x_i). A representative estimating equation is

[ \sum_{i=1}^{n} w(x_i),x_i, \psi!\left( \frac{y_i-x_i^{\mathsf T}\beta} {\sigma,v(x_i)} \right)=0, ]

where (w(x_i)) restricts the magnitude of the predictor contribution and (v(x_i)) can adapt residual scaling to the predictor position. Mallows-type estimators place the leverage weight outside the residual score, whereas Schweppe-type estimators also modify the score’s residual argument.

Such weights depend on a robust description of the predictor distribution. Distances computed from the ordinary sample mean and covariance matrix can themselves be distorted by a cluster of leverage points. Robust multivariate location and scatter estimates supply an alternative geometry for identifying predictor configurations that are distant from the principal data cloud.

Regression equivariance remains a defining structural requirement. If a linear combination of the predictors is added to every response, an equivariant estimator changes its coefficient vector by the corresponding amount. High-breakdown methods are constructed so that resistance to contamination does not depend on an arbitrary choice of response origin or linear predictor parameterization.

High-breakdown estimators

The least median of squares estimator minimizes the median squared residual:

[ \widehat{\beta}_{\mathrm{LMS}}

\operatorname*{arg,min}{\beta} \operatorname{median}{i} \left(y_i-x_i^{\mathsf T}\beta\right)^2. ]

Because fewer than half of the observations cannot generally determine the median residual criterion, the estimator can attain a breakdown point approaching one half. Its objective uses limited information near the center of the residual distribution, which results in relatively low efficiency under Gaussian errors.

The least trimmed squares estimator orders the squared residuals,

[ r_{(1)}^2(\beta)\le \cdots \le r_{(n)}^2(\beta), ]

and minimizes

[ \sum_{i=1}^{h}r_{(i)}^2(\beta) ]

for a designated subset size (h). The parameter (h) determines the relation between breakdown resistance and the fraction of observations contributing directly to the criterion. Unlike a procedure that deletes observations according to a preliminary least-squares fit, trimming is defined simultaneously with coefficient estimation.

S-estimators minimize a robust residual scale subject to a scale equation. Their breakdown point depends on the chosen loss function and calibration constant. MM-estimators begin with a high-breakdown estimate, construct a resistant scale, and define a final (M)-estimate having greater Gaussian efficiency while preserving the initial estimator’s breakdown point under the stated regularity conditions.

Inference and interpretation

Robust regression coefficients retain the interpretation supplied by the underlying conditional model, but changing the loss function can change the population quantity being estimated. Least squares targets a conditional mean when its moment conditions hold. Absolute-deviation regression targets a conditional median under the corresponding specification. Other (M)-estimators target parameters defined by their population estimating equations.

For smooth (M)-estimators, asymptotic covariance commonly has the sandwich form

[ \operatorname{Var}(\widehat{\beta}) \approx \frac{1}{n}A^{-1}BA^{-\mathsf T}, ]

where (A) is the expected derivative of the estimating function and (B) is its variance. This form separates the local curvature of the estimating equation from the variability of individual contributions. Nonsmooth or high-breakdown estimators require results appropriate to their objective functions and may have convergence behavior that differs from that of ordinary smooth estimators.

Robust fitting does not convert arbitrary contamination into verifiable information about the uncontaminated model. When several distinct regression structures each account for a substantial portion of the observations, a high-breakdown estimator can select one structure without establishing that it has a privileged substantive interpretation. The fitted result remains conditional on the regression specification, the contamination model, and the equivariance properties built into the estimator.