Подтвердить что ты не робот

Эффективная память массив массивных numpy в Python

Мне нужно отсортировать ОЧЕНЬ большой геномный набор данных, используя numpy. Я имею массив из 2.6 миллиардов поплавков, размеры = (868940742, 3), который занимает примерно 20 ГБ памяти на моей машине после загрузки и просто сидит там. У меня есть ранний планшет "MacBook Pro" на раннем этапе 2015 года с 16 ГБ оперативной памяти, твердотельный HD 500 ГБ и процессор Intel i7 с тактовой частотой 3,1 ГГц. Просто загрузите массив переполнения в виртуальную память, но не до такой степени, что моя машина страдает, или я должен остановить все остальное, что я делаю.

Я строю этот ОЧЕНЬ большой массив шаг за шагом из 22 меньших подкаров (N, 2).

Функция FUN_1 генерирует 2 новых массива (N, 1), используя каждый из 22 подмассивов, которые я называю sub_arr.

Первый вывод FUN_1 генерируется путем интерполяции значений из sub_arr[:,0] в массиве b = array([X, F(X)]), а второй вывод создается путем размещения sub_arr[:, 0] в бункеры с использованием массива r = array([X, BIN(X)]). Я вызываю эти выходы b_arr и rate_arr соответственно. Функция возвращает 3-кортеж (N, 1) массивов:

import numpy as np

def FUN_1(sub_arr):
    """interpolate b values and rates based on position in sub_arr"""

    b = np.load(bfile)
    r = np.load(rfile)

    b_arr = np.interp(sub_arr[:,0], b[:,0], b[:,1])
    rate_arr = np.searchsorted(r[:,0], sub_arr[:,0])  # HUGE efficiency gain over np.digitize...

    return r[rate_r, 1], b_arr, sub_arr[:,1] 

Я вызываю функцию 22 раза в цикле for и заполняю предварительно выделенный массив нулей full_arr = numpy.zeros([868940742, 3]) со значениями:

full_arr[:,0], full_arr[:,1], full_arr[:,2] = FUN_1

Что касается экономии памяти на этом этапе, я думаю, что это лучшее, что я могу сделать, но я открыт для предложений. В любом случае, я не сталкиваюсь с проблемами до этого момента, и это занимает всего около 2 минут.

Вот процедура сортировки (есть две последовательные сортировки)

for idx in range(2):
    sort_idx = numpy.argsort(full_arr[:,idx])
    full_arr = full_arr[sort_idx]
    # ...
    # <additional processing, return small (1000, 3) array of stats>

Теперь этот вид работал, хотя и медленно (занимает около 10 минут). Тем не менее, я недавно начал использовать более крупную таблицу с более высоким разрешением [X, F(X)] для шага интерполяции выше в FUN_1, который возвращает b_arr, и теперь SORT действительно замедляется, хотя все остальное остается неизменным.

Интересно, что я даже не сортирую на интерполированных значениях на шаге, где сортировка теперь отстает. Вот некоторые фрагменты различных файлов интерполяции: меньший по размеру на 30% меньше в каждом случае и гораздо более равномерный по значениям во втором столбце; более медленное имеет более высокое разрешение и множество других уникальных значений, поэтому результаты интерполяции, вероятно, более уникальны, но я не уверен, что это должно иметь какой-то эффект...?

больший, медленный файл:

17399307    99.4
17493652    98.8
17570460    98.2
17575180    97.6
17577127    97
17578255    96.4
17580576    95.8
17583028    95.2
17583699    94.6
17584172    94

меньший, более равномерный регулярный файл:

1       24  
1001    24  
2001    24  
3001    24  
4001    24  
5001    24
6001    24
7001    24

Я не уверен, что может вызвать эту проблему, и меня бы интересовали любые предложения или просто общий ввод о сортировке в этом типе памяти, ограничивающий регистр!

4b9b3361

Ответ 1

В настоящий момент каждый вызов np.argsort генерирует массив (868940742, 1) индексов int64, который сам по себе займет около 7 ГБ. Кроме того, когда вы используете эти индексы для сортировки столбцов full_arr, вы генерируете еще один массив массивов (868940742, 1), так как причудливая индексация всегда возвращает скорее копию чем вид.

Одним из довольно очевидных улучшений было бы сортировать full_arr на месте с помощью метода .sort(). К сожалению, .sort() не позволяет вам напрямую указывать строку или столбец для сортировки. Однако вы можете указать поле для сортировки для структурированного массива. Поэтому вы можете заставить сортировку inplace по одному из трех столбцов, получив view в свой массив как структурированный массив с тремя полями с плавающей запятой, затем сортировка по одному из этих полей:

full_arr.view('f8, f8, f8').sort(order=['f0'], axis=0)

В этом случае я сортирую full_arr на месте 0-го поля, которое соответствует первому столбцу. Обратите внимание, что я предположил, что есть три столбца float64 ('f8') - вы должны изменить это соответственно, если ваш тип dtype отличается. Это также требует, чтобы ваш массив был смежным и в основном формате, т.е. full_arr.flags.C_CONTIGUOUS == True.

