Doppia somma efficiente di prodotti

Question

Considerare due ndarrays di lunghezza n, arr1 e arr2. Sto calcolando la seguente somma dei prodotti, e di farlo num_runs volte al benchmark:

import numpy as np 
import time 

num_runs = 1000 
n = 100 

arr1 = np.random.rand(n) 
arr2 = np.random.rand(n) 

start_comp = time.clock() 
for r in xrange(num_runs): 
    sum_prods = np.sum([arr1[i]*arr2[j] for i in xrange(n) 
         for j in xrange(i+1, n)]) 

print "total time for comprehension = ", time.clock() - start_comp 

start_loop = time.clock() 
for r in xrange(num_runs): 
    sum_prod = 0.0 
    for i in xrange(n): 
     for j in xrange(i+1, n): 
      sum_prod += arr1[i]*arr2[j] 

print "total time for loop = ", time.clock() - start_loop 

L'uscita è

total time for comprehension = 3.23097066953 
total time for comprehension = 3.9045544426 

in modo da utilizzare di lista appare più veloce.

C'è un'implementazione molto più efficiente, utilizzando le routine NumPy forse, per calcolare una tale somma di prodotti?

Answer 1

ridisporre l'operazione in un algoritmo O (n) runtime invece di O (n^2), e sfruttare NumPy per i prodotti e le somme:

# arr1_weights[i] is the sum of all terms arr1[i] gets multiplied by in the 
# original version 
arr1_weights = arr2[::-1].cumsum()[::-1] - arr2 

sum_prods = arr1.dot(arr1_weights) 

temporizzazione mostra che questo è circa 200 volte più veloce della comprensione della lista per n == 100.

In [21]: %%timeit 
    ....: np.sum([arr1[i] * arr2[j] for i in range(n) for j in range(i+1, n)]) 
    ....: 
100 loops, best of 3: 5.13 ms per loop 

In [22]: %%timeit 
    ....: arr1_weights = arr2[::-1].cumsum()[::-1] - arr2 
    ....: sum_prods = arr1.dot(arr1_weights) 
    ....: 
10000 loops, best of 3: 22.8 µs per loop 

Answer 2

È possibile utilizzare il seguente trasmissione trucco:

a = np.sum(np.triu(arr1[:,None]*arr2[None,:],1)) 
b = np.sum([arr1[i]*arr2[j] for i in xrange(n) for j in xrange(i+1, n)]) 
print a == b # True 

Fondamentalmente, io sto pagando il prezzo di calcolare il prodotto di tutti gli elementi a due a due in arr1 e arr2 per sfruttare la velocità di trasmissione NumPy/vettorializzazione viene eseguito molto più velocemente nel codice di basso livello.

e tempi:

%timeit np.sum(np.triu(arr1[:,None]*arr2[None,:],1)) 
10000 loops, best of 3: 55.9 µs per loop 

%timeit np.sum([arr1[i]*arr2[j] for i in xrange(n) for j in xrange(i+1, n)]) 
1000 loops, best of 3: 1.45 ms per loop 

Answer 3

Un modo vettorizzati: np.sum(np.triu(np.multiply.outer(arr1,arr2),1)).

per un 30x miglioramento:

In [9]: %timeit np.sum(np.triu(np.multiply.outer(arr1,arr2),1)) 
1000 loops, best of 3: 272 µs per loop 

In [10]: %timeit np.sum([arr1[i]*arr2[j] for i in range(n) 
         for j in range(i+1, n)] 
100 loops, best of 3: 7.9 ms per loop 

In [11]: allclose(np.sum(np.triu(np.multiply.outer(arr1,arr2),1)), 
np.sum(np.triu(np.multiply.outer(arr1,arr2),1))) 
Out[11]: True 

Un altro approch veloce è quello di utilizzare numba:

from numba import jit 
@jit 
def t(arr1,arr2): 
    s=0 
    for i in range(n): 
     for j in range(i+1,n): 
      s+= arr1[i]*arr2[j] 
    return s 

per un 10x nuovo fattore:

In [12]: %timeit t(arr1,arr2) 
10000 loops, best of 3: 21.1 µs per loop 

E utilizzando @ user2357112 risposta minima ,

@jit 
def t2357112(arr1,arr2): 
    s=0 
    c=0 
    for i in range(n-2,-1,-1): 
     c += arr2[i+1] 
     s += arr1[i]*c 
    return s 

per

In [13]: %timeit t2357112(arr1,arr2) 
100000 loops, best of 3: 2.33 µs per loop 

, solo facendo le operazioni necessarie.