Tensor Multiplikation mit numpy Tensordot

Question

Ich habe einen Tensor U aus n Matrizen der Dimension (d, k) und einer Matrix V der Dimension (k, n) zusammengesetzt.

Ich möchte sie multiplizieren, so dass das Ergebnis eine Matrix der Dimension zurückkehrt (d, n), in der Spalte j der Matrix zwischen j von U und der Spalte j von V. das Ergebnis der Matrixmultiplikation ist

eine Möglichkeit, dies zu erreichen ist:

for j in range(n): 
    res[:,j] = U[:,:,j] * V[:,j] 

ich frage mich, ob es einen schnelleren Ansatz ist numpy Bibliothek. Insbesondere denke ich an die np.tensordot() Funktion.

Dieses kleine Snippet erlaubt mir, eine einzelne Matrix mit einem Skalar zu multiplizieren, aber die offensichtliche Verallgemeinerung zu einem Vektor bringt nicht zurück, was ich erhofft hatte.

a = np.array(range(1, 17)) 
a.shape = (4,4) 
b = np.array((1,2,3,4,5,6,7)) 
r1 = np.tensordot(b,a, axes=0) 

Irgendwelche Vorschläge?

Answer 1

Es gibt ein paar Möglichkeiten, wie Sie dies tun können. Das erste, was in den Sinn kommt, ist np.einsum:

# some fake data 
gen = np.random.RandomState(0) 
ni, nj, nk = 10, 20, 100 
U = gen.randn(ni, nj, nk) 
V = gen.randn(nj, nk) 

res1 = np.zeros((ni, nk)) 
for k in range(nk): 
    res1[:,k] = U[:,:,k].dot(V[:,k]) 

res2 = np.einsum('ijk,jk->ik', U, V) 

print(np.allclose(res1, res2)) 
# True 

np.einsum verwendet Einstein notation Tensorverjüngung auszudrücken. In dem obigen Ausdruck 'ijk,jk->ik' sind i, j und k Indizes, die den verschiedenen Dimensionen U und V entsprechen. Jede durch Komma getrennte Gruppierung entspricht einem der Operanden, die an np.einsum übergeben wurden (in diesem Fall hat U die Dimensionen ijk und V hat die Dimensionen jk). Der '->ik' Teil gibt die Abmessungen des Ausgabearrays an. Alle Dimensionen mit Indizes, die nicht in der Ausgabezeichenfolge enthalten sind, werden summiert.

np.einsum ist unglaublich nützlich für die Durchführung komplexer Tensor Kontraktionen, aber es kann eine Weile dauern, um Ihren Kopf vollständig umschließen, wie es funktioniert. Sie sollten sich die Beispiele in der Dokumentation (oben verlinkt) ansehen.

---

Einige andere Optionen:

1. Element weise Multiplikation mit broadcasting, durch Summierung gefolgt:

   res3 = (U * V[None, ...]).sum(1) 
   

2. inner1d mit einer Last von der Umsetzung:

   from numpy.core.umath_tests import inner1d 
   
   res4 = inner1d(U.transpose(0, 2, 1), V.T) 
   

Einige Benchmarks:

In [1]: ni, nj, nk = 100, 200, 1000 

In [2]: %%timeit U = gen.randn(ni, nj, nk); V = gen.randn(nj, nk) 
    ....: np.einsum('ijk,jk->ik', U, V) 
    ....: 
10 loops, best of 3: 23.4 ms per loop 

In [3]: %%timeit U = gen.randn(ni, nj, nk); V = gen.randn(nj, nk) 
(U * V[None, ...]).sum(1) 
    ....: 
10 loops, best of 3: 59.7 ms per loop 

In [4]: %%timeit U = gen.randn(ni, nj, nk); V = gen.randn(nj, nk) 
inner1d(U.transpose(0, 2, 1), V.T) 
    ....: 
10 loops, best of 3: 45.9 ms per loop