¿Hay un numpy incorporado para rechazar valores atípicos de una lista?

100

¿Existe un numpy incorporado para hacer algo como lo siguiente? Es decir, tome una lista dy devuelva una lista filtered_dcon los elementos periféricos eliminados en función de una distribución supuesta de los puntos en 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]

Digo 'algo así como' porque la función podría permitir distribuciones variables (poisson, gaussian, etc.) y umbrales de valores atípicos variables dentro de esas distribuciones (como la mque he usado aquí).

aaren
fuente
Relacionado: ¿scipy.stats puede identificar y enmascarar valores atípicos obvios? , aunque esa pregunta parece tratar con situaciones más complejas. Para la tarea simple que describió, un paquete externo parece ser excesivo.
Sven Marnach
Estaba pensando que, dada la cantidad de archivos incorporados en la biblioteca principal, era extraño que no hubiera nada para hacer esto. Parece algo bastante común hacer con datos sin procesar y ruidosos.
aaren

Respuestas:

103

Este método es casi idéntico al suyo, solo que más numpyst (también funciona solo en matrices numpy):

def reject_outliers(data, m=2):
    return data[abs(data - np.mean(data)) < m * np.std(data)]
eumiro
fuente
3
Ese método funciona suficientemente bien si mes lo suficientemente grande (por ejemplo m=6), pero para valores pequeños de meste sufre de la media la varianza no es estimadores robustos.
Benjamin Bannier
30
Sin embargo, eso no es realmente una queja sobre el método, sino una queja sobre la vaga noción de un 'valor atípico'
Eelco Hoogendoorn
¿cómo eliges una m?
john ktejik
1
No he conseguido que esto funcione. Sigo recibiendo un error de retorno de datos [abs (data - np.mean (data)) <m * np.std (data)] TypeError: solo las matrices escalares enteras se pueden convertir a un índice escalar O simplemente congela mi programa
john ktejik
@johnktejik data arg debe ser una matriz numerosa.
Sander van Leeuwen
181

Algo importante cuando se trata de valores atípicos es que se debe intentar utilizar estimadores lo más robustos posible. La media de una distribución estará sesgada por valores atípicos pero, por ejemplo, la mediana será mucho menor.

Basándose en la respuesta de 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]

Aquí he reemplazado la media con la mediana más robusta y la desviación estándar con la mediana de la distancia absoluta a la mediana. Luego escalé las distancias por su (nuevamente) valor mediano para que mesté en una escala relativa razonable.

Tenga en cuenta que para que la data[s<m]sintaxis funcione, datadebe ser una matriz numpy.

Benjamin Bannier
fuente
5
itl.nist.gov/div898/handbook/eda/section3/eda35h.htm este es básicamente el puntaje Z modificado al que se hace referencia aquí, pero con un umbral diferente. Si mis matemáticas son correctas, recomiendan una m de 3.5 / .6745 ~= 5.189(multiplican spor .6745 y especifican una mde 3.5 ... también toman abs(s)). ¿Alguien puede explicar la elección de m? ¿O es algo que identificará a partir de su conjunto de datos en particular?
Charlie G
2
@BenjaminBannier: ¿Puede proporcionar alguna explicación concreta para elegir un valor en mlugar de declaraciones esponjosas como "interacción de pureza y eficiencia"?
stackoverflowuser2010
1
@ stackoverflowuser2010: Como dije, esto depende de sus requisitos específicos, es decir, qué tan limpios necesitamos para señalizar la muestra (falsos positivos), o cuántas mediciones de señal podemos permitirnos desechar para mantener la señal limpia (falsos negativos) . En cuanto a una evaluación de ejemplo específica para un caso de uso determinado, consulte, por ejemplo, desy.de/~blist/notes/whyeffpur.ps.gz .
Benjamin Bannier
2
Recibo el siguiente error cuando llamo a la función con una lista de flotantes:TypeError: only integer scalar arrays can be converted to a scalar index
Vasilis
2
@Charlie, si observa la figura itl.nist.gov/div898/handbook/eda/section3/eda356.htm#MAD , verá que cuando se trata de una distribución normal (que de hecho no es el caso, necesitaría la puntuaciones z modificadas) con SD = 1, tiene MAD ~ 0,68, lo que explica el factor de escala. Por tanto, la elección de m = 3,5 implica que desea descartar el 0,05% de los datos.
Fato39
13

La respuesta de Benjamin Bannier produce un traspaso cuando la mediana de las distancias desde la mediana es 0, por lo que encontré esta versión modificada un poco más útil para los casos que se muestran en el siguiente ejemplo.

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]

Ejemplo:

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

Da:

[[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)
Yigal
fuente
9

Sobre la base de Benjamin, usando pandas.Seriesy reemplazando MAD con 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]

Por ejemplo, si establece iq_range=0.6, los percentiles del rango intercuartílico se convertirían en:, 0.20 <--> 0.80por lo que se incluirán más valores atípicos.

ankostis
fuente
4

Una alternativa es hacer una estimación robusta de la desviación estándar (asumiendo estadísticas gaussianas). Buscando calculadoras en línea, veo que el percentil del 90% corresponde a 1.2815σ y el 95% es 1.645σ ( http://vassarstats.net/tabs.html?#z )

Como un simple ejemplo:

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)

La salida que obtengo es:

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

Que está cerca del valor esperado de 2.

Si queremos eliminar puntos por encima / por debajo de 5 desviaciones estándar (con 1000 puntos esperaríamos 1 valor> 3 desviaciones estándar):

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))

