Diferencia en el rendimiento entre numpy y matlab.

Estoy calculando el algoritmo de propagación hacia atrás para un autoencoder disperso. Lo he implementado en python usando numpy y en matlab . El código es casi el mismo, pero el rendimiento es muy diferente. El tiempo que toma matlab para completar la tarea es de 0.252454 segundos, mientras que el número 0.973672151566 es casi cuatro veces más. Llamaré a este código varias veces después en un problema de minimización, por lo que esta diferencia conlleva varios minutos de retraso entre las implementaciones. ¿Es este un comportamiento normal? ¿Cómo podría mejorar el rendimiento en numpy?

Implementación Numpy:

Sparse.rho es un parámetro de ajuste, sparse.nodes es el número de nodos en la capa oculta (25), sparse.input (64) el número de nodos en la capa de entrada, theta1 y theta2 son las matrices de peso para el primero y segunda capa, respectivamente, con dimensiones 25×64 y 64×25, m es igual a 10000, rhoest tiene una dimensión de (25,), x tiene una dimensión de 10000×64, a3 10000×64 y a2 10000×25.

UPDATE : He introducido cambios en el código siguiendo algunas de las ideas de las respuestas. El rendimiento ahora es numpy: 0.65 vs matlab: 0.25.

 partial_j1 = np.zeros(sparse.theta1.shape) partial_j2 = np.zeros(sparse.theta2.shape) partial_b1 = np.zeros(sparse.b1.shape) partial_b2 = np.zeros(sparse.b2.shape) t = time.time() delta3t = (-(x-a3)*a3*(1-a3)).T for i in range(m): delta3 = delta3t[:,i:(i+1)] sum1 = np.dot(sparse.theta2.T,delta3) delta2 = ( sum1 + sum2 ) * a2[i:(i+1),:].T* (1 - a2[i:(i+1),:].T) partial_j1 += np.dot(delta2, a1[i:(i+1),:]) partial_j2 += np.dot(delta3, a2[i:(i+1),:]) partial_b1 += delta2 partial_b2 += delta3 print "Backprop time:", time.time() -t 

Implementación de Matlab:

 tic for i = 1:m delta3 = -(data(i,:)-a3(i,:)).*a3(i,:).*(1 - a3(i,:)); delta3 = delta3.'; sum1 = W2.'*delta3; sum2 = beta*(-sparsityParam./rhoest + (1 - sparsityParam) ./ (1.0 - rhoest) ); delta2 = ( sum1 + sum2 ) .* a2(i,:).' .* (1 - a2(i,:).'); W1grad = W1grad + delta2* a1(i,:); W2grad = W2grad + delta3* a2(i,:); b1grad = b1grad + delta2; b2grad = b2grad + delta3; end toc 

Sería erróneo decir que “Matlab siempre es más rápido que NumPy” o viceversa. A menudo su rendimiento es comparable. Cuando use NumPy, para obtener un buen rendimiento debe tener en cuenta que la velocidad de NumPy proviene de las funciones subyacentes escritas en C / C ++ / Fortran. Se desempeña bien cuando se aplican esas funciones a matrices completas. En general, obtiene un peor rendimiento cuando llama a esas funciones NumPy en arreglos más pequeños o escalares en un bucle de Python.

¿Qué pasa con un bucle de Python que preguntas? Cada iteración a través del bucle de Python es una llamada a un método next . Cada uso de [] indexación es una llamada a un método __getitem__ . Cada += es una llamada a __iadd__ . Cada búsqueda de atributos de puntos (como en like np.dot ) implica llamadas a funciones. Esas llamadas de función se sumn a un obstáculo significativo para la velocidad. Estos ganchos le dan a Python un poder expresivo, por ejemplo, la indexación de cadenas significa algo diferente a la indexación de dictados. La misma syntax, diferentes significados. La magia se logra dando a los objetos diferentes métodos __getitem__ .

Pero ese poder expresivo tiene un costo en velocidad. Entonces, cuando no necesite toda esa expresividad dinámica, para obtener un mejor rendimiento, intente limitarse a las llamadas de función NumPy en arreglos completos.

