2016-09-09 10 views
3

나는 ndarray에있는 transpose이 matlab의 permute 함수와 동등한 것으로 알고 있지만 간단히 작동하지 않는 특정 유스 케이스가 있음을 이해합니다.순열과 함께 파이썬에서 방송

A는 형상 NxNxM 및 B의 3 차원 텐서
C = @bsxfun(@times, permute(A,[4,2,5,1,3]), permute(B, [1,6,2,7,3,4,5]) 

모양 NxNxMxPxP의 5D 텐서이다 매트랩 I는 다음이있다. 위의 함수는 루프 된 크로네 커 생성물을 벡터화하기위한 것입니다. 나는 Matlab이 A와 B 둘 다에 대해 2 개의 싱글 튼 디멘션을 추가하고 있다고 가정하고 있는데, 이는 왜 그것들을 재정렬 할 수 있는지에 대한 것이다. 이 코드를 파이썬 에 이식하는 방법을 찾고 있습니다.하지만이 추가 차원을 추가 할 수있는 능력이 없다고 생각합니다. . 나는 this 성공적으로 추가 치수를 추가하지만 방송은 동일한 matlab의 bsxfun을 작동하지 않는 것을 발견했습니다.

A = A[...,None,None] 
B = B[...,None,None] 
C = transpose(A,[3,1,4,0,2])*transpose(B,[0,5,1,6,2,3,4]) 

나는 다음과 같은 오류 얻을 : 나는 명백한 번역을 시도 (예 나는이 ndarray 's 및 기능 NumPy와 사용하고 있습니다)

return transpose(axes) 
ValueError: axes don't match array 

내 첫번째 추측은을하는 것입니다 A와 B에 reshape을 추가하면 싱글 톤 차원에 추가됩니까?

지금 다음 오류 얻을 :

mults = transpose(rho_in,[3,1,4,0,2])*transpose(proj,[0,5,1,6,2,3,4]) 
ValueError: operands could not be broadcast together with shapes (1,9,1,9,8) (9,1,9,1,8,40,40) 

편집 : 올바르게 파이썬이 matlab에 곱셈을 방송에 대한 단일 차원을 추가하는 방법에 대한 덜하지만 더로 내 질문에 개정.

답변

6

MATLAB과 numpy의 큰 차이점은 전자는 배열에 대해 열 - 주 형식을 사용하고 후자는 행 주요을 사용한다는 것입니다. 결론적으로 암시 적 싱글 톤 차원은 다르게 처리됩니다.

특히 MATLAB은 후행 단일 차원을 명시 적으로 무시합니다. rand(3,3,1,1,1,1,1)은 실제로는 3x3 행렬입니다. 이 줄을 따라 선도 이 일치하는 경우 을 사용하여 NxNxMxPxP과 호환되는 이 암시 적으로 일치하면 두 배열에서 작동 할 수 있습니다.

Numpy, 반면에 allows implicit singletons up front. permute 배열은 모양이 (1,9,1,9,8) 인 모양이 인 크기가 일치하는 모양 (예 : (40,40,9,1,9,1,8))이고 결과는 모양이 (40,40,9,9,9,9,8)이어야합니다.

더미 예 : 당신이 뭘하려는 건지 아마 numpy.einsum을 사용하여 수행 할 수

>>> import numpy as np 
>>> (np.random.rand(40,40,9,1,9,1,8)+np.random.rand(1,9,1,9,8)).shape 
(40, 40, 9, 9, 9, 9, 8) 

참고. 나는 그것을 자세히 살펴볼 것을 제안한다. 무슨 뜻인지의 예 : 나는 당신이 원하는 것을 수집 질문에서이 작업을 수행 할 수 있습니다 요소를 가지고 A[1:N,1:N,1:M]B[1:N,1:N,1:M,1:P,1:P]

C[i1,i2,i3,i4,i5,i6,i7] = A[i2,i4,i5]*B[i1,i3,i5,i6,i7] 

이 (특정 인덱스 순서가 다를 수 있습니다) C[1:N,1:N,1:N,1:N,1:M,1:P,1:P] 있도록 새로운 배열을 구성.이 맞다면, 당신은 참으로 numpy.einsum()를 사용할 수 있습니다

>>> a = np.random.rand(3,3,2) 
>>> b = np.random.rand(3,3,2,4,4) 
>>> np.einsum('ijk,lmkno->limjkno',a,b).shape 
(3, 3, 3, 3, 2, 4, 4) 
두 가지가 있지만, 주목해야한다

. 첫째, 위 작업은 매우 메모리 집약적 일 것이므로 벡터화 (일반적으로 메모리 요구를 희생시키면서 CPU 시간을 얻는 경우)에서 예상해야합니다. 두 번째로, 코드를 이식 할 때 데이터 모델을 재 배열하는 것을 진지하게 고려해야합니다. 방송이 두 언어로 다르게 작동하는 이유는 열 - 주요/행 - 주요 차이와 복잡하게 연관되어 있습니다. 이것은 또한 A(:,i2,i3) 메모리의 연속 블록에 해당하기 때문에 A(i1,i2,:)하지 않는 동안 MATLAB에서 당신이 먼저 주요 인덱스 작업을해야한다는 것을 의미한다. 반대로 numpy의 경우 A[i1,i2,:]은 인접하지만 A[:,i2,i3]은 그렇지 않습니다.

이러한 고려 사항은 벡터화 작업 바람직 MATLAB에서 인덱스를 선도 NumPy와의 인덱스를 후행와 함께 작동하도록 데이터의 물류를 설정해야하는 것이 좋습니다. 당신은 여전히하지만 당신의 크기는 우리가 코드의 두 버전이 최적의 설정을 사용하는 것으로 가정 적어도 경우, MATLAB에 비해 다른 (아마도 역방향) 순서로해야한다, 작업 자체를 수행 할 numpy.einsum를 사용할 수 있습니다.

+0

이 주제 (특히 'einsum'의 사용)를 확장하기 위해 어떤 종류의 직접적인 메시지로 연락 할 수 있습니까? 나는 나의 미래의 질문이 여기에 관련이 있기에 너무 구체적 일 것 같다. –

+0

@Mike 불행히도 Stack Overflow는 오프 사이트 노력을 분명히 낙담합니다 :) 주제에 대한 우려를 볼 수 있지만 여기에 의견을 남기거나 (매우 밀접하게 관련되어 있고 간단하다면) 새로운 질문을하십시오 (보이는 경우 충분히 큰 질문을 할만 큼) 또는 [MATLAB 회의실] (http://chat.stackoverflow.com/rooms/81987/chatlab-and-talktave) 또는 [python one] (http : //chat.stackoverflow.com/rooms/6/python). –

+0

아마도 permute (A ...)의'(2,1,3)'때문에'np.einsum ('ijk, lmkno-> ljmikno', A, B)이어야할까요? – Divakar

2

당신의 MATLAB 코드를 보면, 당신은 -

C = bsxfun(@times, permute(A,[4,2,5,1,3]), permute(B, [1,6,2,7,3,4,5]) 

그래서 본질적으로 -

B : 1 , 6 , 2 , 7 , 3 , 4 , 5 
A : 4 , 2 , 5 , 1 , 3 

을 이제 우리는 더 높은 사람에서 단일 차원을 빌려했다 MATLAB으로는, 그 이유는 전부 그 에서 가져 오기의 문제는 B에 대한 6, 7을 흐리게하고 A에 대한 45을 어두워집니다. NumPy와에서

, 우리는 np.newaxis/None에 명시 적으로 사람에 가져. 따라서, NumPy와, 우리는과 같이 넣을 수 - N 새로운 축을 나타냅니다

B : 1 , N , 2 , N , 3 , 4 , 5 
A : N , 2 , N , 1 , 3 , N , N 

을. 우리는 A이 정렬을 위해 다른 차원을 앞으로 밀어 넣으려면 끝에 새로운 축을 넣어야한다는 것을 알아 두십시오. 대충, 이것은 MATLAB에서 기본적으로 발생합니다.

크기를 적절하게 만드는 것처럼 보이기 때문에 B이 쉽게 보이며 적절한 위치에 새 축을 추가하면됩니다 (B[:,None,:,None,:,:,:]).

그런 A을 만드는 것은 간단하지 않습니다. AN's을 무시하고, 우리는 것 - A : 2 , 1 , 3합니다. 따라서 시작점은 치수를 바꾸는 것이고 무시 된 두 개의 새 축을 추가합니다 (A.transpose(1,0,2)[None,:,None,:,:,None,None]).

는 지금까지, 우리는이 -

B (new): B[:,None,:,None,:,:,:] 
A (new): A.transpose(1,0,2)[None,:,None,:,:,None,None] 

NumPy와, 우리가 선도하는 새로운 축과 뒤 비 - 싱글 희미을 건너 뛸 수 있습니다.그래서 우리는과 같이 단순화 할 수 -

C = A.transpose(1,0,2)[:,None,...,None,None]*B[:,None,:,None] 

런타임 테스트 안드라스의 게시물에 해당하는 의미 @ 내가 믿는

-

B (new): B[:,None,:,None] 
A (new): A.transpose(1,0,2)[:,None,...,None,None] 

최종 출력이 두 확장 버전 사이의 곱셈 것 np.einsum 구현은 다음과 같습니다. np.einsum('ijk,lmkno->ljmikno',A,B).

In [24]: A = np.random.randint(0,9,(10,10,10)) 
    ...: B = np.random.randint(0,9,(10,10,10,10,10)) 
    ...: 

In [25]: C1 = np.einsum('ijk,lmkno->ljmikno',A,B) 

In [26]: C2 = A.transpose(1,0,2)[:,None,...,None,None]*B[:,None,:,None] 

In [27]: np.allclose(C1,C2) 
Out[27]: True 

In [28]: %timeit np.einsum('ijk,lmkno->ljmikno',A,B) 
10 loops, best of 3: 102 ms per loop 

In [29]: %timeit A.transpose(1,0,2)[:,None,...,None,None]*B[:,None,:,None] 
10 loops, best of 3: 78.4 ms per loop 

In [30]: A = np.random.randint(0,9,(15,15,15)) 
    ...: B = np.random.randint(0,9,(15,15,15,15,15)) 
    ...: 

In [31]: %timeit np.einsum('ijk,lmkno->ljmikno',A,B) 
1 loop, best of 3: 1.76 s per loop 

In [32]: %timeit A.transpose(1,0,2)[:,None,...,None,None]*B[:,None,:,None] 
1 loop, best of 3: 1.36 s per loop 
+0

나는 그것이 실제로'einsum '과 무슨 의미인지를 확인할 수 있습니다. :) –