質問:
非常に大きな行ベースのファイルを効率的にサブセット化するにはどうすればよいですか?
Konrad Rudolph
2017-06-05 17:49:25 UTC
view on stackexchange narkive permalink

これは最近繰り返し発生しています。非常に大きなテキストファイル(数GiBのオーダー)があり、約10,000行の行ベースのサブセット化を実行する必要があります。特定のシナリオ(たとえば、BAMファイルをランダムにサンプリングするための samtools view -s )の解決策はありますが、私のユースケースがこれらのカテゴリに当てはまらない場合があります。

残念ながらナイーブ sed ベースのソリューションは非常に遅い:

  time sed -n -f <(awk -vOFS = '' '{print $ 0、 "p"} 'line_numbers.txt)input_file > selected_lines.txt  

ここで、 line_numbers.txt は、1行に1つの行番号を含むファイルです。

これを10,000行実行するのを忘れてください。すでにわずか1000で停止しています。

これを高速化するにはどうすればよいですか。理想的には、入力ファイルのサイズに合わせてスケーリングし、実行時間はほぼ一定です。サブセット化する行?

この質問は、特定のバイオインフォマティクスファイル形式ではなく、一般的なケースに対処する必要があるため、ここよりもスタックオーバーフローに属していると思います。投票して閉鎖しようとしましたが、どういうわけか、より適切なサイトとしてstackoverflowを選択する選択肢がありません。
@bli私は同意しません。私たちは、バイオインフォマティシャンの仕事の過程で発生する可能性のあるあらゆる質問を受け入れるべきだと思います。多くのバイオインフォマティクスの質問は、単純なテキスト解析操作(たとえば、異なるシーケンス形式間の翻訳を検討してください)に要約できますが、それでもそれらはトピックにあるべきだと感じています。また、特定の移行パスが設定されていて、まだそれらがない場合を除いて、別のサイトに移行することはできません。
@bli私も同意しません。StackOverflowの質問はほとんどの場合言語に関連付けられているためです(またはその逆で、何にも関連付けられておらず、コードソリューションは必要ありません)。対照的に、私は解決策に興味がありますが、テクノロジーについては気にしません。
この質問にもう少しバイオインフォマティクスのコンテキスト/ストーリーを含めると役に立ちます。ユースケースではないものの例を示しましたが、そうではありません。行*番号*に基づいてサブセット化する必要がある状況について頭を悩ませることはできません。いずれかの列の値に基づいてサブセット化するのが一般的です。
awkのみの場合: `awk'BEGIN {while((getline <" line_num.txt ")> 0)l [$ 1] = 1} NR in l'input_file`。
SOに関する関連記事、awkを使用して、[ファイルを行番号と列番号でサブセット化](https://stackoverflow.com/q/40842008/680068)。これは頻繁に出てくることに同意し、この投稿はバイオインフォマティクスのトピックになっていると思います。
トップアンサーのベンチマーク結果を見るといいでしょう。
ファイル全体をメモリに読み込むことを気にしないのであれば、Pythonではこれをすばやく行う必要があります `operator.itemgetter(* [int(line)for line in open(line_file)])(open(input_file).readlines())`
@Chris_Randsジェネレーターで `itertools.islice`を使用してみませんか?これにより、ファイル全体を一度にロードする必要がなくなり、オーバーヘッドが最小限に抑えられます。
@KonradRudolphは、ジェネレーターが消費されるときに相対的な行番号を使用して、行番号の各連続範囲を個別に「islice」するような意味ですか?
@Chris_Randsくそー、私はRを使いすぎています。連続していないインデックス範囲を `islice`に渡すことができると思いました。
四 答え:
#1
+7
Konrad Rudolph
2017-06-05 17:49:25 UTC
view on stackexchange narkive permalink

結局、(サンプル行番号をソートした後)次の候補行を追跡するだけでパフォーマンスの問題が修正され、残りの速度低下のほとんどは実際にファイルを読み取るオーバーヘッドが原因であると思われるため、それほど多くはありません改善するために。

sed でこれを行う方法がわからず、 awk でも簡単ではないため、ここにPerlスクリプトを示します。 :

 #!/ usr / bin / env perluse strict; use warnings; my $ file = $ ARGV [0]; my $ lines_file = $ ARGV [1]; open my $ lines_fh、 '<'、$ lines_file or die "Cannot read file $ lines_file"; chomp(my @lines = < $ lines_fh>); close $ lines_fh; @lines = sort {$ a < = > $ b} @lines; open my $ fh、 '<'、$ file or die "Cannot read file $ file"; my $ line = 1; my $ next_line = 0; while(< $ fh>){last if $ next_line ==スカラー@lines; if($ line ++ == $ lines [$ next_line]){$ next_line ++;印刷; }} close $ fh;  

Rパッケージ用に C ++で同様の関数を実装しました。これはわずかに長いだけです。 Perlスクリプト。テストファイルのPerlスクリプトよりも約3倍高速です。

代わりにハッシュを使用した場合、これははるかに高速ではないでしょうか? `$ line {$ _} ++ for @lines`そして、` while`ループで: `print if $ line {$。}`。
@terdonいいえ。ハッシュルックアップ時間は平均して一定ですが、連続する配列での直接インデックスルックアップよりも遅くなります(さらに時折インデックスが増加します)。理論についてはこれだけです。私もそれをテストしました、そしてそれは実際に成り立ちます。 (そして、予想通り、実際の違いはほとんど無視できます。)
まあ、TIL。おかげで、私は常にハッシュが定義上より速いと思っていました(生物学者にコードを書かせるとそれが得られます:)。
#2
+5
gringer
2017-06-06 00:54:57 UTC
view on stackexchange narkive permalink

行のリストを格納するためにハッシュセットを使用する場合、Perlはこれでかなり高速になるはずです。このような構造は、フィールド値に基づくサブセット化にも機能します。比較は、「$」ではなくフィールドと行われます。

 #!/ usr / bin / perluse strict; use警告; my $ lines_file = $ ARGV [0]; my%include_lines =(); open my $ lines_fh、 '<'、$ lines_file or die "Cannot read file $ lines_file"; while(< $ lines_fh>){chomp; $ include_lines {$ _} = 1;} close $ lines_fh; while(<>){if($ include_lines {$。}){# "$。" -現在のファイル印刷の行番号。 }}  

このSOの回答によると、「$」に注意してください。演算子は厳密には現在の行番号ではなく、さまざまなファイル操作やその他の設定の影響を受ける可能性があります。

編集:ハッシュセットを並べ替えられたリストと比較して、回答の速度に関するコメントを見ただけです。 $ lines [$ next_line] ビットは私には少し奇妙に感じます。ソートされたリストで shift または pop を使用して、次の行をフェッチしてみましたか:

 #!/ usr / bin / perluse strict;警告を使用; my $ lines_file = $ ARGV [0]; open my $ lines_fh、 '<'、$ lines_file or die "Cannot read file $ lines_file"; chomp(my @lines = < $ lines_fh>); close $ lines_fh ; @lines = sort {$ a < = > $ b} @lines; my $ next_line = shift(@lines); while(<>){if($。== $ next_line){$ next_line = shift(@lines) ;印刷;最後の場合(!@lines); }}  

shift pop に変更し(一致するように並べ替え順序を逆にします)、Konradの元のコードで4.8秒、5.8秒、4.1秒の時間を取得しました。ハッシュコードとポップコードはそれぞれ、 / usr / share / dict / british-english-insane から10,000行を25回フェッチします(入力ファイルを / tmp codeにコピーした後) >)。これらはすべて私の観点からは同じです。コマンドを実行するよりも、コマンドを入力するのに時間がかかるほどの速さです。 pop の代わりに shift を使用しても、時間が目立って変わることはないようです。

