Есть ли преимущество в скорости анализа или использовании памяти при использовании HDF5 для хранения больших массивов (вместо плоских двоичных файлов)?

Я обрабатываю большие трехмерные массивы, которые мне часто приходится нарезать различными способами для проведения разнообразного анализа данных. Типичный «куб» может иметь размер ~ 100 ГБ (и, вероятно, в будущем он станет больше).

Кажется, что типичный рекомендуемый формат файла для больших наборов данных в python - использовать HDF5 (либо h5py, либо pytables). Мой вопрос: есть ли какое-либо преимущество в скорости или использовании памяти при использовании HDF5 для хранения и анализа этих кубов по сравнению с их хранением в простых плоских двоичных файлах? Подходит ли HDF5 для табличных данных, а не для больших массивов, с которыми я работаю? Я вижу, что HDF5 может обеспечить хорошее сжатие, но меня больше интересует скорость обработки и проблема переполнения памяти.

Я часто хочу проанализировать только одно большое подмножество куба. Один из недостатков pytables и h5py заключается в том, что, когда я беру часть массива, я всегда получаю обратно массив numpy, используя память. Однако, если я нарежу массивную memmap плоского двоичного файла, я могу получить представление, в котором данные хранятся на диске. Таким образом, мне кажется, что мне легче анализировать определенные сектора моих данных, не переполняя память.

Я исследовал как pytables, так и h5py, и пока не видел пользы от них для моей цели.

Преимущества HDF5: организация, гибкость, совместимость.

Некоторые из основных преимуществ HDF5 - это его иерархическая структура (аналогичная папкам / файлам), необязательные произвольные метаданные, хранящиеся с каждым элементом, и его гибкость (например, сжатие). Эта организационная структура и хранилище метаданных могут показаться банальными, но на практике они очень полезны.

Еще одно преимущество HDF состоит в том, что наборы данных могут быть фиксированного или гибкого размера. Следовательно, легко добавлять данные в большой набор данных без необходимости создавать новую копию.

Кроме того, HDF5 - это стандартизированный формат с библиотеками, доступными практически для любого языка, поэтому совместное использование ваших данных на диске между, скажем, Matlab, Fortran, R, C и Python с HDF очень просто. (Честно говоря, это не так уж сложно и с большим двоичным массивом, если вы знаете упорядочение C и F и знаете форму, dtype и т. Д. Хранимого массива.)

Преимущества HDF для большого массива: более быстрый ввод-вывод произвольного слоя

Так же, как TL / DR: для 3D-массива размером ~ 8 ГБ чтение "полного" среза по любой оси заняло ~ 20 секунд с фрагментированным набором данных HDF5 и от 0,3 секунды (в лучшем случае) до более трех часов (в худшем случае) для запомненный массив тех же данных.

Помимо перечисленных выше вещей, есть еще одно большое преимущество «фрагментированного» * ​​формата данных на диске, такого как HDF5: чтение произвольного фрагмента (акцент на произвольном) обычно выполняется намного быстрее, так как данные на диске более непрерывны. в среднем.

*(HDF5 не обязательно должен быть форматом данных с фрагментами. Он поддерживает фрагменты, но не требует этого. На самом деле h5py, если я правильно помню , по умолчанию для создания набора данных не используется фрагмент.)

По сути, ваша скорость чтения с диска в лучшем случае и скорость чтения с диска в худшем случае для данного фрагмента набора данных будут довольно близки с фрагментированным набором данных HDF (при условии, что вы выбрали разумный размер фрагмента или позволили библиотеке выбрать один за вас). С простым двоичным массивом лучший случай быстрее, но худший - намного хуже.

Одно предостережение: если у вас SSD, вы, вероятно, не заметите огромной разницы в скорости чтения / записи. Однако с обычным жестким диском последовательное чтение происходит намного быстрее, чем случайное. (то есть у обычного жесткого диска долгое seekвремя.) HDF по-прежнему имеет преимущество перед SSD, но это больше связано с другими функциями (например, метаданными, организацией и т. д.), чем с чистой скоростью.

Во-первых, чтобы устранить путаницу, доступ к h5pyнабору данных возвращает объект, который ведет себя довольно аналогично массиву numpy, но не загружает данные в память, пока они не будут разрезаны. (Подобно memmap, но не идентично.) Для получения дополнительной информации ознакомьтесь с h5pyвведением .

Нарезка набора данных загрузит подмножество данных в память, но, предположительно, вы хотите что-то с ними сделать, и в этот момент они вам все равно понадобятся в памяти.

Если вы действительно хотите выполнять вычисления вне ядра, вы можете легко использовать для табличных данных pandasили pytables. Это возможно h5py(лучше для больших массивов ND), но вам нужно опуститься на более низкий уровень и самостоятельно выполнить итерацию.

