IT박스

numpy : max () 및 min () 동시 함수

itboxs 2020. 9. 8. 07:46
반응형

numpy : max () 및 min () 동시 함수


numpy.amax ()는 최대 값을 배열하고 찾을 numpy.amin () 않는 최소 값 동일. 최대 값과 최소값을 모두 찾으려면 두 함수를 모두 호출해야합니다.이 경우 (매우 큰) 배열을 두 번 전달해야합니다.

numpy API에 데이터를 한 번만 통과하여 최대와 최소를 모두 찾는 함수가 있습니까?


numpy API에 데이터를 한 번만 통과하여 최대와 최소를 모두 찾는 함수가 있습니까?

아니요.이 글을 쓰는 시점에는 그러한 기능이 없습니다. (그리고 예, 그러한 함수 있다면 그 성능은 대규모 배열에서 호출 하고 연속적으로 수행하는 것보다 훨씬 더 좋을 것 입니다.)numpy.amin()numpy.amax()


나는 배열을 두 번 전달하는 것이 문제라고 생각하지 않습니다. 다음 의사 코드를 고려하십시오.

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

여기에는 루프가 1 개만 있지만 여전히 2 개의 검사가 있습니다. (각각 1 개의 체크가있는 2 개의 루프를 갖는 대신). 실제로 저장하는 유일한 것은 루프 1 개의 오버 헤드입니다. 배열이 정말로 큰 경우 실제 루프의 작업 부하에 비해 오버 헤드가 적습니다. (이것은 모두 C로 구현되었으므로 루프는 어쨌든 다소 자유 롭습니다).


편집 나에게 찬성하고 믿음을 주신 4 분에게 죄송합니다. 확실히 최적화 할 수 있습니다.

다음을 통해 파이썬 모듈로 컴파일 할 수있는 포트란 코드가 있습니다 f2py(아마도 Cython전문가가이를 최적화 된 C 버전과 비교할 수 있습니다 ...).

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

다음을 통해 컴파일하십시오.

f2py -m untitled -c fortran_code.f90

이제 우리는 그것을 테스트 할 수있는 위치에 있습니다.

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'

결과는 저에게 약간 놀랍습니다.

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

나는 그것을 완전히 이해하지 못한다고 말해야한다. 그냥 np.minminmax1비교 minmax2는 여전히지는 싸움이므로 기억 문제가 아닙니다 ...

참고 -크기를 한 배 10**a늘리고 10**a(문제 크기를 일정하게 유지하는) 한 요소만큼 반복을 줄이면 성능이 변경되지만 메모리 성능과 함수 호출 오버 헤드간에 상호 작용이 있음을 보여주는 일관된 방식은 아닙니다. 파이썬. min포트란 의 간단한 구현을 비교해도 numpy를 약 2 배 능가합니다.


유용한 경우 numpy.ptp 라는 (최대-최소)를 찾는 함수 가 있습니다.

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

하지만 한 번의 순회로 최소값과 최대 값을 모두 찾을 수있는 방법은 없다고 생각합니다.

편집 : ptp는 후드 아래에서 최소 및 최대 호출


LLVM을 사용하는 NumPy 인식 동적 Python 컴파일러 인 Numba를 사용할 수 있습니다 . 그 결과 구현은 매우 간단하고 명확합니다.

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()))

또한 Numpy의 min() & max()구현 보다 빠릅니다 . 그리고 단일 C / Fortran 코드 줄을 작성할 필요가 없습니다.

항상 아키텍처, 데이터, 패키지 버전에 따라 다르므로 자체 성능 테스트를 수행하십시오.


이것은 오래된 실이지만 어쨌든 누군가 이것을 다시 본다면 ...

최소값과 최대 값을 동시에 찾을 때 비교 횟수를 줄일 수 있습니다. 그것이 수레라면 (내가 생각하기에) 이것은 계산 복잡성은 아니지만 시간을 절약 할 수 있습니다.

(Python 코드) 대신 :

_max = ar[0]
_min=  ar[0]
for ii in xrange(len(ar)):
    if _max > ar[ii]: _max = ar[ii]
    if _min < ar[ii]: _min = ar[ii]

