Iglewicz and Hoaglin outlier test with modified z-scores - What should I do if the MAD becomes 0?

Question

I'm a programmer with a small statistics background and I need to find outliers in a small list of integers and floats.

After some search on google I found the Iglewicz and Hoaglin outlier test which creates a modified z-score M_i for every value in the list and check it against an threshold (normally 3.5).

$$M_{i} = \frac{0.6745(x_{i} - \tilde{x})} {\mbox{MAD}}$$

I wrote a litte python script to test it. At first it worked great, but after a few tests I spotted an error.

If you try to find outliers (with my script) in an list with many identically values and one outlier e.g. data = [10, 10, 10, 10, 10, 10, 10, 100] the MAD(median absolute deviation) becomes 0 and this leads my to my question: "What should I do if the MAD becomes 0?".

My first idea was to set the MAD to ∞, but this causes the script to find no outliers.

My second idea was to add very small offsets to the values to make them unique e.g. data = [10.0, 10.00000001, 10.00000002, 10.00000003, 10.00000004, 10.00000004, 10.00000005, 100]. This way the MAD can't become 0 and my script is able to detect the outlier 100.

Does somebody have better ideas?

Am I doing something wrong?

Answer 1

Three facts will help you here.

• What you discovered is called the exact fit property. If a proportion $\alpha > 0.5$ of the observations in your sample have the same value, the mad of your sample will be 0.
• This is not a property of the mad in particular, but of all robust estimators of scale. More precisely: any robust estimator of scale with a breakdown point of $0< \alpha < 0.5$ will have an exact fit property at the level of $1-\alpha$ (see section 3 of Croux et al., 2006, [0], for example).
• Your first proposals amount to replacing the value of $M_i$ by arbitrary numbers in case of exact fit (setting $M_i=0$ in the former and $M_i=O(1/\sigma)$ --where $\sigma$ is the amount by which you perturb the data-- in the latter).

Your proposed solution to the problem (point 3) is not the correct one.

In fact, the correct solution to your problem is much simpler. Keep the MAD, keep the outliers rejection rule. All you need to do is to adopt the convention $0/0:=0$ in the computation of the outliers detection rule. This convention has no impact outside of exact fit cases. Then you can use the rule regardless of whether the MAD is strictly positive or not.

This is because:

In an exact fit situation whereby half or more of the data is tied at an arbitrary value $x$, all observations in your sample that are different from $x$ are severe outliers.

In such a situation, all observations in your sample that are different from $x$ are, after all, infinitely divergent from the pattern of the bulk of the data. Then, adopting the $0/0:=0$ will assign the correct outlyingness score both to those observation equal to $x$ ($M_i=0$) and those different from $x$ ($M_i=\infty$).

The reason you can use this convention is because the exact fit property is bijective:

Mad = 0 $\iff$ more than half of your sample are tied to the same value.

• Algorithms for projection-pursuit robust principal component analysis. (2006). Croux, C. Filzmoser, P. and Oliveira, M. R.

Answer 2

1. A practical suggestion.

Change this part of the code

    if mad == 0:
        mad = 9223372036854775807 # maxint

to

    if mad == 0:
        mad = 2.2250738585072014e-308 # sys.float_info.min

It does the trick. Division by this number blows up the Iglewicz-Hoaglin test statistic – exactly as desired. That is, marking strongly deviant observations as outliers.

---

2. Previous practical suggestion.

What you could do, is check if it works with the closely related definition of mean absolute error (MAE):

$$ \text{MAE} = \frac{1}{n} \sum_{i=1}^n |x_i - \text{median}(x)|, $$

with $e_i = x_i - \text{median}(x)$ the errors (better: residuals, or, deviations).

IBM uses this variant:

$$ M_{i} = \frac{x_{i} - \text{median}(x)} { 1.253314 \cdot \text{MAE} } $$

for the if MAD == 0 case.

---

3. What is going on here? (From a programming perspective)

Consider the two cases:

1. $0/0$,
2. $x/0$ for $x \neq 0$.

Scientific programming languages R, Matlab and Julia have the following behavior:

1. 0/0 returns NaN.
2. 90/0 returns Inf.

Python, on the other hand, throws a ZeroDivisionError in both cases.

Practical suggestion one circumvents both cases for both flavors of zero-division handling.