Python/NumPy: wdrożenie sumę biegu (ale nie całkiem)

Question

Podano dwie tablice o równej długości, jeden Data Holding, jednego gospodarstwa wyniki ale początkowo ustawiony na zero, np:

a = numpy.array([1, 0, 0, 1, 0, 1, 0, 0, 1, 1]) 
b = numpy.array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0]) 

ja jak obliczyć sumę wszystkich możliwych podzbiorów trzech sąsiednich elementów w. Jeżeli suma wynosi 0 lub 1, trzy odpowiadające im elementy w b pozostaną niezmienione; Tylko wtedy, gdy suma przekracza 1 są trzy odpowiadające im elementy b ustawiony na 1, tak aby po obliczenia b się

array([0, 0, 0, 1, 1, 1, 0, 1, 1, 1]) 

prosty pętla osiągnąć:

for x in range(len(a)-2): 
    if a[x:x+3].sum() > 1: 
     b[x:x+3] = 1 

Następnie b ma pożądana forma.

Muszę to zrobić dla dużej ilości danych, więc szybkość jest problemem. Czy w NumPy jest szybszy sposób przeprowadzenia powyższej operacji?

(Rozumiem, że jest to podobne do splotu, ale nie do końca takie samo).

Answer 1

Można zacząć od splotu, wybrać wartości, które przekraczają 1, i wreszcie użyć „rozszerzenie”:

b = numpy.convolve(a, [1, 1, 1], mode="same") > 1 
b = b | numpy.r_[0, b[:-1]] | numpy.r_[b[1:], 0] 

Od tego uniknąć pętli Pythona, powinno być szybsze niż swoim podejściu, ale nie robiłem synchronizacji.

Alternatywą jest użycie drugiego splot rozszerzają:

kernel = [1, 1, 1] 
b = numpy.convolve(a, kernel, mode="same") > 1 
b = numpy.convolve(b, kernel, mode="same") > 0 

Jeśli masz scipy dostępny, kolejna opcja dla dylatacji jest

b = numpy.convolve(a, [1, 1, 1], mode="same") > 1 
b = scipy.ndimage.morphology.binary_dilation(b) 

Edit: Robiąc some timings, Zauważyłem, że to rozwiązanie wydaje się najszybsze dla dużych macierzy:

b = numpy.convolve(a, kernel) > 1 
b[:-1] |= b[1:] # Shift and "smearing" to the *left* (smearing with b[1:] |= b[:-1] does not work) 
b[:-1] |= b[1:] # … and again! 
b = b[:-2] 

Dla tablicy miliona wpisów, było ponad 200 razy szybsze niż oryginalne podejście na moim komputerze. Jak zauważył EOL w komentarzach, to rozwiązanie można uznać za nieco kruche, ponieważ zależy to od szczegółów implementacji NumPy.

Answer 2

można obliczyć „splatania” sum w sposób efektywny z:

>>> a0 = a[:-2] 
>>> a1 = a[1:-1] 
>>> a2 = a[2:] 
>>> a_large_sum = a0 + a1 + a2 > 1 

Uaktualnienie b Następnie można zrobić skutecznie, pisząc coś, co oznacza „co najmniej jeden z trzech sąsiadujących a_large_sum wartości True” : najpierw przedłużyć Ci a_large_sum tablica z powrotem do tej samej liczby elementów jak a (w prawo, w lewo i na prawo, a potem w lewo):

>>> a_large_sum_0 = np.hstack([a_large_sum, [False, False]]) 
>>> a_large_sum_1 = np.hstack([[False], a_large_sum, [False]]) 
>>> a_large_sum_2 = np.hstack([[False, False], a_large_sum]) 

wtedy uzyskać b w sposób efektywny:

>>> b = a_large_sum_0 | a_large_sum_1 | a_large_sum_2 

To daje wynik, który uzyskania, ale w bardzo efektywny sposób, poprzez lewarowanie z NumPy wewnętrznych szybkich pętli.

PS: Takie podejście jest zasadniczo taka sama jak pierwszego rozwiązania Svena, ale jest o wiele bardziej niż pieszy eleganckiego kodu Svena; jest jednak tak szybki. Drugie rozwiązanie Svena (podwójne convolve()) jest jeszcze bardziej eleganckie, i jest dwa razy szybsze.

Answer 3

Możesz również chcieć rzucić okiem na NumPy's stride_tricks. Korzystanie z ustawień taktowania Sven (patrz ogniwem w odpowiedzi Svena), stwierdziliśmy, że dla (bardzo) dużych tablic, jest to także szybki sposób zrobić to, co chcesz (czyli z definicji a):

shape = (len(a)-2,3) 
strides = a.strides+a.strides 
a_strided = numpy.lib.stride_tricks.as_strided(a, shape=shape, strides=strides) 
b = np.r_[numpy.sum(a_strided, axis=-1) > 1, False, False] 
b[2:] |= b[1:-1] | b[:-2] 

Po edycji (patrz komentarze poniżej) nie jest to już najszybszy sposób.

Spowoduje to utworzenie ukośnego widoku na oryginalnej tablicy. Dane w a nie są kopiowane, ale są po prostu oglądane w nowy sposób. Chcemy zasadniczo stworzyć nową tablicę, w której ostatni indeks zawiera podparary, które chcemy zsumować (tj. Trzy elementy, które chcemy zsumować). W ten sposób możemy łatwo podsumować na końcu za pomocą ostatniego polecenia.

Ostatnim elementem tego nowego kształtu musi więc być 3 i pierwszym elementem będzie długość starego a minus 2 (ponieważ możemy tylko sumują się do -2 nd element).

Lista kroków zawiera kroki, w bajtach, które nowa tablica a_strided musi wykonać, aby przejść do następnego elementu w każdym z wymiarów kształtu. Jeśli ustawisz te wartości równe, oznacza to, że a_strided[0,1] i a_strided[1,0] będą zarówno a[1], co jest dokładnie tym, czego chcemy. W normalnej tablicy tak by nie było (pierwszym krokiem byłoby "wymiar pierwszego wymiaru razy długość-tablicy-pierwszego wymiaru (= kształt [0])"), ale w tym przypadku możemy dobrze z niego korzystaj.

Nie jestem pewien, czy wyjaśniłem to wszystko bardzo dobrze, ale po prostu wydrukuj a_strided, a zobaczysz, jaki jest wynik i jak łatwo to zrobi operacja.