numpy : fonction pour max() et min() simultanés

Question

numpy.amax() trouvera la valeur maximale dans un tableau, et numpy.amin() fait la même chose pour la valeur min. Si je veux trouver à la fois le maximum et le minimum, je dois appeler les deux fonctions, ce qui nécessite de passer deux fois sur le (très grand) tableau, ce qui semble lent.

Existe-t-il une fonction dans l'API numpy qui trouve à la fois le maximum et le minimum en un seul passage dans les données ?

Answer 1

Existe-t-il une fonction dans l'API numpy qui trouve à la fois le maximum et le minimum en un seul passage dans les données ?

Non. Au moment où nous écrivons ces lignes, cette fonction n'existe pas. (Et oui, s'il y a étaient une telle fonction, ses performances seraient de manière significative mieux que d'appeler numpy.amin() y numpy.amax() successivement sur un grand tableau).

Answer 2

Je ne pense pas que passer deux fois sur le tableau soit un problème. Considérez le pseudo-code suivant :

minval = array[0]
maxval = array[0]
for i in array:
    if i < minval:
       minval = i
    if i > maxval:
       maxval = i

Bien qu'il n'y ait qu'une seule boucle ici, il y a toujours deux contrôles. (Au lieu d'avoir 2 boucles avec 1 contrôle chacune). La seule chose que vous économisez vraiment est l'overhead d'une boucle. Si les tableaux sont vraiment grands comme vous le dites, cette surcharge est faible par rapport à la charge de travail de la boucle réelle. (Notez que tout ceci est implémenté en C, donc les boucles sont plus ou moins gratuites de toute façon).

---

EDIT Désolé pour les 4 d'entre vous qui ont voté et ont eu confiance en moi. Vous pouvez définitivement optimiser ceci.

Voici du code fortran qui peut être compilé dans un module python via f2py (peut-être un Cython Un gourou peut venir et comparer avec une version C optimisée...) :

subroutine minmax1(a,n,amin,amax)
  implicit none
  !f2py intent(hidden) :: n
  !f2py intent(out) :: amin,amax
  !f2py intent(in) :: a
  integer n
  real a(n),amin,amax
  integer i

  amin = a(1)
  amax = a(1)
  do i=2, n
     if(a(i) > amax)then
        amax = a(i)
     elseif(a(i) < amin) then
        amin = a(i)
     endif
  enddo
end subroutine minmax1

subroutine minmax2(a,n,amin,amax)
  implicit none
  !f2py intent(hidden) :: n
  !f2py intent(out) :: amin,amax
  !f2py intent(in) :: a
  integer n
  real a(n),amin,amax
  amin = minval(a)
  amax = maxval(a)
end subroutine minmax2

Compilez-le via :

f2py -m untitled -c fortran_code.f90

Et maintenant nous sommes dans un endroit où nous pouvons le tester :

import timeit

size = 100000
repeat = 10000

print timeit.timeit(
    'np.min(a); np.max(a)',
    setup='import numpy as np; a = np.arange(%d, dtype=np.float32)' % size,
    number=repeat), " # numpy min/max"

print timeit.timeit(
    'untitled.minmax1(a)',
    setup='import numpy as np; import untitled; a = np.arange(%d, dtype=np.float32)' % size,
    number=repeat), '# minmax1'

print timeit.timeit(
    'untitled.minmax2(a)',
    setup='import numpy as np; import untitled; a = np.arange(%d, dtype=np.float32)' % size,
    number=repeat), '# minmax2'

Les résultats sont un peu stupéfiants pour moi :

8.61869883537 # numpy min/max
1.60417699814 # minmax1
2.30169081688 # minmax2

Je dois dire que je ne le comprends pas complètement. Comparer juste np.min contre minmax1 et minmax2 est toujours une bataille perdue, donc ce n'est pas seulement un problème de mémoire ...

notes -- Augmentation de la taille par un facteur de 10**a et en diminuant la répétition par un facteur de 10**a (en maintenant la taille du problème constante) modifie effectivement les performances, mais pas de manière apparemment cohérente, ce qui montre qu'il existe une certaine interaction entre les performances de la mémoire et le surcoût des appels de fonction dans Python. Même en comparant un simple min L'implémentation en fortran bat celle de numpy par un facteur d'environ 2 ...

Answer 3

Vous pourriez utiliser Numba qui est un compilateur dynamique Python tenant compte de NumPy et utilisant LLVM. L'implémentation qui en résulte est assez simple et claire :

import numpy
import numba

