A Note on High Breakdown

Regression Estimators

by R.W. Oldford

Technical Report No. 39 August 1983

Massachusetts Institute of Technology

Center for Computational Research in Economics and Management Science Alfred P. Sloan School of Management

Abstract

High breakdown, without other measures of estimator resistance, is an inadequate goal for regression estimators. This is shown by constructing an easily computed regression estimator with 50% breakdown. The estimator is essentially least squares.

Acknowledgements

Special thanks go to David Belsley, David Donoho, Peter Rousseeuw and Alexander Samarov for their helpful comments on an earlier draft and to Karen Martel for typing this paper.

1. Introduction Suppose we have the regression model

y_{i} = x_{i}^{T} β + e_{i} i = 1, ..., n

where β ∈ ℜ^{p} is the parameter vector and (x_{i}^{T} , y_{i}) ∈ ℜ^{p+1} is the ith observation vector. Let D_{n} = {(x_{i}^{T} , y_{i}) : i = 1, ..., n} denote the sample of n observations and b(D_{n}) be any estimator of β based on D_{n}. Further let D_{n,m} ⊂ ℜ^{p+1} be any set of cardinality n such that D_{n,m} ∩ D_{n} has cardinality m and let D_{n,m} be the set containing all such sets for fixed m and n.

The finite sample replacement breakdown (Donoho and Huber (1983)) is defined to be ϵ^{∗} = m/n where m is the smallest integer such that

sup_{D_{n,m}∈D_{n,m}} ∥ b(D_{n,m}) − b(D_{n}) ∥= ∞  (1.1)

and ∥ · ∥ denotes the Euclidean norm.

Recently [Siegel (1982), Donoho and Huber (1983), Rousseeuw (1982)] there has been increasing interest in the construction of regression estimators which have high breakdown. It is generally felt that the resulting estimates will provide good starting values for more efficient robust estimation procedures [e.g., Andrews (1974)]. Alternatively, such estimates may be used for exploratory data analytic purposes.

The purpose of the present note is to demonstrate the inadequacy of the definition (1.1) of breakdown for regression estimators. We do this by proposing a new estimator which has breakdown of 1/2, but which is clearly unsuitable for data analysis. In addition, this estimator is equivariant to non-singular transformations of the explanatory variables, highly efficient at the standard Gaussian model and trivially computed. Other high breakdown estimators, repeated medians (RM) [Siegel (1979, 1982)], and least median square (LMS) [Rousseeuw (1982)], do not possess all of those attributes.

The estimator exploits a peculiarity inherent in the definition of breakdown (1.1), namely, that the supremum must be infinite. Since the parameters being estimated are slope and/or intercept parameters, most estimators will give infinite estimates only if some of the |y_{i}|’s → ∞. The proposed estimator prevents this occurrence.

2. The Estimator Let y′ be the median of the y_{i}’s and MAD be their median absolute deviation, median (|y_{i}−y′|).

An estimator b^{∗}(D_{n}) may be constructed as a weighted least squares estimate with weights w_{i} given by

