macでインフォマティクス

macでインフォマティクス

NGS関連のインフォマティクス情報についてまとめています。

SPAdesアセンブラ

 ref.1

 人体から海洋までほとんどの環境のバクテリアは研究所でクローン化できないため、既存のNGS(Next Generation Sequencing)技術を使用してシーケンスを決定することはできない。これは、Human Microbiome Project(HMP)(論文より Gill et al、2006)から抗生物質発見(Li and Vederas、2009)までのさまざまなプロジェクトの主要なボトルネックを表している。例えば、HMPの重要な問題は、バクテリアが互いにどのように相互作用するかである。これらの相互作用は、しばしば、他のバクテリアとのコミュニケーションまたはそれらを殺すために産生される様々なペプチドによって行われる。しかしながら、質量分析(このような研究のための重要な技術)は、かなり完全なプロテオームの知識を必要とするため、ヒトマイクロバイオームのpeptidomics研究は今や限られている。一方、新しいペプチド抗生物質の研究は、Non-Ribosomal Peptide Syntetases(NRPS)をコードする遺伝子のDNA配列決定(Sieber and Marahiel、2005)により大きく恩恵を受けるが、既存のメタゲノミクスアプローチは、通常このような例外的に長い遺伝子(60,000ヌクレオチド)をアセンブリすることができない。

 HMPと新しい抗生物質の発見は、Single-Cell Sequencing(SCS)によって革命を起こす多くのプロジェクトの2つの例にすぎない。 SCSの実験的(Ishoey et al、2008; Navin et al、2011; Islam et al、2011)および計算(Chitsaz et al、2011)の両局面における最近の改良により、1細胞からバクテリアゲノムを決定する可能性が開かれた。特に、Chitsaz et al(2011)は、SCSで生物の代謝を推定するのに十分な多数の遺伝子を捕捉できることを示した。多くのアプリケーション(プロテオミクスや抗生物質の発見を含む)では、大多数の遺伝子を捕捉することは、完全なゲノムと同じくらい有用である。

 現在、MDA(Multiple Displacement Amplification)はシングルセル増幅の主要技術である(Dean et al、2001)。しかし、MDAは極端な増幅バイアス(異なる領域間でカバレッジの差が桁違いに大きい)をもたらし、キメラリードとキメラペアリードが生じ、これが次のアセンブリを複雑にする。既存のアセンブラはこれらの問題に対処できないことから、Rodrigue et al. (2009) は、SCSが直面している課題は実験的ではなく計算的なものになっていると指摘した。最近の論文(Marcy et al、2007; Woyke et al、2010; Youssef et al、2011; Blainey et al、2011; Grindberg et al、2011)は、SCSにおけるフラグメントアセンブリの課題を示している。

 Chitsaz et al (2011)はE+V-SC assemblerを導入し、modified EULER-SR とmodified velvet を組み合わせ、SCSデータのフラグメントアセンブリの大幅な改善を達成した。しかし、Chitsaz et alの共著者である我々(この論文の著者ら)(G.T.およびP.A.P.)は、 (2011)では、SCSの可能性を最大限に活用するために、(既存のツールを変更するだけでなく)アルゴリズム設計を変更する必要があることに気付いた。

 著者らは、SPAdesアセンブラを発表し、SCSと標準(マルチセル)バクテリアの両方のデータセットのために、数多くの新しいアルゴリズムソリューションを導入し、最先端のアセンブラを改良した。フラグメントアセンブリは、そのk-mersのセットから文字列を再構成する問題として、しばしば抽象化される。この抽象化は、多くのフラグメントアセンブリアルゴリズムの基礎となるアセンブリへのde Bruijnアプローチに自然につながる。しかしながら、NGSデータのより高度な抽象化は、文字列中の距離≒dのk-mersのペア(以下、k-bimersと呼ぶ)から文字列を再構築する問題を考慮する。残念ながら、前者の抽象化には単純なアルゴリズムがある一方で、後者の抽象化の分析は、主にde Bruijn graphの後処理ヒューリスティックになっている。ペアリード(以下、bireads2と呼ぶ)の解析のための多くのヒューリスティックが提案されているが(Pevzner et al、2001; Pevzner and Tang、2001; Zerbino and Birney、2008; Butler et al、2008; Simpson et al、2009 ; Chaisson et al、2009; Li et al、2010)、bireads2の利用はアセンブリのステージの中でおそらく最も不十分な段階であると考えられる。 Medvedev et al(2011a)は、後者の抽象化に適した新しい手法であり、標準的なDe Bruijnグラフより重要な利点を持つPaired de Bruijnグラフ(PDBG)を導入した。しかし、PDBGは理論的には実用的ではなくむしろ理論的に導入されており、主にリード間を固定距離とした非現実的なケースを対象としていた。

 IduryとWaterman(1995)がフラグメントアセンブリのためにde Bruijn graphを導入したとき、多くの人々はサンガーのリードのエラー率が高いため、このアプローチを非実用的と見なしていた(ref.3)。Pevzner et al (2001)は、このボトルネックを解消するためにエラー訂正手順を導入した。同様に、PDBGsは、NGSに特徴的なbiread(and thus k-bimer)距離の変動のために実用的でないように見えるかもしれない。 SPAdesは、k-bimersの正確な距離を明らかにするk-bimer調整を導入し、PDBGに触発されたペアのアセンブリグラフを導入することによって、このボトルネックに対処する。特に、SPAdesはペアリードを利用することができる。 E + V-SCは、リードを使用したが、キメラリードペアのレベルの上昇によって生じるミスアセンブリを回避するために、ペアエンドのpairingを無視した。

 読者はPevzner et al(2004 pubmed)で紹介されたA-Bruijn graphの概念に精通していると仮定する 。この論文のDe Bruijn graph、PDBG、およびその他のgraphは、A-Bruijn graphの特別なケースである。 SPAdesは、最初のde Bruijnグラフを作成するためだけにk-mersを使用し、後でそれらについて「忘れる」という意味で普遍的なA-Bruijnアセンブラである。それ以降のステージでは、k-mersによってラベル付けする必要がないgraphに対してのみgraph理論操作を実行する。操作はgraphトポロジーカバレッジ、およびシーケンス長に基づいているが、シーケンスそのものには基づいていない。最後の段階で、コンセンサスDNA配列が復元される。著者らは、同じフレームワーク内でA-Bruijn graphのいくつかのバリエーション(例えば、Bruijn graphの対と多重化)を実装し、これらのgraphが有用であることが証明されている他のアプリケーションに適用する汎用アセンブラを設計した(Bandeira et al、 2007、2008; Pham and Pevzner、2010)。

 セクション2では、SPAdesのステージの概要を説明する。セクション3では、de Bruijn graph、Multisized de Bruijn graph、およびPDBGを定義する。セクション4-6は、SPAdesのさまざまな段階をカバーしている。セクション7では、シングルセルおよび培養された大腸菌(E.coli)データセットについて、SPAdesおよび他のアセンブラベンチマークする。セクション8では、アセンブリgraphの構築とデータ構造に関する追加の詳細を示す。第9章では、ペアアセンブリgraphを構築する詳細な例を示す。セクション10では、ユニバーサルアセンブラの概念についてさらに説明する。

