sparNAはウィルスゲノムのアセンブリツール。ウィルスゲノムはRNA ploplymeraseのエラー率の高さなどの要因でhetero genesityが非常に高いため、特別な仕分け方をしない限りpopulation genomeやmeta genomeのデータセットに近い状態でシーケンス解析が行われる。また、ゲノムサイズが小さいために1万カバレッジ~10万カバレッジオーバーでシーケンスされることが多く、一般のk-merアセンブリではエラーと真の多様性を区別することが難しい(エラーのk-mer配列ですら何百回も出てくる可能性がある)。
sparNAはこのような混じり物で、割合に大きな違いがあるウィルスゲノムを個別にアセンブルして、結果をビジュアルで分析できるツールである。内部ではspadesを動かしているので、やっていることはspadesの総当たりに近いが、複数coverageと複数k-merを全自動で総あたりでアセンブルし、アノテーションもかけてくれるのは大変便利である。
インストール
依存
brewとpipで一括インストールできる。
brew install python3 spades bowtie2 samtools vcftools bcftools seqtk pip3 install argh numpy biopython khmer pandas plotly
本体 (Github)
ヘルプ
./sparna.py help
ラン
アダプタートリミングなしでアセンブル(あらかじめトリミングしておく)
./sparna.py --lca --fwd-fq paired1.fq --rev-fq paired2.fq --norm-c-list 1,5,20,100 --norm-k-list 21,25,31 --out-prefix output --threads 12
- --fwd-fq FWD_FQ
- --rev-fq REV_FQ
- --norm-c-list NORM_C_LIST
- --norm-k-list NORM_K_LIST
- --out-prefix OUT_PREFIX
- -t THREADS(4)
指定したk-merでアセンブルを実行する。またカバレッジを複数指定することで、それをターゲットとして独立にアセンブルを行うらしい。
アダプタートリミングしてからアセンブル。カバレッジのノーマライズは実行しない。
./sparna.py --lca --qual-trim --no-norm --fwd-fq paired1.fq --rev-fq paired2.fq --norm-c-list 1,5,20,100 --norm-k-list 21,25,31 --out-prefix tet --threads 12
- --qual-trim False
- --no-norm False
デフォルトでは Khmerのnormalize-by-median.pyでカバレッジがノーマライズされてアセンブルされるが、--no-normをつけるとノーマライズはスキップされそのままアセンブルが実行される。
結果はhtmlで出力される。
ツールの挙動を確かめるため、6kbのウイルスのシミュレーションデータ(シーケンスエラーがあるだけでほぼ純粋なゲノム)をアセンブルしてみた。
指定したk-merとカバレッジの組み合わせ全てでアセンブリが行われる。上の図で色が違うのはカバレッジまたはk-merサイズの異なるcontigである。例えば濃い緑のcontigはカバレッジ100指定、k-mer=21でアセンブルされたcontigである。カバレッジ指定が異なるcontigもだいたい同じ位置にplotがあるが(緑、黄土色、灰色)、k-merサイズが異なるとGCとcontigサイズが違う位置にcontigがプロットされている(水色と濃い青色)。
領域を指定するとそこが拡大される。
右上のauto scaleボタンかダブルクリックで元に戻る。
blastでhitすればアノテーションも付く。
右下のexport をクリック。
plotlyを使ったcontigの分析リンクが開く(plotly)。
Xボタンを押して消すことで特定のcontigのみ表示できる。
deepシーケンスしたウィルスゲノムなら、k-merのヒストグラムから推定ゲノムサイズを予測し、カバレッジとk-merを振ってランすることになる。
推定カバレッジが10万くらいだとどのような結果が得られるのかテストしてみる。
引用
GitHub - sdwfrost/sparNA: Consensus assembly of deep sequenced RNA virus populations