MATLAB Central Newsreader - how to make the following case vectorized?

how to make the following case vectorized?

James Anderson — Wed, 04 Jun 2008 19:34:03 -0400

Hi,

I have a matrix (say 10000 rows and 50 columns). I want to
find the rank for each entry within the row it belongs to. I
know how to do this by using a loop, suppose A is the matrix

for k = 1:10000;
  [sortx,ind] = sort(A(k,:));
  rankt(ind) = 1:50;
  Rank(k,:) = rankt;
end

The Rank contains what I need. However, this is going to
take a long time, especially when the number of rows get
even larger. I found the function sort can sort the whole
matrix along one direction, however, is there any vectorized
way to get the Rank very fast without doing any looping?

Thanks,

James

Re: how to make the following case vectorized?

roberson@ibd.nrc-cnrc.gc.ca (Walter Roberson) — Wed, 04 Jun 2008 20:02:08 -0400

In article <g26qnb$4nt$1@fred.mathworks.com>,
James Anderson <janderson_net@yahoo.com> wrote:

>I have a matrix (say 10000 rows and 50 columns). I want to
>find the rank for each entry within the row it belongs to. I
>know how to do this by using a loop, suppose A is the matrix

>for k = 1:10000;
> [sortx,ind] = sort(A(k,:));
> rankt(ind) = 1:50;
> Rank(k,:) = rankt;
>end

>The Rank contains what I need. However, this is going to
>take a long time, especially when the number of rows get
>even larger. I found the function sort can sort the whole
>matrix along one direction, however, is there any vectorized
>way to get the Rank very fast without doing any looping?

You can do it as a two-step process:

[ind,ind] = sort(A,2);
[Rank,Rank] = sort(ind,2);

(The repeating of the variables on the left hand side is a signal
to Matlab that the value of the first return parameter is to be
discarded.)
--
"And that's the way it is." -- Walter Cronkite

Re: how to make the following case vectorized?

rych — Thu, 05 Jun 2008 11:52:06 -0400

On Jun 4, 9:02 pm, rober...@ibd.nrc-cnrc.gc.ca (Walter Roberson)
wrote:
>
> [ind,ind] = sort(A,2);
> [Rank,Rank] = sort(ind,2);
>
> (The repeating of the variables on the left hand side is a signal
> to Matlab that the value of the first return parameter is to be
> discarded.)

Where is this feature documented? I mean, the discarding of the first
output parameter?

Re: how to make the following case vectorized?

roberson@ibd.nrc-cnrc.gc.ca (Walter Roberson) — Thu, 05 Jun 2008 16:23:50 -0400

In article <8ffc93c0-2ef6-4875-8d0d-e61b723c9e56@i76g2000hsf.googlegroups.com>,
rych <rychphd@gmail.com> wrote:
>On Jun 4, 9:02 pm, rober...@ibd.nrc-cnrc.gc.ca (Walter Roberson)
>wrote:

>> [Rank,Rank] = sort(ind,2);

>> (The repeating of the variables on the left hand side is a signal
>> to Matlab that the value of the first return parameter is to be
>> discarded.)

>Where is this feature documented? I mean, the discarding of the first
>output parameter?

Good question. I poked around a bit but did not find an explicit
statement in the documentation.

I could imagine the statement instead being interpreted as
assigning the value according to the -first- output argument
with that name and discarding the remaining ones, but I know from
past experience that when a name is repeated, it is the -last-
output value that has effect.
--
"History is a pile of debris" -- Laurie Anderson

Re: how to make the following case vectorized?

Roger Stafford — Thu, 05 Jun 2008 19:33:01 -0400

