Is there a better vectorization technique?

3 次查看(过去 30 天)
I am trying to see if there are other ways of coding this code sample more efficiently. Here, y is an 1xM matrix, (say, 1x1000), and z is an NxM matrix, (say, 5x1000).
mean(ones(N,1)*y.^3 .* z,2)
This code works fine, but I worry of N increases a lot, that the:
ones(N,1)*y.^3
might get too wasteful and make everything slow down.
Thoughts?

采纳的回答

John D'Errico
John D'Errico 2012-8-29
Its not THAT terrible for a matrix that small. However, you can gain from the use of bsxfun in many instances. Here, the matrices are simply too small to really gain anything.
>> N = 5;M =1000;
>> y = rand(1,M);
>> z = rand(N,M);
>> mean(ones(N,1)*y.^3 .* z,2)
ans =
0.12412
0.11669
0.12102
0.11976
0.12196
>> mean(bsxfun(@times,y.^3,z),2)
ans =
0.12412
0.11669
0.12102
0.11976
0.12196
>> z*y.'.^3/M
ans =
0.12412
0.11669
0.12102
0.11976
0.12196
As you can see, all three solutions return the same result. All are equally valid.
Now I'll compare the times required.
>> timeit(@() mean(ones(N,1)*y.^3 .* z,2))
ans =
0.00023018
>> timeit(@() mean(bsxfun(@times,y.^3,z),2))
ans =
0.00026829
>> timeit(@() z*y.'.^3/M)
ans =
0.00016594
As I said, you don't gain much. In fact, bsxfun does not gain at all, and is a bit slower. But you can gain a bit, if you re-write the expression into the third form I've posed. Not much, but a bit.
  2 个评论
Kalamaya
Kalamaya 2012-8-29
I thought about using bsxfun. How is the timing if the N was very large, say, 2000 or something like that? (timeit is your own function?)

请先登录,再进行评论。

更多回答(1 个)

Jan
Jan 2012-8-29
编辑:Jan 2012-8-30
The power() operator is more expensive than the matrix multiplication. Therefore an explicit multiplication saves time:
M = 1000;
N = 5;
y = rand(1, M);
z = rand(N, M);
tic; for i=1:100, a = mean(ones(N,1) * y .^ 3 .* z, 2); end; toc
% Elapsed time is 0.036668 seconds.
tic; for i=1:100, a = z * y.' .^ 3 / M; end; toc
% Elapsed time is 0.026818 seconds.
tic; for i=1:100, a = z * (y .* y .* y)' / M; end; toc
% Elapsed time is 0.001327 seconds.
[EDITED] If the resulting array is very large, a multiplication is faster than a division, but the result can differ due to rounding:
a = z * (y .* y .* y)' * (1 / M);
For the small [5x1] array in the example this does not matter.

类别

Help CenterFile Exchange 中查找有关 Logical 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by