Lo que da:

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

No tengo idea de qué enfoque es el más eficiente / robusto

Chris
fuente
3

Me gustaría proporcionar dos métodos en esta respuesta, una solución basada en la "puntuación z" y una solución basada en "IQR".

El código proporcionado en esta respuesta funciona tanto en una numpymatriz de atenuación única como en una matriz múltiple numpy.

Primero importemos algunos módulos.

import collections
import numpy as np
import scipy.stats as stat
from scipy.stats import iqr

método basado en puntuación z

Este método probará si el número cae fuera de las tres desviaciones estándar. Según esta regla, si el valor es atípico, el método devolverá verdadero, si no, devolverá falso.

def sd_outlier(x, axis = None, bar = 3, side = 'both'):
    assert side in ['gt', 'lt', 'both'], 'Side should be `gt`, `lt` or `both`.'

    d_z = stat.zscore(x, axis = axis)

    if side == 'gt':
        return d_z > bar
    elif side == 'lt':
        return d_z < -bar
    elif side == 'both':
        return np.abs(d_z) > bar

Método basado en IQR

Este método probará si el valor es menor q1 - 1.5 * iqro mayor que q3 + 1.5 * iqr, que es similar al método de trazado de SPSS.

def q1(x, axis = None):
    return np.percentile(x, 25, axis = axis)

def q3(x, axis = None):
    return np.percentile(x, 75, axis = axis)

def iqr_outlier(x, axis = None, bar = 1.5, side = 'both'):
    assert side in ['gt', 'lt', 'both'], 'Side should be `gt`, `lt` or `both`.'

    d_iqr = iqr(x, axis = axis)
    d_q1 = q1(x, axis = axis)
    d_q3 = q3(x, axis = axis)
    iqr_distance = np.multiply(d_iqr, bar)

    stat_shape = list(x.shape)

    if isinstance(axis, collections.Iterable):
        for single_axis in axis:
            stat_shape[single_axis] = 1
    else:
        stat_shape[axis] = 1

    if side in ['gt', 'both']:
        upper_range = d_q3 + iqr_distance
        upper_outlier = np.greater(x - upper_range.reshape(stat_shape), 0)
    if side in ['lt', 'both']:
        lower_range = d_q1 - iqr_distance
        lower_outlier = np.less(x - lower_range.reshape(stat_shape), 0)

    if side == 'gt':
        return upper_outlier
    if side == 'lt':
        return lower_outlier
    if side == 'both':
        return np.logical_or(upper_outlier, lower_outlier)

Por último, si desea filtrar los valores atípicos, use un numpyselector.

Que tengas un buen día.

Pérdidas Don
fuente
3

Tenga en cuenta que todos los métodos anteriores fallan cuando su desviación estándar se vuelve muy grande debido a grandes valores atípicos.

( Simalar como el cálculo promedio falla y debería calcular la mediana. Sin embargo, el promedio es "más propenso a errores como el stdDv". )

Puede intentar aplicar iterativamente su algoritmo o filtrar usando el rango intercuartílico: (aquí, "factor" se relaciona con un rango * sigma, pero solo cuando sus datos siguen una distribución gaussiana)

import numpy as np

def sortoutOutliers(dataIn,factor):
    quant3, quant1 = np.percentile(dataIn, [75 ,25])
    iqr = quant3 - quant1
    iqrSigma = iqr/1.34896
    medData = np.median(dataIn)
    dataOut = [ x for x in dataIn if ( (x > medData - factor* iqrSigma) and (x < medData + factor* iqrSigma) ) ] 
    return(dataOut)
K. Foe
fuente
Lo siento, pasé por alto que ya hay una sugerencia de IQR arriba. ¿Debo dejar esta respuesta de todos modos debido a un código más corto o eliminarlo?
K. Foe
1

Quería hacer algo similar, excepto establecer el número en NaN en lugar de eliminarlo de los datos, ya que si lo elimina, cambia la longitud, lo que puede estropear el trazado (es decir, si solo está eliminando valores atípicos de una columna en una tabla , pero necesita que permanezca igual que las otras columnas para poder trazarlas entre sí).

Para hacerlo, utilicé las funciones de enmascaramiento de numpy :

def reject_outliers(data, m=2):
    stdev = np.std(data)
    mean = np.mean(data)
    maskMin = mean - stdev * m
    maskMax = mean + stdev * m
    mask = np.ma.masked_outside(data, maskMin, maskMax)
    print('Masking values outside of {} and {}'.format(maskMin, maskMax))
    return mask
Alex S
fuente
También puede np.cliplos a los valores mínimos y máximos permitidos para mantener las dimensiones.
Andi R
0

si desea obtener la posición del índice de los valores atípicos idx_list, la devolverá.

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

data_points = np.array([8, 10, 35, 17, 73, 77])  
print(reject_outliers(data_points))

after rejection: [ 8 10 35 17], index positions of outliers: [4 5]
Caner Erden
fuente
0

Para un conjunto de imágenes (cada imagen tiene 3 dimensiones), donde quería rechazar valores atípicos para cada píxel que usé:

mean = np.mean(imgs, axis=0)
std = np.std(imgs, axis=0)
mask = np.greater(0.5 * std + 1, np.abs(imgs - mean))
masked = np.multiply(imgs, mask)

Entonces es posible calcular la media:

masked_mean = np.divide(np.sum(masked, axis=0), np.sum(mask, axis=0))

(Lo uso para la resta de fondo)

ron653
fuente