Тем не менее, будущее многоплановых вычислений вне ядра - за Blaze. Взгляните на него, если действительно хотите пойти по этому маршруту.

Дело "без частей"

Прежде всего, рассмотрим записанный на диск трехмерный C-упорядоченный массив (я смоделирую его, вызвав arr.ravel()и распечатав результат, чтобы сделать вещи более заметными):

In [1]: import numpy as np

In [2]: arr = np.arange(4*6*6).reshape(4,6,6)

In [3]: arr
Out[3]:
array([[[  0,   1,   2,   3,   4,   5],
        [  6,   7,   8,   9,  10,  11],
        [ 12,  13,  14,  15,  16,  17],
        [ 18,  19,  20,  21,  22,  23],
        [ 24,  25,  26,  27,  28,  29],
        [ 30,  31,  32,  33,  34,  35]],

       [[ 36,  37,  38,  39,  40,  41],
        [ 42,  43,  44,  45,  46,  47],
        [ 48,  49,  50,  51,  52,  53],
        [ 54,  55,  56,  57,  58,  59],
        [ 60,  61,  62,  63,  64,  65],
        [ 66,  67,  68,  69,  70,  71]],

       [[ 72,  73,  74,  75,  76,  77],
        [ 78,  79,  80,  81,  82,  83],
        [ 84,  85,  86,  87,  88,  89],
        [ 90,  91,  92,  93,  94,  95],
        [ 96,  97,  98,  99, 100, 101],
        [102, 103, 104, 105, 106, 107]],

       [[108, 109, 110, 111, 112, 113],
        [114, 115, 116, 117, 118, 119],
        [120, 121, 122, 123, 124, 125],
        [126, 127, 128, 129, 130, 131],
        [132, 133, 134, 135, 136, 137],
        [138, 139, 140, 141, 142, 143]]])

Значения будут последовательно сохраняться на диске, как показано в строке 4 ниже. (Давайте пока проигнорируем детали файловой системы и фрагментацию.)

In [4]: arr.ravel(order='C')
Out[4]:
array([  0,   1,   2,   3,   4,   5,   6,   7,   8,   9,  10,  11,  12,
        13,  14,  15,  16,  17,  18,  19,  20,  21,  22,  23,  24,  25,
        26,  27,  28,  29,  30,  31,  32,  33,  34,  35,  36,  37,  38,
        39,  40,  41,  42,  43,  44,  45,  46,  47,  48,  49,  50,  51,
        52,  53,  54,  55,  56,  57,  58,  59,  60,  61,  62,  63,  64,
        65,  66,  67,  68,  69,  70,  71,  72,  73,  74,  75,  76,  77,
        78,  79,  80,  81,  82,  83,  84,  85,  86,  87,  88,  89,  90,
        91,  92,  93,  94,  95,  96,  97,  98,  99, 100, 101, 102, 103,
       104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116,
       117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129,
       130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143])

В лучшем случае возьмем срез по первой оси. Обратите внимание, что это только первые 36 значений массива. Это будет очень быстрое чтение! (один поиск, один чтение)

In [5]: arr[0,:,:]
Out[5]:
array([[ 0,  1,  2,  3,  4,  5],
       [ 6,  7,  8,  9, 10, 11],
       [12, 13, 14, 15, 16, 17],
       [18, 19, 20, 21, 22, 23],
       [24, 25, 26, 27, 28, 29],
       [30, 31, 32, 33, 34, 35]])

Точно так же следующий срез по первой оси будет просто следующими 36 значениями. Чтобы прочитать весь срез по этой оси, нам нужна только одна seekоперация. Если все, что мы собираемся читать, - это различные срезы вдоль этой оси, то это идеальная файловая структура.

Однако давайте рассмотрим наихудший сценарий: срез по последней оси.

In [6]: arr[:,:,0]
Out[6]:
array([[  0,   6,  12,  18,  24,  30],
       [ 36,  42,  48,  54,  60,  66],
       [ 72,  78,  84,  90,  96, 102],
       [108, 114, 120, 126, 132, 138]])

Чтобы прочитать этот фрагмент, нам нужно 36 поисков и 36 считываний, поскольку все значения разделены на диске. Ни один из них не смежный!

Это может показаться незначительным, но по мере того, как мы добираемся до все больших и больших массивов, количество и размер seekопераций быстро растет. Для большого (~ 10 Гб) 3D-массива, хранящегося таким образом и считываемого через memmap, чтение полного среза вдоль «худшей» оси может легко занять десятки минут, даже с современным оборудованием. В то же время срез по лучшей оси может занять меньше секунды. Для простоты я показываю только «полные» срезы вдоль одной оси, но то же самое происходит с произвольными срезами любого подмножества данных.