w_{i} = \{ \begin{matrix}1 & if |y_{i} − y′| ≤ c · MAD \\ 0 & otherwise, \end{matrix}  (2.1)

where c is some constant satisfying 1 ≤ c < ∞. The weights ignore the x-data entirely. Even so, we may prove that it has breakdown 1/2 (as n → ∞).

We make the following assumptions on the elements of D_{n}:

\begin{matrix}A1. & y_{i} ∈ ℜ, & |y_{i}| < ∞ & i = 1, ..., n \\ A2. & x_{i} ∈ ℜ^{p}, & ∥ x_{i} ∥< ∞ & i = 1, ..., n \end{matrix}

A3. x_{1}, ..., x_{n} are in general position, that is no (p + 1) points lie on a (p − 1) dimensional linear manifold.

And we introduce the following notation: let

(i) ∥ A ∥ denote the spectral norm of a matrix A,

(ii) J_{m} ⊂ \{1, ..., n\} of cardinality m and \overline{J}_{m} its complement in this set be defined implicitly by

D_{n,m} = {z_{i} : z_{i}^{T} = (u_{i}^{T} , v_{i}) ∀ i ∈ J_{m} and z_{i}^{T} = (x_{i}^{T} , y_{i})
∀ i ∈ \overline{J}_{m}\} ⊂ ℜ^{p+1},

(iii) W (D) = diag(w_{1}(D), ..., w_{n}(D)) where w_{i}(D)’s are weights as in (2.1) based on a data set D,

(iv) X(D) and y(D) be the usual X matrix and y-vector constructed from a data set D,

(v) λ_{j}(A) the jth largest eigenvalue of A,

(vi) I(D) = {i : w_{i}(D) = 1} ⊂ {1, ..., n},

(vii) [a] denote the largest integer less than or equal to a.

We may now prove the following result.

Theorem: Given A1-A3, the weighted least squares estimator b^{∗}(D_{n}) with weights given by (2.1) has finite sample replacement breakdown of

ϵ^{∗} ≥ \frac{n − [n/2] − p}{n}.  (2.2)

Proof: By A1-A3, ∥ b^{∗}(D_{n}) ∥< ∞ and we need only examine ∥ b^{∗}(D_{n,m}) ∥ in (1.1). By a property of the spectral norm [Wilkinson (1965)]

∥ b^{∗}(D_{n,m}) ∥ ≤ ∥ (X^{T} (D_{n,m})W (D_{n,m})X(D_{n,m}))^{−1}X(D_{n,m})^{T} W (D_{n,m}) ∥
· ∥ W (D_{n,m})y(D_{n,m}) ∥  (2.3)

and ∥ (X^{T} (D_{n,m})W (D_{n,m})X(D_{n,m}))^{−1}X(D_{n,m})^{T} W (D_{n,m}) ∥^{2}

= λ_{p}^{−1} ( ∑_{i∈J_{m}} w_{i}(D_{n,m})u_{i}u_{i}^{T} + ∑_{i∈\overline{J}_{m}} w_{i}(D_{n,m})x_{i}x_{i}^{T} )
≤ λ_{p}^{−1} ( ∑_{i∈\overline{J}_{m}} w_{i}(D_{n,m})x_{i}x_{i}^{T} )
= λ_{p}^{−1} (  ∑_{i∈\overline{J}_{m}∩I(D_{n,m})} x_{i}x_{i}^{T} ) .  (2.4)

The inequality may be found for example in Lawson and Hanson (1974). By A3, (2.4) is finite provided the cardinality of \overline{J}_{m} ∩ I(D_{n,m}) is greater than p. The cardinality of \overline{J}_{m} is n − m. Let that of I(D_{n,m}) be k. The cardinality of their intersection is guaranteed to be greater than p if (n − m) + k − n > p or m < k − p. Since k ≥ n − [n/2] the last quantity is certainly finite if m < n − [n/2] − p. Also,

∥ W (D_{n,m})y(D_{n,m}) ∥^{2}= ∑_{i∈J_{m}} w_{i}(D_{n,m})v_{i}^{2} + ∑_{i∈\overline{J}_{m}} w_{i}(D_{n,m})y_{i}^{2}

which is always finite provided m < [\frac{n}{2}]. Therefore

∥ b^{∗}(D_{n,m}) ∥ < ∞ ∀D_{n,m}

for all m < n − [n/2] − p and the theorem is proved.

Since ϵ^{∗} is bounded from above by \frac{1}{2}, as n → ∞, ϵ^{∗} → \frac{1}{2}. Note also that the theorem could be proven with the inclusion of an intercept term.

3. Discussion In addition to having a high breakdown value, the estimator is equivariant to non-singular

transformations and easily computed. For c sufficiently large in (2.1), the estimator is essentially least squares and therefore highly efficient at the usual Gaussian Model. Further, c can be chosen large enough so that, for practical purposes, the estimator is equivariant to location transformations of the regression parameters.

However, it is clearly not a resistant regression estimator. While sup ∥ b^{∗}(D_{n}) − b(D_{n,m}) ∥ is bounded for m < n − [n/2] − p, it may be quite large for m = 1 (finite-sample-influence or sensitivity function). Figure 1 demonstrates that this can also be the case for the LMS estimator. Here a single point essentially determines the line. Rousseeuw (1982) gives an example with p = 2, where RM is substantially affected by 40% contamination but LMS is not. Clearly these estimators differ on other resistance properties.

Breakdown describes a worst-possible-case scenario. Because of this it is often more easily assessed than an estimator’s sensitivity curve (e.g., LMS). However, without assessment of other resistance properties, it may be misleading. High breakdown implies bounded sensitivity but the bounds may be high enough to be ineffectual.

At least two alternative definitions of breakdown are possible. The first would replace ∞ in (1.1) with some constant k. Breakdown would now be a function ϵ^{∗}(k) and more difficult to assess. Moreover, breakdown would no longer be invariant to a linear transformation of β.

The second possibility is closer in spirit to the original finite sample version of Andrews et al. (1971). It consists of fixing the sets D_{n} and D_{n,m} to be investigated. In such cases ∞ is usually replaced by k < ∞.

Donoho and Rousseeuw (1983; personal communication) have suggested that the exact fit property (EFP) of a regression estimator be investigated. Basically the set D_{n} is chosen to be such that all (x_{i}^{T} , y_{i}) lie exactly on a plane. D_{n,m} is as before and b(D_{n}) is said to have the exact fit property ϵ = m/n if m is the largest integer such that sup_{D_{n,m}∈D_{n,m}} ∥ b(D_{n})−b(D_{n,m}) ∥= 0. Siegel (1979) and Rousseeuw (1983) proved that the RM and LMS respectively have EFP of ϵ = 1/2.

Other possibilities for D_{n} and D_{n,m} should be investigated. It may be the case that estimators which are reasonable on other grounds “break down” only for pairs (D_{n}, D_{n,m}) which rarely occur in practice. It would be helpful to have breakdown type assessments for commonly occurring (D_{n}, D_{n,m}) and exploratory techniques for recognizing those infrequent pairs (D_{n}, D_{n,m}) which are dangerous to the estimator being used.

Figure 1. The LMS estimate for a particular D_{n,m}.

References Andrews, D.F., P.J. Bickel, F.R. Hampel, P.J. Huber, W.H. Rogers, and J.W. Tukey (1972), Robust

Estimates of Location: Survey and Advances, Princeton University Press: Princeton, New Jersey.

Andrews, D.F. (1974), “A Robust Method for Multiple Linear Regression”, Technometrics 16, pp. 523-531.

Donoho, D.L. and P.J. Huber (1983), “The Notion of Breakdown Point”, E.L. Lehmann Festschriftt, Bickel, P.J., K. Doksum and J.L. Hodges, Jr. (eds.), Wadsworth Press.

Lawson, C.R. and R.J. Hanson (1974), Solving Least Squares Problems, Prentice-Hall: Englewood Cliffs, New Jersey.

Rousseeuw, P.J. (1982), “Least Median of Squares Regression”, Technical Report, Centrum voor

Statistiek en Operationeel Onderzoek, Vrije Universiteit Brussel. Siegel, A.F. (1979), “The Repeated Median Algorithm”, unpublished working paper, University

of Wisconsin. Siegel, A.F. (1982), “Robust Regression using Repeated Medians”, Biometrika, 69, pp. 242-244. Wilkinson, J.H. (1965), The Algebraic Eigenvalue Problem, Oxford University Press, London.