Tengo una matriz numpy 2-D como x = array([[ 1., 5.],[ 3., 4.]]), tengo que comparar cada fila con todas las demás filas de la matriz y crear una nueva matriz de valores mínimos de ambas filas y tomar la suma de la fila mínima y guardarla en una nueva matriz Finalmente obtendré una matriz simétrica.

Por ejemplo: comparo array [1,5] consigo mismo. La nueva matriz 2-D es array([[ 1., 5.],[ 1., 5.]]), creo una matriz mínima a lo largo del eje = 0, es decir [1., 5.] luego tomo la suma de la matriz que será 6. De manera similar, repito la operación para todas las filas y Termino con una matriz 2 * 2 array([[ 6, 5.],[ 5, 7.]]).

import numpy as np
x=np.array([[1,5],[3,4]])
y=np.zeros((len(x),len(x)))
for i in range(len(x)):
    array_a=x[i]
    for j in range(len(x)):
       array_b=x[j]
       array_c=np.array([array_a,array_b])
       min_array=np.min(array_c,axis=0)
       array_sum=np.sum(min_array)
       y[i,j]=array_sum

Mi matriz 2-D es muy grande y realizar las operaciones mencionadas anteriormente lleva mucho tiempo. Soy nuevo en Python, por lo que cualquier sugerencia para mejorar el rendimiento será realmente útil.

0
Vinay 15 feb. 2018 a las 09:42

2 respuestas

La mejor respuesta

La forma obvia de acelerar su código sería hacer todos los bucles en numpy. Tuve una primera solución (f2 en el código a continuación), que generaría una matriz que contenía todas las combinaciones que deben compararse y luego reduciría esa matriz al resultado final realizando las np.min y { Comandos {X2}}. Desafortunadamente, ese método consume bastante memoria y, por lo tanto, se vuelve lento cuando las matrices son grandes, porque la matriz intermedia es NxNx2xN para una matriz de entrada NxN.

Sin embargo, encontré una solución diferente que usa un bucle for (f3 a continuación) y parece ser razonablemente rápida. La aceleración al original publicado por el OP es aproximadamente 4 veces para una matriz de 1000x1000. Aquí los códigos con algunas pruebas:

import numpy as np
import timeit

def f(x):
    y = np.zeros_like(x)
    for i in range(x.shape[0]):
        a = x[i]
        for j in range(x.shape[1]):
            b = x[j]
            y[i,j] = np.sum(np.min([a,b], axis=0))
    return y

def f2(x):
    y = np.empty((x.shape[0],1,2,x.shape[0]))
    y[:,0,0,:] = x[:,:]
    y = np.repeat(y, x.shape[0],axis=1)
    y[:,:,1,:] = x[:,:]
    return np.sum(np.min(y,axis=2),axis=2)

def f3(x):
    y = np.empty_like(x)
    for i in range(x.shape[1]):
        y[:,i] = np.sum(np.minimum(x[i,:],x[:,:]),axis=1)
    return y

##some testing that the functions work
x = np.array([[1,5],[3,4]])
a=f(x)
b=f2(x)
c=f3(x)
print(np.all(a==b))
print(np.all(a==c))

x = np.array([[1,7,5],[2,3,8],[5,2,4]])
a=f(x)
b=f2(x)
c=f3(x)
print(np.all(a==b))
print(np.all(a==c))

x = np.random.randint(0,10,(100,100))
a=f(x)
b=f2(x)
c=f3(x)
print(np.all(a==b))
print(np.all(a==c))

##some speed testing:
print('-'*50)
print("speed test small")
x = np.random.randint(0,100,(100,100))
print("original")
print(min(timeit.Timer(
    'f(x)',
    setup = 'from __main__ import f,x',
).repeat(3,10)))

print("using np.repeat")
print(min(timeit.Timer(
    'f2(x)',
    setup = 'from __main__ import f2,x',
).repeat(3,10)))

print("one for loop")
print(min(timeit.Timer(
    'f3(x)',
    setup = 'from __main__ import f3,x',
).repeat(3,10)))

print('-'*50)
print("speed test big")
x = np.random.randint(0,100,(1000,1000))
print("original")
print(min(timeit.Timer(
    'f(x)',
    setup = 'from __main__ import f,x',
).repeat(3,1)))

print("one for loop")
print(min(timeit.Timer(
    'f3(x)',
    setup = 'from __main__ import f3,x',
).repeat(3,1)))

Y aquí la salida:

True
True
True
True
True
True
--------------------------------------------------
speed test small
original
1.3070102719939314
using np.repeat
0.15176948899170384
one for loop
0.029766165011096746
--------------------------------------------------
speed test big
original
17.505746565002482
one for loop
4.437685210024938

En otras palabras, f2 es bastante rápido para matrices que no agotan su memoria, pero especialmente para matrices grandes, f3 es lo más rápido que pude encontrar.

EDITAR :

Inspirado por la respuesta de @ Aguy y esta publicación, aquí todavía hay una modificación que solo calcula el triángulo inferior de la matriz y luego copia el resultados al triángulo superior:

def f4(x):
    y = np.empty_like(x)
    for i in range(x.shape[1]):
        y[i:,i] = np.sum(np.minimum(x[i,:],x[i:,:]),axis=1)
    i_upper = np.triu_indices(x.shape[1],1)
    y[i_upper] = y.T[i_upper]
    return y

La prueba de velocidad para la matriz 1000x1000 ahora da

speed test big
original
18.71281115297461
one for loop over lower triangle
2.0939957330119796

EDITAR 2 :

Aquí todavía hay una versión que usa numba para acelerar. Según esta publicación, es mejor escribir los bucles explícitamente en este caso:

import numba as nb
@nb.jit(nopython=True)
def f_nb(x):
    res = np.empty_like(x)
    for j in range(res.shape[1]):
        for i in range(j,res.shape[0]):
            res[j,i] = res[i,j] = np.sum(np.minimum(x[i,:], x[j,:]))

    return res

Y las pruebas de velocidad relevantes dan:

  • 0.015975199989043176 para una matriz de 100 x 100
  • 0.37946902704425156 para una matriz de 1000 x 1000
  • 467.06363476096885 para una matriz de 10000 x 10000

La prueba de velocidad de 10000x10000 para f4 no parecía querer terminar en absoluto, así que la dejé fuera. Si sus matrices se vuelven mucho más grandes que eso, es posible que tenga problemas de memoria, ¿consideró esto?

0
Thomas Kühn 16 feb. 2018 a las 13:31

La mejora obvia para ahorrar aproximadamente la mitad del tiempo es ejecutar solo en índices i> = j. Por elegancia y algo de ahorro también puede usar menos variables.

import numpy as np
import time

x=np.random.randint(0, 10, (500, 500))
y=np.zeros((len(x),len(x)))

# OP version
t0 = time.time()
for i in range(len(x)):
    array_a=x[i]
    for j in range(len(x)):
       array_b=x[j]
       array_c=np.array([array_a,array_b])
       min_array=np.min(array_c,axis=0)
       array_sum=np.sum(min_array)
       y[i,j]=array_sum
print(time.time() - t0)

z=np.zeros((len(x),len(x)))

# modified version
t0 = time.time()
for i in range(len(x)):
    for j in range(i, len(x)):
        z[i, j]=np.sum(np.min([x[i], x[j]], axis=0))
        z[j, i] = z[i, j]
print(time.time() - t0)

# verify that the result are the same
print(np.all(z == y))

Los resultados en mi máquina:

4.2974278926849365
2.746302604675293
True
1
Aguy 15 feb. 2018 a las 10:03