セクション2

以下では、SCSで特に面倒な問題を扱うSPAdesの4つの段階について説明する。

(1)ステージ1(アセンブリgraph構築)はすべてのNGSアセンブラによって扱われ、しばしばde Bruijn graph簡略化(EULER / Velvetにおけるバルジ/バブル除去)と呼ばれる。この論文では、複数のk-merを使ったde Bruijn graphを使用し、新しいバルジ/チップ除去アルゴリズムを実装し、キメラを検出して除去し、biread情報を距離ヒストグラムに集約し、実行されたgraph操作をバックトラックすることを可能にするアセンブリグラフ構築の新しいアプローチを提案する。
(2)ステージ2(k-bimer調整)は、アセンブリgraphの距離ヒストグラムとパスのジョイント解析を使用して、ゲノム内のk-mer間の正確な距離推定値(アセンブリgraph内のエッジ)を導き出す。
(3)ステージ3は、PDBGアプローチからインスパイアされたペアのアセンブリgraphを構築する。
(4)ステージ4(コンティグ構築)は、サンガーシーケンス時代に十分に研究された(Ewing et al、1998)。 NGSプロジェクトは通常、高いカバレッジを特徴とするため、NGSアセンブラはかなり正確なコンティグを生成する(ただしSCSの精度は低下する)。 SPAdesは、コンティグのDNA配列を構築し、グラフの簡略化(論文 8.6節参照)を目的としてコンティグへリードをマッピングする。
これまでの研究では、さまざまなアセンブラをエラー訂正ツールと組み合わせることで性能が向上することが実証されている(Pevzner et al、2010; Kelley et al、2010; Ilie et al、2010; Gnerre et al、2011)。しかし、 (Quake [Kelley et al、2010]などの)ツールは、暗黙的にほぼ均一なカバレッジを仮定しているため、シングルセルのデータでは機能しない。 Chitsaz et al。(2011)はVelvet-SCとEULERのエラー訂正(Chaisson and Pevzner、2008)を結びつけ、E + V-SCツールとなった。この論文では、SPAdesは、アセンブリ前にエラー修正と品質トリミングのためにHammer(Medvedev et al、2011b)(SCSを対象としています)によるエラー修正を行う。

 

