numpy-discussion

[Numpy-discussion] max argmax combo

[Bill B. <wb...@gm...>]

I find myself often wanting both the max and the argmax of an array.
(And same for the other arg* functions)
Of course I can do argmax first then use fancy indexing to get the max as well.
But the result of argmax isn't really in a format that's readily
usable as an index.
You have to do something like
    a = rand(10,5)
    imax = a.argmax(axis=0)
    vmax = a[(imax, range(5))]

Which isn't terrible, just always takes me a moment to remember the
proper indexing expression.
Would a way to get the argmax and the max at the same time be of
interest to anyone else?  Maybe an extra 'ret' arg on argmax?

   a = rand(10,5)
   imax,vmax = a.argmax(axis=0,retmax=True)

--Bill

Re: [Numpy-discussion] max argmax combo

[Charles R H. <cha...@gm...>]

On 9/18/06, Bill Baxter <wb...@gm...> wrote:
>
> I find myself often wanting both the max and the argmax of an array.
> (And same for the other arg* functions)
> Of course I can do argmax first then use fancy indexing to get the max as
> well.
> But the result of argmax isn't really in a format that's readily
> usable as an index.
> You have to do something like
>     a = rand(10,5)
>     imax = a.argmax(axis=0)
>     vmax = a[(imax, range(5))]
>
> Which isn't terrible, just always takes me a moment to remember the
> proper indexing expression.
> Would a way to get the argmax and the max at the same time be of
> interest to anyone else?  Maybe an extra 'ret' arg on argmax?
>
>    a = rand(10,5)
>    imax,vmax = a.argmax(axis=0,retmax=True)


I don't generally like overloading return values, the function starts to
lose its definition and becomes a bit baroque where simply changing a
keyword value can destroy the viability of the following code. But I can see
the utility of what you want. Hmm,  this problem is not unique to argmax.
Maybe what we need is a general way to extract values, something like

extract(a, imax, axis=0)

to go along with all the single axis functions.

Chuck

Re: [Numpy-discussion] max argmax combo

[Bill B. <wb...@gm...>]

On 9/19/06, Charles R Harris <cha...@gm...> wrote:
> On 9/18/06, Bill Baxter <wb...@gm...> wrote:
> > I find myself often wanting both the max and the argmax of an array.
> > (And same for the other arg* functions)

> > You have to do something like
> >     a = rand(10,5)
> >     imax = a.argmax(axis=0)
> >     vmax = a[(imax, range(5))]
> >
> I don't generally like overloading return values, the function starts to
> lose its definition and becomes a bit baroque where simply changing a
> keyword value can destroy the viability of the following code.

Agreed.  Seems like the only justification is if you get multiple
results from one calculation but only rarely want the extra values.
It doesn't make sense to always return them, but it's also not worth
making a totally different function.


> But I can see the utility of what you want. Hmm,  this problem is not unique to argmax.
> Maybe what we need is a general way to extract values, something like
>
> extract(a, imax, axis=0)
>
> to go along with all the single axis functions.

Yes, I think that would be easier to remember.

It should also work for the axis=None case.
  imax = a.argmax(axis=None)
  v = extract(a, imax, axis=None)



--Bill

Re: [Numpy-discussion] max argmax combo

[Charles R H. <cha...@gm...>]

On 9/18/06, Bill Baxter <wb...@gm...> wrote:
>
> On 9/19/06, Charles R Harris <cha...@gm...> wrote:
> > On 9/18/06, Bill Baxter <wb...@gm...> wrote:
> > > I find myself often wanting both the max and the argmax of an array.
> > > (And same for the other arg* functions)
>
> > > You have to do something like
> > >     a = rand(10,5)
> > >     imax = a.argmax(axis=0)
> > >     vmax = a[(imax, range(5))]
> > >
> > I don't generally like overloading return values, the function starts to
> > lose its definition and becomes a bit baroque where simply changing a
> > keyword value can destroy the viability of the following code.
>
> Agreed.  Seems like the only justification is if you get multiple
> results from one calculation but only rarely want the extra values.
> It doesn't make sense to always return them, but it's also not worth
> making a totally different function.
>
>
> > But I can see the utility of what you want. Hmm,  this problem is not
> unique to argmax.
> > Maybe what we need is a general way to extract values, something like
> >
> > extract(a, imax, axis=0)
> >
> > to go along with all the single axis functions.
>
> Yes, I think that would be easier to remember.
>
> It should also work for the axis=None case.
>   imax = a.argmax(axis=None)
>   v = extract(a, imax, axis=None)


