8000 ENH: minmax function · Issue #9836 · numpy/numpy · GitHub
[go: up one dir, main page]
Skip to content

ENH: minmax function #9836

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom
Open
jakirkham opened this issue Oct 7, 2017 · 32 comments
Open

ENH: minmax function #9836

jakirkham opened this issue Oct 7, 2017 · 32 comments
Labels
01 - Enhancement 15 - Discussion

Comments

@jakirkham
Copy link
Contributor

Would be really handy to have a function that computes and returns both the minimum and maximum of an array. This should be substantially faster than calling min and then max on the same array.

FWIW here's an SO question from some time ago that has gotten a fair bit of attention.
@rbkn99
Copy link
rbkn99 commented Oct 19, 2017

Can I try to implement this function as my first contribution?

@eric-wieser
Copy link
Member

Tricky to make this a ufunc, because max is currently np.maximum.reduce, but a reduction can't compute two values at once (unless we implement #8773)

@jakirkham
Copy link
Contributor Author

Sorry if I'm being a bit slow here, but wouldn't it be possible to just return an array containing two values? Admittedly I've not dug deep into ufuncs internal workings. So am probably missing something.

@eric-wieser
Copy link
Member

The issue is that max(a), being a shorthand for ufunc.reduce, is implemented as f(a[0], f(a[1], f(a[2], ...))), where f = np.maximum, and acts on scalars. That's not enough state to implement a minmax.

You could of course just implement this with NpyIter - I'm just saying this isn't well suited to a ufunc

@rbkn99
Copy link
rbkn99 commented Oct 20, 2017 •
edited
Loading

So, did I understand right that minmax can be implemented in another, not as a ufunc, way?

@eric-wieser
Copy link
Member
eric-wieser commented Oct 20, 2017 •
edited
Loading

Correct, but it would by far be preferable to adapt ufuncs to make it possible.

You could implement it as a gufunc (like linalg.det), but then you lose the axis argument, and it only works in 1D.

@mhvk
Copy link
Contributor
mhvk commented Oct 21, 2017

This is one of quite a number of requests to combine ufuncs in various ways (here, running two in parallel; often, two in series); I still like the idea of making it possible to programmatically make such combined ufuncs. It should actually not be that hard -- for some thoughts, see #8889 (comment). Now if only I had time to work on it...

@jakirkham
Copy link
Contributor Author

Since this can be done without creating a new type of ufuncs, would a proposal that added minmax without a new class of ufuncs still be accepted?

@mhvk
Copy link
Contributor
mhvk commented Nov 6, 2017

@jakirkham - rereading this, I think @eric-wieser is correct that it is in fact not easy to implement minmax. Just as an example, even the calculation of ptp is done by first calculating max and min and then subtracting them (in core/src/multiarray/calculation.c). My own sense would be that effort is better spent on making a general way to combine ufunc, as while it is inefficient to calculate min and max separately, my guess would be that in most cases this is not a very big part of a larger calculation (e.g., that involves any trig function or also needs a median).

@jakirkham
Copy link
Contributor Author

