bsxfun implementation in matrix multiplication

Question

As always trying to learn more from you, I was hoping I could receive some help with the following code.

I need to accomplish the following:

1) I have a vector:

x = [1 2 3 4 5 6 7 8 9 10 11 12]

2) and a matrix:

A =[11    14    1
    5     8    18
    10    8    19
    13    20   16]

I need to be able to multiply each value from x with every value of A, this means:

new_matrix = [1* A
              2* A
              3* A
               ...
              12* A]

This will give me this new_matrix of size (12*m x n) assuming A (mxn). And in this case (12*4x3)

How can I do this using bsxfun from matlab? and, would this method be faster than a for-loop?

Regarding my for-loop, I need some help here as well... I am not able to storage each "new_matrix" as the loop runs :(

for i=x
new_matrix = A.*x(i)
end

Thanks in advance!!

EDIT: After the solutions where given

First solution

clear all
clc
x=1:0.1:50;
A = rand(1000,1000);
tic
val = bsxfun(@times,A,permute(x,[3 1 2]));
out = reshape(permute(val,[1 3 2]),size(val,1)*size(val,3),[]);
toc

Output:

Elapsed time is 7.597939 seconds.

Second solution

clear all
clc
x=1:0.1:50;
A = rand(1000,1000);
tic
Ps = kron(x.',A);
toc

Output:

Elapsed time is 48.445417 seconds.

Answer 1

Send x to the third dimension, so that singleton expansion would come into effect when bsxfun is used for multiplication with A, extending the product result to the third dimension. Then, perform the bsxfun multiplication -

val = bsxfun(@times,A,permute(x,[3 1 2])) 

Now, val is a 3D matrix and the desired output is expected to be a 2D matrix concatenated along the columns through the third dimension. This is achieved below -

out = reshape(permute(val,[1 3 2]),size(val,1)*size(val,3),[])

Hope that made sense! Spread the bsxfun word around! woo!! :)

Answer 2

The kron function does exactly that:

kron(x.',A)

Answer 3

Here is my benchmark of the methods mentioned so far, along with a few additions of my own:

function [t,v] = testMatMult()
    % data
    %{
    x = [1 2 3 4 5 6 7 8 9 10 11 12];
    A = [11 14 1; 5 8 18; 10 8 19; 13 20 16];
    %}
    x = 1:50;
    A = randi(100, [1000,1000]);

    % functions to test
    fcns = {
        @() func1_repmat(A,x)
        @() func2_bsxfun_3rd_dim(A,x)
        @() func2_forloop_3rd_dim(A,x)
        @() func3_kron(A,x)
        @() func4_forloop_matrix(A,x)
        @() func5_forloop_cell(A,x)
        @() func6_arrayfun(A,x)
    };

    % timeit
    t = cellfun(@timeit, fcns, 'UniformOutput',true);

    % check results
    v = cellfun(@feval, fcns, 'UniformOutput',false);
    isequal(v{:})
    %for i=2:numel(v), assert(norm(v{1}-v{2}) < 1e-9), end
end

% Amro
function B = func1_repmat(A,x)
    B = repmat(x, size(A,1), 1);
    B = bsxfun(@times, B(:), repmat(A,numel(x),1));
end

% Divakar
function B = func2_bsxfun_3rd_dim(A,x)
    B = bsxfun(@times, A, permute(x, [3 1 2]));
    B = reshape(permute(B, [1 3 2]), [], size(A,2));
end

% Vissenbot
function B = func2_forloop_3rd_dim(A,x)
    B = zeros([size(A) numel(x)], 'like',A);
    for i=1:numel(x)
        B(:,:,i) = x(i) .* A;
    end
    B = reshape(permute(B, [1 3 2]), [], size(A,2));
end

% Luis Mendo
function B = func3_kron(A,x)
    B = kron(x(:), A);
end

% SergioHaram & TheMinion
function B = func4_forloop_matrix(A,x)
    [m,n] = size(A);
    p = numel(x);
    B = zeros(m*p,n, 'like',A);
    for i=1:numel(x)
        B((i-1)*m+1:i*m,:) = x(i) .* A;
    end
end

% Amro
function B = func5_forloop_cell(A,x)
    B = cell(numel(x),1);
    for i=1:numel(x)
        B{i} = x(i) .* A;
    end
    B = cell2mat(B);
    %B = vertcat(B{:});
end

% Amro
function B = func6_arrayfun(A,x)
    B = cell2mat(arrayfun(@(xx) xx.*A, x(:), 'UniformOutput',false));
end

The results on my machine:

>> t
t =
    0.1650    %# repmat (Amro)
    0.2915    %# bsxfun in the 3rd dimension (Divakar)
    0.4200    %# for-loop in the 3rd dim (Vissenbot)
    0.1284    %# kron (Luis Mendo)
    0.2997    %# for-loop with indexing (SergioHaram & TheMinion)
    0.5160    %# for-loop with cell array (Amro)
    0.4854    %# arrayfun (Amro)

(Those timings can slightly change between different runs, but this should give us an idea how the methods compare)

Note that some of these methods are going to cause out-of-memory errors for larger inputs (for example my solution based on repmat can easily run out of memory). Others will get significantly slower for larger sizes but won't error due to exhausted memory (the kron solution for instance).

I think that the bsxfun method func2_bsxfun_3rd_dim or the straightforward for-loop func4_forloop_matrix (thanks to MATLAB JIT) are the best solutions in this case.

Of course you can change the above benchmark parameters (size of x and A) and draw your own conclusions :)

Answer 4

Just to add an alternative, you maybe can use cellfun to achieve what you want. Here's an example (slightly modified from yours):

x = randi(2, 5, 3)-1;
a = randi(3,3);
%// bsxfun 3D (As implemented in the accepted solution)
val = bsxfun(@and, a, permute(x', [3 1 2])); %//'
out = reshape(permute(val,[1 3 2]),size(val,1)*size(val,3),[]);
%// cellfun (My solution)
val2 = cellfun(@(z) bsxfun(@and, a, z), num2cell(x, 2), 'UniformOutput', false);
out2 = cell2mat(val2); % or use cat(3, val2{:}) to get a 3D matrix equivalent to val and then permute/reshape like for out
%// compare
disp(nnz(out ~= out2));

Both give the same exact result.

For more infos and tricks using cellfun, see: http://matlabgeeks.com/tips-tutorials/computation-using-cellfun/

And also this: https://stackoverflow.com/a/1746422/1121352

Answer 5

If your vector x is of lenght = 12 and your matrix of size 3x4, I don't think that using one or the other would change much in term of time. If you are working with higher size matrix and vector, now that might become an issue.

So first of all, we want to multiply a vector with a matrix. In the for-loop method, that would give something like that :

s = size(A);
new_matrix(s(1),s(2),numel(x)) = zeros;   %This is for pre-allocating. If you have a big vector or matrix, this will help a lot time efficiently.

for i = 1:numel(x)
    new_matrix(:,:,i)= A.*x(i)
end

This will give you 3D matrix, with each 3rd dimension being a result of your multiplication. If this is not what you are looking for, I'll be adding another solution which might be more time efficient with bigger matrixes and vectors.