It shouldn't be too difficult to jig something up given all the example
code. I can do that, but I would like more input first. The questions I have
are these.

1) Should it be done?
2) Should it be a method? (functions being somewhat deprecated)
3) What name should it have?

I think Travis will have to weigh in on this. IIRC, he felt that the number
of methods was getting out of hand.

--Bill


Chuck

Re: [Numpy-discussion] max argmax combo

[Travis O. <oli...@ee...>]

Charles R Harris wrote:

>
>
> On 9/18/06, *Bill Baxter* <wb...@gm... 
> <mailto:wb...@gm...>> wrote:
>
>     On 9/19/06, Charles R Harris <cha...@gm...
>     <mailto:cha...@gm...>> wrote:
>     > On 9/18/06, Bill Baxter <wb...@gm...
>     <mailto:wb...@gm...>> wrote:
>     > > I find myself often wanting both the max and the argmax of an
>     array.
>     > > (And same for the other arg* functions)
>
>     > > You have to do something like
>     > >     a = rand(10,5)
>     > >     imax = a.argmax(axis=0)
>     > >     vmax = a[(imax, range(5))]
>     > >
>     > I don't generally like overloading return values, the function
>     starts to
>     > lose its definition and becomes a bit baroque where simply
>     changing a
>     > keyword value can destroy the viability of the following code.
>
>     Agreed.  Seems like the only justification is if you get multiple
>     results from one calculation but only rarely want the extra values.
>     It doesn't make sense to always return them, but it's also not worth
>     making a totally different function.
>
>
>     > But I can see the utility of what you want. Hmm,  this problem
>     is not unique to argmax.
>     > Maybe what we need is a general way to extract values, something
>     like
>     >
>     > extract(a, imax, axis=0)
>     >
>     > to go along with all the single axis functions.
>
>     Yes, I think that would be easier to remember.
>
>     It should also work for the axis=None case.
>       imax = a.argmax(axis=None)
>       v = extract(a, imax, axis=None)
>
>
> It shouldn't be too difficult to jig something up given all the 
> example code. I can do that, but I would like more input first. The 
> questions I have are these.
>
> 1) Should it be done?
> 2) Should it be a method? (functions being somewhat deprecated)
> 3) What name should it have?
>
> I think Travis will have to weigh in on this. IIRC, he felt that the 
> number of methods was getting out of hand.


I can support adding a *function* that does both.   It can't be named 
extract (that already exists).   There should be one for all the 
"arg"-like functions.

If somebody doesn't add it before 1.0 final, it can wait for 1.0.1

-Travis

Re: [Numpy-discussion] max argmax combo

[Charles R H. <cha...@gm...>]

On 9/18/06, Charles R Harris <cha...@gm...> wrote:
>
>
>
> On 9/18/06, Bill Baxter <wb...@gm...> wrote:
> >
> > On 9/19/06, Charles R Harris <cha...@gm...> wrote:
> > > On 9/18/06, Bill Baxter <wb...@gm...> wrote:
> > > > I find myself often wanting both the max and the argmax of an array.
> > > > (And same for the other arg* functions)
> >
> > > > You have to do something like
> > > >     a = rand(10,5)
> > > >     imax = a.argmax(axis=0)
> > > >     vmax = a[(imax, range(5))]
> > > >
> > > I don't generally like overloading return values, the function starts
> > to
> > > lose its definition and becomes a bit baroque where simply changing a
> > > keyword value can destroy the viability of the following code.
> >
> > Agreed.  Seems like the only justification is if you get multiple
> > results from one calculation but only rarely want the extra values.
> > It doesn't make sense to always return them, but it's also not worth
> > making a totally different function.
> >
> >
> > > But I can see the utility of what you want. Hmm,  this problem is not
> > unique to argmax.
> > > Maybe what we need is a general way to extract values, something like
> > >
> > > extract(a, imax, axis=0)
> > >
> > > to go along with all the single axis functions.
> >
> > Yes, I think that would be easier to remember.
> >
> > It should also work for the axis=None case.
> >   imax = a.argmax(axis=None)
> >   v = extract(a, imax, axis=None)
>
>
> It shouldn't be too difficult to jig something up given all the example
> code. I can do that, but I would like more input first. The questions I have
> are these.
>
> 1) Should it be done?
> 2) Should it be a method? (functions being somewhat deprecated)
> 3) What name should it have?
>
> I think Travis will have to weigh in on this. IIRC, he felt that the
> number of methods was getting out of hand.
>
> --Bill
>
>
> Chuck
>

