17

問題の説明

Python/numpy でモンテカルロ粒子シミュレーター (ブラウン運動と光子放出) を作成中。シミュレーション出力 (>>10GB) をファイルに保存し、2 番目のステップでデータを処理する必要があります。Windows と Linux の両方との互換性は重要です。

粒子の数 ( n_particles) は 10 ～ 100 です。時間ステップ数 ( time_size) は ~10^9 です。

シミュレーションには 3 つのステップがあります (以下のコードはオールイン RAM バージョン用です)。

  1. rate 配列をシミュレート (および保存) しemissionます (ほぼ 0 の要素を多数含む):

    • 形状 ( n_particlesx time_size)、float32、サイズ80GB

  2. 配列を計算countsします (以前に計算されたレートを持つポアソン過程からのランダム値):

    • 形状 ( n_particlesx time_size)、uint8、サイズ20GB

      counts = np.random.poisson(lam=emission).astype(np.uint8)

  3. カウントのタイムスタンプ (またはインデックス) を見つけます。カウントはほとんどの場合 0 であるため、タイムスタンプ配列は RAM に収まります。

    # Loop across the particles
timestamps = [np.nonzero(c) for c in counts]

ステップ 1 を 1 回実行してから、ステップ 2 ～ 3 を何度も (~100 回) 繰り返します。将来的には、計算する前に前処理emission(適用cumsumまたは他の関数)が必要になるかもしれませんcounts。

質問

私は動作中のメモリ内実装を持っており、(はるかに) 長いシミュレーションに拡張できるコア外バージョンを実装するための最良のアプローチは何かを理解しようとしています。

存在してほしいもの

配列をファイルに保存する必要があり、シミュレーションに 1 つのファイルを使用したいと考えています。また、シミュレーション パラメーター (スカラー) のディクショナリを格納および呼び出す "簡単な" 方法も必要です。

理想的には、チャンクを事前に割り当てて埋めることができる、ファイルでバックアップされた numpy 配列が必要です。次に、numpy 配列メソッド ( max、cumsum、...) を透過的に動作させchunksize、各反復でロードする配列の量を指定するキーワードのみを必要とします。

さらに良いことに、キャッシュと RAM の間ではなく、RAM とハード ドライブの間で動作するNumexprが必要です。

実用的なオプションは何ですか

最初のオプションとして、pyTables の実験を開始しましたが、その複雑さと抽象化に満足していません (numpy とは大きく異なります)。さらに、私の現在のソリューション（以下を参照）は醜く、あまり効率的ではありません。

だから私が答えを求める私の選択肢は

  1. 必要な機能を備えた numpy 配列を実装する (どのように?)

  2. よりスマートな方法でpytableを使用する（異なるデータ構造/メソッド）

  3. h5py、blaze、pandas などの別のライブラリを使用します (これまで試したことはありません)。

暫定的な解決策 (pyTables)

シミュレーション パラメーターを'/parameters'グループに保存します。各パラメーターは numpy 配列スカラーに変換されます。詳細なソリューションですが、機能します。

データをチャンクで生成し、新しいチャンクごとに追加する必要があるため (ただし、最終的なサイズはわかっています)、emission拡張可能な配列 ( ) として保存します。EArray貯蓄countsはもっと問題です。pytable 配列のように保存すると、「counts >= 2」のようなクエリを実行するのが難しくなります。したがって、カウントを複数のテーブル (パーティクルごとに 1 つ) [UGLY] として保存し、.get_where_list('counts >= 2'). これがスペース効率に優れているかどうかはわかりません。単一の配列を使用する代わりにこれらすべてのテーブルを生成すると、HDF5 ファイルが大幅に破壊されます。さらに、奇妙なことに、これらのテーブルを作成するには、カスタム dtype を作成する必要があります (標準の numpy dtype の場合でも):

    dt = np.dtype([('counts', 'u1')])        
    for ip in xrange(n_particles):
        name = "particle_%d" % ip
        data_file.create_table(
                    group, name, description=dt, chunkshape=chunksize,
                    expectedrows=time_size,
                    title='Binned timetrace of emitted ph (bin = t_step)'
                        ' - particle_%d' % particle)

各粒子カウント「テーブル」には異なる名前 ( name = "particle_%d" % ip) があり、簡単に反復できるようにそれらを Python リストに入れる必要があります。

編集: この質問の結果は、PyBroMoと呼ばれるブラウン運動シミュレーターです。

4

5 に答える 5

3

Dask.array は、PyTables や h5py などのディスク上の配列でmax、 、 などのチャンク操作を実行できます。cumsum

import h5py
d = h5py.File('myfile.hdf5')['/data']
import dask.array as da
x = da.from_array(d, chunks=(1000, 1000))