オリジナルのSPAdesは以下のようなフローでアセンブリを実行する。

f:id:kazumaxneo:20180614103714j:plain

プレゼンテーションPDFより転載。 de Brujin graph構築にはショートリード情報が用いられ、ロングリード情報はエラーの多寡に関わらず使われない。

 

公式ページ

http://cab.spbu.ru/software/spades/

3.12マニュアル

http://cab.spbu.ru/files/release3.12.0/manual.html 

 

New Frontiers of Genome Assembly with SPAdes 3.1

質疑応答がちょっと面白いです。

ppt1

f:id:kazumaxneo:20180613212251j:plain

ppt2

f:id:kazumaxneo:20180613212329j:plain

 

インストール

mac os 10.13でテストした。

 バイナリをダウンロードする。Linuxは上のマニュアル参照(バイナリが用意されている)。

curl http://cab.spbu.ru/files/release3.12.0/SPAdes-3.12.0-Darwin.tar.gz -o SPAdes-3.12.0-Darwin.tar.gz
tar -zxf SPAdes-3.12.0-Darwin.tar.gz
cd SPAdes-3.12.0-Darwin/bin/

./spades.py 

$ ./spades.py 

SPAdes genome assembler v3.12.0

 

Usage: ./spades.py [options] -o <output_dir>

 

Basic options:

-o <output_dir> directory to store all the resulting files (required)

--sc this flag is required for MDA (single-cell) data

--meta this flag is required for metagenomic sample data

--rna this flag is required for RNA-Seq data 

--plasmid runs plasmidSPAdes pipeline for plasmid detection 

--iontorrent this flag is required for IonTorrent data

--test runs SPAdes on toy dataset

-h/--help prints this usage message

-v/--version prints version

 

Input data:

--12 <filename> file with interlaced forward and reverse paired-end reads

-1 <filename> file with forward paired-end reads

-2 <filename> file with reverse paired-end reads

-s <filename> file with unpaired reads

--merged <filename> file with merged forward and reverse paired-end reads

--pe<#>-12 <filename> file with interlaced reads for paired-end library number <#> (<#> = 1,2,...,9)

--pe<#>-1 <filename> file with forward reads for paired-end library number <#> (<#> = 1,2,...,9)

--pe<#>-2 <filename> file with reverse reads for paired-end library number <#> (<#> = 1,2,...,9)

--pe<#>-s <filename> file with unpaired reads for paired-end library number <#> (<#> = 1,2,...,9)

--pe<#>-m <filename> file with merged reads for paired-end library number <#> (<#> = 1,2,...,9)

--pe<#>-<or> orientation of reads for paired-end library number <#> (<#> = 1,2,...,9; <or> = fr, rf, ff)

--s<#> <filename> file with unpaired reads for single reads library number <#> (<#> = 1,2,...,9)

--mp<#>-12 <filename> file with interlaced reads for mate-pair library number <#> (<#> = 1,2,..,9)

--mp<#>-1 <filename> file with forward reads for mate-pair library number <#> (<#> = 1,2,..,9)

--mp<#>-2 <filename> file with reverse reads for mate-pair library number <#> (<#> = 1,2,..,9)

--mp<#>-s <filename> file with unpaired reads for mate-pair library number <#> (<#> = 1,2,..,9)

--mp<#>-<or> orientation of reads for mate-pair library number <#> (<#> = 1,2,..,9; <or> = fr, rf, ff)

--hqmp<#>-12 <filename> file with interlaced reads for high-quality mate-pair library number <#> (<#> = 1,2,..,9)

--hqmp<#>-1 <filename> file with forward reads for high-quality mate-pair library number <#> (<#> = 1,2,..,9)

