color spaceのアライナー CUSHAW3 - macでインフォマティクス

　次世代シークエンシング（NGS）技術の出現と急速な進歩により、ショートリードアラインメントアルゴリズムの開発には相当量の研究努力が払われてきた。様々なショートリードアライナーが開発されており、機能性の点でさらに2つの世代に分類することができる。第1世代アライナーは、通常、非常に短いリード（通常100 bps以下）で設計され最適化されている。これらのアライナーは、通常、リードがゲノムからの非常に小さな偏差を有すると仮定し、したがって典型的にはミスマッチのみを許容する。いくつかのアライナーがギャップをサポートしているにもかかわらず、ギャップの最大許容数は速度のためにかなり制限されている（通常は1つのギャップ）。第1世代アライナの例としては、RMAP [論文より　ref.1]、MAQ [ref.2]、BFAST [ref.3]、Bowtie [ref.4]、BWA [ref.5]、CUSHAW [ref.6]およびSOAP3 [ref.7]が挙げられる。

　NGSの進歩により、最も広く使用されているilluminaシーケンシングのリードは100bpを超えて着実に増加している。しかしながら、これらのより長いショートリードは、通常、より高いシーケンスエラー率を犠牲にしている。他方、より長いリードはより真の挿入または欠失（indels）を有する傾向がある。そのため第1世代のアライナーは、アライメント品質、速度、またはその両面でより長いショートリードをアライメントするには非効率的となり、第2世代のショートリードアライナーの開発が動機づけられた。

　最近、BWA-SW [ref.8]、GASSST [ref.9]、Bowtie2 [ref.10]、CUSHAW2 [ref.11]、GEM [ref.12]、SeqAlto [ref.13]、SOAP3-dp [ref.14]、 BWA-MEM [ref.15]などいくつかの第二世代アライナーが開発されている（論文執筆時点）。これらのアライナーはすべてシード・アンド・エクステンション・パラダイムに基づいて設計されている。このパラダイムでは、まず、ゲノム上のシード、すなわち短い非ギャップ/ギャップ付きアラインメントを同定し、次に動的プログラミングを用いてアライメントを残りのリードに拡張することによって、リードをアライメントさせる。検索スペースをさらに減らすために、アラインメントの拡張に制約やフィルトレーションが行われることがよくある。様々なアライナによって異なるシードポリシーを使用できる。 BWA-SWは可変長ギャップシードを使用し、Bowtie2は固定長の非ギャップシード（不正確なマッチ）を抽出する。 GUSSTとSeqAltoの両方とも、固定長の正確なマッチ（k-merはk塩基のサブストリングである）シードを、CUSHAW2とBWA-MEMはそれぞれ可変長の最大完全一致（MEM）シードとスーパーMEMシードを同定する。 SOAP3は、Graphics Processing Unit（GPU）コンピューティングに基づくアライナであり、Bowtie2と同様のシードアプローチを採用している。GEMは、フィルタベースの近似文字列マッチングアプローチを採用して、適切なピジョンホールのようなルールで関連する候補マッチを抽出する。（一部略）。これらのアライナーは、高速に大部分のショートリードをマッピングできるが、シミュレートされたものであっても、すべてのリードをヒトゲノムなどの大きなゲノムに正確にアライメントさせることは困難である[ref.11] [ref.16]。したがって、アライメント品質をさらに向上させるために、新しいショートリードアライナーを設計することは非常に重要である。

　この論文では、base spaceとclor spaceの両方のシーケンスに対して、オープンソースの高感度で正確なショートリーダーアライナであるCUSHAW3を紹介する。CUSHAW3では、アライメント品質を向上させる３つのアプローチ; MEM seeds、exact-match k-mer seeds、variable-length seeds derived from local alignments、をハイブリッドシーディングとしてパイプラインに組み込んだ。さらに、正確なペアエンド（PE）アラインメントを容易にするためにweighted seed pairing heuristic、PE alignment pair ranking、read mate rescuing、の3つの技術が提案されている。ハイブリッドシーディングの概念は、ショートリードアライメントのためのいくつかの他の実装においてすでに暗示されている。 1つの例はIlluminaシーケンシングのためのアライナーStampy [ref.17]で、最初はBWA（不正確なマッチシードに基づく）とのアライメントを行い、完全マッチk -merを使った別のシード・アンド・エクステンド・ベースのアプローチでマッピングされていないリードを処理する。別の例は、SSAHA（固定長k-merシード）[ref.18]、BWA、BWA-SWからのアラインメントアプローチを組み込んだイオントレント配列決定のためのアライナーであるTMAPである。

　CUSHAW3のパフォーマンスは、シングルエンド（SE）およびPEを使い、シミュレートされたショートリードとリアルデータのショートリードの両方をヒトゲノムにマッピングすることによって評価された。ベーススペースアライメントのために、アライナーはNovoalign、CUSHAW2、BWA-MEM、Bowtie2およびGEMとも比較された。実験結果は、CUSHAW3が、SEおよびPEアライメントの両方について、CUSHAW2、BWA-MEM、Bowtie2およびGEMより一貫して優れていることを明らかにした。さらにcolor spaceのアラインメントは、CUSHAW3がSHRiMP2 [ref.19]やBFASTと比較して高速かつ最高のアライメント品質を提供する。

インストール

公式サイト

http://cushaw3.sourceforge.net/homepage.htm#latest

解凍してMakefileに従いビルドする。

make

$ ./cushaw3

cushaw3 (v3.0.3) is a fast gapped long-read aligner based on Burrows-Wheeler transform

Usage: cushaw3 <command> [options]

Command:

index build BWT and FM-index

align perform base-space read alignments (e.g. Illumina, 454, Ion Torrent and PacBio)

calign perform color-space read alignments (ABI SOLiD)

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

ラン

残念ながらcsfastaとqualを直接扱うことはできないので、まずfastqに変換しておく。ここではSocWikiでコードが貼られているSOLiD2std.plを使う（リンク）。

perl SOLiD2std.pl -seq input.csfasta -qual input.qual -o output.fastq

1、リファレンスのindexの作成。

cushaw3 index -c reference.fa

-p <string> (prefix of the index, default = fasta name)
-i <int> (interval for reduced suffix array [3], 2^#INT)
-c (build color-space index)

２、SOLiDのcolor space（fastqに変換済み）のマッピング。

先ほど作ったoutput.fastqを使う。

cushaw3 calign -r reference.fa -f output.fastq -t 12 -o output.sam

-r <string> (the file name base for the reference genome)
-f <string> file1 [file2] (single-end sequence files in FASTA/FASTQ format)
-b <string> file1 [file2] (single-end sequence files in BAM format)
-s <string> file1 [file2] (single-end sequence files in SAM format)
-q <string> file1_1 file1_2 [file2_1 file2_2] (paired-end sequence files in FASTA/FASTQ format)
-mode <int> (paired-end/mate-paired reads, default = 1)
- 0 means paired-end reads
- 1 means mate-paired reads