Por lo tanto, eliminar el for-loop; use ecuaciones “vectorizadas” cuando sea posible. Por ejemplo, en lugar de

 for i in range(m): delta3 = -(x[i,:]-a3[i,:])*a3[i,:]* (1 - a3[i,:]) 

puede calcular delta3 para cada i todo a la vez:

 delta3 = -(x-a3)*a3*(1-a3) 

Mientras que en el for-loop delta3 es un vector, usar la ecuación vectorizada delta3 es una matriz.


Algunos de los cálculos en el for-loop no dependen de i y, por lo tanto, deben levantarse fuera del bucle. Por ejemplo, sum2 parece una constante:

 sum2 = sparse.beta*(-float(sparse.rho)/rhoest + float(1.0 - sparse.rho) / (1.0 - rhoest) ) 

Aquí hay un ejemplo ejecutable con una implementación alternativa ( alt ) de su código ( orig ).

Mi punto de referencia de timeit muestra una mejora de velocidad de 6.8x :

 In [52]: %timeit orig() 1 loops, best of 3: 495 ms per loop In [53]: %timeit alt() 10 loops, best of 3: 72.6 ms per loop 

 import numpy as np class Bunch(object): """ http://code.activestate.com/recipes/52308 """ def __init__(self, **kwds): self.__dict__.update(kwds) m, n, p = 10 ** 4, 64, 25 sparse = Bunch( theta1=np.random.random((p, n)), theta2=np.random.random((n, p)), b1=np.random.random((p, 1)), b2=np.random.random((n, 1)), ) x = np.random.random((m, n)) a3 = np.random.random((m, n)) a2 = np.random.random((m, p)) a1 = np.random.random((m, n)) sum2 = np.random.random((p, )) sum2 = sum2[:, np.newaxis] def orig(): partial_j1 = np.zeros(sparse.theta1.shape) partial_j2 = np.zeros(sparse.theta2.shape) partial_b1 = np.zeros(sparse.b1.shape) partial_b2 = np.zeros(sparse.b2.shape) delta3t = (-(x - a3) * a3 * (1 - a3)).T for i in range(m): delta3 = delta3t[:, i:(i + 1)] sum1 = np.dot(sparse.theta2.T, delta3) delta2 = (sum1 + sum2) * a2[i:(i + 1), :].T * (1 - a2[i:(i + 1), :].T) partial_j1 += np.dot(delta2, a1[i:(i + 1), :]) partial_j2 += np.dot(delta3, a2[i:(i + 1), :]) partial_b1 += delta2 partial_b2 += delta3 # delta3: (64, 1) # sum1: (25, 1) # delta2: (25, 1) # a1[i:(i+1),:]: (1, 64) # partial_j1: (25, 64) # partial_j2: (64, 25) # partial_b1: (25, 1) # partial_b2: (64, 1) # a2[i:(i+1),:]: (1, 25) return partial_j1, partial_j2, partial_b1, partial_b2 def alt(): delta3 = (-(x - a3) * a3 * (1 - a3)).T sum1 = np.dot(sparse.theta2.T, delta3) delta2 = (sum1 + sum2) * a2.T * (1 - a2.T) # delta3: (64, 10000) # sum1: (25, 10000) # delta2: (25, 10000) # a1: (10000, 64) # a2: (10000, 25) partial_j1 = np.dot(delta2, a1) partial_j2 = np.dot(delta3, a2) partial_b1 = delta2.sum(axis=1) partial_b2 = delta3.sum(axis=1) return partial_j1, partial_j2, partial_b1, partial_b2 answer = orig() result = alt() for a, r in zip(answer, result): try: assert np.allclose(np.squeeze(a), r) except AssertionError: print(a.shape) print(r.shape) raise 

Consejo: Note que dejé en los comentarios la forma de todos los arreglos intermedios. Saber la forma de las matrices me ayudó a comprender qué hacía su código. La forma de los arreglos puede ayudarlo a guiarlo hacia las funciones NumPy correctas para usar. O al menos, prestar atención a las formas puede ayudarlo a saber si una operación es sensata. Por ejemplo, cuando computas

 np.dot(A, B) 