--hqmp<#>-2 <filename> file with reverse reads for high-quality mate-pair library number <#> (<#> = 1,2,..,9)

--hqmp<#>-s <filename> file with unpaired reads for high-quality mate-pair library number <#> (<#> = 1,2,..,9)

--hqmp<#>-<or> orientation of reads for high-quality mate-pair library number <#> (<#> = 1,2,..,9; <or> = fr, rf, ff)

--nxmate<#>-1 <filename> file with forward reads for Lucigen NxMate library number <#> (<#> = 1,2,..,9)

--nxmate<#>-2 <filename> file with reverse reads for Lucigen NxMate library number <#> (<#> = 1,2,..,9)

--sanger <filename> file with Sanger reads

--pacbio <filename> file with PacBio reads

--nanopore <filename> file with Nanopore reads

--tslr <filename> file with TSLR-contigs

--trusted-contigs <filename> file with trusted contigs

--untrusted-contigs <filename> file with untrusted contigs

 

Pipeline options:

--only-error-correction runs only read error correction (without assembling)

--only-assembler runs only assembling (without read error correction)

--careful tries to reduce number of mismatches and short indels

--continue continue run from the last available check-point

--restart-from <cp> restart run with updated options and from the specified check-point ('ec', 'as', 'k<int>', 'mc', 'last')

--disable-gzip-output forces error correction not to compress the corrected reads

--disable-rr disables repeat resolution stage of assembling

 

Advanced options:

--dataset <filename> file with dataset description in YAML format

-t/--threads <int> number of threads

[default: 16]

-m/--memory <int> RAM limit for SPAdes in Gb (terminates if exceeded)

[default: 250]

--tmp-dir <dirname> directory for temporary files

[default: <output_dir>/tmp]

-k <int,int,...> comma-separated list of k-mer sizes (must be odd and

less than 128) [default: 'auto']

--cov-cutoff <float> coverage cutoff value (a positive float number, or 'auto', or 'off') [default: 'off']

--phred-offset <33 or 64> PHRED quality offset in the input reads (33 or 64)

[default: auto-detect]

ほかにmetaspades.py、rnaspades.py、dipspades.py、plasmidspades.py、/truspades.pyがある。以下のような使い分けになる。

  • metaSPAdes – a pipeline for metagenomic data sets (see metaSPAdes options).
  • plasmidSPAdes – a pipeline for extracting and assembling plasmids from WGS data sets (see plasmidSPAdes options).
  • rnaSPAdes – a de novo transcriptome assembler from RNA-Seq data (see rnaSPAdes manual).
  • truSPAdes – a module for TruSeq barcode assembly (see truSPAdes manual).
  • dipSPAdes – a module for assembling highly polymorphic diploid genomes (deprecated, see dipSPAdes manual).

spadesはbrewでも導入できる。

brew install spades

SPAdesは現在でも開発が続けられており、バージョンを変えると末端部位やカバレッジの低い領域の結果が揺れ動くことがあります。アップグレードする時は注意してください。

 

 

ラン

ペアエンドリードfastqのアセンブリ。非圧縮fastqにも対応。

