Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
45,822 characters · 21 sections · 46 citation commands
On (in)consistency of M-estimators under contamination
In the early 1960s, there was an increasing awareness that standard estimators for normal models may not fare well under deviations from normality. huber1964 proposed maximum likelihood-type ({M}) estimators for location and scale and found that they are more robust to such deviations than traditional estimators. Since then, robustness means that an estimate is only distorted in a bounded way when adding arbitrary observations to the sample hampel1971general. While statistical inference theory was developed for data with infinitesimal contamination HeritierRonchetti1994, little theory is available for cases with more contamination.
We analyze popular M-estimators for location and scale when a fixed proportion of observations is contaminated. We find that the Huber estimator for location is not consistent, whereas the Tukey estimator for location has better consistency properties. To conduct inference, the uncontaminated scale has to be estimated. Unfortunately, robust scale estimators such as the interquartile range (IQR) and the median absolute deviation (MAD) are inconsistent under contamination. Thus, valid statistical inference based on {M}-estimation seems to require detailed modelling of the contamination, which is typically not done in practice.
Huber's article initiated the development of a large body of robust estimators, starting with the location-scale problem and later expanding to regression and many other statistical problems. An attractive feature of Huber's work was the idea of developing estimators that are not strongly based in statistical models. However, when thinking about inference, statistical models are useful. The lack of models is perhaps the reason that inference theory remains incomplete for the contaminated case. A major problem in the inference theory is that we need consistent estimators for both location and scale. Location and scale are entwined. Without knowing one, it is hard to get to know the other. huber1964 was clearly aware of this.
Huber's $\epsilon$-contamination model mixes a normal distribution with an $\epsilon$ proportion of a contamination distribution. This model has become the standard contamination model. It is well suited for analyzing the consequences of infinitesimal contamination hampel1986robust, HeritierRonchetti1994, huber2009; see avella2015robust, yu2017robust, loh2024theoretical for recent reviews and HorowitzManski1995 for a discussion of identification. But, when $\epsilon$ is not vanishing, standard scale estimators are inconsistent and asymptotic inference on the location depends on the contamination, which is usually unknown. Thus, inference for {M}-estimators is fraught under contamination.
Another type of robustness analysis uses the breakdown point hampel1971general. The finite sample breakdown point asks how many arbitrary observations can be added to a sample without distorting the estimator unboundedly donohohuber1983,huber1984finite. Location M-estimators typically have a high breakdown point. In practice, it is common to use {M}-estimators with a high breakdown point, which also maintain high efficiency under the normal model Yohai1987, coakley1993bounded, yu2017robust.
Our approach is to consider a model where a proportion of the observations is contaminated. This matches the breakdown point idea. We consider two classes of {M}-estimators. The first class has an unbounded objective function. It includes the median and the Huber estimator. We find that these estimators are bounded in probability when more than half of the data is uncontaminated, which attests their robustness. Yet, these estimators are typically inconsistent under contamination.
The second class has a bounded objective function and includes the Tukey estimator. We derive a lower bound on the proportion of contamination ensuring that these estimators are bounded in probability and consistent under contamination. However, nuisance parameter free inference requires a consistent scale estimator. We show that robust scale estimators such as the IQR and MAD are inconsistent under contamination.
These results for {M}-estimators are disappointing. In the discussion of the results, we point to the least trimmed squares (LTS) estimator rousseeuw1984least which fares better.
We first revisit {M}-estimation and classical approaches to the robustness problem. Next, we introduce a new contamination model and analyse boundedness, consistency and asymptotic normality properties. This leads to new insights into the behaviour of {M}-estimators under contamination. The results are illustrated in a simulation study.
Suppose we have scalar observations $y_i$ for $i = 1, \dots, n$ and we want to estimate a location-scale specification $y_i = \mu + \sigma \varepsilon_i$ where $\mu\in \mathbb{R}$ and $\sigma > 0$. A model for the errors $\varepsilon_i$ will be specified below.
The case of a known scale, $\sigma=1$ say, is the starting point for the {M}-estimation theory of huber1964. The {M}-estimator $\hat\mu$ for location $\mu$ minimizes the objective function $ \sum_{i=1}^n \rho \big( y_i - \mu \big) . $ We focus on $\rho$ functions that have a constant hence bounded derivative in the tails, as described in the following assumption. Note, this excludes the sample average with $\rho(x) = x^2$, which is not robust.
We distinguish between two types of estimators. Non-redescending estimators have $\psi_* > 0$ and include the median with $\rho(x) = |x|$ and the Huber estimator with $\rho(x) = (x^2/2) \mathbf{1}_{(|x| \leq c)} + c(|x| - c/2)\mathbf{1}_{(|x| > c)}$. Redescending estimators have $\psi_* = 0$ as for the Tukey estimator with $\rho(x) = (c^2/6)\{1 - (1 - x^2/c^2)^3\mathbf{1}_{(|x| \leq c)}\}$. The parameter $c > 0$ is a tuning parameter that the researcher needs to set in advance. In Section (ref), we discuss standard choices that are used in practice.
In practice the scale $\sigma$ will be unknown. Suppose a preliminary estimator $\hat\sigma$ is available. Inserting this estimator in the above objective function for the location case gives a new objective function
Note, the median computed from $\rho(x)=|x|$ is scale equivariant and is computable without knowing the scale.
It is important that the scale estimator is also robust, because otherwise the location estimator can still easily be affected by `outliers'. Standard robust scale estimators include the IQR and the MAD huber2009. Letting $\hat{q}_p$ denote the $p$ quantile of the full sample, we write the estimators as
The tuning parameters, $k$, must be set in advance by the researcher. Some standard choices are discussed in Section (ref). We note that the combination of the median with the MAD is considered a simultaneous {M}-estimate of location and scale huber2009.
The Huber $\epsilon$-contamination model huber1964 has become the standard way of thinking about and modelling contamination. In the $\epsilon$-contaminated normal model, the errors $\varepsilon_i$ are drawn from a mixture $F = (1-\epsilon)\Phi + \epsilon H$ with $\epsilon \in [0,1)$, $\Phi$ denotes the standard normal distribution and $H$ an unknown contamination distribution. Theory often requires symmetric contamination.
In this model, the asymptotic distribution of {M}-estima\-tors depends crucially on two unknown nuisance parameters: $\epsilon$ and $H$. To get valid inference, the user will have to model the contamination. However, in practice, it is common to use the asymptotic theory under the normal model even when there is contamination.
In the infinitesimal approach the contamination proportion vanishes. This leads to the influence function theory, which describes the effect of an infinitesimal contamination of a point mass distribution on the estimator hampel1971general. A condition for (local) robustness is that the influence function of an estimator is bounded.
In Huber's $\epsilon$-contamination model with vanishing $\epsilon$, it is common to conduct inference as if the distribution is normal without any contamination. This eliminates the nuisance parameter problem mentioned in Section (ref). Then, the nuisance parameters $c$ and $k$ are chosen from a perspective of efficiency in the uncontaminated normal model. Usual choices are $c = 1.345$ for the Huber estimator and $c=4.685$ for the Tukey estimator, which give $95\%$ efficiency relative to the average MaronnaMartinYohai2006. The motivation is a trade-off between efficiency and robustness. The sample average is fully efficient but can be arbitrarily distorted by one `outlier'. Thus, a less efficient, but more robust estimator is preferred. Tuning parameters that ensure consistency under normality are $k = 1/(q^\Phi_{3/4}-q^\Phi_{1/4})$ for IQR and $k = 1/q^\Phi_{1/4}$ for MAD, where $q^\Phi_p$ is the standard normal $p$-quantile MaronnaMartinYohai2006.
The infinitesimal approach works well when the contamination proportion $\epsilon$ is very small. In empirical practice, this local property is not always a realistic assumption. Therefore, researchers are also interested in global robustness properties such as the breakdown point.
The finite sample breakdown point with $\epsilon$-contamina\-tion is defined as follows donohohuber1983. Let $X$ be a sample of size $h$ of `good' observations. Adjoin $n-h$ arbitrary values $Y$ to get a sample $X\cup Y$ of size $n$. The sample $X\cup Y$ contains a fraction $\epsilon=(n-h)/n$ of arbitrary values. The breakdown point for an estimator $T$ is then $\epsilon^*(X,T)=\inf\{\epsilon:\sup_{Y}|T(X\cup Y)-T(X)|=\infty\}$. As $X$ is finite, so is $T(X)$ and breakdown happens if $Y$ diverges.
huber1984finite analyzed the finite sample breakdown of {M}-estimators. For non-redescending estimators satisfying Assumption (ref), the breakdown point is $1/2$, while for redescending estimators, breakdown depends on the setup and the breakdown point is in general less than $1/2$.
The infinitesimal approach and the breakdown point approach are not compatible. The contamination proportion is vanishing for the former and non-vanishing for the latter. If one applies the breakdown point approach, then it is not a good idea to set the tuning parameters as if there is no contamination. For instance, it is not guaranteed that the M-estimator is consistent under contamination. Likewise, inference and efficiency considerations based on the normal model are not necessarily relevant.
Although the two approaches are not compatible, it has become standard practice to combine them. There is a focus on the use and development of robust estimators with high efficiency under the normal model and a high breakdown point Yohai1987, coakley1993bounded, GerviniYohai2002. Moreover, it is typically advised not to use estimators that have a high breakdown point, but low asymptotic efficiency under the normal model yu2017robust.
Robustness and efficiency are often seen as unrelated problems. This may not be optimal. For example, even when an estimator has high efficiency under the normal model and a high breakdown point, then it is still not clear how that estimator behaves under contamination. Two estimators that satisfy both conditions, such as the Huber and Tukey estimators, might behave differently under contamination, and it unclear for a practitioner which estimator to choose. Yet, a high breakdown point estimator with a low asymptotic efficiency under the normal model may be preferable under contamination, see Section (ref).
The drawbacks mentioned in Section (ref) stem from the fact that there are no clear asymptotic results for M-estimators under contamination. Therefore, we analyse asymptotic properties such as boundedness, consistency and asymptotic normality of M-estimators under contamination. To this end, we introduce a simple contamination model inspired by the breakdown point theory in Section (ref). We assume that a sample of observations $y_i$, with $i = 1, \dots, n$, has $h \leq n$ `good' observations and $n-h$ `outlier' observations. Let $\zeta_n$ be the set of $h$ elements from $(1, \dots, n)$ with indices of the `good' observations, while $\zeta_n^c$ denotes its complement, consisting of $n-h$ indices of `outliers'. Throughout Section (ref), we impose different assumptions on `good' and `outlying' observations.
We assume that a sample of observations $y_i$, with $i = 1, \dots, n$, has $h \leq n$ `good' observations and $n-h$ `outliers'. The models are indexed by $n$, so that $h \rightarrow \infty$ as $n \rightarrow \infty$ and the proportion of `good' observations satisfies $h/n \rightarrow \lambda$ with $1/2 < \lambda \leq 1$. Let $\zeta_n$ be a deterministic sequence of $h$-sets from $(1, \dots, n)$ of indices of `good' observations. The observations satisfy $y_i = \mu_\circ + \sigma_\circ \varepsilon_i$ for fixed $\mu_\circ, \sigma_\circ$. Throughout, we specify sufficient conditions for errors $\varepsilon_i$ for `good' $i\in\zeta_n$ and `outliers' $j\in \zeta_n^c$.
We show that {M}-estimators are bounded in probability in the asymptotic setup of Section (ref). The result requires a convergent scale estimator $\hat\sigma$. Ideally, the scale estimator should be consistent, but Assumption (ref)$(i)$ allows general limits expressed in terms of a consistency factor $\varsigma$. Indeed, in Section (ref), we show that standard robust scale estimators are typically inconsistent under contamination.
Assumption (ref)$(iii)$ is a uniform law of large numbers. When the scale is known, then a standard law of large numbers suffices. Note, we do not need to make any assumptions on the `outliers'.
In the redescending case, Theorem (ref)$(ii)$ gives a lower bound to the proportion of `outliers' that is greater than a half and which depends on the $\rho$ function and the distribution of the `good' observations. This exactly matches the finding in donohohuber1983,huber1984finite regarding the finite sample breakdown points of these estimators. We explore this connection.
The finite sample addition breakdown point for redescending {M}-estimators is $\epsilon^* = \lceil A \rceil/(h +\lceil A \rceil)$, where $A=h - \sum_{i\in\zeta_n} \rho\{(y_i - \hat\mu)/(\varsigma\sigma_\circ)\}/\rho_*$ huber1984finite. If $\hat\mu \overset{p}{\to} \mu_\circ$, then $A/h$ and $1 - h^{-1} \sum_{i\in\zeta_n} \rho(\varepsilon_i/\varsigma)/\rho_*$ have the same limit $1 -\tilde{\rho}_{\varsigma}/\rho_*$, due to Assumption (ref)$(iii)$. As a result, $1 -\epsilon^*\to 1/(2-\tilde{\rho}_{\varsigma}/\rho_*)$, as derived in Theorem (ref)$(ii)$.
Redescending M-estimators can tolerate more `outliers' when the scale is over-estimated as the lower bound to $\lambda$ in Theorem (ref)$(ii)$ is decreasing in $\varsigma$. However, when the scale is underestimated, then the redescending {M}-estimator can become unbounded. Indeed, as $\varsigma \rightarrow 0$, then $\tilde{\rho}_{\varsigma}/\rho_* \rightarrow 1$ so that $1/(2 - \tilde{\rho}_{\varsigma}/\rho_*) \rightarrow 1$. In Section (ref), we show that the IQR and MAD typically overestimate in the presence of contamination. huber1984finite also noted this in his finite sample breakdown point analysis. In Section (ref), we explore this in simulations.
We derive an oracle property for redescending {M}-estimators, such as the Tukey estimator. That is, the M-estimators for the contamination model are close to the minimizers of $\sum_{i\in\zeta_n} \rho \{(y_i-\mu)/(\sigma_\circ\varsigma)\}$ in which the set of `good' observations is known. This is helpful, as the latter problem has been studied in the literature.
In Assumption (ref), part $(i)$ requires that the smallest `outlier' errors diverges. This holds, for example, if the `good' errors are $\mathsf{N}(0,1)$, while `outliers' are outside the range of the `good' errors. Part $(ii)$ requires that $\rho$ is Lipschitz, which holds for the Tukey estimator. Part $(iii)$ requires that the infeasible {M}-estimator is close to $\mu_\circ$.
We can get an asymptotic distribution for the $M$-estimator under some additional assumptions allowing asymmetric contamination.
Sufficient conditions for Assumption (ref) are given in jureckova1996robust, Chapter 5.3. In particular, if the `good' errors are independent standard normal and $\rho$ is twice differentiable, then $V_{1} = \int_{\mathbb{R}} \{\dot\rho (x)\}^2 d\Phi(x)/\{\int_{\mathbb{R}} \ddot\rho(x) d\Phi(x)\}^2$.
Theorem (ref) gives a basis for discussing asymptotic efficiency within the class of robust estimators with the oracle property that $\hat\mu-\hat\mu_{\zeta_n}$ is small. For example, if the `good' observations are all $\mathsf{N}(0,1)$, then the benchmark is the fully efficient average of the $h$ `good' observations. Redescending {M}-estimators have efficiency $1/V_{\varsigma}$. In particular, the Tukey estimator with known scale and $c=4.685$ achieves 95% efficiency. Thus, the Tukey estimator achieves high efficiency in both the normal model without contamination MaronnaMartinYohai2006 and in the normal model with contamination. This does, however, require a consistent estimator of the scale $\sigma_{\circ}$.
Non-redescending {M}-estimators, such as the median and the Huber estimator, are inconsistent in the presence of asymmetric contamination. We illustrate this using a class of data generating processes matching the setup for finite sample breakdown point considerations.
Note, Assumption (ref)$(i)$-$(iii)$ implies Assumption (ref), so that the boundedness Theorem (ref) applies. Part $(v)$ provides a particular asymmetric contamination.
Theorems (ref)$(ii)$ and (ref) depend on the scale estimator $\hat\sigma$ via the consistency factor $\varsigma$. Thus, it is important to find robust scale estimators that are bounded and consistent under contamination. We first show that the IQR and MAD defined in ((ref)) are bounded in probability.
Theorem (ref) shows that the IQR is bounded when the asymptotic proportion of `good' observations $\lambda > 3/4$, while the MAD requires $\lambda > 1/2$. This aligns with their finite sample breakdown points of $1/4$ for the IQR and $1/2$ for the MAD huber2009.
Next, we show that the two scale estimators are, in general, inconsistent under contamination. We normalize the estimators to be consistent for normal data, see Sections (ref), (ref), but consider `good' errors with a general distribution and allow contamination.
We define the consistency factors
where $d$ solves $\mathsf{F}(c + d) - \mathsf{F}(c - d) = 1/(2\lambda)$ for a $c$ solving $\mathsf{F}(c) = (1/2 - \varrho)/\lambda$, where $q^{\mathsf{F}}_p$ and $q^\Phi_p$ are the $\mathsf{F}$ and normal $p$ quantiles. When the `good' errors are normal and $\lambda=1$, so that $\varrho = 0$, then $\varsigma_{\textsc{IQR}}=\varsigma_{\textsc{MAD}}=1$. In the Appendix, we show that when Assumption (ref) holds with `good' normal errors and $\lambda<1$, then both consistency factors are greater than unity, see Remarks (ref), (ref).
The scale estimators $\hat\sigma_{\textsc{IQR}}$ and $\hat\sigma_{\textsc{MAD}}$ are, in general, inconsistent under contamination. For the non-redescen\-ding median and Huber estimator, boundedness is not affected by the inconsistency of the scale (Theorem (ref)). However, these estimators are inconsistent under asymmetric contamination (Theorem (ref)). For redescending estimators, such as the Tukey estimator, the situation is more subtle. The boundedness result can tolerate more contamination when the scale is overestimated (Theorem (ref)), which it tends to be (Remarks (ref), (ref)). This in line with the breakdown point analysis in huber1984finite. Further, inconsistency of the scale estimator does not affect the consistency of the location estimator. However, it does result in nuisance parameters when conducting inference (Theorem (ref)). In summary, conducting valid inference with {M}-estimators seems to require detailed modelling of the contamination.
We study finite sample properties of four {M}-estimators of location (average, median, Huber, Tukey). We consider bias, the effect of scale estimation and boundedness.
In Section (ref), we used the average, applied infeasibly to the set $\zeta_n$ of $h$ `good' observations, as a benchmark for efficiency. In practice, the set $\zeta_n$ is unknown and must be estimated. The Least Trimmed Squares (LTS) estimator rousseeuw1984least does this and it is asymptotically efficient.
The LTS estimator requires that the user specifies there are $h$ `good' observations and $n-h$ `outliers'. The LTS estimator minimizes the residual sum of squares over all $h$-subsets of the data. In the location-scale case, we can proceed as follows. Let $\zeta$ be a $h$-subset of $(1,\dots,n)$ with indices of the `good' observations. This is estimated by $$ \hat\zeta = \mathrm{argmin}\,_{\zeta} \sum_{i \in \zeta} (y_i - \bar{y}_{\zeta})^2 , $$ while the estimators of location and scale are
The LTS estimator is maximum likelihood in a contamination model satisfying Assumption (ref), where `good' observations are i.i.d.\ normal and `outlier' errors are more extreme than the `good' errors as in Assumption (ref)$(i)$, but otherwise unspecified berenguer2023model. This LTS model matches the setup of the finite sample breakdown point, see Section (ref). Indeed, the LTS estimator has a breakdown point of 50% and this result extends to a regression version of the estimator RousseeuwLeroy1987. When applied to a clean normal sample, then the LTS estimator applied with $h\approx n/2$ has an efficiency of about 7% relative to the full sample average (Butler1982; RousseeuwLeroy1987). However, under the LTS model, the standardized estimator $h^{1/2}(\hat\mu_{LTS}-\mu_\circ)/\hat\sigma_{LTS}$ has the oracle property that it equals the standardized average with large probability, so that it is asymptotically standard normal under the LTS model berenguer2023model. This result generalizes to a regression model that allows leverage effects berenguer2023model,berenguer2023leverage.
In practice, the user will have to estimate the proportion $\lambda$ of `good' observations by estimating $h/n$. A standard method to estimate $\lambda$ is the index plot method RousseeuwLeroy1987. A detailed analysis is likely to show that this method is not consistent and therefore will not result in a fully efficient LTS estimator. Instead, we estimate $\lambda$ using the cumulant based method of berenguer2023model. This estimator is consistent, but it remains work in progress to show that the resulting LTS estimator has the oracle property.
The location estimators are as follows. At first, we consider the scale to be known. The sample average is used as a non-robust benchmark. Neither the average, the median nor the LTS estimator require initial scale estimation. The Huber and Tukey estimators require initial scale estimation and are used with tuning parameters $c = 1.345$ and $c = 4.685$, respectively, to achieve $95\%$ efficiency. The LTS estimator requires a user choice of $h$ to set trimming.
At first, we let the scale $\sigma_{\circ}$ be known and choose $h$ as 80% of $n$ in the LTS estimation. Subsequently, we estimate $\sigma_{\circ}$ with the MAD and we estimate $h$ using a cumulant based normality test berenguer2023model.
We consider six data generating processes (DGPs). For all DGPs we use $\mu_{\circ} = 0$ and $\sigma_{\circ} = 1$, consider sample sizes $n = 25, 100, 400$ and use $10^5$ repetitions.
DGP1--2 have no contamination. DGP1 has independent $\mathsf{N}(0,1)$ errors and DGP2 has independent $t(3)$ errors.
DGP3--6 have contamination: DGP3--5 have $h/n = \lambda = 0.8$; DGP6 has $h=n-1$. The `good' errors are independent $\mathsf{N}(0,1)$ for DGP3,4,6 and independent $t(3)$ for DGP5. The `outlier' errors satisfy $\varepsilon_j = \max_{i \in \zeta_n} \varepsilon_i + \xi$, $j \in \zeta_n^c$. We have $\xi = 1$ for DGP3 and $\xi = 3$ for DGP4--6. DGP6 mimics infinitesimal $\epsilon$-contamination.
Tables (ref)--(ref) report bias for different estimators when scale and trimming are, respectively, known and estimated. This only affects the Huber, Tukey and LTS estimators. As $10^5$ repetitions are used, the maximal Monte Carlo standard error was below $0.0005$.
DGP1 has no contamination and normal errors so that all methods perform well. In this case, the MAD is consistent for $\sigma_{\circ}$. Thus, even with unknown scale, we see good performance of the Huber and Tukey estimators. berenguer2023model give simulation evidence that the cumulant based estimator for $h$ is consistent when the `good' errors are normal, so that the LTS estimator also performs well when estimating $h$.
DGP2 has no contamination and $t(3)$ errors. For this reason, all estimators are consistent. When $\sigma_{\circ}$ is unknown, the MAD overestimates, which does not negatively affect the consistency properties of the Huber and Tukey estimator as there is no contamination. The LTS estimator performs very well, even for unknown $h$.
DGP3--DGP4 have normal `good' errors and extreme `outliers'. The bias of the average grows as the sample size grows because the `outliers' are larger for larger sample sizes. Moreover, the bias is increasing in the `outlier' parameter $\xi$. The bias of the median and the Huber estimator remain constant when $n$ and $\xi$ grow, which shows their robustness. They are, however, inconsistent (Theorem (ref)). When $\sigma_{\circ}$ is unknown, the MAD overestimates (Theorem (ref) and Remark (ref)). Overestimation of the scale results in less downweighting of `outliers', hence more bias. The Tukey estimator is consistent so that the bias decreases as the sample size grows (Theorem (ref)). When $\sigma_{\circ}$ is known, the Tukey estimator does not always penalize the `outliers' in DGP3, so the bias is large in small samples. When $\sigma_{\circ}$ is unknown, the performance of the Tukey estimator worsens and convergence is very slow for DGP3. Extra simulations showed biases of $0.458$ and $0.386$ for $n = 1600, 6400$, respectively. The LTS estimator performs very well, even when $h$ is unknown.
DGP5 has $t(3)$ `good' errors and extreme 'outliers'. As in DGP3--4 the average, median and Huber estimators are biased. The Tukey estimator performs very well. This happens, because the `outliers' are generated by the maximal `good' errors and tend to be larger than the `outliers' in DGP4. These `outliers' are down-weighted, which results in a good performance of the Tukey estimator. With unknown scale, the MAD overestimates, which results in slower convergence than in the known scale scenario. The LTS estimator performs very well, and even better when $h$ is unknown. This happens because $h$ is estimated using a cumulant based normality test, while the `good' observations are $t(3)$, resulting in a beneficial underestimation.
DGP6 mimics the situation of infinitesimal contamination. We see that all estimators are consistent as the bias decreases to zero as sample size increases. This happens, because the influence of the `outlier' diminishes as the samples size increases. We do, however, see that the average is quite biased for smaller sample sizes. In particular, if $\xi \rightarrow \infty$, then the bias of the average would diverge, while the bias of the other estimators would be unaffected.
We explore how the bias of four {M}-estimators depends on the scale estimation. We use $\hat\sigma=\varsigma\sigma_\circ$ and vary the consistency factor $\varsigma$. Figure (ref) has two graphs showing the bias of the estimator for varying $\varsigma$ and $n=100$, $1000$.
We generate data from DGP4, albeit with $h/n = 0.6$. The Huber and Tukey estimators have the same tuning parameters, tuned for $\sigma_{\circ} = 1$, as before.
The average and median are biased, but constant in $\varsigma$ due to their scale invariance. For large $n$, the median approaches the $5/6$ normal quantile, which is 0.97. The average diverges as $0.4 (\max_{i\in\zeta_n}\varepsilon_i+3)$. We have $\max_{i\in\zeta_n}\varepsilon_i/\{2\log(0.6n)\}^{1/2}\overset{p}{\to} 1$, see berenguer2023leverage, Appendix B.
Consider the Tukey estimator. When $\varsigma = 1$, the estimator has zero bias. When $\varsigma > 1$, the bias is zero up until $\varsigma = 1.2\ (n = 100)$ or $\varsigma = 1.3\ (n = 1000)$ and converges towards the bias of the average. When $\varsigma$ is large, the `outliers' are not down-weighted so that the Tukey estimator behaves as the average. Nevertheless, it is consistent since the `outliers' diverge and therefore for any fixed $\varsigma > 1$ the `outliers' are, at some point, down-weighted. When $\varsigma < 1$, the bias is zero down to $\varsigma = 0.65\ (n=100)$ or $\varsigma = 0.55\ (n = 1000)$. Below those values, the Tukey estimator hones in on the outliers so that its bias diverges as $\max_{i\in\zeta_n}\varepsilon_i+3$. This and the dependency on $\xi$ shows unboundedness of the estimator. The red vertical line in Figure (ref) indicates the smallest value of $\varsigma$ for which boundedness holds according to Theorem (ref)$(ii)$. The simulations indicate that the first part of the bias curve converges to a function with a step at the red line as expected.
Consider the Huber estimator. When $\varsigma \rightarrow 0$, it behaves as the median as the width of the central quadratic part of the objective function shrinks, so that the linear parts dominate. When $\varsigma \rightarrow \infty$, it behaves as the average as the quadratic part of the objective function now dominates.
Robustness is the property that an estimator changes in a bounded way when adding contaminated observations to a sample. This is traditionally measured in terms of the finite sample breakdown point. {M}-estimators are known to have a high finite sample breakdown point for location-scale models donohohuber1983,huber1984finite. The results on asymptotic boundedness confirm this for location estimators (Theorem (ref)) and for scale estimators (Theorem (ref)). These results do not generalize to regression M-estimators he1990tail.
Robust estimators that are consistent and have simple inferential properties under contamination are particularly desirable. Theorems (ref), (ref) showed that redescending {M}-estimators are consistent under contamination, while non-redescending {M}-estimators are not. Even so, simple inference based on redescending M-estimators requires a scale estimator that is consistent under contamination. This fails for robust scale estimators such as the IQR and MAD (Theorem (ref)). It remains an open problem to find robust {M}-estimators with simple inference when a proportion of the observations are contaminated.
By simulation, we found that the LTS estimator performs better than various {M}-estimators, making it an attractive alternative. LTS estimators have a high breakdown point in location-scale models and in regression RousseeuwLeroy1987. LTS estimators for location and scale, and more generally for regression coefficients, have the oracle property that inference is the same as for the infeasible least squares estimator on the `good' observations only berenguer2023model,berenguer2023leverage. Inference is, therefore, nuisance parameter free and it achieves full efficiency in the contaminated model (Section (ref)). These results do, however, require that the proportion of `good' observations is known. One could say that the problem of estimating the scale has been replaced by the problem of estimating the proportion of the `good' observations. This appears to be an easier problem to solve theoretically. Indeed, an estimator proposed in berenguer2023model appeared successful in the simulations.