vectorize array python performance numpy for-loop numpy-ufunc

python - array - numpy ufuncs velocidad vs velocidad de bucle



numpy vectorize (1)

Este es el comportamiento normal y esperado. Simplemente está demasiado simplificado para aplicar la afirmación "evitar para bucles con numpy" . Si se trata de bucles internos, es (casi) siempre cierto. Pero en el caso de los bucles externos (como en su caso) hay muchas más excepciones. Especialmente si la alternativa es usar transmisión porque esto acelera su operación al usar mucha más memoria.

Solo para agregar un poco de fondo a esa declaración de "evitar los bucles con numpy" :

Las matrices NumPy se almacenan como matrices contiguas con tipos de c . ¡Python int no es lo mismo que C int ! Por lo tanto, cada vez que itere sobre cada elemento de una matriz, debe conectar el elemento de la matriz, convertirlo a Python int y luego hacer lo que quiera con él y, finalmente, puede necesitar convertirlo a entero CA (llamado boxeo y desempaquetar el valor). Por ejemplo, desea sum los elementos en una matriz usando Python:

import numpy as np arr = np.arange(1000) %%timeit acc = 0 for item in arr: acc += item # 1000 loops, best of 3: 478 µs per loop

Es mejor que uses numpy:

%timeit np.sum(arr) # 10000 loops, best of 3: 24.2 µs per loop

Incluso si presionas el lazo en el código Python C estás muy lejos del rendimiento numpy:

%timeit sum(arr) # 1000 loops, best of 3: 387 µs per loop

Puede haber excepciones a esta regla, pero estas serán realmente dispersas. Al menos mientras haya una funcionalidad numpy equivalente. Entonces, si iteraras sobre elementos individuales, deberías usar numpy.

A veces, un ciclo simple de python es suficiente. No se anuncia ampliamente, pero las funciones numpy tienen una gran sobrecarga en comparación con las funciones de Python. Por ejemplo, considere una matriz de 3 elementos:

arr = np.arange(3) %timeit np.sum(arr) %timeit sum(arr)

¿Cuál será más rápido?

Solución: la función Python funciona mejor que la solución numpy:

# 10000 loops, best of 3: 21.9 µs per loop <- numpy # 100000 loops, best of 3: 6.27 µs per loop <- python

Pero, ¿qué tiene esto que ver con tu ejemplo? No mucho de hecho porque siempre usa funciones numpy en arreglos (no elementos individuales y ni siquiera algunos elementos) por lo que su bucle interno ya usa las funciones optimizadas. Es por eso que ambos realizan más o menos lo mismo (+/- un factor 10 con muy pocos elementos a un factor 2 en aproximadamente 500 elementos). Pero eso no es realmente el ciclo de arriba, es la llamada de la función de arriba!

Su solución de bucle

Usando line-profiler y una resolution = 100 :

def fun_func(tim, prec, values): for i, ti in enumerate(tim): values[i] = np.sum(np.sin(prec * ti)) %lprun -f fun_func fun_func(tim, prec, values) Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 101 752 7.4 5.7 for i, ti in enumerate(tim): 3 100 12449 124.5 94.3 values[i] = np.sum(np.sin(prec * ti))

El 95% se gasta dentro del bucle, incluso he dividido el cuerpo del bucle en varias partes para verificar esto:

def fun_func(tim, prec, values): for i, ti in enumerate(tim): x = prec * ti x = np.sin(x) x = np.sum(x) values[i] = x %lprun -f fun_func fun_func(tim, prec, values) Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 101 609 6.0 3.5 for i, ti in enumerate(tim): 3 100 4521 45.2 26.3 x = prec * ti 4 100 4646 46.5 27.0 x = np.sin(x) 5 100 6731 67.3 39.1 x = np.sum(x) 6 100 714 7.1 4.1 values[i] = x

El momento en que los consumidores son np.multiply , np.sin , np.sum aquí, como puede verificar fácilmente al comparar su tiempo por llamada con su sobrecarga:

arr = np.ones(1, float) %timeit np.sum(arr) # 10000 loops, best of 3: 22.6 µs per loop