`shift` /` pop`は配列を変更しますが、結果としてメモリが移動すると非常に遅くなる可能性があります(ただし、Perlがそれを行うかどうかはわかりません)。インデックス付きアクセスよりも確実に高速であってはなりません。正直なところ、コードは少なくなります。
申し訳ありません。私はこれらをテストせずにすばやく書き上げました。これらのバグを修正し(12行目で%を$に変更)、18行目の先頭に「my」を配置しました。
#3
+4
bli
2017-06-05 19:17:42 UTC
view on stackexchange narkive permalink

いくつかの関連する質問が他のサイトに表示され、興味深い解決策が考えられます。ここで報告します:

空でない行の約1%をサンプリングするには:

  awk'BEGIN {srand()}!/ ^ $ / {if(rand()< = .01)print $ 0} 'input_file  

(from https:// stackoverflow .com / a / 692321/1878788

1000行のランダム行を選択するには:

  shuf -n 1000 input_file  

https://stackoverflow.com/a/15065490/1878788、および https://unix.stackexchange.com/a/108604/55127から)

編集:行のリストを使用したPythonソリューション

行インデックスのセットを使用し、セットのメンバーシップをテストして行を選択します:

 #!/ usr / bin / env python3import syswith open(sys.argv [2]、 "r")as line_numbers_file:line_indices = set(int(line)-1 for line in line_numbers_file)with open(sys.argv [1]、 "r" )as input_file:print(*(line.strip()for(idx、line)in enumerate(input_file)if idx in line_indices)、sep = "\ n") 

numpyブール配列を itertools.compress と一緒に使用する:

 #!/ usr / bin / env python3import sysfrom itertools import compressfrom numpy import zeroswith open( sys.argv [2]、 "r")as line_numbers_file:line_indices = [int(line)-1 for line in line_numbers_file] selector = zeros(max(line_indices)+ 1、dtype = bool)selector [line_indices] = 1with open (sys.argv [1]、 "r")as input_file:print(*(line.strip()for line in compress(input_file、selector))、sep = "\ n") 

15774756samレコードと10000個の事前生成されたランダム行番号のリストを含むファイルでいくつかのテストを行いました。

KonradRudolph( https:// bioinformatics)によって提案されたperlスクリプト。 stackexchange.com/a/454/292)は約5.3秒で実行されます。

セットメンバーシップテストのPythonソリューションは約4.45秒で実行されます。

圧縮ベースのソリューション約3.4秒で実行されます。 反復回数はブール配列の長さに依存するため、これは必要な最大の行番号によって大きく異なる可能性があると思います。ここでの最大の行番号は15773768で、合計行数に比べてかなり多いです。

Python3.6で試してみました。 Python 2.7の方が少し速いかもしれないと思いますが、テストはしていません。

問題は、ランダムなサンプルだけでなく、ファイルから特定の行のセットを選択することだと思います。
@terdonは正しいです。最近遭遇したユーザーのケースでは、特定のランダムでない行が必要でした。それにもかかわらず、これは良い補足です。
@KonradRudolph `samtools view -s`についておっしゃっていたので、残りの質問はランダムに選択することを考えていて、行番号のリストはランダムに生成されたと思いました。
@bli私の意図は、これを反例として、問題をどのように解決できないかを示すことでした。
`int(line.strip())`は `int(line)`に簡略化できます。 `int`はとにかく先頭/末尾の空白を無視します
@Chris_Rands私はそれを知りませんでした。ありがとう。
#4
+1
Alex Reynolds
2017-06-10 01:27:11 UTC
view on stackexchange narkive permalink

Githubにある subset というコマンドライン(C ++ 14)ツールを作成しました: https://github.com/alexpreynolds/subset

これは、適度にメモリ効率が高く、高速である必要があります。 subset ツールは、入力行をテーブルに格納しませんが、代わりにファイルを1回ストリーミングし、入力ファイルの4または8kのバッファーチャンクを格納します(OSによって異なります)。

行番号を配列に格納しますが、整数あたり8バイト* 100kは800kBであり、その使用例ではメモリはそれほど多くありません。

O(nlogn)ソートがあります行番号配列にペナルティがありますが、このリストはクエリファイルよりもはるかに小さく、整数の並べ替えはかなり最適化されているため、ヒットは小さいはずです。

行番号リストがすでに並べ替えられている場合は、並べ替えをスキップするオプションを追加できます。それが役立つかどうか教えてください。

フィルタリングステップでは、行番号配列と入力ファイルを直線的に調べ、インデックスが一致する行を印刷し、残りをスキップします。

実際、クエリする行番号がなくなると、 subset は入力ファイルの解析の早い段階で終了します。したがって、この機能は、非常に大きなクエリファイルのフィルタリングを高速化する場合に特に役立ちます。 (たとえば、クエリファイルに100万行があり、対象の最後の行番号が12345の場合、ファイルの残りの部分を読み取る理由はありません。)

次のように取得、ビルド、インストールできます。だから:

  $ git clone https://github.com/alexpreynolds/subset.git$ cdサブセット$ make $ cpサブセット/ usr / local / bin  

バイナリがパスに含まれると、それを使用する方法がいくつかあります。

たとえば、開始インデックスと長さの値を指定できます。以下は、33行目から始まる7行を取得します(0インデックス値として32):

  $サブセット--prefix-with-indices-s 32 -n 7 -i query.txt > answer.txt  

または、それぞれが別々の行にある行番号を含むテキストファイルを指定することもできます。以下は、 line-numbers.txt というファイルを読み込み、それを使用して query.txt をフィルタリングします。

  $サブセット--prefix -with-indices -l line-numbers.txt -i query.txt > answer.txt  

line-numbers.txt のインデックスは正である必要があり、 0-インデックス付きの整数。 subset が番号のリストを並べ替えるので、番号のリストを並べ替える必要はありません。これは、入力/クエリファイルの効率的なシングルパスを実行できるようにするためです。

-prefix-with-indices を省略して、デバッグプレフィックスを省略できます。これは、結果のサニティチェックを実行できるようにするためのものです。

test / makefile テストは、2種類のフィルタリングのオプションと使用法を示しています。



このQ&Aは英語から自動的に翻訳されました。オリジナルのコンテンツはstackexchangeで入手できます。これは、配布されているcc by-sa 3.0ライセンスに感謝します。
Loading...