Creating and managing 3D matrices in Matlab

Question

I am desperate to avoid the for loop in Matlab, but I cannot figure out how to do this. Here's the situation:

I have two matrices mxn A and B and two vectors v and w length d . I want the outer one to multiply A and v to get the matrix mxnxd , where the entry is (i,j,k) A_(i,j) * v_k , and similarly for B and w .

Then I want to add the resulting mxnxd matrices, and then take mean along the last dimension to return the mxn matrix.

I'm sure I can handle the last part, but the first part is completely stuck to me. I tried using bsxfun no avail. Does anyone know an effective way to do this? Many thanks!

EDIT: This revision comes after three wonderful answers below. gnovice has the best answer to the question I asked, no doubt. However, the question I wanted to ask includes the squaring of each record before taking the average value. I forgot to mention this part initially. Given this annoyance, both other answers work well, but the clever trick of doing algebra before coding doesn't help this time. Thanks for the help, everyone!

Answer 1

EDIT:

Although the issue in the question has been updated, an algebraic approach can still be used to simplify the questions. You still don't need to worry about 3D matrices. Your result will be as follows:

 output = mean(v.^2).*A.^2 + 2.*mean(v.*w).*A.*B + mean(w.^2).*B.^2; 

If your matrices and vectors are large, this solution will give you much better performance due to less memory required compared to solutions using BSXFUN or REPMAT .

Explanation:

Assuming M is the m-by-n-by-d matrix that you get as a result, before you take the average of the third dimension, this is what the span of the third dimension will cover:

 M(i,j,:) = A(i,j).*v + B(i,j).*w; 

In other words, the vector v scaled by A(i,j) plus the vector w scaled by B(i,j) . And this is what you get when you apply elemental squaring:

 M(i,j,:).^2 = (A(i,j).*v + B(i,j).*w).^2; = (A(i,j).*v).^2 + ... 2.*A(i,j).*B(i,j).*v.*w + ... (B(i,j).*w).^2; 

Now, when you take the average of the third dimension, the result for each element of output(i,j) will be as follows:

 output(i,j) = mean(M(i,j,:).^2); = mean((A(i,j).*v).^2 + ... 2.*A(i,j).*B(i,j).*v.*w + ... (B(i,j).*w).^2); = sum((A(i,j).*v).^2 + ... 2.*A(i,j).*B(i,j).*v.*w + ... (B(i,j).*w).^2)/d; = sum((A(i,j).*v).^2)/d + ... sum(2.*A(i,j).*B(i,j).*v.*w)/d + ... sum((B(i,j).*w).^2)/d; = A(i,j).^2.*mean(v.^2) + ... 2.*A(i,j).*B(i,j).*mean(v.*w) + ... B(i,j).^2.*mean(w.^2); 

Answer 2

Try changing the vectors v and w to 1 x 1 xd :

  mean (bsxfun(@times, A, reshape(v, 1, 1, [])) ... + bsxfun(@times, B, reshape(w, 1, 1, [])), 3) 

Here I use [] in the reshape argument to tell him to fill out this dimension depending on the product of all other dimensions and the total number of elements in the vector.

Answer 3

Use repmat to split the matrix in the third dimension.

 A = 1 2 3 4 5 6 >> repmat(A, [1 1 10]) ans(:,:,1) = 1 2 3 4 5 6 ans(:,:,2) = 1 2 3 4 5 6 

and etc.

Answer 4

You still do not need to resort to explicit loops or an indirect loop using bsxfun et al. for your updated requirements. You can achieve what you want with a simple vector solution, as follows

 output = reshape(mean((v(:)*A(:)'+w(:)*B(:)').^2),size(A)); 

Since the OP only says that v and w are vectors of length d , the above solution should work for both row vectors and columns. If they are known as column vectors, v(:) can be replaced by v and similarly for w .

---

You can check if this matches Lambdageek's answer (modified to round conditions) as follows

 outputLG = mean ((bsxfun(@times, A, reshape(v, 1, 1, [])) ... + bsxfun(@times, B, reshape(w, 1, 1, []))).^2, 3); isequal(output,outputLG) ans = 1