[解決済み] 大容量アレイストレージにHDF5を使用することで、解析速度やメモリ使用量に利点がありますか（フラットバイナリファイルの代わりに）？

2022-10-25 16:53:01

質問

私は大規模な 3D 配列を処理しており、さまざまなデータ分析を行うために、しばしばさまざまな方法でスライスする必要があります。典型的なキューブは 100GB ほどになります (将来はもっと大きくなる予定です)。

Pythonで大規模なデータセットのために推奨される典型的なファイル形式は、HDF5（h5pyまたはpytablesのいずれか）を使用することであるようです。私の質問は、単純なフラットバイナリファイルにキューブを格納するよりも、HDF5 を使用して格納および分析することに、速度またはメモリ使用の利点があるかどうかということです。HDF5は、私が扱っているような大きな配列ではなく、表形式のデータに適しているのでしょうか? HDF5 が優れた圧縮を提供できることは知っていますが、処理速度とメモリのオーバーフローに対処することにもっと関心があります。

私は頻繁に、キューブの 1 つの大きなサブセットのみを分析したいと思います。pytables と h5py の両方の欠点は、私が配列のスライスを取るとき、私は常に numpy 配列を返し、メモリを使用することだと思われます。しかし、フラットバイナリファイルのnumpy memmapをスライスすると、データをディスク上に保持するビューを得ることができます。そのため、メモリをオーバーランすることなく、データの特定のセクターをより簡単に分析することができるようです。

私はpytablesとh5pyの両方を調査しましたが、私の目的のために今のところどちらの利点も見出せませんでした。

どのように解決するのですか？

HDF5の利点。組織化、柔軟性、相互運用性

HDF5の主な利点は、階層構造（フォルダ/ファイルに似ている）、各項目で保存される任意のメタデータ、および柔軟性（圧縮など）です。この組織構造とメタデータの保存は些細なことに聞こえるかもしれませんが、実際にはとても便利なものなのです。

HDFのもう一つの利点は、データセットが固定サイズの または のどちらかを選択できることです。したがって、全く新しいコピーを作成することなく、大きなデータセットにデータを追加することが簡単にできます。

さらに、HDF5は標準化されたフォーマットで、ほとんどすべての言語で利用可能なライブラリが用意されています。(公平に見て、C対Fの順序を意識し、格納された配列の形状、dtypeなどを知っている限り、大きなバイナリ配列でもそれほど難しくありません)。

大きな配列のためのHDFの利点。任意のスライスのより速いI/O

TL/DRと同じように 8GBの3D配列の場合、HDF5データセットをチャンクして、任意の軸に沿ったquot;full"のスライスを読み込むのに20秒かかり、0.3秒（ベストケース）かかっています。 3 時間以上 (最悪の場合) であった。

任意のスライス (任意であることを強調) の読み取りは、ディスク上のデータが平均してより連続的であるため、一般的にはるかに高速になります。

* (HDF5 はチャンクされたデータフォーマットである必要はありません。チャンキングをサポートしていますが、必須ではありません。実際、データセットを作成する際のデフォルトは h5py でのデータセット作成のデフォルトは、私の記憶が正しければ、チャンクしないことです)。

基本的に、データセットの特定のスライスに対する最良の場合のディスク読み取り速度と最悪の場合のディスク読み取り速度は、チャンクされたHDFデータセットでかなり近くなります（あなたが妥当なチャンクサイズを選択したか、ライブラリがあなたに代わって選択したと仮定しています）。単純なバイナリ配列の場合、ベストケースはより速くなりますが、ワーストケースは はるかに より悪いです。

1 つの注意点として、SSD を使用している場合、おそらく読み取り/書き込み速度に大きな違いは感じないでしょう。しかし、通常のハードドライブでは、シーケンシャル読み取りはランダム読み取りよりもはるかに高速です。 (つまり、通常のハードディスク・ドライブには、長い seek の時間がかかります)。HDF は SSD 上でもまだ優位性がありますが、それは生の速度よりも他の機能 (たとえば、メタデータ、整理など) によるものです。

まず最初に、混乱を避けるために、HDF にアクセスすることは h5py datasetはnumpyの配列にかなり似た振る舞いをするオブジェクトを返しますが、スライスされるまではデータをメモリにロードしません。 (memmapに似ていますが、同一ではありません。) h5py 導入をクリックすると詳細が表示されます。

データセットをスライスすると、データのサブセットがメモリにロードされますが、おそらくあなたはそれを使って何かをしたいのでしょう、その時にはとにかくメモリにあるデータが必要になります。

もしアウトオブコアの計算をしたいのであれば、表形式のデータに対して、かなり簡単に pandas あるいは pytables . で可能です。 h5py でも可能ですが (大きな N-D 配列ではより適切)、タッチダウンして自分で反復を処理する必要があります。

