Czy istnieje numpy wbudowane do odrzucania wartości odstających z listy


101

Czy jest wbudowana funkcja numpy, która wykonuje coś podobnego do następującego? Oznacza to, że weź listę di zwróć listę filtered_dz usuniętymi elementami zewnętrznymi na podstawie pewnego zakładanego rozkładu punktów w 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]

Mówię „coś podobnego”, ponieważ funkcja może zezwalać na różne rozkłady (poissona, gaussa itp.) I różne progi wartości odstających w tych rozkładach (tak jak mużyłem tutaj).


Powiązane: Czy scipy.stats może zidentyfikować i zamaskować oczywiste wartości odstające? , chociaż to pytanie wydaje się dotyczyć bardziej złożonych sytuacji. W przypadku prostego zadania, które opisałeś, zewnętrzny pakiet wydaje się być przesadą.
Sven Marnach

Myślałem, że biorąc pod uwagę liczbę wbudowanych w głównej bibliotece numpy, było dziwne, że nie ma nic do zrobienia. Wydaje się, że jest to dość powszechne w przypadku surowych, zaszumionych danych.
aaren

Odpowiedzi:


104

Ta metoda jest prawie identyczna z twoją, tylko bardziej numpyst (działa również tylko na tablicach numpy):

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

3
Ta metoda działa wystarczająco dobrze, jeśli mjest dostatecznie duża (np. m=6), Ale dla małych wartości mta obarczona jest średnią wariancją, która nie jest solidnymi estymatorami.
Benjamin Bannier

30
nie jest to jednak skarga dotycząca metody, ale skarga na niejasne pojęcie „wartości odstającej”
Eelco Hoogendoorn

jak wybrać m?
john ktejik

1
Nie udało mi się to. Ciągle otrzymuję błąd zwracany dane [abs (dane - np.mean (dane)) <m * np.std (dane)] TypeError: tylko całkowite tablice skalarne można przekonwertować na indeks skalarny LUB po prostu zawiesza mój program
john ktejik

@johnktejik argument danych musi być tablicą numpy.
Sander van Leeuwen

181

W przypadku wartości odstających ważne jest, aby starać się używać estymatorów tak solidnych, jak to tylko możliwe. Średnia rozkładu będzie obciążona wartościami odstającymi, ale np. Mediana będzie znacznie mniejsza.

Opierając się na odpowiedzi 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]

Tutaj zastąpiłem średnią bardziej solidną medianą, a odchylenie standardowe medianą bezwzględnej odległości od mediany. Następnie przeskalowałem odległości przez ich (ponownie) medianę, tak aby mbyła na rozsądnej skali względnej.

Zauważ, że aby data[s<m]składnia działała, datamusi być tablicą numpy.


5
itl.nist.gov/div898/handbook/eda/section3/eda35h.htm jest to zasadniczo zmodyfikowany punkt Z-score, o którym mowa w tym miejscu, ale z innym progiem. Jeśli moja matematyka jest prawidłowa, zalecają m 3.5 / .6745 ~= 5.189(mnożą sprzez 0,6745 i określają m3,5 ... również biorą abs(s)). Czy ktoś może wyjaśnić wybór m? A może jest to coś, co zidentyfikujesz w swoim konkretnym zbiorze danych?
Charlie G

2
@BenjaminBannier: Czy możesz podać konkretne wyjaśnienie wyboru wartości mzamiast puszystych stwierdzeń, takich jak „wzajemne oddziaływanie czystości i wydajności”?
stackoverflowuser2010

1
@ stackoverflowuser2010: Tak jak powiedziałem, zależy to od twoich konkretnych wymagań, tj. jak czysto musimy sygnalizować próbkę (fałszywie dodatnie) lub ile pomiarów sygnału możemy sobie pozwolić na wyrzucenie, aby utrzymać czysty sygnał (fałszywe negatywy) . Jeśli chodzi o konkretną przykładową ocenę dla określonego przypadku użycia, patrz np . Desy.de/~blist/notes/whyeffpur.ps.gz .
Benjamin Bannier

2
Otrzymuję następujący błąd, kiedy wywołuję funkcję z listą pływaków:TypeError: only integer scalar arrays can be converted to a scalar index
Vasilis

2
@Charlie, jeśli spojrzysz na rysunek itl.nist.gov/div898/handbook/eda/section3/eda356.htm#MAD , zobaczysz, że mając do czynienia z rozkładem normalnym (co w rzeczywistości nie ma miejsca, potrzebujesz zmodyfikowane wyniki z) z SD = 1, masz MAD ~ 0,68, co wyjaśnia współczynnik skalowania. Wybór m = 3,5 oznacza zatem, że chcesz odrzucić 0,05% danych.
Fato39