I note that argsort also produces indexes that are hard to use in the ndim>1
case. The two problems aren't quite equivalent because argsort maintains
ndim while argmax reduces ndim by one, but it would be nice if there was
something that would work for both. Using imax[...,newaxis] would make the
shapes compatible for input but then the output would need to lose the
newaxis on return. Hmmm. These are both examples of an indirect indexing
problem where the indices on the specified axis are a function of the
indices on the other axis. Using the other indices in the argmax case yields
a scalar index while in the argsort case it yields a 1d array that can be
used for fancy indexing. I'm just floating around here trying to think of a
consistent way to regard these things. Ummm, and I think this could work. In
fact, the arrays could be even deeper and fancy indexing on the specified
axis would produce what was essentially an array of arrays. This is sort of
the photo-negative version of take.

Apropos the overloaded return types, I think that that is precisely what
makes it tricky to use functions that may return either copies or views. The
results should really be marked read only because otherwise unexpected
results can arise.

Chuck

Re: [Numpy-discussion] max argmax combo

[Francesc A. <fa...@ca...>]

A Dimarts 19 Setembre 2006 07:18, Charles R Harris va escriure:
> I note that argsort also produces indexes that are hard to use in the
> ndim>1 case.

Perhaps it is worth to mention here that I've always liked to have a sort()=
=20
and argsort() functionality merged in one shot function (or method). Having=
=20
both separated implies two sorts, while if I'd have such a combo available,=
=20
the resulting operation would take quite less time.

One can easily stablish kind of a pattern for situations where this could=20
happen: in all the 'arg*()' functions. Perhaps I'm wrong in my count, but=20
there appear to be only three of such functions in NumPy, namely, argmax,=20
argmin and argsort. Adding three additional 'combos' doesn't seem a lot to =
my=20
mind, but it can be just 'too much' for more common sense minds.

Cheers,

=2D-=20
>0,0<   Francesc Altet =C2=A0 =C2=A0 http://www.carabos.com/
V   V   C=C3=A1rabos Coop. V. =C2=A0=C2=A0Enjoy Data
 "-"

Re: [Numpy-discussion] max argmax combo

[Charles R H. <cha...@gm...>]

On 9/19/06, Francesc Altet <fa...@ca...> wrote:
>
> A Dimarts 19 Setembre 2006 07:18, Charles R Harris va escriure:
> > I note that argsort also produces indexes that are hard to use in the
> > ndim>1 case.
>
> Perhaps it is worth to mention here that I've always liked to have a
> sort()
> and argsort() functionality merged in one shot function (or method).
> Having
> both separated implies two sorts, while if I'd have such a combo
> available,
> the resulting operation would take quite less time.


Do you want both the indexes and index sorted array returned, or just sort
the array using indexes, i.e., something like

a.sort(kind="quicksort", method="indirect")

IIRC, the old numeric argsort actually did something like this under the
covers but discarded the sorted array.

One can easily stablish kind of a pattern for situations where this could
> happen: in all the 'arg*()' functions. Perhaps I'm wrong in my count, but
> there appear to be only three of such functions in NumPy, namely, argmax,
> argmin and argsort. Adding three additional 'combos' doesn't seem a lot to
> my
> mind, but it can be just 'too much' for more common sense minds.


Chuck

Re: [Numpy-discussion] max argmax combo

[Francesc A. <fa...@ca...>]

