Schleife über (oder Vektorisieren) variabler Länge in Matrizen Theano

Question

Ich habe eine Liste von Matrizen L, wobei jedes Element eine Mx*n Matrix (x eine Variable ist, n ist eine Konstante).

Ich möchte die Summe von M'*M für alle Positionen in L berechnen (M' ist die transponierte M) wie die folgende Python Code tut:

for M in L: 
    res += np.dot(M.T, M) 

Eigentlich möchte ich dies in Theano implementieren (die doesn 't unterstützen mehrdimensionale Arrays mit variabler Länge, und ich möchte nicht alle Matrizen auf die gleiche Größe auffüllen, da dies zu viel Speicherplatz verschwenden wird (einige der Matrizen können sehr groß sein).

Gibt es einen besseren Weg, dies zu tun?

bearbeiten:

L vor der Kompilation Theano bekannt.

bearbeiten:

erhielt zwei ausgezeichnete Antworten von @DanielRenshaw und @Divakar, emotional schwierig zu wählen zu akzeptieren.

Answer 1

Sie können die Eingabearrays nur entlang der ersten Achse auffüllen, die alle x summieren. Somit würden wir mit einem großen (X,n) Array enden, wobei X =x1+x2+x3+..... Dies kann transponiert werden und sein Skalarprodukt mit seinem Selbst wäre die gewünschte Ausgabe der Form (n,n). All dies wird mit einer reinen vektorisierten Lösung erreicht, die das leistungsfähige Punktprodukt nutzt.Somit würde die Umsetzung sein -

# Concatenate along axis=0 
Lcat = np.concatenate(L,axis=0) 

# Perform dot product of the transposed version with self 
out = Lcat.T.dot(Lcat) 

Runtime-Tests und überprüfen Ausgabe -

In [116]: def vectoized_approach(L): 
    ...: Lcat = np.concatenate(L,axis=0) 
    ...: return Lcat.T.dot(Lcat) 
    ...: 
    ...: def original_app(L): 
    ...: n = L[0].shape[1] 
    ...: res = np.zeros((n,n)) 
    ...: for M in L: 
    ...:  res += np.dot(M.T, M) 
    ...: return res 
    ...: 

In [117]: # Input 
    ...: L = [np.random.rand(np.random.randint(1,9),5)for iter in range(1000)] 

In [118]: np.allclose(vectoized_approach(L),original_app(L)) 
Out[118]: True 

In [119]: %timeit original_app(L) 
100 loops, best of 3: 3.84 ms per loop 

In [120]: %timeit vectoized_approach(L) 
1000 loops, best of 3: 632 µs per loop 

Answer 2

Vorausgesetzt, dass die Anzahl der Matrizen bekannt ist, bevor die Theano-Kompilierung stattfinden muss, kann man einfach reguläre Python-Listen von Theano-Matrizen verwenden.

Hier ist ein vollständiges Beispiel, das den Unterschied zwischen numpy und Theano-Versionen zeigt.

Dieser Code wurde aktualisiert, um einen Vergleich mit dem vektorisierten Ansatz von @ Divakar zu ermöglichen, der bessere Ergebnisse liefert. Zwei vektorisierte Ansätze sind für Theano möglich, eine, bei der Theano die Verkettung vornimmt, und eine, bei der numpy die Verkettung durchführt, deren Ergebnis dann an Theano übergeben wird.

import timeit 
import numpy as np 
import theano 
import theano.tensor as tt 


def compile_theano_version1(number_of_matrices, n, dtype): 
    assert number_of_matrices > 0 
    assert n > 0 
    L = [tt.matrix() for _ in xrange(number_of_matrices)] 
    res = tt.zeros(n, dtype=dtype) 
    for M in L: 
     res += tt.dot(M.T, M) 
    return theano.function(L, res) 


def compile_theano_version2(number_of_matrices): 
    assert number_of_matrices > 0 
    L = [tt.matrix() for _ in xrange(number_of_matrices)] 
    concatenated_L = tt.concatenate(L, axis=0) 
    res = tt.dot(concatenated_L.T, concatenated_L) 
    return theano.function(L, res) 


def compile_theano_version3(): 
    concatenated_L = tt.matrix() 
    res = tt.dot(concatenated_L.T, concatenated_L) 
    return theano.function([concatenated_L], res) 


def numpy_version1(*L): 
    assert len(L) > 0 
    n = L[0].shape[1] 
    res = np.zeros((n, n), dtype=L[0].dtype) 
    for M in L: 
     res += np.dot(M.T, M) 
    return res 


def numpy_version2(*L): 
    concatenated_L = np.concatenate(L, axis=0) 
    return np.dot(concatenated_L.T, concatenated_L) 


def main(): 
    iteration_count = 100 
    number_of_matrices = 20 
    n = 300 
    min_x = 400 
    dtype = 'float64' 
    theano_version1 = compile_theano_version1(number_of_matrices, n, dtype) 
    theano_version2 = compile_theano_version2(number_of_matrices) 
    theano_version3 = compile_theano_version3() 
    L = [np.random.standard_normal(size=(x, n)).astype(dtype) 
     for x in range(min_x, number_of_matrices + min_x)] 

    start = timeit.default_timer() 
    numpy_res1 = np.sum(numpy_version1(*L) 
         for _ in xrange(iteration_count)) 
    print 'numpy_version1', timeit.default_timer() - start 

    start = timeit.default_timer() 
    numpy_res2 = np.sum(numpy_version2(*L) 
         for _ in xrange(iteration_count)) 
    print 'numpy_version2', timeit.default_timer() - start 

    start = timeit.default_timer() 
    theano_res1 = np.sum(theano_version1(*L) 
         for _ in xrange(iteration_count)) 
    print 'theano_version1', timeit.default_timer() - start 

    start = timeit.default_timer() 
    theano_res2 = np.sum(theano_version2(*L) 
         for _ in xrange(iteration_count)) 
    print 'theano_version2', timeit.default_timer() - start 

    start = timeit.default_timer() 
    theano_res3 = np.sum(theano_version3(np.concatenate(L, axis=0)) 
         for _ in xrange(iteration_count)) 
    print 'theano_version3', timeit.default_timer() - start 

    assert np.allclose(numpy_res1, numpy_res2) 
    assert np.allclose(numpy_res2, theano_res1) 
    assert np.allclose(theano_res1, theano_res2) 
    assert np.allclose(theano_res2, theano_res3) 


main() 

Wenn diese Drucke laufen (so etwas wie)

numpy_version1 1.47830819649 
numpy_version2 1.77405482179 
theano_version1 1.3603150303 
theano_version2 1.81665318145 
theano_version3 1.86912039489 

Das behauptet Pass zeigt, dass die Theano und numpy Versionen beide das gleiche Ergebnis zu hohen Grad an Genauigkeit berechnen. Diese Genauigkeit wird deutlich reduziert, wenn float32 anstelle von float64 verwendet wird.

Die Timing-Ergebnisse zeigen, dass der vektorisierte Ansatz möglicherweise nicht vorzuziehen ist, sondern von den Matrixgrößen abhängt. Im obigen Beispiel sind die Matrizen groß und der Nicht-Verkettungsansatz ist schneller, aber wenn die Parameter n und min_x in der main-Funktion geändert werden, um viel kleiner zu sein, dann ist der Verkettungsansatz schneller. Andere Ergebnisse können beim Ausführen auf einer GPU enthalten sein (nur Theano-Versionen).

Answer 3

Neben @ DanielRenshaw Antwort, wenn wir die Anzahl der Matrizen bis 1000 zu erhöhen, die compile_theano_version1 Funktion RuntimeError: maximum recursion depth exceeded ergeben und compile_theano_version2 scheint ewig dauern, um zu kompilieren.

Es hierfür ein Fix ist durch typed_list mit:

def compile_theano_version4(number_of_matrices, n): 
    import theano.typed_list 
    L = theano.typed_list.TypedListType(tt.TensorType(theano.config.floatX, broadcastable=(None, None)))() 
    res, _ = theano.scan(fn=lambda i: tt.dot(L[i].T, L[i]), 
         sequences=[theano.tensor.arange(number_of_matrices, dtype='int64')]) 
    return theano.function([L], res.sum(axis=0)) 

Außerdem habe ich den Datentyp aller relevanten Variablen zu float32 gesetzt und lief auf GPU @ DanielRenshaw Skript, es stellte sich heraus, dass @ Vorschlag des Divakar (theano_version3) ist in diesem Fall am effizientesten. Obwohl @DanielRenshaw sagte, ist die Verwendung einer riesigen Matrix nicht immer eine gute Übung.

Die folgenden sind Einstellungen und Ausgänge auf meinem Rechner.

iteration_count = 100 
number_of_matrices = 200 
n = 300 
min_x = 20 
dtype = 'float32' 
theano.config.floatX = dtype 


numpy_version1 5.30542397499 
numpy_version2 3.96656394005 
theano_version1 5.26742005348 
theano_version2 1.76983904839 
theano_version3 1.03577589989 
theano_version4 5.58366179466