y A.shape = (n, m) y B.shape = (m, p) , entonces np.dot(A, B) será una matriz de forma (n, p) .


Puede ayudar a construir los arreglos en orden C_CONTIGUOUS (al menos, si se usa np.dot ). Puede haber tanto como una aceleración de 3x al hacerlo:

A continuación, x es lo mismo que xf excepto que x es C_CONTIGUOUS y xf es F_CONTIGUOUS, y la misma relación para y y yf .

 import numpy as np m, n, p = 10 ** 4, 64, 25 x = np.random.random((n, m)) xf = np.asarray(x, order='F') y = np.random.random((m, n)) yf = np.asarray(y, order='F') assert np.allclose(x, xf) assert np.allclose(y, yf) assert np.allclose(np.dot(x, y), np.dot(xf, y)) assert np.allclose(np.dot(x, y), np.dot(xf, yf)) 

%timeit puntos de referencia de %timeit muestran la diferencia de velocidad:

 In [50]: %timeit np.dot(x, y) 100 loops, best of 3: 12.9 ms per loop In [51]: %timeit np.dot(xf, y) 10 loops, best of 3: 27.7 ms per loop In [56]: %timeit np.dot(x, yf) 10 loops, best of 3: 21.8 ms per loop In [53]: %timeit np.dot(xf, yf) 10 loops, best of 3: 33.3 ms per loop 

Respecto al benchmarking en Python:

Puede ser engañoso usar la diferencia en pares de llamadas time.time() para comparar la velocidad del código en Python. Necesitas repetir la medida muchas veces. Es mejor deshabilitar el recolector de basura automático. También es importante medir grandes períodos de tiempo (como un mínimo de 10 segundos de repeticiones) para evitar errores debidos a una resolución deficiente en el temporizador del reloj y para reducir la importancia de la time.time llamadas de time.time . En lugar de escribir todo ese código, Python le proporciona el módulo timeit . Básicamente, lo uso para sincronizar las piezas de código, excepto que lo estoy llamando a través de un terminal IPython para mayor comodidad.

No estoy seguro de si esto está afectando sus puntos de referencia, pero tenga en cuenta que podría hacer una diferencia. En la pregunta a la que me time.time , de acuerdo con time.time dos piezas de código diferían en un factor de 1.7x, mientras que los puntos de referencia que usaban timeit mostraban que las piezas de código corrían en cantidades de tiempo esencialmente idénticas.

Comenzaría con operaciones in situ para evitar asignar nuevas matrices cada vez:

 partial_j1 += np.dot(delta2, a1[i,:].reshape(1,a1.shape[1])) partial_j2 += np.dot(delta3, a2[i,:].reshape(1,a2.shape[1])) partial_b1 += delta2 partial_b2 += delta3 

Puede reemplazar esta expresión:

 a1[i,:].reshape(1,a1.shape[1]) 

De una forma más sencilla y rápida (gracias a Bi Rico ):

 a1[i:i+1] 

Además, esta línea:

 sum2 = sparse.beta*(-float(sparse.rho)/rhoest + float(1.0 - sparse.rho) / (1.0 - rhoest)) 

parece ser el mismo en cada bucle, no es necesario volver a calcularlo.

Y, una optimización probablemente menor, puede reemplazar todas las apariciones de x[i,:] con x[i] .

Finalmente, si puede permitirse asignar la memoria m veces más, puede seguir la sugerencia de unutbu y vectorizar el bucle:

 for m in range(m): delta3 = -(x[i]-a3[i])*a3[i]* (1 - a3[i]) 

con:

 delta3 = -(x-a3)*a3*(1-a3) 

Y siempre puede usar Numba y ganar en velocidad significativamente sin vectorizar (y sin usar más memoria).

La diferencia en el rendimiento entre numpy y matlab siempre me ha frustrado. A menudo, al final se reducen a las bibliotecas de paquetes subyacentes. Por lo que sé, matlab usa el paquete completo de atlas por defecto, mientras que numpy usa una luz de paquete. Matlab reconoce que las personas no se preocupan por el espacio y el volumen, mientras que Numpy reconoce que las personas sí lo hacen. Pregunta similar con una buena respuesta.