pętli nad (lub vectorize) Zmienna długość w matryce Theano

Question

mam listę macierzy L, w którym każdy element M jest x*n macierzy (x jest zmienna, n jest stała).

Chcę obliczyć sumę M'*M dla wszystkich elementów w L (M' jest transpozycją M) za pomocą następującego kodu Pythona robi:

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

Właściwie chcę zaimplementować to w Theano (co robi 't obsługują wielowymiarowe tablice o zmiennej długości) i nie chcę podkładać wszystkich macierzy do tego samego rozmiaru, ponieważ to zmarnuje zbyt dużo miejsca (niektóre macierze mogą być bardzo duże).

Czy jest lepszy sposób to zrobić?

Edit:

L jest znany przed kompilacją Theano.

Edit:

otrzymała dwie doskonałe odpowiedzi od @DanielRenshaw i @Divakar, emocjonalnie trudno wybrać jedną, aby zaakceptować.

Answer 1

Można po prostu podkładać tablice wejściowe wzdłuż pierwszej osi, która jest sumą x. W ten sposób uzyskalibyśmy wysoką tablicę (X,n), gdzie X =x1+x2+x3+..... To może być transponowane, a jego produkt z kropką sam w sobie będzie pożądaną formą kształtu (n,n). Wszystko to osiąga się dzięki czystemu wektoryzowanemu rozwiązaniu wykorzystującemu potężny produkt dotowy.Zatem realizacja byłaby -

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

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

badania środowiska i sprawdzić wyjście -

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

Biorąc pod uwagę, że liczba macierzy jest znana, zanim kompilacja Theano musi się odbyć, można po prostu użyć zwykłych list Pythona macierzy Theano.

Oto pełny przykład pokazujący różnicę między wersjami numpy i Theano.

Ten kod został zaktualizowany, aby uwzględnić porównanie z wektoryzowanym podejściem @ Divakar, które działa lepiej. Możliwe są dwa wektoryzowane podejścia do Theano, w którym Theano wykonuje konkatenację, i takie, w którym numpy wykonuje konkatenację, której wynik jest następnie przekazywany Theano.

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() 

Kiedy uruchomić ten nadrukami (coś podobnego)

numpy_version1 1.47830819649 
numpy_version2 1.77405482179 
theano_version1 1.3603150303 
theano_version2 1.81665318145 
theano_version3 1.86912039489 

twierdzi przepustkę, pokazując, że Theano i NumPy wersje zarówno obliczyć ten sam wynik dużą dokładnością. Oczywiście ta dokładność zmniejszy się, jeśli użyjemy float32 zamiast float64.

Wyniki pomiaru czasu pokazują, że podejście wektoryzacji może nie być preferowane, zależy to od rozmiarów matrycy. W powyższym przykładzie macierze są duże, a podejście bez konkatenacji jest szybsze, ale jeśli parametry i min_x są zmienione w funkcji main, aby były znacznie mniejsze, wówczas podejście łączenia jest szybsze. Inne wyniki mogą zachowywać się podczas pracy na GPU (tylko wersje Theano).

Answer 3

Oprócz użytkownika @ DanielRenshaw odpowiedź, jeśli wzrost liczby matryc do 1000, funkcja compile_theano_version1 przyniesie RuntimeError: maximum recursion depth exceeded , i compile_theano_version2 wydaje się trwać wiecznie, aby skompilować.

Jest to poprawka do tego za pomocą typed_list:

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)) 

Ponadto ustawić typ danych wszystkich istotnych zmiennych float32 i pobiegł @ skryptu DanielRenshaw dotyczącego GPU, okazało się, że @ sugestią Divakar za (theano_version3) jest najbardziej efektywny w tym przypadku. Chociaż, jak powiedział @DanielRenshaw, używanie ogromnej matrycy może nie zawsze być dobrą praktyką.

Poniżej znajdują się ustawienia i wyjścia na moim komputerze.

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