먼저 배열에서 인접한 두 값을 비교 한 다음 작은 값만 현재 최소값과 비교하고 큰 값은 현재 최대 값과 비교할 수 있습니다.

## for an even-sized array
_max = ar[0]
_min = ar[0]
for ii in xrange(0, len(ar), 2)):  ## iterate over every other value in the array
    f1 = ar[ii]
    f2 = ar[ii+1]
    if (f1 < f2):
        if f1 < _min: _min = f1
        if f2 > _max: _max = f2
    else:
        if f2 < _min: _min = f2
        if f1 > _max: _max = f1

여기에있는 코드는 속도를 위해 C 또는 Fortran 또는 Cython을 사용하는 것이 분명한 Python으로 작성되었지만, 이렇게하면 3/2 * len (ar) 비교를 제공하는 len (ar) / 2 반복을 사용하여 반복 당 3 회 비교를 수행합니다. 반대로 "명백한 방법"을 비교하면 반복 당 두 번의 비교를 수행하여 2 * len (ar) 비교를 수행합니다. 비교 시간을 25 % 절약합니다.

언젠가 누군가가 유용하다고 생각할 것입니다.


일반적으로 한 번에 두 개의 요소를 처리하고 더 작은 요소를 임시 최소값과 더 큰 요소 만 임시 최대 값과 비교하여 최소값 알고리즘에 대한 비교 량을 줄일 수 있습니다. 평균적으로 순진한 접근 방식보다 3/4의 비교 만 필요합니다.

이는 c 또는 fortran (또는 기타 저수준 언어)으로 구현 될 수 있으며 성능 측면에서 거의 타의 추종을 불허합니다. 원칙을 설명하고 매우 빠르고 dtype 독립적 인 구현을 얻기 위해 사용 하고 있습니다.

import numba as nb
import numpy as np

@nb.njit
def minmax(array):
    # Ravel the array and return early if it's empty
    array = array.ravel()
    length = array.size
    if not length:
        return

    # We want to process two elements at once so we need
    # an even sized array, but we preprocess the first and
    # start with the second element, so we want it "odd"
    odd = length % 2
    if not odd:
        length -= 1

    # Initialize min and max with the first item
    minimum = maximum = array[0]

    i = 1
    while i < length:
        # Get the next two items and swap them if necessary
        x = array[i]
        y = array[i+1]
        if x > y:
            x, y = y, x
        # Compare the min with the smaller one and the max
        # with the bigger one
        minimum = min(x, minimum)
        maximum = max(y, maximum)
        i += 2

    # If we had an even sized array we need to compare the
    # one remaining item too.
    if not odd:
        x = array[length]
        minimum = min(x, minimum)
        maximum = max(x, maximum)

    return minimum, maximum

Peque가 제시 한 순진한 접근 방식보다 확실히 빠릅니다 .

arr = np.random.random(3000000)
assert minmax(arr) == minmax_peque(arr)  # warmup and making sure they are identical 
%timeit minmax(arr)            # 100 loops, best of 3: 2.1 ms per loop
%timeit minmax_peque(arr)      # 100 loops, best of 3: 2.75 ms per loop

예상대로 새로운 minmax 구현은 순진한 구현에 걸린 시간의 약 3/4 밖에 걸리지 않습니다 ( 2.1 / 2.75 = 0.7636363636363637).


아무도 numpy.percentile을 언급하지 않았기 때문에 그렇게 할 것이라고 생각했습니다. [0, 100]백분위 수 를 요청 하면 최소 (0 번째 백분위 수)와 최대 (100 번째 백분위 수)의 두 요소 배열이 제공됩니다.

그러나 OP의 목적을 충족시키지 못합니다. 최소 및 최대 개별적으로 빠르지 않습니다. 그 때문에 (더 단단한 문제, 비 극단적 인 백분위 수 있도록 할 몇 가지 기계에 아마 해야 오래 걸립니다).

In [1]: import numpy

In [2]: a = numpy.random.normal(0, 1, 1000000)

In [3]: %%timeit
   ...: lo, hi = numpy.amin(a), numpy.amax(a)
   ...: 
100 loops, best of 3: 4.08 ms per loop