Между прочим, есть несколько форматов файлов, которые используют это преимущество и в основном хранят на диске три копии огромных трехмерных массивов: одну в порядке C, одну в порядке F и одну в промежуточном положении между ними. (Примером этого является формат D3D от Geoprobe, хотя я не уверен, что он где-либо задокументирован.) Кого волнует, если конечный размер файла составляет 4 ТБ, хранение стоит дешево! Самое безумное в этом заключается в том, что, поскольку основной вариант использования - извлечение одного суб-фрагмента в каждом направлении, чтение, которое вы хотите выполнить, выполняется очень и очень быстро. Работает очень хорошо!

Простой случай "по частям"

Скажем, мы храним 2x2x2 «фрагмента» 3D-массива как непрерывные блоки на диске. Другими словами, что-то вроде:

nx, ny, nz = arr.shape
slices = []
for i in range(0, nx, 2):
    for j in range(0, ny, 2):
        for k in range(0, nz, 2):
            slices.append((slice(i, i+2), slice(j, j+2), slice(k, k+2)))

chunked = np.hstack([arr[chunk].ravel() for chunk in slices])

Итак, данные на диске будут выглядеть так chunked:

array([  0,   1,   6,   7,  36,  37,  42,  43,   2,   3,   8,   9,  38,
        39,  44,  45,   4,   5,  10,  11,  40,  41,  46,  47,  12,  13,
        18,  19,  48,  49,  54,  55,  14,  15,  20,  21,  50,  51,  56,
        57,  16,  17,  22,  23,  52,  53,  58,  59,  24,  25,  30,  31,
        60,  61,  66,  67,  26,  27,  32,  33,  62,  63,  68,  69,  28,
        29,  34,  35,  64,  65,  70,  71,  72,  73,  78,  79, 108, 109,
       114, 115,  74,  75,  80,  81, 110, 111, 116, 117,  76,  77,  82,
        83, 112, 113, 118, 119,  84,  85,  90,  91, 120, 121, 126, 127,
        86,  87,  92,  93, 122, 123, 128, 129,  88,  89,  94,  95, 124,
       125, 130, 131,  96,  97, 102, 103, 132, 133, 138, 139,  98,  99,
       104, 105, 134, 135, 140, 141, 100, 101, 106, 107, 136, 137, 142, 143])

И чтобы показать, что это блоки размером 2x2x2 arr, обратите внимание, что это первые 8 значений chunked:

In [9]: arr[:2, :2, :2]
Out[9]:
array([[[ 0,  1],
        [ 6,  7]],

       [[36, 37],
        [42, 43]]])

Чтобы читать в любом срезе вдоль оси, мы читали либо 6, либо 9 смежных фрагментов (вдвое больше данных, чем нам нужно), а затем сохраняли только ту часть, которую хотели. Это наихудший случай: максимум 9 поисков по сравнению с максимум 36 поисками для версии без фрагментов. (Но в лучшем случае все равно 6 поисков против 1 для массива memmapped.) Поскольку последовательные чтения очень быстрые по сравнению с поисками, это значительно сокращает время, необходимое для чтения произвольного подмножества в память. Еще раз, этот эффект усиливается при увеличении массивов.

HDF5 делает еще несколько шагов вперед. Чанки не обязательно должны храниться непрерывно, и они индексируются B-деревом. Более того, они не обязательно должны быть одинакового размера на диске, поэтому сжатие может применяться к каждому фрагменту.

Разделенные массивы с h5py

По умолчанию h5pyне создает на диске фрагментированные файлы HDF (я думаю pytables, наоборот). Однако, если вы укажете chunks=Trueпри создании набора данных, вы получите на диске фрагментированный массив.

В качестве быстрого минимального примера:

import numpy as np
import h5py

data = np.random.random((100, 100, 100))

with h5py.File('test.hdf', 'w') as outfile:
    dset = outfile.create_dataset('a_descriptive_name', data=data, chunks=True)
    dset.attrs['some key'] = 'Did you want some metadata?'

Обратите внимание, что chunks=Trueговорит нам h5pyавтоматически выбрать размер блока. Если вы знаете больше о вашем наиболее распространенном варианте использования, вы можете оптимизировать размер / форму фрагмента, указав кортеж формы (например, (2,2,2)в простом примере выше). Это позволяет сделать чтение по определенной оси более эффективным или оптимизировать чтение / запись определенного размера.

Сравнение производительности ввода-вывода

Чтобы подчеркнуть эту мысль, давайте сравним чтение в срезах из разбитого на блоки набора данных HDF5 и большого (~ 8 ГБ) упорядоченного по Фортрану трехмерного массива, содержащего одни и те же точные данные.

Я очищал все кеши ОС между запусками, поэтому мы наблюдаем «холодную» производительность.

