Jaki jest najbardziej efektywny sposób mapowania funkcji na tablicy numpy? Sposób, w jaki robiłem to w moim bieżącym projekcie, jest następujący:

import numpy as np 

x = np.array([1, 2, 3, 4, 5])

# Obtain array of square of each element in x
squarer = lambda t: t ** 2
squares = np.array([squarer(xi) for xi in x])

Wydaje się jednak, że jest to prawdopodobnie bardzo nieefektywne, ponieważ używam interpretacji listy, aby skonstruować nową tablicę jako listę Pythona przed przekształceniem jej z powrotem w tablicę numpy.

Czy możemy zrobić lepiej?

Testowałem wszystkie sugerowane metody plus np.array(map(f, x))z perfplot(moim małym projektem).

Wiadomość nr 1: Jeśli możesz korzystać z natywnych funkcji Numpy, zrób to.

Jeśli funkcja, którą próbujesz wektoryzować, jest już wektoryzowana (jak w x**2przykładzie w oryginalnym poście), użycie jej jest znacznie szybsze niż cokolwiek innego (zwróć uwagę na skalę logu):

Jeśli faktycznie potrzebujesz wektoryzacji, nie ma większego znaczenia, którego wariantu używasz.

Kod do odtworzenia wykresów:

import numpy as np
import perfplot
import math


def f(x):
    # return math.sqrt(x)
    return np.sqrt(x)


vf = np.vectorize(f)


def array_for(x):
    return np.array([f(xi) for xi in x])


def array_map(x):
    return np.array(list(map(f, x)))


def fromiter(x):
    return np.fromiter((f(xi) for xi in x), x.dtype)


def vectorize(x):
    return np.vectorize(f)(x)


def vectorize_without_init(x):
    return vf(x)


perfplot.show(
    setup=lambda n: np.random.rand(n),
    n_range=[2 ** k for k in range(20)],
    kernels=[f, array_for, array_map, fromiter, vectorize, vectorize_without_init],
    xlabel="len(x)",
)

Co powiesz na używanie numpy.vectorize.

import numpy as np
x = np.array([1, 2, 3, 4, 5])
squarer = lambda t: t ** 2
vfunc = np.vectorize(squarer)
vfunc(x)
# Output : array([ 1,  4,  9, 16, 25])

TL; DR

Jak zauważył @ user2357112 , „bezpośrednia” metoda zastosowania funkcji jest zawsze najszybszym i najprostszym sposobem mapowania funkcji na tablicach Numpy:

import numpy as np
x = np.array([1, 2, 3, 4, 5])
f = lambda x: x ** 2
squares = f(x)

Ogólnie należy unikać np.vectorize, ponieważ nie działa on dobrze i ma (lub miał) wiele problemów . Jeśli masz do czynienia z innymi typami danych, możesz sprawdzić inne metody przedstawione poniżej.

Porównanie metod

Oto kilka prostych testów służących do porównania trzech metod mapowania funkcji, w tym przykładzie z użyciem Pythona 3.6 i NumPy 1.15.4. Po pierwsze, funkcje konfiguracji do testowania:

import timeit
import numpy as np

f = lambda x: x ** 2
vf = np.vectorize(f)

def test_array(x, n):
    t = timeit.timeit(
        'np.array([f(xi) for xi in x])',
        'from __main__ import np, x, f', number=n)
    print('array: {0:.3f}'.format(t))

def test_fromiter(x, n):
    t = timeit.timeit(
        'np.fromiter((f(xi) for xi in x), x.dtype, count=len(x))',
        'from __main__ import np, x, f', number=n)
    print('fromiter: {0:.3f}'.format(t))

def test_direct(x, n):
    t = timeit.timeit(
        'f(x)',
        'from __main__ import x, f', number=n)
    print('direct: {0:.3f}'.format(t))

def test_vectorized(x, n):
    t = timeit.timeit(
        'vf(x)',
        'from __main__ import x, vf', number=n)
    print('vectorized: {0:.3f}'.format(t))