Кредит за этот метод должен пойти Джо Кингтону за его ответ здесь.


Хотя для этого требуется меньше памяти, сортировка структурированного массива по полю, к сожалению, намного медленнее по сравнению с использованием np.argsort для создания массива индексов, как вы упомянули в комментариях ниже (см. этот предыдущий вопрос). Если вы используете np.argsort для получения набора индексов для сортировки, вы можете увидеть умеренное усиление производительности, используя np.take, а не прямую индексацию, чтобы получить отсортированный массив:

 %%timeit -n 1 -r 100 x = np.random.randn(10000, 2); idx = x[:, 0].argsort()
x[idx]
# 1 loops, best of 100: 148 µs per loop

 %%timeit -n 1 -r 100 x = np.random.randn(10000, 2); idx = x[:, 0].argsort()
np.take(x, idx, axis=0)
# 1 loops, best of 100: 42.9 µs per loop

Однако я бы не ожидал увидеть какую-либо разницу в использовании памяти, так как оба метода будут генерировать копию.


Относительно вашего вопроса о том, почему сортировка второго массива выполняется быстрее - да, вы должны ожидать, что любой разумный алгоритм сортировки будет быстрее, когда в массиве будет меньше уникальных значений, потому что в среднем меньше работы для этого. Предположим, что у меня есть случайная последовательность цифр от 1 до 10:

5  1  4  8  10  2  6  9  7  3

Есть 10!= 3628800 возможные способы расположения этих цифр, но только один, в котором они находятся в порядке возрастания. Теперь предположим, что есть только 5 уникальных цифр:

4  4  3  2  3  1  2  5  1  5

Теперь есть 2⁵ = 32 способа упорядочения этих цифр в порядке возрастания, так как я могу поменять любую пару одинаковых цифр в отсортированном векторе, не нарушая порядок.

По умолчанию np.ndarray.sort() использует Quicksort. Вариант qsort работает путем рекурсивного выбора элемента "pivot" в массиве, а затем переупорядочивает массив таким образом, чтобы все элементы перед ним помещается меньше значения поворота, а после него располагаются все элементы, превышающие значение поворота. Значения, равные оси, уже отсортированы. Имея меньше уникальных значений означает, что, в среднем, больше значений будет равно значению поворота на любой развертки, и, следовательно, меньше свип необходимы, чтобы полностью сортировать массив.

Например:

%%timeit -n 1 -r 100 x = np.random.random_integers(0, 10, 100000)
x.sort()
# 1 loops, best of 100: 2.3 ms per loop

%%timeit -n 1 -r 100 x = np.random.random_integers(0, 1000, 100000)
x.sort()
# 1 loops, best of 100: 4.62 ms per loop

В этом примере типы двух массивов одинаковы. Если ваш меньший массив имеет меньший размер элемента по сравнению с большим массивом, тогда стоимость его копирования из-за фантазии индексации также будет меньше.

Ответ 2

EDIT: в случае, если кто-то новый для программирования и numpy попадает на этот пост, я хочу указать на важность рассмотрения используемого np.dtype. В моем случае я действительно смог уйти с использованием плавающей запятой с половинной точностью, т.е. np.float16, которая уменьшила объект на 20 ГБ в памяти до 5 ГБ и сделала сортировку более управляемой. Значение по умолчанию, используемое numpy, равно np.float64, что является большой точностью, что вам может и не понадобиться. Ознакомьтесь с doc здесь, где описывается емкость различных типов данных. Спасибо @ali_m за то, что указали это в комментариях.

Я плохо справился с этим вопросом, но я обнаружил некоторые полезные обходные пути, которые, по моему мнению, будут полезны для всех, кому нужно отсортировать массивный массив numpy.

Я создаю очень большой массив numpy из 22 "поддиапазонов" данных генома человека, содержащих элементы [position, value]. В конечном счете, окончательный массив должен быть численно отсортирован "на месте" на основе значений в конкретном столбце и без перетасовки значений внутри строк.

Размеры подматрицы следуют форме:

arr1.shape = (N1, 2)
...
arr22.shape = (N22, 2)

sum([N1..N2]) = 868940742 т.е. существует около 1BN позиций для сортировки.

Сначала я обрабатываю 22 под-массива с помощью функции process_sub_arrs, которая возвращает 3-кортеж 1D-массивов той же длины, что и вход. Я складываю 1D массивы в новый массив (N, 3) и вставляю их в массив np.zeros, инициализированный для полного набора данных:

    full_arr = np.zeros([868940742, 3])
    i, j = 0, 0

    for arr in list(arr1..arr22):  
        # indices (i, j) incremented at each loop based on sub-array size
        j += len(arr)
        full_arr[i:j, :] = np.column_stack( process_sub_arrs(arr) )
        i = j

    return full_arr

EDIT: Поскольку я понял, что мой набор данных может быть представлен полуточными поплавками, теперь я инициализирую full_arr следующим образом: full_arr = np.zeros([868940742, 3], dtype=np.float16), который только 1/4 размера и намного проще сортировать.

Результат - массивный массив размером 20 ГБ:

full_arr.nbytes = 20854577808

