2024/08/26 誤字修正
Kraken2はシークエンシングリードから菌叢解析を行うプログラムです。事前に作成されたRefSeq完全長ゲノムのDBなどを使って、シークエンシングリードの超高速な分類学的分類を実行できます。
kraken2はkraken1と比べてメモリ使用率が大きく減っていますが、それでもメモリ使用量は多い傾向があります。この理由は、高速にプロファイリングするため、krakenがデフォルト設定ではk-mer(minimizer)データベースをメモリにロードしてk-mer一致を探索するためです。例えば2024年初頭の時点のスタンダードデータベースを使うと、60GB近いRAM使用量になります。
どんな菌叢かを調べる目的ならstandardやnrのような巨大なデータベースを使うしかないのですが、研究の目的によっては関心がある分類が既にあって、その分類の株を読んだリードだけ分類できればよいこともあります。そのような場合、カスタムデータベースを作成することでメモリ使用量を大きく減らすことができます。以前、GTDBのカスタムデータベース作成手順を公式の説明に従って簡単に説明しました。ヘルパースクリプトの助けを借りて、該当するゲノムのfastaをダウンロードしてデータベースを作成する流れです。
この記事と同じGTDB tax.ベースのカスタムデータベースを作るなら、GTDBサーバーのgtdb_taxonomy.tsvファイルからその分類群を含む行だけgrepで抽出してデータベースを作成します(*1)。
自分の場合、シアノバクテリア門にだけ関心があったので、GTDB R220に含まれるゲノムのp__cyanobacteriaを含む行だけ取り出してGTDB tax.ベースのkrakenデータベースを作成しました(半日かかった)。ゲノムには2200個ほどゲノムが含まれていましたが、このデータベースを使用すると、kraken2のメモリ使用量は0.8GBほどに留まりました。データベースロード時間も短縮され、データベースロードタイムはrefseq DBが20秒以上かかったのに対して1.5秒程度で完了しました(PCI-E SSD使用)。
kraken2は非常に高速なため、データベースを短時間でメモリに読み込めれば10秒以内に数十万リードのプロファイリングを完了でき、したがって1,000個程度の小さなfastqなら数時間でプロファイリングすることができます。メモリ使用量が少ないので、スレッド使用数を控えめにしてkraken2を並列実行すればさらに短時間に多数のfastqのプロファイリングができるかもしれません。
このように、関心のある分類群が決まっているなら、その分類に限定したデータベースを使うことでkrakenをさらに高速実行できます。一点注意点として、複数の分類にまたがって完全一致する遺伝子やゲノム領域が存在する場合、kraken2の分類精度が低下する可能性もあります。全ての結果を信じる前に、どこかの時点で正しい結果かどうかの検証はするべきかなと思います。
引用
Improved metagenomic analysis with Kraken 2
Derrick E. Wood, Jennifer Lu & Ben Langmead
Genome Biology volume 20, Article number: 257 (2019)
参考
Kraken2の事前ビルド済みのデータべースは、https://benlangmead.github.io/aws-indexes/k2 からダウンロードできます。
Library列からは使用されているゲノム情報をダウンロードでき、Inspect列からはkrakenのinspectコマンドで閲覧できるデータベースの内容を示したテキストをダウンロードできる。このテキストには、データベースのミニマイザーのうち、さまざまな分類群/クレードにマップされたものの数が示されている(kraken2は羊土社から出版されているメタゲノム本の年度更新2年目の1つの項でも説明してます)。
コメント
SRAから少数のリードずつだけダウンロードしてkraken2で分析することで、公開されているシークエンシングデータの中で、どのようなメタデータのサンプルに関心のある菌が多く存在しているのか調査できる。メタゲノムのshallow metagenomic approach(わずか500リードでも主要な菌叢は分かるという論文が出ている)をデジタルで行うので、自分はこれを勝手にdigital shallow metagenomic approachと呼んでいる。
もしランダムサンプリングやフローセルの座標の指定でなく、RNAseqのデジタル正規化のようなことがサーバー側でできれば、存在量の少ない菌のリードを維持しつつ存在量の多い菌のリードを減らしてダウンロードの負担を下げることができるはずだが現状それば難しい(菌がいるかいないかだけわかれば良い時の話)。
*1
記事はR89となっているが古いのでR220のほうを使いたい。R220ならGTDBのサーバーからR220のgtdb_taxonomy.tsvをダウンロードして進める(ゲノムはヘルパースクリプトでダウンロードされる)。時間がかかるので、最初に少数のゲノムだけでテストすることを推奨。
データベースロード時間が短いなら、kraken2 folkのように複数サンプルを同時に扱えなくても必ずしもランタイムは長くならない。