Kraken 2は、k-merの完全一致を利用したシークエンシングリードの分類学的プロファイリングツールで、メタゲノムやメタアンプリコンの分類や汚染のチエックなどに幅広く使用されている。データベースは自分で作成することもできるが(custom databaseの項を参照)、おそそらくは公開されているスタンダードデータベースがよく使われている(standard databseの項を参照)。このスタンダードデータベースはバクテリア、古細菌、ウイルスドメインのRefSeqの全ゲノム、ヒトゲノム、既知のベクターのコレクションなどが含まれている(構築後のサイズは2022年8月でhash.k2dは57GBくらい)。また、より大きなGTDB互換kraken2データベースも有志の手によって公開されている(構築後のhash.k2dのサイズは308GBくらい)。kraken2の実行時は、このデータベースがRAMにロードされる。スタンダードデータベースではこのデータベースのロードにより50GB以上のメモリを必要とする(v1よりかなり抑えられている)。kraken2 GTDB-D.Bだと6倍ほど多くメモリを使う。
kraken2を大量のサンプルに対して実行する場合、このD.Bのロード時間が一番のボトルネックとなる。例えば自分の計算環境ではシーケンシャルリード6500MB/sくらいのエンタープライズ向けPCI接続SSDにスタンダードデータベースを配置してD.Bのロードを行っているが、それでもおよそ30秒かかる。それに対して、ロード後の100万リードの分類は素晴らしく高速で5秒程度で終わる。SATA接続SSDだとデータベースのロードにはより長い時間がかかり、HDDから読み出すとだといつまで経っても終わらないような感覚になる。サンプル数が増えるにつれて、データベースのロード時間が全ランタイムに対して支配的になる(1,000サンプルではD.Bロードに30,000秒(8.3-h)、分析に5,000秒(1.4-h))*1。。データベースロード時間の問題を抜本的に解決するには、複数サンプルを扱う時もデータベースのロードが1回であることが望ましい (すなわち、1,000サンプルではD.Bロードに30秒、分析に5,000秒)。しかし2022年12月現在では本家kraken2は複数サンプル指定に対応しておらず、ワイルドカードでまとめて指定してもスキップされる。
そこで、ここではマルチサンプルに対応したkraken2のフォークを試す。これを使う事で、データベースのロードは一度のみで大量のfastqの分類学的プロファイリングを行うことができる。
インストール
kraken2.0.9-betaを試した。
git clone https://github.com/daydream-boost/kraken2.git
cd kraken2/
./install_kraken2.sh <path>/<to>/pre-build_kraken_database
実行方法
各サンプルのシークエンシングリードを順番に指定する。フルパスでの指定には対応していないので相対パスで指定する。拡張子はfqやfastqなどに対応している。圧縮している場合、解凍してから使用する(本家kraken2はgzipped fastqに対応している)。サンプルが多いときはワイルドカードと組み合わせて指定する。
mkdir out_dir
kraken2 --paired --threads 20 -db kraken_db --confidence 0.5 --output out_dir/ cleandata/sample1_clean_R1.fq cleandata/sample1_clean_R2.fq ... cleandata/samplen_clean_R1.fq cleandata/samplen_clean_R2.fq
#ワイルドカード指定
kraken2 --paired --threads 20 -db kraken_db --confidence 0.5 --output out_dir/ fastq_dir/sample*.fastq
目的の分類群が存在するかしないか調査のため、ゼロカウントも全て出力。
mkdir kraken2_taxonomic_profiling
kraken2 --use-names --report-zero-counts --paired --threads 20 -db <path>/<to>/kraken2-full-database/ --output kraken2_taxonomic_profiling/ --report report.txt fastq_dir/SRR*_{1,2}.fastq
kraken2_taxonomic_profiling/
カレントパスにはreport.txtが出力される。
引用
The second version of the Kraken taxonomic sequence classification system
https://github.com/DerrickWood/kraken2
参考
https://github.com/DerrickWood/kraken2/issues/87
*1
最新のEPYC計算機構成でRAM diskを試す事も思いつくが、いずれにせよRAM diskは共用の計算機環境では使えない