13

Odpowiedź Benjamina Banniera daje wynik pass-through, gdy mediana odległości od mediany wynosi 0, więc uważam, że ta zmodyfikowana wersja jest nieco bardziej pomocna w przypadkach, jak podano w poniższym przykładzie.

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]

Przykład:

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

Daje:

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

9

Opierając się na Benjamin's, używając pandas.Seriesi zastępując MAD przez 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]

Na przykład, jeśli ustawisz iq_range=0.6, percentyle z międzykwartylowego zasięgu staną: 0.20 <--> 0.80tak więcej odstających zostaną uwzględnione.


4

Alternatywą jest dokonanie solidnego oszacowania odchylenia standardowego (przy założeniu statystyki Gaussa). Patrząc na kalkulatory online, widzę, że 90% percentyl odpowiada 1,2815σ, a 95% to 1,645σ ( http://vassarstats.net/tabs.html?#z )

Jako prosty przykład:

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)

Wynik, który otrzymuję, to:

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

Co jest zbliżone do oczekiwanej wartości 2.

Jeśli chcemy usunąć punkty powyżej / poniżej 5 odchyleń standardowych (przy 1000 punktów oczekiwalibyśmy 1 wartość> 3 odchylenia standardowe):

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

Co daje:

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

Nie mam pojęcia, które podejście jest bardziej wydajne / solidne


3

W tej odpowiedzi chciałbym podać dwie metody, rozwiązanie oparte na „z score” i rozwiązanie oparte na „IQR”.

Kod podany w tej odpowiedzi działa zarówno na pojedynczej numpytablicy dim, jak i na wielu numpytablicach.

Najpierw zaimportujmy niektóre moduły.

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

Metoda oparta na punktacji z

Ta metoda sprawdzi, czy liczba wykracza poza trzy odchylenia standardowe. Na podstawie tej reguły, jeśli wartość jest odstająca, metoda zwróci true, jeśli nie, zwróci false.

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

Metoda oparta na IQR

Ta metoda sprawdzi, czy wartość jest mniejsza q1 - 1.5 * iqrlub większa niż q3 + 1.5 * iqr, co jest podobne do metody wykresu 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)

Wreszcie, jeśli chcesz odfiltrować wartości odstające, użyj numpyselektora.

Miłego dnia.


3

Weź pod uwagę, że wszystkie powyższe metody zawodzą, gdy odchylenie standardowe staje się bardzo duże z powodu dużych wartości odstających.

( Podobnie jak średnia kalkulacja zawodzi i powinna raczej obliczyć medianę. Chociaż średnia jest „bardziej podatna na taki błąd jak stdDv” ).

Możesz spróbować iteracyjnie zastosować swój algorytm lub filtrować za pomocą zakresu międzykwartylowego: (tutaj „współczynnik” odnosi się do zakresu * sigma, ale tylko wtedy, gdy dane są zgodne z rozkładem Gaussa)

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)

Przepraszam, przeoczyłem, że jest już sugestia dotycząca IQR powyżej. Czy mimo to powinienem zostawić tę odpowiedź ze względu na krótszy kod, czy go usunąć?
K. Foe

1

Chciałem zrobić coś podobnego, z wyjątkiem ustawienia liczby na NaN zamiast usuwania jej z danych, ponieważ jeśli ją usuniesz, zmienisz długość, która może zepsuć wykres (tj. Jeśli usuwasz wartości odstające tylko z jednej kolumny w tabeli , ale potrzebujesz, aby pozostała taka sama jak inne kolumny, aby można było wykreślić je względem siebie).

Aby to zrobić, użyłem funkcji maskujących 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

Możesz też np. Przyciąć je do minimalnych i maksymalnych dozwolonych wartości, aby zachować wymiary.
Andi R

0

jeśli chcesz uzyskać pozycję indeksu wartości odstających idx_list, zwróci ją.

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]

0

Dla zestawu obrazów (każdy obraz ma 3 wymiary), w którym chciałem odrzucić wartości odstające dla każdego piksela, którego użyłem:

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)

Wtedy można obliczyć średnią:

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

(Używam go do odejmowania w tle)

Korzystając z naszej strony potwierdzasz, że przeczytałeś(-aś) i rozumiesz nasze zasady używania plików cookie i zasady ochrony prywatności.
Licensed under cc by-sa 3.0 with attribution required.