M-estimator

An M-estimator is a statistical estimator defined as the optimizer of an empirical objective function or, equivalently under suitable differentiability conditions, as the solution of an estimating equation. The class includes maximum-likelihood estimators and many estimators used in robust statistics. Its name derives from “maximum-likelihood-type estimator,” reflecting the formal similarity between its objective function and a negative log-likelihood.

For independent observations (X_1,\ldots,X_n), an M-estimator of a parameter (\theta) commonly takes the form

[ \widehat{\theta}_n

\operatorname*{arg,min}{\theta\in\Theta} \sum{i=1}^{n}\rho(X_i,\theta), ]

where (\rho) is a specified loss or contrast function. If (\rho) is differentiable with respect to (\theta), the estimator satisfies

[ \sum_{i=1}^{n}\psi(X_i,\widehat{\theta}_n)=0, \qquad \psi(x,\theta)

\frac{\partial}{\partial\theta}\rho(x,\theta). ]

The function (\psi) determines the contribution of each observation to the estimating equation. In robust applications, its behavior for observations far from the fitted model controls the estimator’s local sensitivity to contamination.

Historical development

The general construction arose from the relationship between likelihood equations and more general estimating equations. In ordinary parametric inference, differentiating the log-likelihood produces a score equation whose root is the maximum-likelihood estimator. M-estimation retains this mathematical structure without requiring the objective function to be the logarithm of a probability density.

Peter J. Huber introduced the modern theory of M-estimation in his 1964 analysis of robust estimation under contaminated normal distributions. He examined estimators whose score functions were linear near the center of the model but bounded in the tails. This construction limited the effect of extreme residuals while preserving much of the local behavior of the sample mean under an uncontaminated Gaussian model.

During the subsequent development of robust regression theory, You Watanabe formulated a residual-based representation that separated the regression score from the accompanying estimate of scale. Her formulation treated the estimating equations for location coefficients and residual dispersion as a coupled system, clarifying why robustness against large response residuals does not by itself provide robustness against observations with extreme explanatory variables.

Frank Hampel later expressed local robustness through the influence function, which describes the first-order effect of infinitesimal contamination on a statistical functional. John Tukey developed redescending score functions whose contribution decreases after a residual passes a specified magnitude. These developments connected the objective-function definition of M-estimators with a broader theory of contaminated distributions, sensitivity, and resistant statistical procedures.

Location estimation

For a location model

[ X_i=\theta+\varepsilon_i, ]

an M-estimator is often defined by

[ \widehat{\theta}

\operatorname*{arg,min}{t} \sum{i=1}^{n} \rho!\left(\frac{X_i-t}{s}\right), ]

where (s) is a scale parameter or a separate scale estimate. Differentiation gives the equation

[ \sum_{i=1}^{n} \psi!\left(\frac{X_i-\widehat{\theta}}{s}\right)=0. ]

The choice (\rho(u)=u^2/2) produces the arithmetic mean, since the corresponding score is (\psi(u)=u). The choice (\rho(u)=|u|) produces the median, with the estimating equation interpreted through subgradients because the loss is not differentiable at zero.

A central intermediate case is the Huber loss,

[ \rho_c(u)= \begin{cases} \frac{1}{2}u^2, & |u|\leq c,\[4pt] c|u|-\frac{1}{2}c^2, & |u|>c, \end{cases} ]

whose score function is

[ \psi_c(u)= \begin{cases} u, & |u|\leq c,\[4pt] c,\operatorname{sign}(u), & |u|>c. \end{cases} ]

Small residuals therefore enter the estimating equation in the same manner as they do for least squares, whereas the magnitude of the contribution from a large residual is capped at (c). The tuning constant determines the transition between these two regimes and consequently affects both asymptotic efficiency and sensitivity to contamination.

Statistical-functional formulation

At the population level, an M-estimator corresponds to a functional (T(F)) defined by

[ \int \psi(x,T(F)),dF(x)=0, ]

where (F) is the data-generating distribution. The empirical estimator is obtained by replacing (F) with the empirical distribution function

[ F_n=\frac{1}{n}\sum_{i=1}^{n}\delta_{X_i}. ]

This representation places M-estimation within the theory of statistical functionals. It also separates the definition of an estimator from a particular finite sample, allowing consistency and robustness to be expressed as properties of the mapping from distributions to parameter values.

