2024/09/03 追記, 9/５追記
微生物パンゲノム解析は、原核生物ゲノム中の遺伝子の有無を同定する。しかし、現在のツールでは、配列の多様性が高い種や、属や科のような分類学上の上位の種を解析する場合には限界がある。Roary ILP Bacterial core Annotation Pipeline (RIBAP)は、コア遺伝子を同定するためにRoaryが予測した遺伝子クラスターを改良する整数線形計画法（integer linear programming approach）を用いている。RIBAPは、クラミジア、クレブシエラ、ブルセラ、エンテロコッカスゲノムの複雑さと多様性をうまく処理し、属レベルで包括的なコア遺伝子を同定する他の確立された最新のパンゲノムツールを凌駕している。RIBAPは、github.com/hoelzer-lab/ribapおよびzenodo.org/doi/10.5281/zenodo.10890871で自由に利用できるNextflowパイプラインである。
（省略）
パンゲノミクスにおける課題の一つは、複数のゲノムを解析することで得られる大量のデータを扱うことである。これには多大な計算資源とバイオインフォマティクスの専門知識が必要となる。加えて、パンゲノミクス研究ではしばしば多様な集団を分析することになるが、これらの集団を正確に定義しサンプリングすることは困難である。このような状況では、コアゲノムに属する遺伝子間の配列類似性が低い場合、例えば、属レベルで多様な種のパンゲノムを計算する場合に、特に問題が生じる。このような場合、遺伝子をコアゲノムに正しく割り当てることが難しくなり、誤ってアクセサリーゲノムの独立したグループになってしまう可能性がある。このように、配列の類似性だけに基づいて相同性を定義することは、特に種や属の境界を越えてゲノムを比較する場合、真のコアゲノムを過小評価することが多い。確立されたパンゲノム解析ツールは、通常、種レベル（at the species level）でのパンゲノムを計算することに重点を置いており、この進化的なレベルで評価している。しかし、本著者らの経験では、入力ゲノムの構成や配列の類似性（進化的な関連性を反映する）は、計算パンゲノムツールにとって課題となりうる。特に、生物種レベルを超えると、デフォルトの配列類似度のしきい値が高すぎて、コア遺伝子セットの過小評価につながる可能性があり、一方、しきい値が低すぎると、偽陽性の割り当てが多くなる可能性がある。既存のバイオインフォマティクスのパイプラインを拡張するだけでは、急速に増大する多様なゲノムデータセットの可能性を完全に実現することはできない。したがって、計算パンゲノミクスの分野を発展させるためには、質的に異なる新しい計算手法とパラダイムが必要である。
他のパンゲノムツールはオルソログとパラログを区別するためにある程度の遺伝子近傍情報を用いるが、配列の同一性とは無関係にオルソロジーそのものを推測するために遺伝子近傍情報を用いることはない。例えば、Roaryは保存遺伝子の近傍性の情報を使ってパラログを含む相同グループを真のオルソログを含むグループに分割する。しかし、そもそも配列の類似性が十分に高くないために相同グループが形成されなかったとしたらどうだろうか？ここでは、高い配列類似性に基づいた遺伝子クラスターの正確な初期計算と、これらのクラスターをより大きな遺伝子ファミリーグループにまとめるための、より厳密でないスキャフォールアプローチを組み合わせた新しい方法を提案する。このスキャフォールディングを行うために、整数線形計画法（ILP）を採用している。ILPは、制約条件の下で、特定のターゲットを最適化することによって、選択肢から最良の解を見つけるために用いられる数学的手法である。パンゲノミクスの文脈では、ILPは遺伝子が時間の経過とともに移動、重複、消失した可能性が最も高い方法を理解するためのツールと考えることができる。言い換えれば、本著者らはシンテニー情報（遺伝子の順序と向きの保存）をILPの定式化を通してモデル化し、高い配列類似性を超えて遺伝子クラスターを拡張するために遺伝子の順序の保存を最適化することを可能にする。
インストール
nextflow version 23.10.1.5891とlocal, dockerプロファイルでテストした（ubuntu22.04, 5995WX, 768GBメモリ）。
nextflow pull hoelzer-lab/ribap
> nextflow info hoelzer-lab/ribap
#v1.0.3を選択、プロファイルはlocal & docker
nextflow run hoelzer-lab/ribap -r 1.0.3 --fasta "$HOME/.nextflow/assets/hoelzer-lab/ribap/data/*.fasta" -profile local,docker
dockerではなくプロファイルはlocal & conda
#conda
nextflow run hoelzer-lab/ribap -r 1.0.3 --fasta "$HOME/.nextflow/assets/hoelzer-lab/ribap/data/*.fasta" -profile local,conda
#mamba
nextflow run hoelzer-lab/ribap -r 1.0.3 --fasta "$HOME/.nextflow/assets/hoelzer-lab/ribap/data/*.fasta" -profile local,mamba
ゲノムは--list --fasta genomes.csvでリストとしても指定できる。
出力例
workが作業ディレクトリ（ゲノムを増やして再解析する時に有利になる）、resultが最終出力のディレクトリとなる。遺伝子ファミリーごとのCDSアノテーション、アラインメント、遺伝子ツリー、コア遺伝子グループが含まれる。
ribap_roary95_summary.txt
upset plot
同定されたRIBAPグループを要約し、ユーザーが探索できるインタラクティブなHTMLテーブルも出力される。
テストデータでは5株の遺伝子を比較している。右端の5列が株を表し、その株の列で保存されていない遺伝子はNA表記となる。
HTMLテーブルでは、関心のある遺伝子ファミリーについて調べることができるようになっている。左端の緑色マークをクリックすると展開され、Roaryクラスターのヒートマップや遺伝子ツリーがplotされている。
30S ribosomal protein S12の場合
適当なhypoの場合
ヒートマップは異なる配列同一性閾値でのRoaryクラスターを表す。Roaryの60 %閾値では、5株のhypothetical proteinをコードするORF1つのクラスターにクラスタリングされ、より厳しい閾値（80%以上）では一番下の株のこのORＦは別のクラスターに分類されている。樹形でも、この配列だけrootから分岐して長い枝長が観察されている。
遺伝子の近接性もクラスタリングに関わっているので、roaryの60％閾値で別クラスターになっていても共クラスター化されているグループもある。
同定されたコア遺伝子全てから、IQ-TREEツリーを使った系統推定を実行する（1000 bootstrap replicates未満の値は使用できない）。
nextflow run hoelzer-lab/ribap -r 1.0.3 --fasta '*.fasta' --tree --bootstrap 1000 --outdir ribap
レポジトリより
- ゲノムFASTAとそれに対応するGenBankファイルをCSV形式で提供する場合は、--listと--reference パラメータを使用する（"--list --fasta genomes.csv --reference refs.csv"）。
-
RIBAPは現在のところ、数百、数千の入力ゲノムを扱うことは想定していない。
-
ILPとその結果を保存するためにかなりのディスク容量が必要になる。そのため、RIBAPはデフォルトで中間ILPファイルとその結果を保存しない。必要であれば、--keepILPsを使って保存できる。Chlamydia psittaciゲノムの数を変え、--keepILPsオプションを使用した場合と使用しなかった場合の実行時間とディスク容量の使用量がレポジトリの表に示されている（--keepILPsフラグを立てて32ゲノムを分析すると作業ディレクトリのworkのサイズが80GB以上となっている。このフラグ無しだと2.2GB）。
-
ローカルおよびHPCでの実行には--chunks 8が、デフォルトで使用されている。特にHPC上では、--chunksの値を大きくする（より並列にILPを解く）ことで実行時間を改善可能。
コメント
30ゲノム使って実行すると約1日半かかった。profile docker,localで実行したが最後のwebのステップでエラーを出して終了した。-resumeを付けてprofile mamba,localでそのままやり直すと、一部のステップの計算がやり直しになりさらに半日かかったが、今度はエラーなくランできた。
引用
RIBAP: a comprehensive bacterial core genome annotation pipeline for pangenome calculation beyond the species level
Kevin Lamkiewicz, Lisa-Marie Barf, Konrad Sachse & Martin Hölzer
Genome Biology volume 25, Article number: 170 (2024)
