Matlab - How to vectorize an indexing for loop?

Author

How to vectorize an indexing for loop?

Arthur G

2008-01-30, 9:35 am

Hi,

I'm trying to figure out if there's an easy and efficient way to
convert a particular kind of for loop into some sort of vectorized
statement. The for loop creates a column vector from a 2D matrix,
selecting one element per row.

For example, let's say the variable index says which element I want
from each row. Normally, I'd write something like this:
for i=1:nColumns
ySubset(i) = y(i, index(i));
end

If I try to write this as a single statement with indexing, e.g.
ySubset1 = y(:, index);
I incorrectly take entire columns from y.

I did realize that the diagonal of ySubset1 is what I want, so I
*could* do this:
ySubset2 = diag(y(:,index));

But that seems inefficient, because it could create a huge temporary
array if y is large.

Any ideas out there?
--Arthur

Randy Poe

2008-01-30, 9:35 am

On Jan 30, 8:30 am, Arthur G <gorramfr...@gmail.com> wrote:
> Hi,
>
> I'm trying to figure out if there's an easy and efficient way to
> convert a particular kind of for loop into some sort of vectorized
> statement. The for loop creates a column vector from a 2D matrix,
> selecting one element per row.
>
> For example, let's say the variable index says which element I want
> from each row. Normally, I'd write something like this:
> for i=1:nColumns
> ySubset(i) = y(i, index(i));
> end

This is where IND2SUB and SUB2IND come in handy.

Any Matlab array can be indexed in two ways: using
subscripts (row, column, page, etc), or considering it
as one long contiguous sequence in memory and using
an index into that sequence (which is how it's actually
stored). That is, if y is a 10 x 10 array and you want
to get the (5,5) element, you could ask for either y(5,5)
or y(45).

SUB2IND and IND2SUB convert between these
indexing schemes. They take as first argument the array
size of the array you are trying to index.

ySubset = y( sub2ind(size(y), 1:nColumns, index) );

And if you wanted to set these elements:

y( sub2ind(size(y), 1:nColumns, index) ) = ySubset;

- Randy

Peter Boettcher

2008-01-30, 9:35 am

Arthur G <gorramfreak@gmail.com> writes:

> Hi,
>
> I'm trying to figure out if there's an easy and efficient way to
> convert a particular kind of for loop into some sort of vectorized
> statement. The for loop creates a column vector from a 2D matrix,
> selecting one element per row.
>
> For example, let's say the variable index says which element I want
> from each row. Normally, I'd write something like this:
> for i=1:nColumns
> ySubset(i) = y(i, index(i));
> end
>
> If I try to write this as a single statement with indexing, e.g.
> ySubset1 = y(:, index);
> I incorrectly take entire columns from y.

ySubset = y(sub2ind(size(y), i, index));

-Peter

Arthur G

2008-01-30, 8:20 pm

On Jan 30, 8:57=A0am, Peter Boettcher <boettc...@ll.mit.edu> wrote:
> Arthur G <gorramfr...@gmail.com> writes:
>
>
>
>
> ySubset =3D y(sub2ind(size(y), i, index));
>
> -Peter

For anyone interested, I did some benchmarking, and for matrices with
few rows (say, less than 100), the diag approach is faster than
sub2ind. For larger matrices, sub2ind is faster. However, both are
always slower than just using a for loop. For larger matrices, it's
maybe 50 times faster (on my computer) to use the for loop than to use
sub2ind.

Thanks folks!
--Arthur
--Arthur

John D'Errico

2008-01-30, 11:14 pm

Arthur G <gorramfreak@gmail.com> wrote in message <0fe9c00d-52fa-
471c-914e-2f06cbf144f7@e10g2000prf.googlegroups.com>...
> For anyone interested, I did some benchmarking, and for matrices with
> few rows (say, less than 100), the diag approach is faster than
> sub2ind. For larger matrices, sub2ind is faster. However, both are
> always slower than just using a for loop. For larger matrices, it's
> maybe 50 times faster (on my computer) to use the for loop than to use
> sub2ind.

Perhaps your code was written poorly then.

A = rand(1000);
rowindex = (1:1000)';
colindex = ceil(rand(1000,1)*1000);

tic,Y0 = A(rowindex+(colindex-1)*1000);toc
Elapsed time is 0.000607 seconds.

tic,Y1 = zeros(1000,1);for i=1:1000,Y1(i)=A(rowindex(i),colindex(i)
);end,toc
Elapsed time is 0.007446 seconds.

tic,Y2 = A(sub2ind(size(A),rowindex,colindex));to
c
Elapsed time is 0.001511 seconds.

I'll stick with the alternatives to loops where I can.

John

Jesse Lai

2008-01-31, 10:18 am

John D'Errico wrote:
>
> Perhaps your code was written poorly then.
>
> A = rand(1000);
> rowindex = (1:1000)';
> colindex = ceil(rand(1000,1)*1000);
>
> tic,Y0 = A(rowindex+(colindex-1)*1000);toc
> Elapsed time is 0.000607 seconds.
>
> tic,Y1 = zeros(1000,1);for i=1:1000,Y1(i)=A(rowindex(i),colindex(i)
);end,toc
> Elapsed time is 0.007446 seconds.
>
> tic,Y2 = A(sub2ind(size(A),rowindex,colindex));to
c
> Elapsed time is 0.001511 seconds.
>
> I'll stick with the alternatives to loops where I can.
>
> John

John,

I tried this code on my machine, and indeed the for loop is faster. I
typically run code inside a for loop multiple times if the execution
time is less than a second just to get a better feel for the result.

If you run it in cell mode, then the for loop is much slower because JIT
is not applied (I ran in 2007a). But if you run the code with JIT
enabled, the for loop method is faster.

And I also got the for loop to be faster in the single case by about ten
times when JIT was enabled (0.000045 sec [for loop] vs. 0.000472 sec
[ind2sub]).

clc
A = rand(1000);
rowindex = (1:1000)';
colindex = ceil(rand(1000,1)*1000);

tic,Y0 = A(rowindex+(colindex-1)*1000);toc

Y1 = zeros(1000,1);
tic
for index =1:1000;
for i=1:1000
Y1(i)=A(rowindex(i),colindex(i));
end
end;
toc

tic
for index=1:1000;
Y2 = A(sub2ind(size(A),rowindex,colindex));
end;
toc

--------------------------------

Elapsed time is 0.000157 seconds.
Elapsed time is 0.041587 seconds.
Elapsed time is 0.220166 seconds.

Sponsored Links