Для каждого типа файла мы протестируем чтение в «полном» x-срезе по первой оси и «полном» z-срезе по последней оси. Для массива memmapped, упорядоченного по Фортрану, срез «x» является наихудшим случаем, а срез «z» - наилучшим случаем.

Используемый код находится в сути (включая создание hdfфайла). Я не могу легко поделиться данными, которые здесь используются, но вы можете смоделировать их массивом нулей той же формы ( 621, 4991, 2600)и типа np.uint8.

В chunked_hdf.pyвыглядит следующим образом :

import sys
import h5py

def main():
    data = read()

    if sys.argv[1] == 'x':
        x_slice(data)
    elif sys.argv[1] == 'z':
        z_slice(data)

def read():
    f = h5py.File('/tmp/test.hdf5', 'r')
    return f['seismic_volume']

def z_slice(data):
    return data[:,:,0]

def x_slice(data):
    return data[0,:,:]

main()

memmapped_array.pyаналогичен, но имеет немного большую сложность, чтобы гарантировать, что срезы действительно загружаются в память (по умолчанию memmappedбудет возвращен другой массив, который не будет сравнением яблок с яблоками).

import numpy as np
import sys

def main():
    data = read()

    if sys.argv[1] == 'x':
        x_slice(data)
    elif sys.argv[1] == 'z':
        z_slice(data)

def read():
    big_binary_filename = '/data/nankai/data/Volumes/kumdep01_flipY.3dv.vol'
    shape = 621, 4991, 2600
    header_len = 3072

    data = np.memmap(filename=big_binary_filename, mode='r', offset=header_len,
                     order='F', shape=shape, dtype=np.uint8)
    return data

def z_slice(data):
    dat = np.empty(data.shape[:2], dtype=data.dtype)
    dat[:] = data[:,:,0]
    return dat

def x_slice(data):
    dat = np.empty(data.shape[1:], dtype=data.dtype)
    dat[:] = data[0,:,:]
    return dat

main()

Давайте сначала посмотрим на производительность HDF:

jofer at cornbread in ~ 
$ sudo ./clear_cache.sh

jofer at cornbread in ~ 
$ time python chunked_hdf.py z
python chunked_hdf.py z  0.64s user 0.28s system 3% cpu 23.800 total

jofer at cornbread in ~ 
$ sudo ./clear_cache.sh

jofer at cornbread in ~ 
$ time python chunked_hdf.py x
python chunked_hdf.py x  0.12s user 0.30s system 1% cpu 21.856 total

«Полный» x-срез и «полный» z-срез занимают примерно одинаковое количество времени (~ 20 секунд). Учитывая, что это массив на 8 ГБ, это неплохо. Большую часть времени

И если мы сравним это с временем отображения запомненного массива (он упорядочен в Фортране: «z-срез» - лучший случай, «x-срез» - худший случай):

jofer at cornbread in ~ 
$ sudo ./clear_cache.sh

jofer at cornbread in ~ 
$ time python memmapped_array.py z
python memmapped_array.py z  0.07s user 0.04s system 28% cpu 0.385 total

jofer at cornbread in ~ 
$ sudo ./clear_cache.sh

jofer at cornbread in ~ 
$ time python memmapped_array.py x
python memmapped_array.py x  2.46s user 37.24s system 0% cpu 3:35:26.85 total

Да, вы правильно прочитали. 0,3 секунды для одного направления и ~ 3,5 часа для другого.

Время нарезки в направлении «x» намного больше, чем время, необходимое для загрузки всего массива 8 ГБ в память и выбора нужного фрагмента! (Опять же, это массив, упорядоченный по Фортрану. Противоположная синхронизация среза x / z будет иметь место для массива, упорядоченного по C.)

Однако, если мы всегда хотим взять срез в лучшем случае, большой двоичный массив на диске очень хорош. (~ 0,3 сек!)

При использовании массива memmapped вы застряли с этим несоответствием ввода-вывода (или, возможно, анизотропия - лучший термин). Однако для разбитого на части набора данных HDF вы можете выбрать размер фрагмента так, чтобы доступ был либо равным, либо оптимизированным для конкретного варианта использования. Это дает вам больше гибкости.

В итоге

Надеюсь, это поможет прояснить одну часть вашего вопроса, во всяком случае. HDF5 имеет много других преимуществ по сравнению с "сырыми" мемокартами, но у меня нет возможности раскрыть их здесь все. Сжатие может ускорить некоторые вещи (данные, с которыми я работаю, не сильно выигрывают от сжатия, поэтому я редко использую его), а кэширование на уровне ОС часто лучше работает с файлами HDF5, чем с "необработанными" мемокартами. Помимо этого, HDF5 - действительно фантастический контейнерный формат. Он дает вам большую гибкость в управлении вашими данными и может использоваться более или менее на любом языке программирования.

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