A Dimarts 19 Setembre 2006 19:21, Charles R Harris va escriure:
> On 9/19/06, Francesc Altet <fa...@ca...> wrote:
> > A Dimarts 19 Setembre 2006 07:18, Charles R Harris va escriure:
> > > I note that argsort also produces indexes that are hard to use in the
> > > ndim>1 case.
> >
> > Perhaps it is worth to mention here that I've always liked to have a
> > sort()
> > and argsort() functionality merged in one shot function (or method).
> > Having
> > both separated implies two sorts, while if I'd have such a combo
> > available,
> > the resulting operation would take quite less time.
>
> Do you want both the indexes and index sorted array returned, or just sort
> the array using indexes, i.e., something like
>
> a.sort(kind=3D"quicksort", method=3D"indirect")

Uh, I don't understand what do you mean by "sort the array using indexes",=
=20
sorry.

What I'd like is to get the indexes and the sorted array returned. Or, anot=
her=20
option could be getting the sorted indexes returned and a sort of the array=
=20
done in-place (which is surely faster).=20

Regards,

=2D-=20
>0,0<   Francesc Altet =C2=A0 =C2=A0 http://www.carabos.com/
V   V   C=C3=A1rabos Coop. V. =C2=A0=C2=A0Enjoy Data
 "-"

Re: [Numpy-discussion] max argmax combo

[Charles R H. <cha...@gm...>]

On 9/19/06, Francesc Altet <fa...@ca...> wrote:
>
> A Dimarts 19 Setembre 2006 19:21, Charles R Harris va escriure:
> > On 9/19/06, Francesc Altet <fa...@ca...> wrote:
> > > A Dimarts 19 Setembre 2006 07:18, Charles R Harris va escriure:
> > > > I note that argsort also produces indexes that are hard to use in
> the
> > > > ndim>1 case.
> > >
> > > Perhaps it is worth to mention here that I've always liked to have a
> > > sort()
> > > and argsort() functionality merged in one shot function (or method).
> > > Having
> > > both separated implies two sorts, while if I'd have such a combo
> > > available,
> > > the resulting operation would take quite less time.
> >
> > Do you want both the indexes and index sorted array returned, or just
> sort
> > the array using indexes, i.e., something like
> >
> > a.sort(kind="quicksort", method="indirect")
>
> Uh, I don't understand what do you mean by "sort the array using indexes",
> sorry.
>
> What I'd like is to get the indexes and the sorted array returned. Or,
> another
> option could be getting the sorted indexes returned and a sort of the
> array
> done in-place (which is surely faster).


Don't know about the faster part, it depends on the architecture and the
size of the elements being swapped. It is the fact that the new quicksort
does 1/2 as many swaps together with the type specific code that makes it
faster than the original numpy version. The overhead in the current indirect
sorts comes mostly from the indirect references to array elements in the
innermost loop and it is quite possible that doing swaps on both the indexes
and the array would be a bigger hit than the current argsort followed by a
suitable take type function.

Lexsort is another function that would benefit from what we are discussing.

Chuck

Re: [Numpy-discussion] max argmax combo

[Bill B. <wb...@gm...>]

On 9/20/06, Francesc Altet <fa...@ca...> wrote:
> A Dimarts 19 Setembre 2006 19:21, Charles R Harris va escriure:
> >
> > Do you want both the indexes and index sorted array returned, or just sort
> > the array using indexes, i.e., something like
> >
> > a.sort(kind="quicksort", method="indirect")
>
> Uh, I don't understand what do you mean by "sort the array using indexes",
> sorry.
>

I think he meant do an argsort first, then use fancy indexing to get
the sorted array.
For a 1-d array that's just

  ind = A.argsort()
  Asorted = A[ind]

That should be O(N lg N + N), aka O(N lg N)
For A >1-d, you need an indexing expression that's a little more
complicated, hence the discussion about making an "extract" function
for that purpose.  Then you could say:

  ind = A.argsort(axis=d)
  Asorted = A.extract(ind,axis=d)

and it should just work no matter what 'd' is.

--bb

Re: [Numpy-discussion] max argmax combo

[Charles R H. <cha...@gm...>]