Как отметил @ali_m в своем подробном сообщении, моя ранняя процедура была неэффективной:

sort_idx = np.argsort(full_arr[:,idx])
full_arr = full_arr[sort_idx]

массив sort_idx, который составляет 33% от размера full_arr, зависает и уничтожает память после сортировки full_arr. Этот тип, предположительно, генерирует копию full_arr из-за "фантазии" индексации, что потенциально увеличивает память до 233% того, что уже используется для хранения массивного массива! Это медленный шаг, который длится около десяти минут и в значительной степени опирается на виртуальную память.

Я не уверен, что "причудливый" вид делает постоянную копию. Наблюдая за использованием памяти на моей машине, кажется, что full_arr = full_arr[sort_idx] удаляет ссылку на несортированный оригинал, потому что примерно через 1 секунду все, что осталось, это память, используемая отсортированным массивом и индексом, даже если есть временная копия,

Более компактное использование argsort() для сохранения памяти - это следующее:

    full_arr = full_arr[full_arr[:,idx].argsort()]

Это все еще вызывает всплеск во время назначения, где и временный индексный массив, и временная копия, но память почти мгновенно освобождается снова.

@ali_m указал на хороший трюк (зачисленный Джо Кингтону) за создание де-факто структурированного массива с view на full_arr. Преимущество состоит в том, что они могут быть отсортированы "на месте", поддерживая стабильный порядок строк:

full_arr.view('f8, f8, f8').sort(order=['f0'], axis=0)

Представления отлично подходят для выполнения математических операций массива, но для сортировки он слишком неэффективен даже для одного подматрица из моего набора данных. В общем, структурированные массивы просто не очень хорошо масштабируются, даже если они обладают действительно полезными свойствами. Если кто-нибудь знает, почему это было бы интересно узнать.

Одним из хороших вариантов минимизации потребления памяти и повышения производительности с помощью очень больших массивов является создание конвейера небольших простых функций. Функции очищают локальные переменные после их завершения, поэтому, если промежуточные структуры данных наращивают и подрывают память, это может быть хорошим решением.

Это эскиз конвейера, который я использовал для ускорения массива сортировки массива:

def process_sub_arrs(arr):
    """process a sub-array and return a 3-tuple of 1D values arrays"""

    return values1, values2, values3

def build_arr():
    """build the initial array by joining processed sub-arrays"""

    full_arr = np.zeros([868940742, 3])
    i, j = 0, 0

    for arr in list(arr1..arr22):  
        # indices (i, j) incremented at each loop based on sub-array size
        j += len(arr)
        full_arr[i:j, :] = np.column_stack( process_sub_arrs(arr) )
        i = j

    return full_arr

def sort_arr():
    """return full_arr and sort_idx"""

    full_arr = build_arr()
    sort_idx = np.argsort(full_arr[:, index])

    return full_arr[sort_idx]

def get_sorted_arr():
    """call through nested functions to return the sorted array"""

    sorted_arr = sort_arr()
    <process sorted_arr>

    return statistics

стек вызовов: get_sorted_arr → sort_arr → build_arr → process_sub_arrs

Как только каждая внутренняя функция завершена, get_sorted_arr(), наконец, просто удерживает отсортированный массив, а затем возвращает небольшой массив статистики.

EDIT: Здесь также стоит отметить, что даже если вы можете использовать более компактный dtype для представления своего огромного массива, вы захотите использовать более высокую точность для сводных вычислений. Например, поскольку full_arr.dtype = np.float16, команда np.mean(full_arr[:,idx]) пытается вычислить среднее значение в точке с плавающей запятой с половинной точностью, но это быстро переполняется при суммировании по массивному массиву. Использование np.mean(full_arr[:,idx], dtype=np.float64) предотвратит переполнение.

Сначала я задал этот вопрос, потому что меня озадачило то, что набор данных одинакового размера внезапно начал задыхаться в моей системной памяти, хотя была большая разница в пропорции уникальных значений в новом "медленном" наборе. @ali_m отметил, что, действительно, более простые данные с меньшим количеством уникальных значений легче сортировать:

Вариант qsort Quicksort работает путем рекурсивного выбора 'pivot' в массиве, а затем переупорядочить массив таким образом, чтобы все перед ним помещаются элементы, меньшие, чем опорная величина, и все элементов, превышающих значение поворота, после него. Значения, которые равны стержню, уже отсортированы, поэтому интуитивно, чем меньше уникальных значений в массиве, тем меньше число которые необходимо сделать.

В этой заметке окончательное изменение, которое я решил сделать, чтобы попытаться решить эту проблему, заключалось в том, чтобы заблаговременно объединять новый набор данных, поскольку из этапа интерполяции был излишне высокий уровень десятичной точности. Это в конечном счете имело еще больший эффект, чем другие шаги сохранения памяти, показывая, что сам алгоритм сортировки был лимитирующим фактором в этом случае.

С нетерпением ждем других комментариев или предложений, которые могут быть у кого-либо по этой теме, и я почти наверняка оговорился о некоторых технических проблемах, поэтому я был бы рад услышать ответ: -)