"James Anderson" <janderson_net@yahoo.com> wrote in message <g26qnb
$4nt$1@fred.mathworks.com>...
> Hi,
>
> I have a matrix (say 10000 rows and 50 columns). I want to
> find the rank for each entry within the row it belongs to. I
> know how to do this by using a loop, suppose A is the matrix
>
> for k = 1:10000;
> [sortx,ind] = sort(A(k,:));
> rankt(ind) = 1:50;
> Rank(k,:) = rankt;
> end
>
> The Rank contains what I need. However, this is going to
> take a long time, especially when the number of rows get
> even larger. I found the function sort can sort the whole
> matrix along one direction, however, is there any vectorized
> way to get the Rank very fast without doing any looping?
>
> Thanks,
>
> James
--------------
Here is a way that avoids a second sort, though I don't know if it's any faster
than Walter's. It depends on how large n is. The sort would be of order
m*n*log(n) for large n, as against order m*n for the last two lines below.

[m,n] = size(A);
[ig,ind] = sort(A,2);
[r,Rank] = ndgrid(1:m,1:n);
Rank(r+m*(ind-1)) = Rank;

Roger Stafford

Re: how to make the following case vectorized?

Jos — Thu, 05 Jun 2008 21:18:01 -0400

"Roger Stafford"
<ellieandrogerxyzzy@mindspring.com.invalid> wrote in
message <g29f1d$ra2$1@fred.mathworks.com>...
> "James Anderson" <janderson_net@yahoo.com> wrote in
message <g26qnb
> $4nt$1@fred.mathworks.com>...
> > Hi,
> >
> > I have a matrix (say 10000 rows and 50 columns). I want
to
> > find the rank for each entry within the row it belongs
to. I
> > know how to do this by using a loop, suppose A is the
matrix
> >
> > for k = 1:10000;
> > [sortx,ind] = sort(A(k,:));
> > rankt(ind) = 1:50;
> > Rank(k,:) = rankt;
> > end
> >
> > The Rank contains what I need. However, this is going to
> > take a long time, especially when the number of rows get
> > even larger. I found the function sort can sort the
whole
> > matrix along one direction, however, is there any
vectorized
> > way to get the Rank very fast without doing any
looping?
> >
> > Thanks,
> >
> > James
> --------------
> Here is a way that avoids a second sort, though I don't
know if it's any faster
> than Walter's. It depends on how large n is. The sort
would be of order
> m*n*log(n) for large n, as against order m*n for the last
two lines below.
>
> [m,n] = size(A);
> [ig,ind] = sort(A,2);
> [r,Rank] = ndgrid(1:m,1:n);
> Rank(r+m*(ind-1)) = Rank;
>
> Roger Stafford
>
>

and a solution that avoids ndgrid as well as a seconc call
to sort:

[m,n] = size(A);
[ig,ind] = sort(A,2);
Rank3 = ind + repmat(n*(0:m-1),n,1).' ;
Rank3(Rank3) = repmat(1:n,m,1) ;
Rank3 = reshape(Rank3,[],m).' ;

hth
Jos

Re: how to make the following case vectorized?

Yi Cao — Thu, 05 Jun 2008 22:09:02 -0400

"Jos " <DELjos@jasenDEL.nl> wrote in message <g29l69
$90t$1@fred.mathworks.com>...
> "Roger Stafford"
> <ellieandrogerxyzzy@mindspring.com.invalid> wrote in
> message <g29f1d$ra2$1@fred.mathworks.com>...
> > "James Anderson" <janderson_net@yahoo.com> wrote in
> message <g26qnb
> > $4nt$1@fred.mathworks.com>...
> > > Hi,
> > >
> > > I have a matrix (say 10000 rows and 50 columns). I
want
> to
> > > find the rank for each entry within the row it
belongs
> to. I
> > > know how to do this by using a loop, suppose A is the
> matrix
> > >
> > > for k = 1:10000;
> > > [sortx,ind] = sort(A(k,:));
> > > rankt(ind) = 1:50;
> > > Rank(k,:) = rankt;
> > > end
> > >
> > > The Rank contains what I need. However, this is going
to
> > > take a long time, especially when the number of rows
get
> > > even larger. I found the function sort can sort the
> whole
> > > matrix along one direction, however, is there any
> vectorized
> > > way to get the Rank very fast without doing any
> looping?
> > >
> > > Thanks,
> > >
> > > James
> > --------------
> > Here is a way that avoids a second sort, though I
don't
> know if it's any faster
> > than Walter's. It depends on how large n is. The sort
> would be of order
> > m*n*log(n) for large n, as against order m*n for the
last
> two lines below.
> >
> > [m,n] = size(A);
> > [ig,ind] = sort(A,2);
> > [r,Rank] = ndgrid(1:m,1:n);
> > Rank(r+m*(ind-1)) = Rank;
> >
> > Roger Stafford
> >
> >
>
>
> and a solution that avoids ndgrid as well as a seconc
call
> to sort:
>
> [m,n] = size(A);
> [ig,ind] = sort(A,2);
> Rank3 = ind + repmat(n*(0:m-1),n,1).' ;
> Rank3(Rank3) = repmat(1:n,m,1) ;
> Rank3 = reshape(Rank3,[],m).' ;
>
> hth
> Jos

