Existe-t-il un utilitaire numpy pour rejeter les valeurs aberrantes d'une liste ?

Question

Existe-t-il un buildin numpy pour faire quelque chose comme ce qui suit ? C'est-à-dire, prendre une liste d et retourner une liste filtered_d en éliminant tous les éléments aberrants sur la base d'une distribution supposée des points dans d .

import numpy as np

def reject_outliers(data):
    m = 2
    u = np.mean(data)
    s = np.std(data)
    filtered = [e for e in data if (u - 2 * s < e < u + 2 * s)]
    return filtered

>>> d = [2,4,5,1,6,5,40]
>>> filtered_d = reject_outliers(d)
>>> print filtered_d
[2,4,5,1,6,5]

Je dis "quelque chose comme" parce que la fonction pourrait permettre de varier les distributions (poisson, gaussienne, etc.) et de varier les seuils des valeurs aberrantes au sein de ces distributions (comme la fonction m que j'ai utilisé ici).

Answer 1

Il est important, lorsqu'on traite des valeurs aberrantes, d'essayer d'utiliser des estimateurs aussi robustes que possible. La moyenne d'une distribution sera biaisée par les valeurs aberrantes, mais la médiane, par exemple, le sera beaucoup moins.

En s'appuyant sur la réponse d'eumiro :

def reject_outliers(data, m = 2.):
    d = np.abs(data - np.median(data))
    mdev = np.median(d)
    s = d/mdev if mdev else 0.
    return data[s<m]

Ici, j'ai remplacé la moyenne par la médiane, plus robuste, et l'écart-type par la distance absolue à la médiane. J'ai ensuite mis à l'échelle les distances par leur valeur (encore) médiane de sorte que m est sur une échelle relative raisonnable.

Notez que pour le data[s<m] syntaxe pour travailler, data doit être un tableau numpy.

Answer 2

Cette méthode est presque identique à la vôtre, mais plus numpyst (elle ne fonctionne également que sur les tableaux numpy) :

def reject_outliers(data, m=2):
    return data[abs(data - np.mean(data)) < m * np.std(data)]

Answer 3

La réponse de Benjamin Bannier donne un résultat positif lorsque la médiane des distances à la médiane est égale à 0. J'ai donc trouvé cette version modifiée un peu plus utile pour les cas donnés dans l'exemple ci-dessous.

def reject_outliers_2(data, m=2.):
    d = np.abs(data - np.median(data))
    mdev = np.median(d)
    s = d / (mdev if mdev else 1.)
    return data[s < m]

Exemple :

data_points = np.array([10, 10, 10, 17, 10, 10])
print(reject_outliers(data_points))
print(reject_outliers_2(data_points))

Donne :

[[10, 10, 10, 17, 10, 10]]  # 17 is not filtered
[10, 10, 10, 10, 10]  # 17 is filtered (it's distance, 7, is greater than m)

Answer 4

En s'appuyant sur celui de Benjamin, en utilisant pandas.Series et en remplaçant MAD avec IQR :

def reject_outliers(sr, iq_range=0.5):
    pcnt = (1 - iq_range) / 2
    qlow, median, qhigh = sr.dropna().quantile([pcnt, 0.50, 1-pcnt])
    iqr = qhigh - qlow
    return sr[ (sr - median).abs() <= iqr]

Par exemple, si vous définissez iq_range=0.6 les percentiles de l'écart interquartile deviendraient : 0.20 <--> 0.80 donc plus les valeurs aberrantes seront incluses.

Answer 5

Une alternative consiste à faire une estimation robuste de l'écart-type (en supposant des statistiques gaussiennes). En consultant les calculateurs en ligne, je vois que le percentile de 90% correspond à 1,2815σ et celui de 95% à 1,645σ ( http://vassarstats.net/tabs.html?#z )

Un exemple simple :

import numpy as np

# Create some random numbers
x = np.random.normal(5, 2, 1000)

# Calculate the statistics
print("Mean= ", np.mean(x))
print("Median= ", np.median(x))
print("Max/Min=", x.max(), " ", x.min())
print("StdDev=", np.std(x))
print("90th Percentile", np.percentile(x, 90))

# Add a few large points
x[10] += 1000
x[20] += 2000
x[30] += 1500

# Recalculate the statistics
print()
print("Mean= ", np.mean(x))
print("Median= ", np.median(x))
print("Max/Min=", x.max(), " ", x.min())
print("StdDev=", np.std(x))
print("90th Percentile", np.percentile(x, 90))

# Measure the percentile intervals and then estimate Standard Deviation of the distribution, both from median to the 90th percentile and from the 10th to 90th percentile
p90 = np.percentile(x, 90)
p10 = np.percentile(x, 10)
p50 = np.median(x)
# p50 to p90 is 1.2815 sigma
rSig = (p90-p50)/1.2815
print("Robust Sigma=", rSig)

rSig = (p90-p10)/(2*1.2815)
print("Robust Sigma=", rSig)

Le résultat que j'obtiens est le suivant :

Mean=  4.99760520022
Median=  4.95395274981
Max/Min= 11.1226494654   -2.15388472011
Sigma= 1.976629928
90th Percentile 7.52065379649

Mean=  9.64760520022
Median=  4.95667658782
Max/Min= 2205.43861943   -2.15388472011
Sigma= 88.6263902244
90th Percentile 7.60646688694

Robust Sigma= 2.06772555531
Robust Sigma= 1.99878292462

Ce qui est proche de la valeur attendue de 2.

Si nous voulons supprimer les points au-dessus/au-dessous de 5 écarts types (avec 1000 points, nous nous attendons à une valeur > 3 écarts types) :

y = x[abs(x - p50) < rSig*5]

# Print the statistics again
print("Mean= ", np.mean(y))
print("Median= ", np.median(y))
print("Max/Min=", y.max(), " ", y.min())
print("StdDev=", np.std(y))

Ce qui donne :

Mean=  4.99755359935
Median=  4.95213030447
Max/Min= 11.1226494654   -2.15388472011
StdDev= 1.97692712883

Je n'ai aucune idée de l'approche la plus efficace/robuste.