Reputation: 1144
I was trying to vectorize a certain weighted sum but couldn't figure out how to do it. I have created a simple minimal working example below. I guess the solution involves either bsxfun or reshape and kronecker products but I still have not managed to get it working.
rng(1);
N = 200;
T1 = 5;
T2 = 7;
T3 = 10;
A = rand(N,T1,T2,T3);
w1 = rand(T1,1);
w2 = rand(T2,1);
w3 = rand(T3,1);
B = zeros(N,1);
for i = 1:N
for j1=1:T1
for j2=1:T2
for j3=1:T3
B(i) = B(i) + w1(j1) * w2(j2) * w3(j3) * A(i,j1,j2,j3);
end
end
end
end
A = B;
For the two dimensional case there is a smart answer here.
Upvotes: 3
Views: 109
Reputation: 781
If we're going the route of having functions anyway, and are favoring performance over elegance/brevity, then consider this:
function B = weightReduce(A, varargin)
B = A;
for i = length(varargin):-1:1
N = length(varargin{i});
B = reshape(B, [], N) * varargin{i};
end
end
This is the performance comparison I see:
tic;
for i = 1:10000
W = createWeights(w1,w2,w3);
B = reshape(A, size(A,1), [])*W(:);
end
toc
Elapsed time is 0.920821 seconds.
tic;
for i = 1:10000
B2 = weightReduce(A, w1, w2, w3);
end
toc
Elapsed time is 0.484470 seconds.
Upvotes: 1
Reputation: 22215
This is the logic behind it:
ww1 = repmat (permute (w1, [4, 1, 2, 3]), [N, 1, T2, T3]);
ww2 = repmat (permute (w2, [3, 4, 1, 2]), [N, T1, 1, T3]);
ww3 = repmat (permute (w3, [2, 3, 4, 1]), [N, T1, T2, 1 ]);
B = ww1 .* ww2 .* ww3 .* A;
B = sum (B(:,:), 2)
You can avoid permute
by creating w1
, w2
, and w3
in the appropriate dimension in the first place. Also you can use bsxfun
instead of repmat
as appropriate for extra performance, I'm just showing the logic here and repmat
is easier to follow.
EDIT: Generalised version for arbitrary input dimensions:
Dims = {N, T1, T2, T3}; % add T4, T5, T6, etc as appropriate
Params = cell (1, length (Dims));
Params{1} = rand (Dims{:});
for n = 2 : length (Dims)
DimSubscripts = ones (1, length (Dims)); DimSubscripts(n) = Dims{n};
RepSubscripts = [Dims{:}]; RepSubscripts(n) = 1;
Params{n} = repmat (rand (DimSubscripts), RepSubscripts);
end
B = times (Params{:});
B = sum (B(:,:), 2)
Upvotes: 1
Reputation: 65430
You can use an additional multiplication to modify the w1 * w2'
grid from the previous answer to then multiply by w3
as well. You can then use matrix multiplication again to multiply with a "flattened" version of A
.
W = reshape(w1 * w2.', [], 1) * w3.';
B = reshape(A, size(A, 1), []) * W(:);
You could wrap the creation of weights into it's own function and make this generalizable to N
weights. Since this uses recursion, N
is limited to your current recursion limit (500 by default).
function W = createWeights(W, varargin)
if numel(varargin) > 0
W = createWeights(W(:) * varargin{1}(:).', varargin{2:end});
end
end
And use it with:
W = createWeights(w1, w2, w3);
B = reshape(A, size(A, 1), []) * W(:);
Update
Using part of @CKT's very good suggestion to use kron
, we could modify createWeights
just a little bit.
function W = createWeights(W, varargin)
if numel(varargin) > 0
W = createWeights(kron(varargin{1}, W), varargin{2:end});
end
end
Upvotes: 5
Reputation: 781
Again, you couldn't generalize this that well for N-D unless you made some function to construct the Kronecker product vector, but how about
A = reshape(A, N, []) * kron(w3, kron(w2, w1));
Upvotes: 1