Entonces, tan pronto como la sobrecarga de llamada de la función comulativa sea pequeña en comparación con el tiempo de ejecución de cálculo, tendrá tiempos de ejecución similares. Incluso con 100 elementos, estás bastante cerca del tiempo de sobrecarga. El truco es saber en qué punto se equilibran. Con 1000 elementos, la sobrecarga de llamadas sigue siendo significativa:

%lprun -f fun_func fun_func(tim, prec, values) Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 1001 5864 5.9 2.4 for i, ti in enumerate(tim): 3 1000 42817 42.8 17.2 x = prec * ti 4 1000 119327 119.3 48.0 x = np.sin(x) 5 1000 73313 73.3 29.5 x = np.sum(x) 6 1000 7287 7.3 2.9 values[i] = x

Pero con una resolution = 5000 la sobrecarga es bastante baja en comparación con el tiempo de ejecución:

Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 5001 29412 5.9 0.9 for i, ti in enumerate(tim): 3 5000 388827 77.8 11.6 x = prec * ti 4 5000 2442460 488.5 73.2 x = np.sin(x) 5 5000 441337 88.3 13.2 x = np.sum(x) 6 5000 36187 7.2 1.1 values[i] = x

Cuando gastas 500us en cada llamada np.sin ya no te importan los 20us de gastos generales.

Una palabra de precaución es probable en orden: line_profiler probablemente incluye algunos gastos generales adicionales por línea, tal vez también por llamada de función, por lo que el punto en el que la sobrecarga de llamada de función se vuelve insignificante podría ser menor!

Su solución de transmisión

Empecé por perfilar la primera solución, hagamos lo mismo con la segunda solución:

def fun_func(tim, prec, values): x = tim[:, np.newaxis] x = x * prec x = np.sin(x) x = np.sum(x, axis=1) return x

De nuevo usando line_profiler con resolution=100 :

%lprun -f fun_func fun_func(tim, prec, values) Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 1 27 27.0 0.5 x = tim[:, np.newaxis] 3 1 638 638.0 12.9 x = x * prec 4 1 3963 3963.0 79.9 x = np.sin(x) 5 1 326 326.0 6.6 x = np.sum(x, axis=1) 6 1 4 4.0 0.1 return x

Esto ya excede significativamente el tiempo de sobrecarga y, por lo tanto, terminamos un factor 10 más rápido en comparación con el ciclo.

También hice el perfil para la resolution=1000 :

Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 1 28 28.0 0.0 x = tim[:, np.newaxis] 3 1 17716 17716.0 14.6 x = x * prec 4 1 91174 91174.0 75.3 x = np.sin(x) 5 1 12140 12140.0 10.0 x = np.sum(x, axis=1) 6 1 10 10.0 0.0 return x

y con precision=5000 :

Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def fun_func(tim, prec, values): 2 1 34 34.0 0.0 x = tim[:, np.newaxis] 3 1 333685 333685.0 11.1 x = x * prec 4 1 2391812 2391812.0 79.6 x = np.sin(x) 5 1 280832 280832.0 9.3 x = np.sum(x, axis=1) 6 1 14 14.0 0.0 return x

El tamaño de 1000 es aún más rápido, pero como hemos visto, la sobrecarga de llamadas aún no era insignificante en la solución de bucle. Pero para la resolution = 5000 el tiempo que se pasa en cada paso es casi idéntico (algunos son un poco más lentos, otros más rápidos pero en general bastante similares)

Otro efecto es que la transmisión real cuando haces la multiplicación se vuelve significativa. Incluso con las soluciones numpy muy inteligentes esto todavía incluye algunos cálculos adicionales. Para la resolution=10000 , verá que la multiplicación de difusión comienza a ocupar más "% de tiempo" en relación con la solución de bucle:

Line # Hits Time Per Hit % Time Line Contents ============================================================== 1 def broadcast_solution(tim, prec, values): 2 1 37 37.0 0.0 x = tim[:, np.newaxis] 3 1 1783345 1783345.0 13.9 x = x * prec 4 1 9879333 9879333.0 77.1 x = np.sin(x) 5 1 1153789 1153789.0 9.0 x = np.sum(x, axis=1) 6 1 11 11.0 0.0 return x Line # Hits Time Per Hit % Time Line Contents ============================================================== 8 def loop_solution(tim, prec, values): 9 10001 62502 6.2 0.5 for i, ti in enumerate(tim): 10 10000 1287698 128.8 10.5 x = prec * ti 11 10000 9758633 975.9 79.7 x = np.sin(x) 12 10000 1058995 105.9 8.6 x = np.sum(x) 13 10000 75760 7.6 0.6 values[i] = x

Pero hay otra cosa además del tiempo real: el consumo de memoria. Su solución de bucle requiere memoria O(n) porque siempre procesa n elementos. La solución de radiodifusión, sin embargo, requiere memoria O(n*n) . Probablemente tengas que esperar un tiempo si usas resolution=20000 con tu loop pero solo necesitaría 8bytes/element * 20000 element ~= 160kB pero con la transmisión necesitarás ~3GB . ¡Y esto descuida el factor constante (como matrices temporales también conocidas como matrices intermedias)! ¡Supongamos que va más allá y se quedará sin memoria muy rápido!

Es hora de resumir los puntos de nuevo:

  • Si haces un ciclo de python sobre elementos individuales en una matriz numpy, lo estás haciendo mal.
  • Si recorre subarreglos de una matriz numpy, asegúrese de que la sobrecarga de llamada de función en cada bucle sea insignificante en comparación con el tiempo empleado en la función.
  • Si transmite matrices numpy, asegúrese de no quedarse sin memoria.

Pero el punto más importante sobre la optimización sigue siendo:

  • ¡Solo optimice el código si es demasiado lento! Si es demasiado lento, optimícelo solo después de crear un perfil de su código.

  • No confíe ciegamente en declaraciones simplificadas y nunca más optimice sin crear perfiles.

Un último pensamiento:

Las funciones que requieren un bucle o difusión pueden implementarse fácilmente usando cython , numba o numexpr si no hay una solución ya existente en numpy o scipy .

Por ejemplo, una función numba que combina la eficiencia de memoria de la solución de bucle con la velocidad de la solución de difusión a baja resolutions se vería así:

from numba import njit import math @njit def numba_solution(tim, prec, values): size = tim.size for i in range(size): ti = tim[i] x = 0 for j in range(size): x += math.sin(prec[j] * ti) values[i] = x

Como se señala en los comentarios, numexpr también puede evaluar el cálculo emitido muy rápido y sin requerir memoria O(n*n) :

>>> import numexpr >>> tim_2d = tim[:, np.newaxis] >>> numexpr.evaluate(''sum(sin(tim_2d * prec), axis=1)'')

He leído mucho "evitar los bucles con numpy". Entonces, lo intenté. Estaba usando este código (versión simplificada). Algunos datos auxiliares:

In[1]: import numpy as np resolution = 1000 # this parameter varies tim = np.linspace(-np.pi, np.pi, resolution) prec = np.arange(1, resolution + 1) prec = 2 * prec - 1 values = np.zeros_like(tim)

Mi primera implementación fue con for loop:

In[2]: for i, ti in enumerate(tim): values[i] = np.sum(np.sin(prec * ti))

Luego, me deshice del ciclo explícito y logré esto:

In[3]: values = np.sum(np.sin(tim[:, np.newaxis] * prec), axis=1)

Y esta solución fue más rápida para matrices pequeñas, pero cuando amplié la escala tuve una dependencia tan grande del tiempo:

¿Qué me falta o es un comportamiento normal? Y si no lo es, ¿dónde cavar?

EDITAR : Según los comentarios, aquí hay algo de información adicional. Los tiempos se midieron con %timeit y %%timeit , cada ejecución se realizó en el kernel nuevo. Mi computadora portátil es acer aspire v7-482pg (i7, 8GB) . Estoy usando:

  • Python 3.5.2
  • numpy 1.11.2 + mkl
  • Windows 10