Qual é a maneira mais eficiente de mapear uma função em uma matriz numpy? A maneira como eu faço isso no meu projeto atual é a seguinte:

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

No entanto, isso parece provavelmente muito ineficiente, pois estou usando uma compreensão de lista para construir a nova matriz como uma lista Python antes de convertê-la novamente em uma matriz numpy.

Podemos fazer melhor?

Eu testei todos os métodos sugeridos mais np.array(map(f, x))com perfplot(um pequeno projeto meu).

Mensagem 1: Se você pode usar as funções nativas do numpy, faça isso.

Se a função que você está tentando vetorizar já estiver vetorizada (como o x**2exemplo na postagem original), usar isso é muito mais rápido do que qualquer outra coisa (observe a escala de log):

Se você realmente precisa de vetorização, não importa muito qual variante você usa.

Código para reproduzir as parcelas:

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

Que tal usar 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

Conforme observado por @ user2357112 , um método "direto" de aplicar a função é sempre a maneira mais rápida e simples de mapear uma função sobre matrizes Numpy:

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

Geralmente evite np.vectorize, pois não apresenta um bom desempenho e possui (ou teve) vários problemas . Se você estiver lidando com outros tipos de dados, convém investigar os outros métodos mostrados abaixo.

Comparação de métodos

Aqui estão alguns testes simples para comparar três métodos para mapear uma função, usando este exemplo com Python 3.6 e NumPy 1.15.4. Primeiro, as funções de configuração para teste:

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

Teste com cinco elementos (classificados do mais rápido para o mais lento):

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

Com centenas de elementos:

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

E com milhares de elementos de matriz ou mais:

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

Versões diferentes do Python / NumPy e otimização do compilador terão resultados diferentes, portanto faça um teste semelhante para o seu ambiente.

Existem numexpr , numba e cython , o objetivo desta resposta é levar essas possibilidades em consideração.

Mas primeiro vamos declarar o óbvio: não importa como você mapeie uma função Python em um array numpy, ela permanece uma função Python, o que significa para todas as avaliações:

• O elemento numpy-array deve ser convertido em um objeto Python (por exemplo, a Float).
• todos os cálculos são feitos com objetos Python, o que significa ter a sobrecarga de intérprete, despacho dinâmico e objetos imutáveis.

Portanto, quais máquinas são usadas para fazer um loop na matriz não desempenham um grande papel por causa da sobrecarga mencionada acima - ele permanece muito mais lento do que usar a funcionalidade incorporada do numpy.

Vamos dar uma olhada no seguinte exemplo:

# 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.vectorizeé escolhido como um representante da classe de abordagens da função python puro. Usando perfplot(veja o código no apêndice desta resposta) obtemos os seguintes tempos de execução:

Podemos ver que a abordagem numpy é 10x-100x mais rápida que a versão python pura. A diminuição do desempenho para tamanhos de matriz maiores é provavelmente porque os dados não se ajustam mais ao cache.

Vale mencionar também, que vectorizetambém usa muita memória, e muitas vezes o uso da memória é o gargalo (consulte a pergunta SO relacionada ). Observe também que a documentação da numpy np.vectorizeafirma que é "fornecida principalmente por conveniência, não por desempenho".

Outras ferramentas devem ser usadas, quando o desempenho é desejado, além de escrever uma extensão C a partir do zero, existem as seguintes possibilidades:

Ouve-se com frequência que o desempenho numpy é tão bom quanto ele ganha, porque é puro C sob o capô. No entanto, há muito espaço para melhorias!

A versão numpy vetorizada usa muita memória e acessos à memória adicionais. A biblioteca Numexp tenta agrupar as matrizes numpy e, assim, obter uma melhor utilização do cache:

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

Leva à seguinte comparação:

Não posso explicar tudo no gráfico acima: podemos ver uma sobrecarga maior para a biblioteca numexpr no início, mas como ela utiliza melhor o cache, é cerca de 10 vezes mais rápida para matrizes maiores!

