Outer product of two rectangular matrices

조회 수: 332 (최근 30일)
Dario Cortese
Dario Cortese 2019년 2월 19일
답변: Branco Juan 2021년 12월 7일
If I have a vector r , I can easily calculate its inner product
r=[1 2 3];
inner = r*r'
inner =
14
Same goes for the outer product (please see here for complete definition)
outer=r'*r
outer =
1 2 3
2 4 6
3 6 9
outer, as it should be, has components (where N is the total number of components, here 3). inner, on the other hand has components (where m is the number of rows, here 1).
I want to be able to do this standard thing to rectangular matrices too. The inner product of rectangular matrices is easy enough:
r =
1 2 3
1 1 1
inner=r*r'
inner =
14 6
6 3
inner has components (2x2=4) and this is what I expect from the matrix multiplication of r with its transponse.
Clearly, though, it is not clear how I should calculate the outer product of r with itself, because now the definition of "inner product with transponse" and "outer product with itself" have the same syntax in matlab. Indeed, if I try to repeat what I have done for the vector r, I obtain:
outer=r'*r
outer =
2 3 4
3 5 7
4 7 10
which is not the outer product of r with itself, as it's evident from the fact that it does not have =36, but only components (where n is the number of column). What matlab has interpreted my calculation to be is the inner product of r transponse and r.
How do I obtain the proper outer product, whose components are all the combinations of products between components of r?

답변 (3개)