On 9/19/06, Bill Baxter <wb...@gm...> wrote:
>
> On 9/20/06, Francesc Altet <fa...@ca...> wrote:
> > A Dimarts 19 Setembre 2006 19:21, Charles R Harris va escriure:
> > >
> > > Do you want both the indexes and index sorted array returned, or just
> sort
> > > the array using indexes, i.e., something like
> > >
> > > a.sort(kind="quicksort", method="indirect")
> >
> > Uh, I don't understand what do you mean by "sort the array using
> indexes",
> > sorry.
> >
>
> I think he meant do an argsort first, then use fancy indexing to get
> the sorted array.
> For a 1-d array that's just
>
>   ind = A.argsort()
>   Asorted = A[ind]
>
> That should be O(N lg N + N), aka O(N lg N)
> For A >1-d, you need an indexing expression that's a little more
> complicated, hence the discussion about making an "extract" function
> for that purpose.  Then you could say:
>
>   ind = A.argsort(axis=d)
>   Asorted = A.extract(ind,axis=d)


I'm also thinking of the name argtake.

Chuck

Re: [Numpy-discussion] max argmax combo

[Francesc A. <fa...@ca...>]

A Dimarts 19 Setembre 2006 21:41, Bill Baxter va escriure:
> I think he meant do an argsort first, then use fancy indexing to get
> the sorted array.
> For a 1-d array that's just
>
>   ind =3D A.argsort()
>   Asorted =3D A[ind]
>
> That should be O(N lg N + N), aka O(N lg N)

I see. Thanks.  OTOH, maybe your estimations are right, but the effect of t=
he=20
constants in these O(whatever) estimations can be surprisingly high:

In [3]: from timeit import Timer

In [4]: Timer("b=3Dnumpy.argsort(a);c=3Da[b]", "import numpy;=20
a=3Dnumpy.arange(100000,-1,-1)").repeat(3,100)
Out[4]: [1.6653108596801758, 1.670341968536377, 1.6632120609283447]

In [5]: Timer("b=3Dnumpy.argsort(a);c=3Dnumpy.sort(a)", "import numpy;=20
a=3Dnumpy.arange(100000,-1,-1)").repeat(3,100)
Out[5]: [1.6533238887786865, 1.6272940635681152, 1.6253311634063721]

In [6]: Timer("b=3Dnumpy.argsort(a);a.sort();c=3Da", "import numpy;=20
a=3Dnumpy.arange(100000,-1,-1)").repeat(3,100)
Out[6]: [0.95492100715637207, 0.90312504768371582, 0.90426898002624512]

so, it seems that argsorting first and fancy indexing later on is the most=
=20
expensive procedure for relatively large arrays (100000).

Interestingly, figures above seems to indicate that in-place sort is=20
stunningly fast:

In [7]: Timer("a.sort()","import numpy;=20
a=3Dnumpy.arange(100000,-1,-1)").repeat(3,100)
Out[7]: [0.32840394973754883, 0.2746579647064209, 0.2770991325378418]

and much faster indeed than fancy indexing

In [8]: Timer("b[a]","import numpy;=20
a=3Dnumpy.arange(100000,-1,-1);b=3Da.copy()").repeat(3,100)
Out[8]: [0.79876089096069336, 0.74172186851501465, 0.74209499359130859]

i.e. in-place sort seems 2.7x faster than fancy indexing (at least for thes=
e=20
datasets).

Mmmm, with this, I really ponder if a combo that makes the argsort() and=20
sort() in one shot really makes any sense, at least from the point of view =
of=20
speed:

In [10]: Timer("b=3Dnumpy.argsort(a);a.sort();c=3Da","import numpy;=20
a=3Dnumpy.arange(100000,-1,-1)").repeat(3,100)
Out[10]: [0.98506593704223633, 0.89880609512329102, 0.89982390403747559]

In [11]: Timer("b=3Dnumpy.argsort(a)","import numpy;=20
a=3Dnumpy.arange(100000,-1,-1)").repeat(3,100)
Out[11]: [0.92959284782409668, 0.85385990142822266, 0.87773990631103516]

So, it seems that doing an in-place sort() immediately after an argsort()=20
operation is very efficient (cache effects here?), and would avoid the need=
=20
of the combo function (from the point of view of efficiency, I repeat).

Cheers,

=2D-=20
>0,0<   Francesc Altet =A0 =A0 http://www.carabos.com/
V   V   C=E1rabos Coop. V. =A0=A0Enjoy Data
 "-"