Testowanie z pięcioma elementami (posortowanymi od najszybszego do najwolniejszego):

x = np.array([1, 2, 3, 4, 5])
n = 100000
test_direct(x, n)      # 0.265
test_fromiter(x, n)    # 0.479
test_array(x, n)       # 0.865
test_vectorized(x, n)  # 2.906

Z setkami elementów:

x = np.arange(100)
n = 10000
test_direct(x, n)      # 0.030
test_array(x, n)       # 0.501
test_vectorized(x, n)  # 0.670
test_fromiter(x, n)    # 0.883

I z tysiącami elementów tablicy lub więcej:

x = np.arange(1000)
n = 1000
test_direct(x, n)      # 0.007
test_fromiter(x, n)    # 0.479
test_array(x, n)       # 0.516
test_vectorized(x, n)  # 0.945

Różne wersje Python / NumPy i optymalizacja kompilatora przyniosą różne wyniki, więc wykonaj podobny test dla swojego środowiska.

Wokół jest numexpr , numba i cython , celem tej odpowiedzi jest wzięcie pod uwagę tych możliwości.

Ale najpierw określmy oczywiste: bez względu na to, jak mapujesz funkcję Pythona na tablicę numpy, pozostaje ona funkcją Pythona, co oznacza dla każdej oceny:

• element numpy-array musi zostać przekonwertowany na obiekt Python (np Float ).
• wszystkie obliczenia są wykonywane za pomocą obiektów Python, co oznacza, że ​​mamy narzut interpretera, dynamicznej wysyłki i niezmiennych obiektów.

Tak więc, która maszyna jest używana do przechodzenia przez tablicę, nie odgrywa dużej roli z powodu wspomnianego powyżej narzutu - pozostaje znacznie wolniejsza niż korzystanie z wbudowanej funkcjonalności numpy.

Rzućmy okiem na następujący przykład:

# numpy-functionality
def f(x):
    return x+2*x*x+4*x*x*x

# python-function as ufunc
import numpy as np
vf=np.vectorize(f)
vf.__name__="vf"

np.vectorizejest wybierany jako reprezentant klasy metod czysto pythonowych. Używając perfplot(patrz kod w załączniku do tej odpowiedzi) otrzymujemy następujące czasy działania:

Widzimy, że podejście numpy jest 10x-100x szybsze niż wersja czysto pythonowa. Spadek wydajności w przypadku większych rozmiarów macierzy jest prawdopodobnie spowodowany tym, że dane nie pasują już do pamięci podręcznej.

Warto również wspomnieć, że vectorizerównież zużywa dużo pamięci, więc często użycie pamięci to szyjka butelki (patrz powiązane pytanie SO ). Zauważ też, że dokumentacja tego numpy nanp.vectorize stwierdza, że ​​jest „przede wszystkim dla wygody, a nie dla wydajności”.

Gdy pożądana jest wydajność, należy użyć innych narzędzi, oprócz napisania rozszerzenia C od zera, istnieją następujące możliwości:

Często słyszy się, że numpy-wydajność jest tak dobra, jak to tylko możliwe, ponieważ pod maską jest czystym C. Ale jest jeszcze wiele do zrobienia!

Wektoryzowana wersja numpy wykorzystuje wiele dodatkowej pamięci i dostęp do pamięci. Biblioteka Numexp próbuje kafelkować tablice numpy, a tym samym uzyskać lepsze wykorzystanie pamięci podręcznej:

# less cache misses than numpy-functionality
import numexpr as ne
def ne_f(x):
    return ne.evaluate("x+2*x*x+4*x*x*x")

Prowadzi do następującego porównania:

Nie mogę wyjaśnić wszystkiego na powyższym wykresie: na początku możemy zobaczyć większy narzut dla biblioteki numexpr, ale ponieważ lepiej wykorzystuje pamięć podręczną, dla większych tablic jest około 10 razy szybszy!

