macでインフォマティクス

macでインフォマティクス

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

高速なショートリードとロングリードのアライナ Kart

 

 次世代シーケンシング(NGS)により、生物学者はヌクレオチド分解能でゲノム全体の変異を調べることができる。数多くの画期的な発見に寄与し、DNAの配列決定や集団内の変異の特徴付けに非常に一般的な手法となっている。新しいシークエンシング技術は、1日に数百万/ 10億塩基のオーダでリードを生成できるため、多くのNGSアプリケーションでは非常に高速なアライメントアルゴリズムが必要になる。 BLAST(Altschul et al、1990)またはBLAT(Kent、2002)のような従来の配列アラインメントアプローチは、膨大な量のショートリードを効率的に処理することができない。その結果、近年、NGSのためのアライナーが多く開発されている。index付け方法に応じて、ハッシュテーブルとsuffix array/BWTの2つのカテゴリに分類できる。ハッシュテーブルに基づくアライナは、k-mersのすべてのサブシーケンスを使用して位置を取得する。対照的に、suffix array/BWTベースのアライナは、リードとリファレンスゲノムとの間の最大完全一致(MEM)を見出す。リードアライナの各カテゴリーにはそれぞれ独自のメリットと欠陥がある。しかしながら、suffix array/BWTベースのアライナは、メモリ消費効率が良好でより一般的である。

 ハッシュテーブルに基づくアライナには、CloudBurst(Schatz、2009)、Eland(proprietary)、MAQ(Li et al、2008a)、RMAP(Smith et al、2008)、SeqMap(Jiang and Wong、2008)、SHRiMP NovaAlign(proprietary)、SSAHA(Ning et al、2001)およびSOAPv1(Li et al、2009)、ZOOM(Lin et al、2008)、BFAST(Homer et al、2009))。ほとんどのハッシュテーブルベースのアライナは、基本的に同じシード・アンド・エクステンション戦略に従います(Li and Homer、2010)。この戦略の代表的なアルゴリズムはBLASTである。 BLASTはデータベースシーケンスの各k-merの出現位置をハッシュテーブルに保持し、次にハッシュテーブルを検索することにより、所与のクエリーのk-merを用いて正確な一致を検索して見つける。完全一致は、クエリと参照配列との間のSmith-Watermanアルゴリズムを用いてアラインメントを拡張するためのシードとして使用される。

 カテゴリのほとんどのアライナは、最大完全一致(MEMと呼ばれる)を特定する suffix arrayに依存し、完全一致に基づいてアライメントを作成する。これはシード・アンド・extensionの方法と同様である。 1つの例外は、シーケンスから複数のシードを使用してマッピングされたゲノム位置を決定するシードアンドボトムステップを採用するSubreadアライナである。サフィックス配列を使用する主な利点は、反復サブシーケンスが単一のパス上に折りたたまれているため、1回のみアライメントする必要があることである(Li and Homer、2010)。

 現在のショートリードアライナは、NGSテクノロジによって生成された膨大な量のリードシーケンスをマッピングするためのソリューションを提供するが、一部では高速ではなく、一部では十分正確ではない。さらに、第3世代シーケンシング技術は、非常に長いシーケンスとはるかに高いエラー率というデータ分析の課題をさらに引き起している。例えば、Pac​​Bio RS IIシーケンサーは、平均して5500〜8500bpの長さのリードを生成できるが、単一リードの精度は約87%である。ほとんどのショートリードアライナは、それらのリードシーケンスを処理するのが難しい。

 これらの課題をすべて念頭に置いて、Kartというアラインメントアルゴリズムを開発した。このアルゴリズムは、BWT配列とハッシュテーブルの両方を使用する。Kartは、リードをアライメントしやすい領域とギャップアライメントを必要とする領域とに分け、最終的なアライメントを構成するために各領域を独立にアライメントさせる divide-and-conquer 戦略を採用する。著者らの実験では、ギャップアライメントを必要とする断片の平均サイズは、元のリード長にかかわらず約20である。合成データセットの実験では、Kartは、ほとんどのアライナよりも短い時間でロングリードをアライメントし、エラー率が15%でも信頼性の高いアライメントを生成する。リアルデータセットの実験で、Kartがシーケンシング品質の低いリードを処理できることが示した。

 

 

インストール

Github

https://github.com/hsinnan75/Kart/

git clone https://github.com/hsinnan75/Kart.git
cd Kart/
make

> ./kart

$ ./kart 

kart v2.4.4 (Hsin-Nan Lin & Wen-Lian Hsu)

 

Usage: ./kart -i Index_Prefix -f <ReadFile_A1 ReadFile_B1 ...> [-f2 <ReadFile_A2 ReadFile_B2 ...>] > out.sam

 

Options: -t INT        number of threads [4]

         -f            files with #1 mates reads (format:fa, fq, fq.gz)

         -f2           files with #2 mates reads (format:fa, fq, fq.gz)

         -o            output filename [stdout]

         -m            output multiple alignments

         -g INT        max gaps (indels) [5]

         -p            paired-end reads are interlaced in the same file

         -pacbio       pacbio data

bwa_indexとkartをパスの通ったディレクトリに移動しておく。

 

 

ラン

はじめにリファレンスのindexを作成する。

bwt_index ref.fasta ref

 

シングルエンドのマッピング。シーケンスデータはfasta かfastq。fastqはgz圧縮にも対応している。

kart -i ref -f ReadA_1.fq -o out.sam
  • -t     number of threads [16]
  • -i     index prefix [BWT based (BWA), required]
  • -f     read filename [required, fasta or fastq or fq.gz]
  • -f2   read filename2 [optional, fasta or fastq or fq.gz], f and f2 are files with paired reads
  • -o    alignment output

ペアエンドは-f-f2で指定する。interleaveのfastqを使うなら-fの代わりに-pを使う。

 

複数ファイルある場合、順に記載する。

kart -i ref -f ReadA_1.fq ReadB_1.fq ReadC_1.fq -f2 ReadA_2.fq ReadB_2.fq ReadC_2.fq -o out.sam.gz

  -o の拡張子を.gzにすると自動でgz出力に変更される。

 

 出力をsamtoolsにパイプしてbam出力する。 

kart -i ref -f Read1.fa.gz -f2 Read2.fa.gz | samtools view -bo out.bam

 

ショートリードのアライメント時間をbwa memと同一スレッドで比較すると、3倍ほどkartのほうが高速でした。

 

引用

Kart: a divide-and-conquer algorithm for NGS read alignment

Lin HN, Hsu WL

Bioinformatics. 2017 Aug 1;33(15):2281-2287