Why is it not easy to implement minmax? Are you just meaning as a ufunc? If so, you may be right. However if you mean the algorithm itself, I have to disagree. In fact, gave it a try and wrapped up a little Cython package that does just that ( https://github.com/jakirkham/cyminmax ). Turns out there are well defined bounds on the Big-O of such an algorithm as well. Though would not be surprised if NumPy has all sorts of neat tricks for something like this that I have missed in my naive implementation.

As to the usefulness of minmax, as with all tools their usefulness is dependent on the problem being solved. For me, taking min and max to normalize images is a regular occurrence. Even if there might be a computation that takes somewhat longer, a significant number of calls are made to min and max. At this stage, I find most other computations I'm doing have been well optimized excepting min and max. Having minmax fixes that. Given the large amount of traffic on the SO question linked above, would think that other people agree that minmax is worth solving too.

@eric-wieser
Copy link
Member
eric-wieser commented Nov 10, 2017 •
edited
Loading

Why is it not easy to implement minmax?

It's a little messy to implement minmax(1d_array) for all input dtypes in C. As you show, using cython makes it a lot easier, but right now there's not much precedent for using cython inside numpy.

It's harder to implement minmax(3d_array, axis=(0,2)), which you need to do for this to actually work as an underlying implementation for np.ptp. I suppose you could expose a public minmax python function that falls back on calling min and max separately if it's not a simple axis=None case.

I might prefer the name np.extents to np.minmax, were such a function to be created.

@jakirkham
Copy link
Contributor Author

It's true C is pretty rough when it comes to parameterizing functions over types. My impression is that most C programs leverage m4 to solve this problem. Though that is probably not a reasonable requirement to add to building NumPy. Have long wondered if people might start using Jinja to do the same thing. It looks like the Chromium project explored this enough to write up a developer guide. Though this is a bit off topic.

Certainly Cython helps as well. Though there is still some boilerplate required to glue everything together. This is mainly caused by Cython not supporting arbitrary dimensions for a memory view (e.g. could be 1-D or 3-D).

One option might be to rip the C generated code out from Cython and use that as a base for NumPy. Would need some tidying though. While it is not common to use Cython in NumPy, some random number generation code does use it and it is always possible to ship the C. Perhaps the lower maintenance burden is already a reason for NumPy to start thinking about it. 😉

Have a WIP implementation of axes in PR ( jakirkham/cyminmax#31 ). Though have no idea how well it performs yet. The approach is at least similar to what things 8000 like apply_along_axis. Namely transposing the axes. Also supports out in the same PR. Have not tried keepdims yet. Though think this should be possible with some more changes.

So minmax is the name used by MATLAB and C++. Though C++ oddly calls the traversal function minmax_element. Julia also calls it minmax. Interestingly so does Ruby. I don't think R has something like this in its core library, but there appear to be other packages that provided it and call it minmax. Would think it best to stick with minmax for consistency with users from these other languages that are likely to crossover into Python for numerical use cases and thus use NumPy.

@mhvk
Copy link
Contributor
mhvk commented Nov 10, 2017

@jakirkham - my comment really wasn't that minmax is not useful, just that in my opinion one can get greater reward for effort spent in trying to make combining ufuncs generally easier. But obviously that is just my opinion, and where I will spend the time (if I have some).

In this respect, in your example of normalizing an image, one also needs to subtract min and divide by max-min, so there are four operations, and having a minmax would thus produce at best a 33% improvement (if one ignores that for large arrays subtract and divide are substantially slower than min and max, as many more writes to memory have to be done). Writing a way to combine ufuncs (subtract and divide) would give you a similar improvement.

@jakirkham
Copy link
Contributor Author

Sure. Fusing operations together would be very powerful and widely useful. Though that feels like a different topic to me. Is there already an issue open on this topic? Would be very interested in that discussion. :)

Would add that the minmax algorithm is actually different than just doing min and max as it iterates over pairs of values and uses them to avoid roughly half of a comparison per value checked. Just adding this in case it is not obvious.

@eric-wieser
Copy link
Member
eric-wieser commented Nov 11, 2017 •
edited
Loading

Just adding this in case it is not obvious.

It certainly wasn't to me. This makes a good case for it needing to be a gufunc and not just a reduction ufunc.

Do you mean it does something like this?

max = min = data[0]
for i in range(1,n, 2):
    try:
        a, b = data[i:i+2]
    except ValueError:
        a = b = data[i]
    if b > a:
        a, b = b, a
    if a > max:
        max = a
    if b < min
        min = b

@eric-wieser
Copy link
Member

Though C++ oddly calls the traversal function minmax_element.

I think that if we add one of these to numpy, we should add both:

min, max = np.minimum_maximum(a, b)   # elementwise min max, like np.minimum
min, max = np.minmax(a)   # array-wise min max, like np.min

@jakirkham
Copy link
Contributor Author

Do you mean it does something like this?

Yes, correct.

It's true. It might not have been obvious. Hence why I drew attention to it. ;) Here's the equivalent code in cyminmax. The last if outside the for-loop is just handling the except case that you had above.

I think that if we add one of these to numpy, we should add both

Agreed. The elementwise variant looks nice as well.

@juliantaylor
Copy link
Contributor
juliantaylor commented Nov 11, 2017 •
edited
Loading

to get performance out of min/max on pcs you have to abandon your CS comparison count ideas :)
any half decent cpu can compute minmax without branching via special instructions. Unfortunately compilers tend to be picky when they emit them. For example most of the time you have to put an explicit if ... else ... instead of just an else, despite that looking like it is more operations, on the cpu it are less when the compiler catches on.
Another issue is that this works well for min or max but min and max together often still confuses compiler.

tldr: write it explicitly in compiler intrinsics to get good code, see the min/max implementation in numpy in numpy/core/src/umath/simd.inc.src

Some more hints, you can probably get axis reduction for free via np.lib._ureduce like median does, I also could have sworn we already had a PR for this lying around since years, but that might have also been for sumabs/sumsquared etc not minmax...

@eric-wieser
Copy link
Member

to get performance out of min/max on pcs you have to abandon your CS comparison count ideas

Those ideas at least work for dtype=object :)

@jakirkham
Copy link
Contributor Author

Thanks for the pointers on intrinsics, @juliantaylor. IIUC that only applies for integers (thinking of the packed integer intrinsics), correct? Or are there intrinsics for floats as well?

@jakirkham jakirkham mentioned this issue Aug 17, 2018
Open
@jakirkham
Copy link
Contributor Author

cc @jni (as we were discussing this previously)

@hofingermarkus
Copy link

Since it is related. I just wanted to add that it would also be great to get the according indices at the same time. At least by setting an according input flag. So maybe something like:

vmin = np.min(A,axis=1)
vmin, idmin = np.min(A, ret_ids = True, axis=1)
(vmin, idmin), (vmax, idmax) = np.minmax(A, ret_ids = True, axis=1)

@jni jni mentioned this issue May 1, 2019
Open
@jakirkham
Copy link
Contributor Author

Today one would probably use argmin and argmax to get those indices and then select the values after. My guess is we could add an argminmax that would be similar in nature. This of course hinges on adding minmax in the first place though.

@WarrenWeckesser
Copy link
Member

For what it's worth, I implemented a simple minmax as a gufunc in ufunclab. Currently, it doesn't handle datetime64 or timedelta64 (ran into the Cliffs of Casting on my journey up the NumPy C API learning curve. :)

@TomAugspurger TomAugspurger mentioned this issue May 28, 2020
Closed
@vadimkantorov vadimkantorov mentioned this issue Jul 17, 2020
Closed
@WarrenWeckesser WarrenWeckesser mentioned this issue Feb 24, 2021
Open
@iver56
Copy link
iver56 commented Jan 13, 2024 •
edited
Loading

For what it's worth, I implemented a SIMD-optimized minmax function for numpy arrays here: https://github.com/nomonosound/numpy-minmax

It's mainly optimized for float32, but is also partially optimized for int16

@hprodh
Copy link
hprodh commented Oct 28, 2024 •
edited
Loading

If minmax cannot be a reduction like min and max are, and at some point another function has to be implemented, the gain provided by this function would basically be due to the fact that the number of "iterative scanning" would be reduced to one.

I think that function should not be called minmax unless it was a reduction, and np.extents, as proposed earlier, seem fairly good. The bothering of implementing it and adding it as a new additional function might justify making it even more versatile, or generalized.

Therefore I think a more complete signature in the API could be like :
(vmins, idmins), (vmaxs, idmaxs) = np.extents(A, nmin=1, nmax=1, ret_ids=True, ignore_nans=True, axis=1)
where a number of min and max values could be provided as an argument, and also an option for ignoring nans.

This function could reduce the cost (and the code length) of use cases where np.max, np.min, np.argmax, np.argmin are used sequentially, including their nan-ignoring equivalents.

@jakirkham
Copy link
Contributor Author

Wonder if someone could do very a mild tweak to the array before passing it to a (private) ufunc leveraged for minmax by adding an extra dimension of length 2 with a view onto the original array like so to get the array to be more ufunc friendly

import numpy as np

a = np.random.random((10, 11))
av = np.lib.stride_tricks.as_strided(a,
                                     shape=(2, *a.shape),
                                     strides=(0, *a.strides),
                                     writeable=False)

assert (av[0] == a).all()
assert (av[1] == a).all()

# a_min, a_max = np._minmax(av)

@jakirkham
Copy link
Contributor Author

It is worth noting that np.divmod returns 2 values and is a ufunc

In [1]: import numpy as np

In [2]: np.divmod
Out[2]: <ufunc 'divmod'>

FWIW np.divmod was added in PR: #9063

@mhvk
Copy link
Contributor
mhvk commented Nov 14, 2024

Indeed, (g)ufunc can have multiple outputs. It really is a question of (i) deciding whether one should have a minmax and (ii) implementing it...

@jni
Copy link
Contributor
jni commented Nov 14, 2024

(i) I'm a definite yes
(ii) please. 😂

@mhvk
Copy link
Contributor
mhvk commented Nov 15, 2024

For the implementation, clearly there is already the one provided by @WarrenWeckesser (see #9836 (comment)), so one could start there. But the trick will be to actually get it to good performance with vectorization, etc., because this will not be able to use much of the machinery in numpy/_core/src/umath/fast_loop_macros.h yet one does not want to start from scratch.

@WarrenWeckesser
Copy link
Member

FYI & FWIW: The core calculations that I implemented are in https://github.com/WarrenWeckesser/ufunclab/blob/main/src/minmax/minmax_gufunc.h (which also includes calculations for argmin, argmax, argminmax, min_argmin and max_argmax). All the boilerplate gufunc wrapping is handled by code generation scripts, based on the configuration file https://github.com/WarrenWeckesser/ufunclab/blob/main/src/minmax/define_cxx_gufunc_extmod.py. I implemented just the integer and floating point types, and the object type.

My implementation is not optimized; it is not much more than the serial loop in C++ that you would expect. I have experimented with a SIMD implementation in C++ using xsimd in https://github.com/WarrenWeckesser/experiments/blob/main/c%2B%2B/xsimd/reductions/minmax.hpp (no NumPy, just C++). For NumPy, one would use the Google Highway library. (Adding SIMD implementations to ufunclab is on my aspirational to-do list, but not my immediate to-do list.)

I'd be happy to work on adding minmax to NumPy (assuming there was agreement that it should be added), but at the moment I have too many languishing work-in-progress PRs and projects, so I can't commit to something new (at least not for a few months).

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
01 - Enhancement 15 - Discussion
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
0