X は numpy 配列のように見え、API の多くをコピーします。x に対する操作は、必要に応じてディスクからストリーミングする複数のコアを使用して効率的に実行されるメモリ内操作の DAG を作成します。

da.exp(x).mean(axis=0).compute()

http://dask.pydata.org/en/latest/

conda install dask
or 
pip install dask
于 2015-05-01T19:16:00.607 に答える
2

PyTable ソリューション

Pandas が提供する機能は不要であり、処理は非常に遅くなるため (以下のノートブックを参照)、最善の方法は PyTables または h5py を直接使用することです。これまでのところ、pytables アプローチのみを試しました。

すべてのテストは、このノートブックで実行されました。

pytables データ構造の紹介

参照: PyTables の公式ドキュメント

Pytables では、配列とテーブルの 2 種類の形式で HDF5 ファイルにデータを格納できます。

配列

Array配列には 、 、の 3 種類がありCArrayますEArray。それらはすべて、numpy スライスに似た表記法で (多次元) スライスを格納および取得できます。

# Write data to store (broadcasting works)
array1[:]  = 3

# Read data from store
in_ram_array = array1[:]

一部のユースケースでの最適化のために、作成時にCArrayサイズを選択できる「チャンク」に保存されます。chunk_shape

ArrayCArrayサイズは作成時に固定されます。ただし、作成後に配列をチャンクごとに塗りつぶす/書き込むことができます。逆にメソッドEArrayで拡張できます。.append()

テーブル

はtableまったく異なる獣です。それは基本的に「テーブル」です。1D インデックスのみがあり、各要素は行です。各行内には「列」データ型があり、各列は異なる型を持つことができます。numpy record-arraysに精通している場合、テーブルは基本的に 1D レコード配列であり、各要素は列として多くのフィールドを持ちます。

1 次元または 2 次元の numpy 配列をテーブルに格納できますが、もう少し注意が必要です。行データ型を作成する必要があります。たとえば、1D uint8 numpy 配列を格納するには、次のようにする必要があります。

table_uint8 = np.dtype([('field1', 'u1')])
table_1d = data_file.create_table('/', 'array_1d', description=table_uint8)

では、なぜテーブルを使用するのでしょうか。配列とは異なり、テーブルは効率的にクエリできるためです。たとえば、巨大なディスクベースのテーブルで要素 > 3 を検索したい場合は、次のように実行できます。

index = table_1d.get_where_list('field1 > 3')

シンプルであるだけでなく (ファイル全体をチャンクでスキャンindexしてループで構築する必要がある配列と比較して)、非常に高速でもあります。

シミュレーション パラメータの保存方法

シミュレーション パラメータを格納する最良の方法は、グループ (つまり/parameters) を使用し、各スカラーを numpy 配列に変換して として格納することCArrayです。

emission" "の配列

emission生成され、順次読み取られる最大の配列です。この使用パターンの場合、適切なデータ構造はEArrayです。~50% のゼロ要素を持つ「シミュレートされた」データでは、blosc圧縮 ( level=5) は 2.2x の圧縮率を達成します。2^18 (256k) のチャンク サイズが最小の処理時間であることがわかりました。

counts「 」を格納しています

" " も保存countsすると、ファイル サイズが 10% 増加し、タイムスタンプの計算に 40% の時間がかかります。最終counts的に必要なのはタイムスタンプだけなので、保存すること自体は利点ではありません。

利点は、1 つのコマンド ( ) で完全な時間軸をクエリするため、インデックス (タイムスタンプ) の再構築がより簡単になること.get_where_list('counts >= 1')です。逆に、チャンク処理では、少しトリッキーないくつかのインデックス演算を実行する必要があり、維持するのが負担になる可能性があります。

ただし、コードの複雑さは、両方の場合に必要な他のすべての操作 (並べ替えとマージ) に比べて小さい場合があります。

timestamps「 」を格納しています

タイムスタンプは RAM に蓄積できます。ただし、開始前に配列のサイズがわからないためhstack()、リストに格納されているさまざまなチャンクを「マージ」するために最後の呼び出しが必要です。これによりメモリ要件が 2 倍になるため、RAM が不足する可能性があります。

を使用して、進行中のタイムスタンプをテーブルに保存できます.append()。最後に、 を使用してテーブルをメモリにロードできます.read()。これはオールインメモリ計算よりも 10% 遅いだけですが、「ダブル RAM」の要件を回避できます。さらに、最終的な全ロードを回避し、RAM の使用量を最小限に抑えることができます。

H5Py

H5pyはpytablesよりもはるかに単純なライブラリです。この (主に) 順次処理のユースケースでは、pytables よりも適しているようです。唯一欠けている機能は、'blosc' 圧縮の欠如です。これによりパフォーマンスが大幅に低下する場合は、まだテストが必要です。

于 2014-01-09T00:25:47.847 に答える