2013-06-16 18 views
8

Sto cercando di fare qualcosa di semplice in numpy, e sono sicuro che ci dovrebbe essere un modo semplice per farlo.Numpy: prodotto esterno di n vettori

Fondamentalmente, ho un elenco di n vettori con varie lunghezze. Se v1[i] è la i 'th entrata del primo vettore allora voglio trovare una gamma dimensionale n, A, in modo tale che

A[i,j,k...] = v1[i] v2[j] v3[k] ... 

mio problema è che:

  1. outer vogliono solo due vettore argomenti.

  2. einsum richiede un parametro come "abcd ..." che sembra non necessario.

  3. kron richiede quello che sembra un rimodellamento piuttosto complesso e richiede solo due argomenti.

Vorrei evitare il più possibile la complessità, in modo da evitare l'introduzione di bug. Quindi preferibilmente vorrei un singolo comando.

Finora, il migliore che ho un po 'in mente è:

vs = [v1, v2, v3 ...] 
shape = map(len, vs) 

# specify the orientation of each vector 
newshapes = diag(array(shape)-1)+1 
reshaped = [x.reshape(y) for x,y in zip(vs, newshapes)] 

# direct product 
A = reduce(lambda a,b: a*b, reshaped, 1) 
+0

Il numero di vettori è sconosciuto fino al runtime? – DarenW

+0

@DarenW sì, è corretto. – Lucas

+1

Mi piace questo 'riduci (lambda a, b: a [..., np.newaxis] * b, vs)' ma non sono sicuro che possa essere considerato "un singolo comando". O se ci sono modi più veloci. – jorgeca

risposta

8

si utilizza l'uso seguente codice di una riga:

reduce(np.multiply, np.ix_(*vs)) 

np.ix_() farà la trasmissione esterna, è necessario ridurre, ma puoi passare l'ufunc np.multiply senza la funzione lambda.

Ecco il confronto tra:

import numpy as np 
vs = [np.r_[1,2,3.0],np.r_[4,5.0],np.r_[6,7,8.0]] 
shape = map(len, vs) 

# specify the orientation of each vector 
newshapes = np.diag(np.array(shape)-1)+1 
reshaped = [x.reshape(y) for x,y in zip(vs, newshapes)] 

# direct product 
A = reduce(lambda a,b: a*b, reshaped, 1) 
B = reduce(np.multiply, np.ix_(*vs)) 

np.all(A==B) 

Il reuslt:

True 
+5

+1 Può essere un miglioramento trascurabile, a meno che le liste vettoriali non siano molto lunghe, ma penserei che numpy 'np.multiply.reduce (np-ix_ (vs))' probabilmente funzionerebbe meglio del costrutto Python. – Jaime

1

C'è una linea alternativa di codice:

reduce(np.multiply.outer, vs) 

E 'più trasparente per me che np.ix_(*vs) costruzione e supporta matrici multidimensionali come in this question.

Timings sono le stesse entro una tolleranza:

import numpy as np 
from functools import reduce 

def outer1(*vs): 
    return np.multiply.reduce(np.ix_(*vs)) 
def outer2(*vs): 
    return reduce(np.multiply.outer, vs) 

v1 = np.random.randn(100) 
v2 = np.random.randn(200) 
v3 = np.random.randn(300) 
v4 = np.random.randn(50) 

%timeit outer1(v1, v2, v3, v4) 
# 1 loop, best of 3: 796 ms per loop 

%timeit outer2(v1, v2, v3, v4) 
# 1 loop, best of 3: 795 ms per loop 

np.all(outer1(v1, v2, v3, v4) == outer2(v1, v2, v3, v4)) 
# True 
+0

Sì, preferirei usarlo. – Lucas

Problemi correlati