Có sản phẩm chấm có dấu chấm/scipy không, chỉ tính các mục chéo của kết quả?

Question

Hãy tưởng tượng có mảng 2 NumPy:

> A, A.shape = (n,p) 
> B, B.shape = (p,p) 

Điển hình p là một số nhỏ hơn (p < = 200), trong khi n có thể tùy tiện lớn.

Tôi đang làm như sau:

result = np.diag(A.dot(B).dot(A.T)) 

Như bạn thấy, tôi chỉ giữ lại các mục chéo n, tuy nhiên có một trung gian (n x n) mảng tính từ mà chỉ có các mục chéo được lưu giữ.

Tôi muốn có một hàm như diag_dot(), chỉ tính toán các mục nhập đường chéo của kết quả và không cấp phát bộ nhớ hoàn chỉnh.

Một kết quả sẽ là:

> result = diag_dot(A.dot(B), A.T) 

Có một chức năng premade như thế này và điều này có thể được thực hiện một cách hiệu quả mà không cần phân bổ (n x n) mảng trung gian?

Answer 1

Tôi nghĩ rằng tôi đã nhận nó trên của riêng tôi, nhưng vẫn sẽ chia sẻ các giải pháp:

kể từ khi nhận được chỉ là đường chéo của một phép nhân ma trận

> Z = N.diag(X.dot(Y)) 

tương đương với số tiền cá nhân của sản phẩm vô hướng hàng của X và cột của Y, báo cáo kết quả trước đó là tương đương với:

> Z = (X * Y.T).sum(-1) 

Đối với các biến độc đáo này có nghĩa là:

01.

> result = (A.dot(B) * A).sum(-1) 

Xin vui lòng sửa tôi nếu tôi sai nhưng điều này nên được nó ...

Answer 2

Bạn có thể nhận được hầu hết mọi thứ bạn đã bao giờ mơ ước với numpy.einsum. Cho đến khi bạn bắt đầu nhận được hang của nó, nó về cơ bản có vẻ như voodoo đen ...

>>> a = np.arange(15).reshape(5, 3) 
>>> b = np.arange(9).reshape(3, 3) 

>>> np.diag(np.dot(np.dot(a, b), a.T)) 
array([ 60, 672, 1932, 3840, 6396]) 
>>> np.einsum('ij,ji->i', np.dot(a, b), a.T) 
array([ 60, 672, 1932, 3840, 6396]) 
>>> np.einsum('ij,ij->i', np.dot(a, b), a) 
array([ 60, 672, 1932, 3840, 6396]) 

EDIT Bạn thực sự có thể nhận được toàn bộ sự việc trong một shot duy nhất, thật nực cười ...

>>> np.einsum('ij,jk,ki->i', a, b, a.T) 
array([ 60, 672, 1932, 3840, 6396]) 
>>> np.einsum('ij,jk,ik->i', a, b, a) 
array([ 60, 672, 1932, 3840, 6396]) 

CHỈNH SỬA Bạn không muốn để nó quá nhiều vào chính nó ... Thêm câu trả lời của OP vào câu hỏi của chính nó để so sánh.

n, p = 10000, 200 
a = np.random.rand(n, p) 
b = np.random.rand(p, p) 

In [2]: %timeit np.einsum('ij,jk,ki->i', a, b, a.T) 
1 loops, best of 3: 1.3 s per loop 

In [3]: %timeit np.einsum('ij,ij->i', np.dot(a, b), a) 
10 loops, best of 3: 105 ms per loop 

In [4]: %timeit np.diag(np.dot(np.dot(a, b), a.T)) 
1 loops, best of 3: 5.73 s per loop 

In [5]: %timeit (a.dot(b) * a).sum(-1) 
10 loops, best of 3: 115 ms per loop 

Answer 3

Một câu trả lời cho người đi bộ, mà tránh việc xây dựng các mảng trung gian lớn là:

result=np.empty([n.], dtype=A.dtype) 
for i in xrange(n): 
    result[i]=A[i,:].dot(B).dot(A[i,:])