Innym podejściem jest skompilowanie funkcji przez jit, a tym samym uzyskanie prawdziwego UFunc w czystym C. Oto podejście Numba:

# runtime generated C-function as ufunc
import numba as nb
@nb.vectorize(target="cpu")
def nb_vf(x):
    return x+2*x*x+4*x*x*x

Jest 10 razy szybszy niż oryginalne podejście numpy:

Jednak zadanie jest kłopotliwie równoległe, dlatego moglibyśmy również użyć prangedo równoległego obliczenia pętli:

@nb.njit(parallel=True)
def nb_par_jitf(x):
    y=np.empty(x.shape)
    for i in nb.prange(len(x)):
        y[i]=x[i]+2*x[i]*x[i]+4*x[i]*x[i]*x[i]
    return y

Zgodnie z oczekiwaniami funkcja równoległa jest wolniejsza dla mniejszych wejść, ale szybsza (prawie współczynnik 2) dla większych rozmiarów:

Podczas gdy numba specjalizuje się w optymalizacji operacji za pomocą tablic numpy, Cython jest bardziej ogólnym narzędziem. Bardziej skomplikowane jest wyodrębnienie tej samej wydajności, co w przypadku numba - często jest to zależne od llvm (numba) vs lokalnego kompilatora (gcc / MSVC):

%%cython -c=/openmp -a
import numpy as np
import cython

#single core:
@cython.boundscheck(False) 
@cython.wraparound(False) 
def cy_f(double[::1] x):
    y_out=np.empty(len(x))
    cdef Py_ssize_t i
    cdef double[::1] y=y_out
    for i in range(len(x)):
        y[i] = x[i]+2*x[i]*x[i]+4*x[i]*x[i]*x[i]
    return y_out

#parallel:
from cython.parallel import prange
@cython.boundscheck(False) 
@cython.wraparound(False)  
def cy_par_f(double[::1] x):
    y_out=np.empty(len(x))
    cdef double[::1] y=y_out
    cdef Py_ssize_t i
    cdef Py_ssize_t n = len(x)
    for i in prange(n, nogil=True):
        y[i] = x[i]+2*x[i]*x[i]+4*x[i]*x[i]*x[i]
    return y_out

Cython powoduje nieco wolniejsze funkcje:

Wniosek

Oczywiście testowanie tylko jednej funkcji niczego nie dowodzi. Należy również pamiętać, że dla wybranego przykładu funkcji przepustowość pamięci była szyjką butelki dla rozmiarów większych niż 10 ^ 5 elementów - dlatego mieliśmy taką samą wydajność dla numba, numexpr i cython w tym regionie.

Ostatecznie ostateczna odpowiedź zależy od rodzaju funkcji, sprzętu, dystrybucji Pythona i innych czynników. Na przykład Anaconda-dystrybucji używa Intela VML dla funkcji NumPy i tym samym przewyższa Numba (chyba że korzysta SVML, zobacz ten SO-post ) łatwo za transcendentalne funkcje jak exp, sin, cosi podobne - patrz np następującym SO-post .

Jednak na podstawie tego dochodzenia i z dotychczasowego doświadczenia powiedziałbym, że numba wydaje się najłatwiejszym narzędziem o najlepszym działaniu, o ile nie są zaangażowane żadne funkcje transcendentalne.

Rysowanie czasów pracy za pomocą pakietu perfplot:

import perfplot
perfplot.show(
    setup=lambda n: np.random.rand(n),
    n_range=[2**k for k in range(0,24)],
    kernels=[
        f, 
        vf,
        ne_f, 
        nb_vf, nb_par_jitf,
        cy_f, cy_par_f,
        ],
    logx=True,
    logy=True,
    xlabel='len(x)'
    )

squares = squarer(x)

Operacje arytmetyczne na tablicach są automatycznie stosowane elementarnie, z wydajnymi pętlami na poziomie C, które pozwalają uniknąć narzutu interpretera, który miałby zastosowanie do pętli lub rozumienia na poziomie Pythona.

