如何用numpy表達這個張量計算？

Question

假設x是一個形狀爲(a,)的向量，T是一個形狀爲(b, a, a)的張量。 如果我想計算(x^T)Tx，我可以使用x.dot(w.dot(x).transpose())來完成。

例如：

x = np.array([1.,2.,3.,4.,5.]) 

w = np.array([[[1.,2.,3.,4.,5.], 
       [1.,2.,3.,4.,5.], 
       [1.,2.,3.,4.,5.]], 
       [[1.,2.,3.,4.,5.], 
       [1.,2.,3.,4.,5.], 
       [1.,2.,3.,4.,5.]]]) 

x.dot(w.dot(x).transpose()) 

但是，如果我想分解T成兩個張量P和Q（低階快遞）與形狀(b,a,r)和(b,r,a)和r<<a所以每個矩陣T這是a*a分解到什麼a*r和r*a，這些數據減少了很多。那我該如何計算(x^T)PQx與numpy？

Answer 1

你的例子有問題。

x.shape (5,) 
w.shape (2,3,5) 
x.dot(w.dot(x).transpose()) 
    ValueError: matrices are not aligned 

但使用你的描述：

`x` `(a,)`, `T` `(b,a,a)`; `(x^T)Tx` 

我喜歡思考複雜的產品時使用einsum（愛因斯坦總和）。我覺得你的x'Tx是：

np.einsum('i,kij,j->k', x, T, x) 

牛逼分解爲：P(b,a,r)，Q(b,r,a);

np.einsum('kir,krj->kij', P,Q) == T 

在一起表達式是：

np.einsum('i,kir,krj,j->k', x, P, Q, x) 

einsum不是最佳時的尺寸是大的，由於k,i,j,r組合迭代空間可以是大的。這仍然是考慮問題的有用方法。

我認爲它可以改寫爲3 dots：

P1 = np.einsum('i,kir->kr', x, P) 
Q1 = np.einsum('krj,j->kr', Q, x) 
np.einsum('kr,kr->k', P1, Q1) 

---

樣本計算：

In [629]: a,b,r = 5,3,2 
In [630]: x=np.arange(1.,a+1) 
In [632]: P=np.arange(b*a*r).reshape(b,a,r) 
In [633]: Q=np.arange(b*a*r).reshape(b,r,a) 

In [635]: T=np.einsum('kir,krj->kij',P,Q) 

In [636]: P 
Out[636]: 
array([[[ 0, 1], 
     [ 2, 3], 
     [ 4, 5], 
     [ 6, 7], 
     ... 
     [24, 25], 
     [26, 27], 
     [28, 29]]]) 

In [637]: Q 
Out[637]: 
array([[[ 0, 1, 2, 3, 4], 
     [ 5, 6, 7, 8, 9]], 
     ... 
     [[20, 21, 22, 23, 24], 
     [25, 26, 27, 28, 29]]]) 

In [638]: T 
Out[638]: 
array([[[ 5, 6, 7, 8, 9], 
     [ 15, 20, 25, 30, 35], 
     [ 25, 34, 43, 52, 61], 
     [ 35, 48, 61, 74, 87], 
     ... 
     [1105, 1154, 1203, 1252, 1301], 
     [1195, 1248, 1301, 1354, 1407], 
     [1285, 1342, 1399, 1456, 1513]]]) 

In [639]: T.shape 
Out[639]: (3, 5, 5) 

In [640]: R1=np.einsum('i,kij,j->k',x,T,x) 
... 
In [642]: R1 
Out[642]: array([ 14125., 108625., 293125.]) 

In [643]: R2=np.einsum('i,kir,krj,j->k',x,P,Q,x) 
In [644]: R2 
Out[644]: array([ 14125., 108625., 293125.]) 

In [645]: P1=np.einsum('i,kir->kr',x,P) 
In [646]: Q1=np.einsum('krj,j->kr',Q,x) 
In [647]: R3=np.einsum('kr,kr->k',P1,Q1) 
In [648]: R3 
Out[648]: array([ 14125., 108625., 293125.]) 

In [649]: P1 
Out[649]: 
array([[ 80., 95.], 
     [ 230., 245.], 
     [ 380., 395.]]) 

In [650]: Q1 
Out[650]: 
array([[ 40., 115.], 
     [ 190., 265.], 
     [ 340., 415.]]) 

計算的最後一組可與dot

In [656]: np.dot(x,P) 
Out[656]: 
array([[ 80., 95.], 
     [ 230., 245.], 
     [ 380., 395.]]) 

In [657]: np.dot(Q,x) 
Out[657]: 
array([[ 40., 115.], 
     [ 190., 265.], 
     [ 340., 415.]]) 

In [658]: np.dot(np.dot(x,P),np.dot(Q,x).T) 
Out[658]: 
array([[ 14125., 40375., 66625.], 
     [ 37375., 108625., 179875.], 
     [ 60625., 176875., 293125.]]) 

完成

但我們只想要最後的dot的對角線。更簡單的產品總和是更好的：

In [661]: (P1*Q1).sum(axis=1) 
Out[661]: array([ 14125., 108625., 293125.])