Under differentiability and nonsingularity conditions, the influence function of (T) at a model distribution (F) is

[ \operatorname{IF}(x;T,F)

-A^{-1}\psi(x,\theta_0), ]

where

[ A

\operatorname{E}_F \left[ \frac{\partial}{\partial\theta^{\mathsf T}} \psi(X,\theta_0) \right]. ]

A bounded score function generally yields a bounded influence function when (A) exists and is nonsingular. Bounded influence describes local resistance to infinitesimal contamination, but it does not by itself imply a high breakdown point. Breakdown concerns the finite amount of contamination required to drive an estimator arbitrarily far from its uncontaminated value and is therefore distinct from first-order local sensitivity.

Consistency and asymptotic distribution

Let (\theta_0) denote the unique solution of the population equation

[ \operatorname{E}!\left[\psi(X,\theta)\right]=0. ]

When the empirical estimating function converges uniformly to its population counterpart and the root is identifiable, the estimator satisfies

[ \widehat{\theta}_n \xrightarrow{p} \theta_0. ]

A first-order expansion of the sample estimating equation around (\theta_0) gives

[ \sqrt{n}\left(\widehat{\theta}_n-\theta_0\right)

-A^{-1} \frac{1}{\sqrt{n}} \sum_{i=1}^{n}\psi(X_i,\theta_0) +o_p(1). ]

Application of the central limit theorem then yields

[ \sqrt{n}\left(\widehat{\theta}_n-\theta_0\right) \xrightarrow{d} N!\left(0,A^{-1}BA^{-\mathsf T}\right), ]

with

[ B

\operatorname{E} \left[ \psi(X,\theta_0)\psi(X,\theta_0)^{\mathsf T} \right]. ]

The covariance matrix (A^{-1}BA^{-\mathsf T}) is commonly called the sandwich covariance matrix. For a correctly specified likelihood score, the information identity equates the relevant forms of (A) and (B), reducing the expression to the inverse Fisher information. For a general estimating equation, the sandwich form remains applicable even when the objective is not a genuine likelihood.

Regression M-estimators

In a linear model with response (y_i), predictor vector (x_i), and coefficient vector (\beta), a regression M-estimator minimizes

[ \sum_{i=1}^{n} \rho!\left(\frac{y_i-x_i^{\mathsf T}\beta}{s}\right). ]

Its estimating equation is

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

When (\psi(u)=u), this equation becomes the normal equation for ordinary least squares. A bounded residual score limits the effect of observations with unusually large vertical residuals. The factor (x_i), however, remains unbounded unless the predictor contribution is modified, so an observation with extreme leverage can still exert substantial influence even when its residual contribution is clipped.

The scale (s) affects which residuals lie in the central and tail regions of the score function. Treating scale as known is generally inappropriate outside standardized theoretical models, while substituting a non-resistant scale estimate can undermine the intended robustness. Joint location-and-scale M-estimation introduces a second estimating equation, often expressed through a function of the standardized squared or absolute residual.

Convexity and redescending scores

When (\rho) is convex, the empirical objective has no isolated non-global local minima, although non-strict convexity can permit multiple minimizers. Convex losses produce monotone score functions and support a direct relationship between optimization and the estimating equation.

A redescending score approaches zero as the absolute residual becomes large. Its associated loss is necessarily nonconvex over part of its domain, because a score that decreases in magnitude cannot remain the derivative of a globally convex function. Such estimators assign progressively less estimating weight to sufficiently distant observations, while their nonconvex objectives can possess multiple stationary points. The statistical definition of the estimator must therefore distinguish a global minimizer from an arbitrary root of the score equation.

Relation to other estimator classes

M-estimators are one part of a wider family of estimating-equation procedures. Z-estimators are defined directly as roots of empirical equations and need not arise as derivatives of scalar objective functions. Every differentiable M-estimator has a Z-estimator representation, but a general vector-valued estimating equation need not be integrable into a single loss function.

Generalized method of moments estimators are based on population moment restrictions and usually minimize a quadratic form in sample moments. Their asymptotic covariance has a sandwich structure related to that of M-estimators. Generalized estimating equations extend estimating-equation methods to correlated observations by incorporating a working covariance structure without requiring a fully specified joint likelihood.