spades.py -1 pair1.fq.gz -2 pair2.fq.gz -k auto --careful -t 12 -o output
  • -1     file with forward paired-end reads-1 <filename> file with forward paired-end reads
  • -2     file with reverse paired-end reads
  • -t     number of threads-t/--threads <int> number of threads [default: 16]
  • -k     comma-separated list of k-mer sizes (must be odd and-k <int,int,...> comma-separated list of k-mer sizes (must be odd and less than 128) [default: 'auto']
  • --careful     tries to reduce number of mismatches and short indels--careful tries to reduce number of mismatches and short indels
  • -o   directory to store all the resulting files (required)-o <output_dir> directory to store all the resulting files (required)

インターレースfastqは--12、シングルエンドfastqは-sで指定する。シングルセルシーケンスなどで複数fastqがある場合は、数値をつけて指定する。例えばインターレースfastqを複数指定するなら"--pe1-12 sample1.fq --pe2-12 sample2.fq --pe3-12 sample3.fq ..."。詳細はマニュアルの Specifying multiple librariesの節を参照。

 

ペアエンドリードfastqとONTのロングリードのハイブリッドアセンブリ

spades.py -1 pair1.fq.gz -2 pair2.fq.gz --nanopore long.fa -k auto --careful -t 12 -o output

 

ランが終わるとk-merごとにアセンブリデータフォルダができる。scaffolds.fastaが最終出力となる。correctedにはエラー訂正されたfq.gzも保存されている。

f:id:kazumaxneo:20180613230424j:plain

 

ケースバイケースですが、うまくパラメータを使えばパフォーマンスが上がる可能性があります。気づいたことを書いておきます。

  1. --carefulをつけると、付属のbwaを使い、アセンブリ後にミスマッチを訂正するステップが挿入される。フラグを立てると、マッピングの時間分だけランタイムが増えるが、スモールゲノムではフラグをつけることが推奨されている(defaultではオフ)。ラージゲノムでは非推奨。純粋な培養株のシーケンスを想定して修復するため、rnaSPADdesとmetaSPAdesにはこのオプションはない。
  2. たとえハイクオリティなPacbioのロングリードでも、ロングリードだけでSPAdesを動かすことは非推奨になってます。
  3. --only-error-correctionをつけると、BayesHammerによるエラー訂正だけが実行されます。エラー訂正されたfastqが欲しいときに使えます。
  4. --only-assemblerをつけると、エラー訂正なしでアセンブリが実行されます。非推奨ですが、ほかのエラー訂正ツールで修復したfastqを使う時にも利用できます。
  5. カビくらい大きなゲノムになるとかなりメモリを食います。論文を読んで、ランできそうか確認してください。defaultでは-m --memoryは最大250GBで、越えるとランを停止します。
  6. 軽度のコンタミネーションは、--cov-cutoff <NUM>でカバレッジcut-offをつけることで抑制できます。アセンブリ後にカバレッジGC、taxonomy profileを分析して、異種生物ゲノムのcontigのコンタミカバレッジで区別できそうなら使用して下さい。

 

 

http://quast.bioinf.spbau.ru

 

  • Unicycler

 

引用

ref.1

SPAdes: A New Genome Assembly Algorithm and Its Applications to Single-Cell Sequencing

Anton Bankevich, Sergey Nurk, Dmitry Antipov, Alexey A. Gurevich, Mikhail Dvorkin, Alexander S. Kulikov, Valery M. Lesin, Sergey I. Nikolenko, Son Pham, Andrey D. Prjibelski, Alexey V. Pyshkin, Alexander V. Sirotkin, Nikolay Vyahhi, Glenn Tesler, Max A. Alekseyev,Pavel A. Pevzner, 

J Comput Biol. 2012 May; 19(5): 455–477.

 

ref.2

ExSPAnder: a universal repeat resolver for DNA fragment assembly

Prjibelski AD, Vasilinetc I, Bankevich A, Gurevich A, Krivosheeva T, Nurk S, Pham S, Korobeynikov A, Lapidus A, Pevzner PA.

Bioinformatics. 2014 Jun 15;30(12):i293-301

 

ref.3

Assembling short reads from jumping libraries with large insert sizes

Vasilinetc I, Prjibelski AD, Gurevich A, Korobeynikov A, Pevzner PA.

Bioinformatics. 05 Oct 5;(0):6-8.

 

ref.4

hybridSPAdes: an algorithm for hybrid assembly of short and long reads

Antipov D, Korobeynikov A, McLean JS, Pevzner PA

Bioinformatics. 2016 Apr 1;32(7):1009-15.

 

ref.5

TruSPAdes: barcode assembly of TruSeq synthetic long reads

Bankevich A, Pevzner PA

Nat Methods. 2016 Mar;13(3):248-50.

 

ref.6

metaSPAdes: a new versatile metagenomic assembler

Nurk S, Meleshko D, Korobeynikov A, Pevzner PA

Genome Res. 2017 May;27(5):824-834.

 

ref.7

BayesHammer: Bayesian clustering for error correction in single-cell sequencing

Nikolenko SI, Korobeynikov AI, Alekseyev MA

BMC Genomics. 2013;14 Suppl 1:S7.

 

 

参考

バクテリアのゲノムアセンブリ

Tychus: a whole genome sequencing pipeline for as- sembly, annotation and phylogenetics of bacterial genomes

https://www.biorxiv.org/content/biorxiv/early/2018/03/16/283101.full.pdf

 

アセンブリの進め方。

http://compgenomics2017.biology.gatech.edu/images/a/ad/Genome_Assembly_Final_Pipeline_and_Results.pdf