Jan
Jan 2019년 2월 19일
편집: Jan 2019년 2월 19일
If you have a [N x M] array and a [R x S x T] array, the output product becomes [N x M x R x S x T]. This can be done with nested for loops:
function C = OuterProduct(A, B) % version 1
sizeA = size(A);
sizeB = size(B);
C = zeros([sizeA, sizeB]);
if ndims(A) == 2 && ndims(B) == 2
for b2 = 1:sizeB(2)
for b1 = 1:sizeB(1)
for a2 = 1:sizeA(2)
for a1 = 1:sizeA(1)
C(a1, a2, b1, b2) = A(a1, a2) * B(b1, b2);
end
end
end
end
else
error('Not implemented yet.');
end
end
Is this correct so far? This is a primitive approach and it can be improved with respect to run-time. But this is a proof of concept yet. Now an improvement:
function C = OuterProduct(A, B) % version 2
sizeA = size(A);
sizeB = size(B);
C = zeros([sizeA, sizeB]);
if ndims(A) == 2 && ndims(B) == 2
for b2 = 1:sizeB(2)
for b1 = 1:sizeB(1)
%for a2 = 1:sizeA(2)
% for a1 = 1:sizeA(1)
C(:, :, b1, b2) = A * B(b1, b2);
% end
%end
end
end
else
error('Not implemented yet.');
end
end
A general solution is still ugly:
function C = OuterProduct(A, B) % version 3
C = zeros([size(A), size(B)]);
q = ndims(A) + 1;
index = cell(1, q);
index(1:q-1) = {':'}; % A cell as index vector of dynamic length
for k = 1:numel(B)
index{q} = k; % Linear index for elements of B
C(index{:}) = A * B(k);
end
end
Nicer and faster with bsxfun or the modern automatic expanding and the elementwise product:
function C = OuterProduct(A, B) % version 4
C = A .* reshape(B, [ones(1, ndims(A)), size(B)]);
% Matlab < R2016b:
% C = bsxfun(@times, A, reshape(B, [ones(1, ndims(A)), size(B)]))
end
Note: There is no reason to reshape A also:
AA = reshape(A, [size(A), ones(1, ndims(B))])
size(A)
size(AA) % Has still the same size!
Trailing dimensions of length 1 are ignored in Matlab.
Alternatively (see Matt J's solution):
function C = OuterProduct(A, B) % version 5
C = reshape(A(:) * B(:).', [size(A), size(B)]);
end
I love Matlab for the last two versions.

Matt J
Matt J 2019년 2월 19일
편집: Matt J 2019년 2월 19일
outer = r(:)*r(:).'

Branco Juan
Branco Juan 2021년 12월 7일
Notation: majuscule for Tensors, Matrices and Vectors; minuscule for their elements
In MatLab, the operator * is always the Matrix Product of matrices (tensors), which means it is the Dot Product for Vectors in Euclidean space (the Inner Product <V1;V2>=V1'.M.V2 with M being the identity).
It requires matrix dimensions to agree and suffices all due properties of such (preserving the order of the two tensors being operated, combining the dimensions involved and changing the number of elements)
Being for example V1 a raw vector of order 1 and dimension dim_v, that can be seen as a tensor of order 2 and dimensions (1;dim_v); V2 a column vector of order 1 and dimension dim_v, that can be seen as a tensor of order 2 and dimensions (dim_v;1)
DP = V1*V2 = sum( v1(i)*v2(i) ) for all i=1:dim_v, of order 2 (although it seems like order 1)
So that DP is a tensor (1x1) and the element dp(1,1) = v1(1,1)*v2(1,1) + v1(1,2)*v2(2,1) + ...v1(1,dim_v)*v2(dim_v,1)
So far, standard stuff.
Alternatively, being A a tensor of order 2 and dimensions (dim_1,dim_2); B a tensor of order 2 and dimensions (dim_2,dim_3), more general than your r and r' in the question above, where dim_1 = dim_3.
MP = A*B, of order 2
So that MP is a tensor (dim_1 x dim_3) and the elements mp(i,j) = sum( a(i,k)*b(k,j) ) for all k=1:dim_2
Do not confuse the operator * with the concept of Outer Product (well defined in the reference you've provided)
OP = A X B, of order 4
So that OP is a tensor (dim_1 x dim_2 x dim_2 x dim_3) and the elements op(i,j,k,g) = a(i,j)*b(k,g)
Matlab practitioners use the operator * to compute the Outer Product of Vectors because by chance it gives the same result, as you can see here:
MPV = V2*V1, of order 2
So that MPV is a tensor (dim_v x dim_v) and the elements mpv(i,j) = sum( a(i,k)*b(k,j) ) for all k=1:1;
thus mpv(i,j) = v1(i,1)*v2(1,j)
OPV = V1 X V2, of order 4
So that OPV is a tensor (1 x dim_v x dim_v x 1) and the elements opv(1,i,j,1) = v1(1,i)*v2(j,1) et voilà!!
You'll see the result as a (dim_v x dim_v) Matrix because we tend to ignore the dimensions of 1. But pay attention to the possible meaning of such possitions with respect to the final application.
So, once all the maths behind your problem are clear (I hope so), let's compute Outer Product of Tensors:
OPM = B X A
So that OPM should be a tensor of order 4, (dim_2 x dim_3 x dim_1 x dim_2) and the elements opm(k,g,i,j) = b(k,g)*a(i,j)
Solved via the new Matrix G of order 4 with the elements of B in its first two dimensions and repeated in the next two, and a new Matrix H of order 4 with the elements of A but in its first two dimensions, repeated in the next two, and finally permuted properly.
G should be (dim_2 x dim_3 x dim_1 x dim_2)
H should be (dim_2 x dim_3 x dim_1 x dim_2)
OPM = G.*H the element-wise product, or element-by-element multiplication.
So that OPM should be a tensor (dim_2 x dim_3 x dim_1 x dim_2) and the elements omp(m,n,i,j) = g(m,n)*h(i,j)
To solve your precise example:
A=[1,2,3;1,1,1]; % This is your r, a matrix (2x3)
B=A'; % This is your r', a matrix (3x2)
% you want to solve (r') Outer_Product (r), so Result = B Outer_Product A a
% Matrix (3x2x2x3) with 36 elements in total
% let's create the new matrices
G = repmat(B,[1 1 size(A)]); % this is a matrix (3x2x2x3)
Aux = repmat(A,[1 1 size(B)]); % this is a matrix (2x3x3x2), so we must reposition its elements
H = permute(Aux,[3 4 1 2]); % this is a matrix (3x2x2x3)
% the order [3 4 1 2] is applicable to any size matrices if they are of
% order 2.
Result = G.*H; % That's it
% Faster, without explanation
Result = B.*(permute(repmat(A,[1 1 size(B)]),[3 4 1 2]));
% No need to reshape B because MatLab is smart.
% Elapsed time is min 0.000481 seconds and max 0.001386 seconds
% reshape(B(:) * A(:).', [size(B), size(A)]); takes min 0.000244 seconds
% and max 0.000721 seconds
% B .* reshape(A, [ones(1, ndims(B)), size(A)]); takes min 0.000635 seconds
% and max 0.001464 seconds
% reshape(B, [size(B), ones(1, ndims(A))]); gives a wrong result...

카테고리

Help CenterFile Exchange에서 Matrices and Arrays에 대해 자세히 알아보기

태그

제품


릴리스

R2018b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by