2016-07-28 3 views
0

모양 (8,9)이있는 scipy 스파 스 행렬과 모양 (9, 12, 17)이있는 다른 배열이 있습니다. 나는 (8,9) 행렬이 효과적으로 첫 번째 차원 만 곱한 크기/행렬 (8,12,17)을 얻도록 이러한 값을 곱합니다. 이 작업을 수행하기 위해 크로네 커 제품을 사용해야합니까, 아니면 numpy로 간단한 방법이 있습니까?numpy, scipy에서 다차원 행렬 곱셈을 수행하는 방법

+0

출력이 너무 적습니까? – Divakar

+0

아니요, 출력은 아마 같을 것입니다. – user1827975

+0

그래서'.todense()'로 밀집 포맷으로 변환하고'np.dot/np.tensordot'를 사용합니까? – Divakar

답변

0

다음은 작동시키기위한 몇 가지 방법입니다. 두 번째 테스트가 더 좋고 12 배 빠릅니다.

def multiply_3D_dim_zero_slow(matrix, array): 
shape = array.shape 
final_shape = (matrix.shape[0], array.shape[1], array.shape[2]) 
result = np.zeros(final_shape) 
for i in xrange(shape[1]): 
    for j in xrange(shape[2]): 
     result[:, i, j] = matrix * array[:, i, j] 
return result.reshape(final_shape) 

그리고 더 빠른 버전으로 다차원 배열을 2D 배열로 만듭니다.

def multiply_3D_dim_zero(matrix, array): 
shape = array.shape 
final_shape = (matrix.shape[0], array.shape[1], array.shape[2]) 
array_reshaped = array.reshape(shape[0], shape[1] * shape[2]) 
return (matrix * array_reshaped).reshape(final_shape)ode here 

이것은 첫 번째 차원에서 효과가 있지만 필요한 것은 일반화 할 수 있습니다. hpaulj이 코멘트에 알 수 있듯이

+0

'tensordot'는 이런 종류의 재 형성을하므로'점 '을 호출 할 수 있습니다. 그것은 올바른 순서로 그들을 얻기 위해 swapaxes 수도 있습니다. – hpaulj

1

는이 작업을 수행하는 가장 쉬운 방법은 조밀 한 매트릭스 np.einsum입니다 :

>>> a = np.random.randn(8, 9) 
>>> b = np.random.randn(9, 12, 17) 
>>> c = np.einsum('ij,jkl->ikl', a, b) 
>>> c.shape 
(8, 12, 17) 
2

m1는 2D 희소 행렬 인 경우, m1.A는 조밀 한 배열의 형태입니다. dinmensions는 실질적으로 einsum 표현을 씁니다. 예를 들어

np.einsum('ij,jkl->ikl', m1.A, m2) 

:

In [506]: M = sparse.random(8, 9, 0.1) 
In [507]: A = np.ones((9, 12, 17)) 
In [508]: np.einsum('ij,jkl->ikl', M.A, A).shape 
Out[508]: (8, 12, 17) 
1

@Divakar는 np.tensordot을 권장하고 @hpaulj 및 @Praveen는 np.einsum을 제안했다. 또 다른 방법은 전치되는 축 : 당신이 인용 작은 치수

(a @ b.transpose((2, 0, 1))).transpose((1, 2, 0)) 

np.einsum 및 전위가 빠른 것 같다. 하지만 일단 곱하는 축의 크기를 확대하기 시작하면 np.tensordot이 다른 두 개의 수를 상회합니다.

import numpy as np 

m, n, k, l = 8, 9, 12, 17 
a = np.random.random((m, n)) 
b = np.random.random((n, k, l)) 

%timeit np.tensordot(a, b, axes=([1], [0])) 
# => 10000 loops, best of 3: 22 µs per loop 
%timeit np.einsum("ij,jkl->ikl", a, b) 
# => 100000 loops, best of 3: 10.1 µs per loop 
%timeit (a @ b.transpose((2, 0, 1))).transpose((1, 2, 0)) 
# => 100000 loops, best of 3: 11.1 µs per loop 

m, n, k, l = 8, 900, 12, 17 
a = np.random.random((m, n)) 
b = np.random.random((n, k, l)) 

%timeit np.tensordot(a, b, axes=([1], [0])) 
# => 1000 loops, best of 3: 198 µs per loop 
%timeit np.einsum("ij,jkl->ikl", a, b) 
# => 1000 loops, best of 3: 868 µs per loop 
%timeit (a @ b.transpose((2, 0, 1))).transpose((1, 2, 0)) 
# => 1000 loops, best of 3: 907 µs per loop 

m, n, k, l = 8, 90000, 12, 17 
a = np.random.random((m, n)) 
b = np.random.random((n, k, l)) 

%timeit np.tensordot(a, b, axes=([1], [0])) 
# => 10 loops, best of 3: 21.7 ms per loop 
%timeit np.einsum("ij,jkl->ikl", a, b) 
# => 10 loops, best of 3: 164 ms per loop 
%timeit (a @ b.transpose((2, 0, 1))).transpose((1, 2, 0)) 
# => 10 loops, best of 3: 166 ms per loop 
+2

'@'는'np.matmul'입니다 (새로운 Python이 없다면). docs에 대한'help (np.matmul)'. – hpaulj

관련 문제