Większość funkcji, które chcesz zastosować elementarnie do tablicy NumPy, po prostu będzie działać, choć niektóre mogą wymagać zmian. Na przykład ifnie działa elementarnie. Chcesz przekonwertować je na konstrukcje takie jak numpy.where:

def using_if(x):
    if x < 5:
        return x
    else:
        return x**2

staje się

def using_where(x):
    return numpy.where(x < 5, x, x**2)

Wierzę, że w nowszej wersji (używam 1.13) numpy możesz po prostu wywołać funkcję, przekazując tablicę numpy do funkcji, którą napisałeś dla typu skalarnego, automatycznie zastosuje wywołanie funkcji do każdego elementu nad tablicą numpy i zwróci ci kolejna tablica liczb

>>> import numpy as np
>>> squarer = lambda t: t ** 2
>>> x = np.array([1, 2, 3, 4, 5])
>>> squarer(x)
array([ 1,  4,  9, 16, 25])

W wielu przypadkach numpy.apply_along_axis najlepszym wyborem będzie . Zwiększa wydajność o około 100x w porównaniu z innymi podejściami - i to nie tylko w przypadku trywialnych funkcji testowych, ale także w przypadku bardziej złożonych kompozycji funkcji od numpy i scipy.

Kiedy dodam metodę:

def along_axis(x):
    return np.apply_along_axis(f, 0, x)

do kodu perfplot, otrzymuję następujące wyniki:

Wygląda na to, że nikt nie wspominał o wbudowanej fabrycznej metodzie produkcji ufuncw paczkach numpy: np.frompyfuncktórą przetestowałem ponownie np.vectorizei osiągnąłem lepsze wyniki o około 20 ~ 30%. Oczywiście będzie działał dobrze jak przepisany kod C, a nawet numba(którego nie testowałem), ale może być lepszą alternatywą niżnp.vectorize

f = lambda x, y: x * y
f_arr = np.frompyfunc(f, 2, 1)
vf = np.vectorize(f)
arr = np.linspace(0, 1, 10000)

%timeit f_arr(arr, arr) # 307ms
%timeit vf(arr, arr) # 450ms

Testowałem także większe próbki, a poprawa jest proporcjonalna. Zobacz dokumentację również tutaj

Jak wspomniano w tym poście , wystarczy użyć wyrażeń generatora takich jak:

numpy.fromiter((<some_func>(x) for x in <something>),<dtype>,<size of something>)

Wszystkie powyższe odpowiedzi dobrze się porównują, ale jeśli potrzebujesz użyć niestandardowej funkcji do mapowania, a masz numpy.ndarrayi musisz zachować kształt tablicy.

Porównałem tylko dwa, ale zachowa kształt ndarray. Do porównania użyłem tablicy z 1 milionem wpisów. Tutaj używam funkcji kwadratowej, która jest również wbudowana w numpy i ma świetne zwiększenie wydajności, ponieważ tam, gdzie było coś potrzebne, możesz użyć wybranej funkcji.

import numpy, time
def timeit():
    y = numpy.arange(1000000)
    now = time.time()
    numpy.array([x * x for x in y.reshape(-1)]).reshape(y.shape)        
    print(time.time() - now)
    now = time.time()
    numpy.fromiter((x * x for x in y.reshape(-1)), y.dtype).reshape(y.shape)
    print(time.time() - now)
    now = time.time()
    numpy.square(y)  
    print(time.time() - now)

Wynik

>>> timeit()
1.162431240081787    # list comprehension and then building numpy array
1.0775556564331055   # from numpy.fromiter
0.002948284149169922 # using inbuilt function

tutaj możesz wyraźnie zobaczyć, że numpy.fromiterdziała świetnie, biorąc pod uwagę proste podejście, a jeśli jest dostępna wbudowana funkcja, skorzystaj z niej.

Posługiwać się numpy.fromfunction(function, shape, **kwargs)

Zobacz „ https://docs.scipy.org/doc/numpy/reference/generated/numpy.fromfunction.html ”