しかし、numpy的なアウトオブコア計算の未来はBlazeです。見てみてくださいを見てみてください。

チャンクされないケース

まず最初に、ディスクに書き込まれた3D C-ordered配列を考えてみましょう（ここでは 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個の値だけになります。この軸に沿った完全なスライスを読み取るには、1つの 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 操作の数とサイズは急速に大きくなります。この方法で保存された大規模な (~10Gb) 3D 配列を memmap で読み込む場合、quot; worst" 軸に沿った完全なスライスの読み取りには、最新のハードウェアでも数十分かかることがあります。同時に、最適な軸に沿ったスライスは 1 秒もかかりません。簡単のために、1 つの軸に沿った完全なスライスのみを示していますが、データの任意のサブセットの任意のスライスでもまったく同じことが起こります。

ちなみに、これを利用したファイルフォーマットがいくつかあり、基本的には 巨大な 1 つは C 次、もう 1 つは F 次、そしてもう 1 つはその中間の次です。 (この例は Geoprobe の D3D フォーマットですが、どこにも文書化されていないようです)。最終的なファイルサイズが4TBになろうが関係ない、ストレージは安いのだから! これの面白いところは、主な用途が各方向に1つのサブスライスを抽出することなので、必要な読み取りがとてもとても速いということです。非常によく機能します!

単純な "chunked" の場合

3D配列の2x2x2個のquot;チャンクを連続したブロックとしてディスクに保存するとします。言い換えれば、次のようなものです。

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 個の連続したチャンク (必要なデータの 2 倍) を読み込んで、必要な部分のみを保持することになります。この場合、最悪の場合、最大9シークとなり、非チャンク化バージョンの最大36シークとなります。(シーケンシャルリードはシークに比べて非常に高速であるため、任意のサブセットをメモリに読み込むのにかかる時間を大幅に短縮することができます。繰り返しになりますが、この効果はアレイが大きくなるほど大きくなります。

HDF5では、これをさらに数歩進めています。チャンクは連続的に保存される必要はなく、B-Treeによってインデックスが作成されます。さらに、ディスク上で同じサイズである必要はないため、各チャンクに圧縮を適用することができます。

チャンクされたアレイに `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) のように）。これにより、特定の軸に沿った読み取りをより効率的にしたり、特定のサイズの読み取り/書き込みに対して最適化したりすることができます。

I/O パフォーマンスの比較

ポイントを強調するために、チャンクされた HDF5 データセットからのスライスの読み込みと、同じ正確なデータを含む大きな (~8GB) Fortran 順序の 3D 配列を比較してみましょう。

私はすべての OS キャッシュをクリアしました。そのため、コールドパフォーマンスを見ることができます。

各ファイルタイプについて、最初の軸に沿った "full" x スライスと最後の軸に沿った "full" z スライスの読み取りをテストしてみます。 Fortran の順序付き memmapped 配列では、"x" スライスが最悪のケースであり、"z" スライスが最良のケースとなります。

使用されるコードは gistにある (を作成することを含む）。 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 秒) を要します。これが 8GB のアレイであることを考えると、それほど悪くはないでしょう。ほとんどの場合

また、これを memmapped 配列の時間と比較すると (これは Fortran 順序です。A "z-slice" がベストケースで、"x-slice" がワーストケースです)。

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

はい、その通りです。1つのスライス方向に対して0.3秒、そして～3.5 時間です。

x"方向にスライスする時間は はるか 方向にスライスする時間は、8GB の配列全体をメモリにロードして、必要なスライスを選択する時間よりもはるかに長いのです! (繰り返しますが、これは Fortran の順序付き配列です。 C 言語による順序付き配列では、x/z スライスのタイミングは逆になります)。

しかし、常にベストケース方向に沿ってスライスを取りたいのであれば、ディスク上の大きなバイナリ配列は非常に優れています。(～0.3秒！)

memmapped配列では、このI/Oの不一致（あるいは異方性と言った方がいいかもしれません）に悩まされます。しかし、チャンクされたHDFデータセットでは、アクセスが同等になるか、特定のユースケースに最適化されるようにチャンクサイズを選択することができます。これによって、より多くの柔軟性が得られます。

要約すると

何はともあれ、質問の一部分をクリアするのに役立つことを願っています。 HDF5 には "raw" memmaps に対して他にも多くの利点がありますが、ここでそのすべてについて説明する余裕はありません。また、OSレベルのキャッシュは、"raw" memmapsよりもHDF5ファイルの方がうまく機能することが多いのです。それ以上に、HDF5は本当に素晴らしいコンテナフォーマットです。HDF5はデータ管理に多くの柔軟性を与えてくれますし、多かれ少なかれどんなプログラミング言語からでも使用することができます。

全体として、あなたのユースケースでうまく機能するかどうか試してみてください。きっと驚かれることでしょう。