Colonnes d'effondrement de Numpy selon la liste

Question

En numpy, j'ai un tableau adxn A et une liste L de longueur n décrivant où je veux que chaque colonne de A se retrouve dans la matrice B. L'idée est que la colonne i de la matrice B est la somme de toutes les colonnes de A pour lesquelles la valeur correspondante dans L est i. Je peux le faire avec une boucle:

A = np.arange(15).reshape(3,5) 
L = [0,1,2,1,1] 
n_cols = 3 
B = np.zeros((len(A), n_cols)) 
# assume I know the desired number of columns, 
# which is also one more than the maximum value of L 
for i, a in enumerate(A.T): 
    B[:, L[i]] += a 

Est-il possible de le faire en tableau A trancher (ou d'une certaine manière autrement vectorisé)?

Answer 1

Vous êtes somme réduction/effondrement des colonnes de A à l'aide L pour sélectionner ces colonnes. En outre, vous mettez à jour les colonnes du tableau de sortie en fonction de l'unicité de L elems.

Ainsi, vous pouvez utiliser np.add.reduceat pour une solution vectorisée, comme si -

sidx = L.argsort() 
col_idx, grp_start_idx = np.unique(L[sidx],return_index=True) 
B_out = np.zeros((len(A), n_cols)) 
B_out[:,col_idx] = np.add.reduceat(A[:,sidx],grp_start_idx,axis=1) 

---

Test d'exécution -

In [129]: def org_app(A,n_cols): 
    ...:  B = np.zeros((len(A), n_cols)) 
    ...:  for i, a in enumerate(A.T): 
    ...:   B[:, L[i]] += a 
    ...:  return B 
    ...: 
    ...: def vectorized_app(A,n_cols): 
    ...:  sidx = L.argsort() 
    ...:  col_idx, grp_start_idx = np.unique(L[sidx],return_index=True) 
    ...:  B_out = np.zeros((len(A), n_cols)) 
    ...:  B_out[:,col_idx] = np.add.reduceat(A[:,sidx],grp_start_idx,axis=1) 
    ...:  return B_out 
    ...: 

In [130]: # Setup inputs with an appreciable no. of cols & lesser rows 
    ...: # so as that memory bandwidth to work with huge number of 
    ...: # row elems doesn't become the bottleneck 
    ...: d,n_cols = 10,5000 
    ...: A = np.random.rand(d,n_cols) 
    ...: L = np.random.randint(0,n_cols,(n_cols,)) 
    ...: 

In [131]: np.allclose(org_app(A,n_cols),vectorized_app(A,n_cols)) 
Out[131]: True 

In [132]: %timeit org_app(A,n_cols) 
10 loops, best of 3: 33.3 ms per loop 

In [133]: %timeit vectorized_app(A,n_cols) 
100 loops, best of 3: 1.87 ms per loop 

Comme le nombre de lignes deviennent comparables avec le nombre de Col. Dans A, les exigences de largeur de bande de mémoire haute pour l'approche vectorisée compenserait toute accélération notable de celle-ci.

Answer 2

Cette itération sur `` B` est-elle la même (non testée)?

for I in range(B.shape[1]): 
     B[:, I] = A[:, L==I].sum(axis=1) 

Les boucles de nombre seront moins nombreuses. Mais plus important encore, il peut donner d'autres idées de vectorisation. (Edit) lors du test, cela fonctionne, mais il est beaucoup plus lent.

+ ========

scipy.sparse ne résume de la colonne avec le produit de la matrice avec une matrice d'unités. Pouvons-nous faire la même chose ici? Faire C tableau avec 1s dans les colonnes de droite

def my2(A,L): 
    n_cols = L.shape[0] 
    C = np.zeros((n_cols,n_cols),int) 
    C[np.arange(n_cols), L] = 1 
    return A.dot(C) 

Ceci est 7x plus rapide que votre boucle, et un peu plus rapide que le code @Divakarsreduceat.

==========

In [126]: timeit vectorized_app(A,L) 
The slowest run took 8.16 times longer than the fastest. This could mean that an intermediate result is being cached. 
10000 loops, best of 3: 99.7 µs per loop 
In [127]: timeit val2 = my2(A,L) 
The slowest run took 10.46 times longer than the fastest. This could mean that an intermediate result is being cached. 
10000 loops, best of 3: 81.6 µs per loop 
In [128]: timeit org1(A,n_cols) 
1000 loops, best of 3: 623 µs per loop 
In [129]: d,n_cols 
Out[129]: (10, 100)