In [4]: %%timeit
   ...: lo, hi = numpy.percentile(a, [0, 100])
   ...: 
100 loops, best of 3: 17.2 ms per loop

In [5]: numpy.__version__
Out[5]: '1.14.4'

Numpy의 향후 버전 [0, 100]은 요청 된 경우에만 정상적인 백분위 수 계산을 건너 뛰는 특별한 경우를 넣을 수 있습니다. 인터페이스에 아무것도 추가하지 않고 Numpy에게 한 번의 호출에서 최소 및 최대를 요청하는 방법이 있지만 (허용 된 답변에서 말한 것과는 반대로) 라이브러리의 표준 구현은이 경우를 활용하지 않습니다. 할 보람 있는.


언뜻보기에 트릭을 수행하는 것처럼 보입니다 .numpy.histogram

count, (amin, amax) = numpy.histogram(a, bins=1)

...하지만 해당 함수 소스살펴보면 단순히 호출 a.min()하고 a.max()독립적으로 수행되므로이 ​​질문에서 설명한 성능 문제를 피할 수 없습니다. :-(

마찬가지로 scipy.ndimage.measurements.extrema가능성처럼 보이지만, 단순히 호출 a.min()하고 a.max()독립적으로 수행합니다.


어쨌든 저에게는 노력할만한 가치가 있었으므로 관심이있는 사람을 위해 여기에서 가장 어렵고 덜 우아한 해결책을 제안 할 것입니다. 내 솔루션은 C ++의 한 패스 알고리즘에서 다중 스레드 최소-최대를 구현하고 이것을 사용하여 Python 확장 모듈을 만드는 것입니다. 이 작업에는 Python 및 NumPy C / C ++ API 사용 방법을 배우는 데 약간의 오버 헤드가 필요하며, 여기에서는 코드를 보여주고이 경로를 따르고 자하는 사람을위한 간단한 설명과 참조를 제공합니다.

다중 스레드 최소 / 최대

여기에는 너무 흥미로운 것이 없습니다. 배열은 크기의 덩어리로 나뉩니다 length / workers. 의 각 청크에 대해 최소 / 최대가 계산 future된 다음 전역 최소 / 최대에 대해 스캔됩니다.

    // mt_np.cc
    //
    // multi-threaded min/max algorithm

    #include <algorithm>
    #include <future>
    #include <vector>

    namespace mt_np {

    /*
     * Get {min,max} in interval [begin,end)
     */
    template <typename T> std::pair<T, T> min_max(T *begin, T *end) {
      T min{*begin};
      T max{*begin};
      while (++begin < end) {
        if (*begin < min) {
          min = *begin;
          continue;
        } else if (*begin > max) {
          max = *begin;
        }
      }
      return {min, max};
    }

    /*
     * get {min,max} in interval [begin,end) using #workers for concurrency
     */
    template <typename T>
    std::pair<T, T> min_max_mt(T *begin, T *end, int workers) {
      const long int chunk_size = std::max((end - begin) / workers, 1l);
      std::vector<std::future<std::pair<T, T>>> min_maxes;
      // fire up the workers
      while (begin < end) {
        T *next = std::min(end, begin + chunk_size);
        min_maxes.push_back(std::async(min_max<T>, begin, next));
        begin = next;
      }
      // retrieve the results
      auto min_max_it = min_maxes.begin();
      auto v{min_max_it->get()};
      T min{v.first};
      T max{v.second};
      while (++min_max_it != min_maxes.end()) {
        v = min_max_it->get();
        min = std::min(min, v.first);
        max = std::max(max, v.second);
      }
      return {min, max};
    }
    }; // namespace mt_np

Python 확장 모듈

여기에서 상황이 추악 해지기 시작합니다 ... Python에서 C ++ 코드를 사용하는 한 가지 방법은 확장 모듈을 구현하는 것입니다. 이 모듈은 distutils.core표준 모듈을 사용하여 구축 및 설치할 수 있습니다 . 이에 수반되는 내용에 대한 전체 설명은 Python 문서 https://docs.python.org/3/extending/extending.html 에서 다룹니다 . 참고 : https://docs.python.org/3/extending/index.html#extending-index 를 인용하여 유사한 결과를 얻을 수있는 다른 방법이 있습니다 .

이 가이드는이 버전의 CPython의 일부로 제공되는 확장을 만들기위한 기본 도구 만 다룹니다. Cython, cffi, SWIG 및 Numba와 같은 타사 도구는 Python 용 C 및 C ++ 확장을 만드는 데 더 간단하고 정교한 접근 방식을 제공합니다.

본질적으로이 경로는 실용적 이라기보다는 학문적 일 것입니다. 그 말을 듣고 다음에 내가 한 일은 튜토리얼에 매우 가깝게 모듈 파일을 만드는 것입니다. 이것은 본질적으로 distutils가 여러분의 코드로 무엇을해야하는지 알고 파이썬 모듈을 생성하기위한 상용구입니다. 이 작업을 수행하기 전에 시스템 패키지를 오염시키지 않도록 Python 가상 환경 을 만드는 것이 현명 할 것입니다 ( https://docs.python.org/3/library/venv.html#module-venv 참조 ).

다음은 모듈 파일입니다.

// mt_np_forpy.cc
//
// C++ module implementation for multi-threaded min/max for np

#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION

#include <python3.6/numpy/arrayobject.h>

#include "mt_np.h"

#include <cstdint>
#include <iostream>

using namespace std;

/*
 * check:
 *  shape
 *  stride
 *  data_type
 *  byteorder
 *  alignment
 */
static bool check_array(PyArrayObject *arr) {
  if (PyArray_NDIM(arr) != 1) {
    PyErr_SetString(PyExc_RuntimeError, "Wrong shape, require (1,n)");
    return false;
  }
  if (PyArray_STRIDES(arr)[0] != 8) {
    PyErr_SetString(PyExc_RuntimeError, "Expected stride of 8");
    return false;
  }
  PyArray_Descr *descr = PyArray_DESCR(arr);
  if (descr->type != NPY_LONGLTR && descr->type != NPY_DOUBLELTR) {
    PyErr_SetString(PyExc_RuntimeError, "Wrong type, require l or d");
    return false;
  }
  if (descr->byteorder != '=') {
    PyErr_SetString(PyExc_RuntimeError, "Expected native byteorder");
    return false;
  }
  if (descr->alignment != 8) {
    cerr << "alignment: " << descr->alignment << endl;
    PyErr_SetString(PyExc_RuntimeError, "Require proper alignement");
    return false;
  }
  return true;
}

template <typename T>
static PyObject *mt_np_minmax_dispatch(PyArrayObject *arr) {
  npy_intp size = PyArray_SHAPE(arr)[0];
  T *begin = (T *)PyArray_DATA(arr);
  auto minmax =
      mt_np::min_max_mt(begin, begin + size, thread::hardware_concurrency());
  return Py_BuildValue("(L,L)", minmax.first, minmax.second);
}

static PyObject *mt_np_minmax(PyObject *self, PyObject *args) {
  PyArrayObject *arr;
  if (!PyArg_ParseTuple(args, "O", &arr))
    return NULL;
  if (!check_array(arr))
    return NULL;
  switch (PyArray_DESCR(arr)->type) {
  case NPY_LONGLTR: {
    return mt_np_minmax_dispatch<int64_t>(arr);
  } break;
  case NPY_DOUBLELTR: {
    return mt_np_minmax_dispatch<double>(arr);
  } break;
  default: {
    PyErr_SetString(PyExc_RuntimeError, "Unknown error");
    return NULL;
  }
  }
}

static PyObject *get_concurrency(PyObject *self, PyObject *args) {
  return Py_BuildValue("I", thread::hardware_concurrency());
}

static PyMethodDef mt_np_Methods[] = {
    {"mt_np_minmax", mt_np_minmax, METH_VARARGS, "multi-threaded np min/max"},
    {"get_concurrency", get_concurrency, METH_VARARGS,
     "retrieve thread::hardware_concurrency()"},
    {NULL, NULL, 0, NULL} /* sentinel */
};

static struct PyModuleDef mt_np_module = {PyModuleDef_HEAD_INIT, "mt_np", NULL,
                                          -1, mt_np_Methods};

PyMODINIT_FUNC PyInit_mt_np() { return PyModule_Create(&mt_np_module); }

이 파일에는 NumPy API뿐만 아니라 Python이 많이 사용됩니다. 자세한 내용은 https://docs.python.org/3/c-api/arg.html#c.PyArg_ParseTuple 및 NumPy를 참조하십시오. : https://docs.scipy.org/doc/numpy/reference/c-api.array.html .

모듈 설치

다음으로 할 일은 distutils를 사용하여 모듈을 설치하는 것입니다. 여기에는 설정 파일이 필요합니다.

# setup.py

from distutils.core import setup,Extension

module = Extension('mt_np', sources = ['mt_np_module.cc'])

setup (name = 'mt_np', 
       version = '1.0', 
       description = 'multi-threaded min/max for np arrays',
       ext_modules = [module])

마지막으로 모듈을 설치하려면 python3 setup.py install가상 환경에서 실행 하십시오.

모듈 테스트

마지막으로 C ++ 구현이 실제로 NumPy의 순진한 사용을 능가하는지 테스트 할 수 있습니다. 이를 위해 다음은 간단한 테스트 스크립트입니다.

# timing.py
# compare numpy min/max vs multi-threaded min/max

import numpy as np
import mt_np
import timeit

def normal_min_max(X):
  return (np.min(X),np.max(X))

print(mt_np.get_concurrency())

for ssize in np.logspace(3,8,6):
  size = int(ssize)
  print('********************')
  print('sample size:', size)
  print('********************')
  samples = np.random.normal(0,50,(2,size))
  for sample in samples:
    print('np:', timeit.timeit('normal_min_max(sample)',
                 globals=globals(),number=10))
    print('mt:', timeit.timeit('mt_np.mt_np_minmax(sample)',
                 globals=globals(),number=10))

이 모든 작업을 통해 얻은 결과는 다음과 같습니다.

8  
********************  
sample size: 1000  
********************  
np: 0.00012079699808964506  
mt: 0.002468645994667895  
np: 0.00011947099847020581  
mt: 0.0020772050047526136  
********************  
sample size: 10000  
********************  
np: 0.00024697799381101504  
mt: 0.002037393998762127  
np: 0.0002713389985729009  
mt: 0.0020942929986631498  
********************  
sample size: 100000  
********************  
np: 0.0007130410012905486  
mt: 0.0019842900001094677  
np: 0.0007540129954577424  
mt: 0.0029724110063398257  
********************  
sample size: 1000000  
********************  
np: 0.0094779249993735  
mt: 0.007134920000680722  
np: 0.009129883001151029  
mt: 0.012836456997320056  
********************  
sample size: 10000000  
********************  
np: 0.09471094200125663  
mt: 0.0453535050037317  
np: 0.09436299200024223  
mt: 0.04188535599678289  
********************  
sample size: 100000000  
********************  
np: 0.9537652180006262  
mt: 0.3957935369980987  
np: 0.9624398809974082  
mt: 0.4019058070043684  

이는 스레드의 초기에 나타난 결과보다 훨씬 덜 고무적입니다. 이는 약 3.5 배의 속도 향상을 나타내며 멀티 스레딩을 통합하지 않았습니다. 내가 달성 한 결과는 다소 합리적이며, 스레딩 오버 헤드가 어레이가 매우 커질 때까지 시간을 지배 할 것으로 예상되며,이 시점에서 성능 증가가 std::thread::hardware_concurrencyx 증가 에 가까워지기 시작할 것 입니다.

결론

일부 NumPy 코드에 대한 애플리케이션 특정 최적화의 여지가 있습니다. 특히 멀티 스레딩과 관련하여 보일 것입니다. 노력할 가치가 있는지 여부는 분명하지 않지만 확실히 좋은 운동 (또는 무언가)처럼 보입니다. Cython과 같은 "타사 도구"중 일부를 배우는 것이 시간을 더 잘 활용할 수 있다고 생각하지만 누가 알겠습니까?

참고 URL : https://stackoverflow.com/questions/12200580/numpy-function-for-simultaneous-max-and-min

반응형