@numba.jit
def minmax(x):
    maximum = x[0]
    minimum = x[0]
    for i in x[1:]:
        if i > maximum:
            maximum = i
        elif i < minimum:
            minimum = i
    return (minimum, maximum)

numpy.random.seed(1)
x = numpy.random.rand(1000000)
print(minmax(x) == (x.min(), x.max()))

Elle devrait également être plus rapide que la méthode de Numpy min() & max() mise en œuvre. Et tout cela sans avoir à écrire une seule ligne de code en C/Fortran.

Effectuez vos propres tests de performance, car ils dépendent toujours de votre architecture, de vos données, des versions de vos paquets...

Answer 4

Il existe une fonction pour trouver (max-min) appelée numpy.ptp si cela vous est utile :

>>> import numpy
>>> x = numpy.array([1,2,3,4,5,6])
>>> x.ptp()
5

mais je ne pense pas qu'il y ait un moyen de trouver à la fois min et max avec une seule traversée.

EDIT : ptp appelle juste min et max sous le capot

Answer 5

Juste pour avoir une idée des chiffres auxquels on peut s'attendre, compte tenu des approches suivantes :

import numpy as np

def extrema_np(arr):
    return np.max(arr), np.min(arr)

import numba as nb

@nb.jit(nopython=True)
def extrema_loop_nb(arr):
    n = arr.size
    max_val = min_val = arr[0]
    for i in range(1, n):
        item = arr[i]
        if item > max_val:
            max_val = item
        elif item < min_val:
            min_val = item
    return max_val, min_val

import numba as nb

@nb.jit(nopython=True)
def extrema_while_nb(arr):
    n = arr.size
    odd = n % 2
    if not odd:
        n -= 1
    max_val = min_val = arr[0]
    i = 1
    while i < n:
        x = arr[i]
        y = arr[i + 1]
        if x > y:
            x, y = y, x
        min_val = min(x, min_val)
        max_val = max(y, max_val)
        i += 2
    if not odd:
        x = arr[n]
        min_val = min(x, min_val)
        max_val = max(x, max_val)
    return max_val, min_val

%%cython -c-O3 -c-march=native -a
#cython: language_level=3, boundscheck=False, wraparound=False, initializedcheck=False, cdivision=True, infer_types=True

import numpy as np

cdef void _extrema_loop_cy(
        long[:] arr,
        size_t n,
        long[:] result):
    cdef size_t i
    cdef long item, max_val, min_val
    max_val = arr[0]
    min_val = arr[0]
    for i in range(1, n):
        item = arr[i]
        if item > max_val:
            max_val = item
        elif item < min_val:
            min_val = item
    result[0] = max_val
    result[1] = min_val

def extrema_loop_cy(arr):
    result = np.zeros(2, dtype=arr.dtype)
    _extrema_loop_cy(arr, arr.size, result)
    return result[0], result[1]

%%cython -c-O3 -c-march=native -a
#cython: language_level=3, boundscheck=False, wraparound=False, initializedcheck=False, cdivision=True, infer_types=True

import numpy as np

cdef void _extrema_while_cy(
        long[:] arr,
        size_t n,
        long[:] result):
    cdef size_t i, odd
    cdef long x, y, max_val, min_val
    max_val = arr[0]
    min_val = arr[0]
    odd = n % 2
    if not odd:
        n -= 1
    max_val = min_val = arr[0]
    i = 1
    while i < n:
        x = arr[i]
        y = arr[i + 1]
        if x > y:
            x, y = y, x
        min_val = min(x, min_val)
        max_val = max(y, max_val)
        i += 2
    if not odd:
        x = arr[n]
        min_val = min(x, min_val)
        max_val = max(x, max_val)
    result[0] = max_val
    result[1] = min_val

def extrema_while_cy(arr):
    result = np.zeros(2, dtype=arr.dtype)
    _extrema_while_cy(arr, arr.size, result)
    return result[0], result[1]

(le extrema_loop_*() Les approches sont similaires à ce qui est proposé aquí alors que extrema_while_*() sont basées sur le code de aquí )

Les horaires suivants :

indiquent que le extrema_while_*() sont les plus rapides, avec extrema_while_nb() étant le plus rapide. Dans tous les cas, les extrema_loop_nb() y extrema_loop_cy() surpassent les performances de l'approche NumPy seule (en utilisant np.max() y np.min() séparément).

Enfin, notez qu'aucune de ces options n'est aussi flexible que l'option np.min() / np.max() (en termes de support n-dim, axis paramètre, etc.).

(le code complet est disponible aquí )