reshape can also be avoid:

[m,n]=size(A);
[Rank,ind]=sort(A,2);
Rank(bsxfun(@plus,(1:m)',(ind-1)*m))=repmat(1:n,m,1);

hth
Yi

Re: how to make the following case vectorized?

Sebastiaan — Fri, 06 Jun 2008 08:58:02 -0400

roberson@ibd.nrc-cnrc.gc.ca (Walter Roberson) wrote in
message <g293um$gs8$1@canopus.cc.umanitoba.ca>...
> In article
<8ffc93c0-2ef6-4875-8d0d-e61b723c9e56@i76g2000hsf.googlegroups.com>,
> rych <rychphd@gmail.com> wrote:
> >On Jun 4, 9:02 pm, rober...@ibd.nrc-cnrc.gc.ca (Walter
Roberson)
> >wrote:
>
> >> [Rank,Rank] = sort(ind,2);
>
> >> (The repeating of the variables on the left hand side
is a signal
> >> to Matlab that the value of the first return parameter
is to be
> >> discarded.)
>
> >Where is this feature documented? I mean, the discarding
of the first
> >output parameter?
>
> Good question. I poked around a bit but did not find an
explicit
> statement in the documentation.
>
> I could imagine the statement instead being interpreted as
> assigning the value according to the -first- output argument
> with that name and discarding the remaining ones, but I
know from
> past experience that when a name is repeated, it is the -last-
> output value that has effect.

Be carefull with things like this. It works allright:

void mexFunction(int nlhs, mxArray* plhs[], int nrhs, const
mxArray* prhs[])
{
    double *A, *B;
    if (nlhs!=2)
        mexErrMsgTxt("Two output arguments are required: [a,
b] = tester;");

    plhs[0] = mxCreateDoubleMatrix(1, 1, mxREAL);
    plhs[1] = mxCreateDoubleMatrix(1, 1, mxREAL);

    A = mxGetPr(plhs[0]);
    B = mxGetPr(plhs[1]);

    B[0] = 2; /* first initialize second argument, then */
    A[0] = 1; /* overwrite the first */
    return;
}

>> [a, a] = tester
a =
     2

Going OT here, but do not expect things always to work as
expected. Take this code:
{
    double *A;
    if (nrhs!=1)
        mexErrMsgTxt("One input/output argument is required:
tester2(A);");

    A = mxGetPr(prhs[0]);

    A[0] += 1;

    return;
}

i.e. overwrites the input argument (useful when working with
large arrays - saves the need of reallocating memory when
updating only a small part of a matrix).

At a first glance, it seems to work fine, altough the
variable is not updated in the 'workspace' tab:
>> a = 1;
>> testout2(a)
>> a
a =
     2

Now try this:
>> b = a;
>> testout2(a)
>> a
a =
     3
>> b
b =
     3

Which is wrong, since b was initialized when a==2. To the
end-user, it is unknown how Matlab handles its pointers
internally, and might change in between Matlab versions.

Ergo: only use these kind of constructions when necessary,
and only in a tested part of the code.

Greetings,
Sebastiaan