3D 매트릭스에 2D 행렬에 곱하십시오
-
20-09-2019 - |
문제
내가 있다고 가정 해 AXBXC 행렬 X
그리고 a BXD 행렬 Y
.
각각을 곱할 수있는 비 루프 방법이 있습니까? 씨 AXB 매트릭스 Y
?
해결책
함수를 사용하여 한 줄로이 작업을 수행 할 수 있습니다. Num2Cell 매트릭스를 끊습니다 X
셀 어레이로 셀프 세포를 가로 질러 작동하려면 :
Z = cellfun(@(x) x*Y,num2cell(X,[1 2]),'UniformOutput',false);
결과 Z
a 1 x-c 각 셀에 an이 포함 된 셀 어레이 a-by-d 행렬. 네가 원한다면 Z
an a-by-d-by-c 매트릭스, 당신은 사용할 수 있습니다 고양이 기능:
Z = cat(3,Z{:});
노트: 내 오래된 솔루션이 사용되었습니다 MAT2CELL 대신에 Num2Cell, 간결하지 않은 것 :
[A,B,C] = size(X);
Z = cellfun(@(x) x*Y,mat2cell(X,A,B,ones(1,C)),'UniformOutput',false);
다른 팁
개인적인 취향으로, 나는 내 코드가 가능한 한 간결하고 읽을 수있는 것을 좋아합니다.
'No-loops'요구 사항을 충족하지는 않지만 다음은 다음과 같습니다.
for m = 1:C
Z(:,:,m) = X(:,:,m)*Y;
end
결과가 발생합니다 A x d x c 행렬 지.
물론, 당신은 항상 Z를 사전 할당하여 사용하여 속도를 높일 수 있습니다. Z = zeros(A,D,C);
.
다음은 한 줄 솔루션입니다 (3 차원으로 분할하려면 2 개) :
A = 2;
B = 3;
C = 4;
D = 5;
X = rand(A,B,C);
Y = rand(B,D);
%# calculate result in one big matrix
Z = reshape(reshape(permute(X, [2 1 3]), [A B*C]), [B A*C])' * Y;
%'# split into third dimension
Z = permute(reshape(Z',[D A C]),[2 1 3]);
따라서 지금 : Z(:,:,i)
결과를 포함합니다 X(:,:,i) * Y
설명:
위의 내용은 혼란스러워 보일 수 있지만 아이디어는 간단합니다. 먼저 나는 3 차원의 3 차원으로 시작합니다 X
첫 번째 희미한 점을 따라 수직 연합을 수행하십시오.
XX = cat(1, X(:,:,1), X(:,:,2), ..., X(:,:,C))
... 어려움은 그랬습니다 C
변수이므로 사용하여 해당 표현식을 일반화 할 수 없습니다. 고양이 또는 Vertcat. 다음으로 우리는 이것을 곱합니다 Y
:
ZZ = XX * Y;
마침내 나는 그것을 3 차원으로 다시 나누었다 :
Z(:,:,1) = ZZ(1:2, :);
Z(:,:,2) = ZZ(3:4, :);
Z(:,:,3) = ZZ(5:6, :);
Z(:,:,4) = ZZ(7:8, :);
그래서 당신은 그것을 하나의 행렬 곱셈이 필요하다는 것을 알 수 있지만 재구성 전후에 매트릭스.
나는 가장 효율적인 방법을 위해 똑같은 문제에 접근하고 있습니다. 외부 라이브러리를 사용하지 않는 약 3 가지 접근 방식이 있습니다 (즉, mtimesx):
- 3D 매트릭스의 슬라이스를 통해 루프
- repmat-and-permute 마법사
- Cellfun 곱셈
최근에 세 가지 방법을 모두 비교하여 가장 빠른 방법을 확인했습니다. 내 직관은 (2) 승자가 될 것입니다. 코드는 다음과 같습니다.
% generate data
A = 20;
B = 30;
C = 40;
D = 50;
X = rand(A,B,C);
Y = rand(B,D);
% ------ Approach 1: Loop (via @Zaid)
tic
Z1 = zeros(A,D,C);
for m = 1:C
Z1(:,:,m) = X(:,:,m)*Y;
end
toc
% ------ Approach 2: Reshape+Permute (via @Amro)
tic
Z2 = reshape(reshape(permute(X, [2 1 3]), [A B*C]), [B A*C])' * Y;
Z2 = permute(reshape(Z2',[D A C]),[2 1 3]);
toc
% ------ Approach 3: cellfun (via @gnovice)
tic
Z3 = cellfun(@(x) x*Y,num2cell(X,[1 2]),'UniformOutput',false);
Z3 = cat(3,Z3{:});
toc
세 가지 접근법 모두 동일한 출력 (Phew!)을 생성했지만 놀랍게도 루프는 가장 빠릅니다.
Elapsed time is 0.000418 seconds.
Elapsed time is 0.000887 seconds.
Elapsed time is 0.001841 seconds.
시간은 한 번의 시험마다 다를 수 있으며 때로는 (2) 가장 느리게 나옵니다. 이러한 차이는 더 큰 데이터에서 더욱 극적입니다. 하지만 함께 많이 더 큰 데이터, (3) 비트 (2). 루프 방법이 여전히 가장 좋습니다.
% pretty big data...
A = 200;
B = 300;
C = 400;
D = 500;
Elapsed time is 0.373831 seconds.
Elapsed time is 0.638041 seconds.
Elapsed time is 0.724581 seconds.
% even bigger....
A = 200;
B = 200;
C = 400;
D = 5000;
Elapsed time is 4.314076 seconds.
Elapsed time is 11.553289 seconds.
Elapsed time is 5.233725 seconds.
그러나 루프 방법 ~할 수 있다 루프 차원이 다른 치수보다 훨씬 큰 경우 (2)보다 느리게하십시오.
A = 2;
B = 3;
C = 400000;
D = 5;
Elapsed time is 0.780933 seconds.
Elapsed time is 0.073189 seconds.
Elapsed time is 2.590697 seconds.
따라서 (2)이 (아마도 극단적 인) 사건에서 큰 요인으로 승리합니다. 모든 경우에 최적의 접근 방식은 없을 수도 있지만 루프는 여전히 꽤 좋으며 많은 경우에 가장 좋습니다. 가독성 측면에서도 가장 좋습니다. 루프를 멀리!
아니요. 여러 가지 방법이 있지만 항상 직접적이거나 간접적으로 루프로 나옵니다.
내 호기심을 기쁘게하기 위해, 왜 어쨌든 그것을 원하십니까?
질문에 답하기 위해 그리고 가독성은 다음을 참조하십시오.
- ndmult, Ajuanpi (Juan Pablo Carbajal), 2013, Gnu Gpl
입력
- 2 배열
- 어둑한
예시
nT = 100;
t = 2*pi*linspace (0,1,nT)’;
# 2 experiments measuring 3 signals at nT timestamps
signals = zeros(nT,3,2);
signals(:,:,1) = [sin(2*t) cos(2*t) sin(4*t).^2];
signals(:,:,2) = [sin(2*t+pi/4) cos(2*t+pi/4) sin(4*t+pi/6).^2];
sT(:,:,1) = signals(:,:,1)’;
sT(:,:,2) = signals(:,:,2)’;
G = ndmult (signals,sT,[1 2]);
원천
원래 소스. 인라인 댓글을 추가했습니다.
function M = ndmult (A,B,dim)
dA = dim(1);
dB = dim(2);
# reshape A into 2d
sA = size (A);
nA = length (sA);
perA = [1:(dA-1) (dA+1):(nA-1) nA dA](1:nA);
Ap = permute (A, perA);
Ap = reshape (Ap, prod (sA(perA(1:end-1))), sA(perA(end)));
# reshape B into 2d
sB = size (B);
nB = length (sB);
perB = [dB 1:(dB-1) (dB+1):(nB-1) nB](1:nB);
Bp = permute (B, perB);
Bp = reshape (Bp, sB(perB(1)), prod (sB(perB(2:end))));
# multiply
M = Ap * Bp;
# reshape back to original format
s = [sA(perA(1:end-1)) sB(perB(2:end))];
M = squeeze (reshape (M, s));
endfunction
나는 당신이 사용하는 것이 좋습니다 MMX 도구 상자 Matlab. N- 차원 행렬에 가능한 빨리 곱할 수 있습니다.
의 장점 MMX 이다:
- 그것은이다 쉬운 사용.
- 곱하다 N- 차원 매트릭스 (실제로는 2D 행렬의 배열을 곱할 수 있습니다)
- 그것은 다른 것을 수행합니다 매트릭스 작업 (전환, 2 차 곱하기, chol 분해 등)
- 사용합니다 C 컴파일러 그리고 멀티 스레드 속도를 높이기위한 계산.
이 문제의 경우이 명령 만 작성하면됩니다.
C=mmx('mul',X,Y);
다음은 가능한 모든 방법에 대한 벤치 마크입니다. 자세한 내용은이를 참조하십시오 의문.
1.6571 # FOR-loop
4.3110 # ARRAYFUN
3.3731 # NUM2CELL/FOR-loop/CELL2MAT
2.9820 # NUM2CELL/CELLFUN/CELL2MAT
0.0244 # Loop Unrolling
0.0221 # MMX toolbox <===================
나는 재귀를 생각하지만 그것은 당신이 할 수있는 유일한 다른 비 루프 방법입니다.
루프를 "롤"할 수 있습니다. 즉, 루프에서 발생할 모든 곱셈을 순차적으로 작성합니다.
다음과 같은 문제에 대한 답변을 공유하고 싶습니다.
1) 두 텐서의 텐서 생성물 (모든 원자가);
2) 치수를 따라 두 개의 텐서의 수축.
다음은 첫 번째 및 두 번째 작업에 대한 서브 루틴입니다.
1) 텐서 제품 :
function [C] = tensor(A,B)
C = squeeze( reshape( repmat(A(:), 1, numel(B)).*B(:).' , [size(A),size(B)] ) );
end
2) 수축 : 여기서 A와 B는 각각 다임 션 I 및 J를 따라 수축 될 텐서입니다. 물론 이러한 치수의 길이는 동일해야합니다. 이것에 대한 점검은 없지만 (이것은 코드를 가릴 것입니다) 이것을 제외하고는 잘 작동합니다.
function [C] = tensorcontraction(A,B, i,j)
sa = size(A);
La = length(sa);
ia = 1:La;
ia(i) = [];
ia = [ia i];
sb = size(B);
Lb = length(sb);
ib = 1:Lb;
ib(j) = [];
ib = [j ib];
% making the i-th dimension the last in A
A1 = permute(A, ia);
% making the j-th dimension the first in B
B1 = permute(B, ib);
% making both A and B 2D-matrices to make use of the
% matrix multiplication along the second dimension of A
% and the first dimension of B
A2 = reshape(A1, [],sa(i));
B2 = reshape(B1, sb(j),[]);
% here's the implicit implication that sa(i) == sb(j),
% otherwise - crash
C2 = A2*B2;
% back to the original shape with the exception
% of dimensions along which we've just contracted
sa(i) = [];
sb(j) = [];
C = squeeze( reshape( C2, [sa,sb] ) );
end
비평가?