Outra abordagem é compilar rapidamente a função e, assim, obter um UFunc C puro puro. Esta é a abordagem da 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

É 10 vezes mais rápido que a abordagem numpy original:

No entanto, a tarefa é embaraçosamente paralelelizável, portanto, também poderíamos usar prangepara calcular o loop em paralelo:

@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

Como esperado, a função paralela é mais lenta para entradas menores, mas mais rápida (quase fator 2) para tamanhos maiores:

Enquanto a numba se especializa em otimizar operações com matrizes numpy, o Cython é uma ferramenta mais geral. É mais complicado extrair o mesmo desempenho que o numba - geralmente é o llvm (numba) versus o compilador local (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

O Cython resulta em funções um pouco mais lentas:

Conclusão

Obviamente, testar apenas uma função não prova nada. Também devemos ter em mente que, para o exemplo de função escolhido, a largura de banda da memória era o gargalo para tamanhos maiores que 10 ^ 5 elementos - portanto, tivemos o mesmo desempenho para numba, numexpr e cython nessa região.

No final, a resposta definitiva depende do tipo de função, hardware, distribuição Python e outros fatores. Por exemplo Anaconda-de distribuição usa VML da Intel para funções de numpy e assim Supera numba (a menos que ele usa SVML, consulte este SO-post ) facilmente para funções transcendentais como exp, sin, cose semelhante - ver, por exemplo o seguinte SO-post .

No entanto, a partir desta investigação e da minha experiência até agora, eu afirmaria que o numba parece ser a ferramenta mais fácil com melhor desempenho, desde que nenhuma função transcendental esteja envolvida.

Plotando tempos de execução com perfplot -package :

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)

As operações aritméticas em matrizes são aplicadas automaticamente de maneira elementar, com loops eficientes no nível C que evitam toda a sobrecarga do interpretador que se aplicaria a um loop ou compreensão no nível do Python.

A maioria das funções que você deseja aplicar a uma matriz NumPy elementwise funcionará, embora algumas possam precisar de alterações. Por exemplo, ifnão funciona de maneira elementar. Você deseja convertê-los para usar construções como numpy.where:

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

torna-se

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

Eu acredito que na versão mais recente (eu uso o 1.13) do numpy, você pode simplesmente chamar a função passando a matriz numpy para a função que você escreveu para o tipo escalar, ela aplicará automaticamente a chamada de função a cada elemento na matriz numpy e retornará você outra matriz numpy

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

Em muitos casos, numpy.apply_along_axis será a melhor opção. Aumenta o desempenho em cerca de 100x em comparação com as outras abordagens - e não apenas para funções triviais de teste, mas também para composições de funções mais complexas, como numpy e scipy.

Quando adiciono o método:

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

para o código perfplot, obtenho os seguintes resultados:

Parece que ninguém mencionou um método de fábrica embutido para produzir ufuncem embalagens numpy: np.frompyfuncque eu testei novamente np.vectorizee superei em cerca de 20 a 30%. Obviamente, ele funcionará bem como o código C prescrito ou mesmo numba(que eu não testei), mas pode ser uma alternativa melhor do quenp.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

Também testei amostras maiores e a melhoria é proporcional. Veja a documentação também aqui

Conforme mencionado neste post , basta usar expressões geradoras como estas :

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

Todas as respostas acima se comparam bem, mas se você precisar usar a função personalizada para mapeamento, e tiver numpy.ndarray, e precisar manter a forma da matriz.

Comparei apenas dois, mas ele manterá a forma de ndarray. Eu usei a matriz com 1 milhão de entradas para comparação. Aqui eu uso a função quadrada, que também está embutida em numpy e tem um ótimo desempenho, já que, como havia necessidade de algo, você pode usar a função de sua escolha.

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)

Resultado

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

aqui você pode ver claramente que numpy.fromiterfunciona muito bem considerando uma abordagem simples e, se a função embutida estiver disponível, use-a.

Usar numpy.fromfunction(function, shape, **kwargs)

Consulte " https://docs.scipy.org/doc/numpy/reference/generated/numpy.fromfunction.html "