macでインフォマティクス

macでインフォマティクス

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

エラーコレクションツール karect

 

 De novo assemblyでは、シーケンスエラーを間違って足場に組み込んでしまうと、dead-endのグラフができたり、並列なパスができたり、あるいはキメラのパスができてしまう可能性がある。そのため、アセンブル前にエラーと思われるような低クオリティのリードはトリムする必要があるが、トリムすると低カバレッジ領域が増え、dead-endのコンティグが増えアセンブルのパフォーマンスが大きく落ちてしまう。低頻度のmerを除く方法論も同じで、やはり低カバレッジの領域の組み込みが困難になる。そこで、エラーと思われるような低頻度の配列を頻出する配列であらかじめ置き換えてしまう方法論が必要となる。一般的にエラーコレクションツール(ECツール)はこのような作業を担っている。

 ECツールとして、SPadesのBayesHammerやMaSuRCAのQuORUMのようにde novo  assemblyのツールの内部で動作するツールや、単独で動作するツールがある。ECツールはいくつか発表されているが、使用するデータやECツールによってエラーの検出感度や間違って置き換える偽陽性率には差がある。特に頻出する短い繰り返し配列ではそのようなことが起こりやすい。このエラーの検出感度とエラーの修正精度はトレードオフの関係になりやすく、recallもprecisonも良くてかつF measureも良いような方法論の構築は簡単ではない。karectはエラー補正を行うツールで、ミスマッチの他にindelも含めて補正を行う。原論文ではkarectによりNGA50が10%改善したと主張されている。エラー補正がアセンブルに及ぼす影響を13のECツールで比較した論文中では(ref.1)、NGA50最も改善されていた。ただし、ピークメモリとランタイムは他のECツールより多い傾向にあった。

 

 

 

インストール

Github

https://github.com/aminallam/karect

git clone https://github.com/aminallam/karect.git 
cd karect/
make

./karect #動作確認

 

user$ ./karect 

-Please specify the tool you want to run: (correct-align-eval-merge).

-Run "karect -[correct|align|eval|merge]" to find information about how to run a specific tool.

  1) correct: a tool to correct assembly reads from fasta/fastq files.

  2) align:   a tool to align original assembly reads as pre-processing for evaluation.

  3) eval:    a tool to evaluate assembly reads correction.

  4) merge:   a tool to concatenate fasta/fastq files.

-Example: "./karect -correct"

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

ヘルプはそれぞれのコマンドを打ち表示させる。例えば"./karect -correct"

 

ラン

 

1、Githubからリンクされているテストデータ(ペアードエンドのfastq) をダウンロードして解凍する。

2、エラー補正を実行。

karect -correct -threads=12 -matchtype=hamming -celltype=haploid -inputfile=./frag_1.fastq -inputfile=./frag_2.fastq -kmer=9
  • -correct a tool to correct assembly reads from fasta/fastq files.
  • -inputfile=f Specify an input fasta/fastq file. This option can be repeated for multiple files
  • -threads=i Specify the number of threads [Default=16].
  • -celltype=[haploid|diploid] Specify the cell type. Use "haploid" for bacteria and viruses.
  • -matchtype=[edit|hamming|insdel] Specify the matching type. "hamming" allows substitution errors only. "edit" allows insertions, deletion, and substitutions with equal costs. "insdel" is the same as "edit", but the cost of substitutions is doubled. Use "hamming" for Illumina datasets, and "edit" for 454 datasets.
  • -kmer=i Specify the minimum kmer size (will increase according to '-kmerfactor') [Default=9].

ref.1の論文では、k=9の時が一番パフォーマンスは良くなっている。

 

終わると複数ファイルが出力される。karect_frag~がエラーが修正されたfastq。

f:id:kazumaxneo:20171219194351j:plain

 

 

3、genomeファイルをダウンロードし、それにアライメントすることでエラー補正を評価する。

alignコマンドを実行。

karect -align -threads=12 -matchtype=hamming -inputfile=./frag_1.fastq -inputfile=./frag_2.fastq -refgenomefile=./genome.fasta -alignfile=./align.txt
  • -align a tool to align original assembly reads as pre-processing for evaluation.
  • -refgenomefile=f Specify the file containing the reference genome (to be aligned with).
  • -alignfile=f Specify the output alignment file.

 

 

4、最後にevalコマンドで評価する。

karect -eval -threads=12 -matchtype=hamming -inputfile=./frag_1.fastq -inputfile=./frag_2.fastq -resultfile=./karect_frag_1.fastq -resultfile=./karect_frag_2.fastq -refgenomefile=./genome.fasta -alignfile=./align.txt -evalfile=./eval.txt 
  • -eval a tool to evaluate assembly reads correction.
  • -resultfile=f Specify a result fasta/fastq file (resulting from running "karect -correct"
  • -alignfile=f Specify the alignment file resulted from running "karect -align".
  • -evalfile=f Specify the output evaluation file.

出力

Evaluation Started

Reference Genome Size = 2903081

..........

-------------------------------------------------------------------------------

Whole Read Statistics:

TN) num_correct_stayed_correct      = 382060 [of 1011910] (37.756322 %)

TP) num_wrong_correctly_fixed       = 622666 [of 1011910] (61.533733 %)

FP) num_correct_wrongly_fixed       = 91 [of 1011910] (0.008993 %)

FN) num_wrong_stayed_wrong          = 7093 [of 1011910] (0.700952 %)

*) Recall (Sensitivity)             = 98.873696 %

*) Precision                        = 99.985388 %

*) FScore                           = 99.426434 %

*) Gain                             = 98.859246 %

-------------------------------------------------------------------------------

Bases Statistics:

TN) num_correct_stayed_correct      = 99587505 [of 102202910] (97.440968 %)

TP) num_wrong_correctly_fixed       = 2601317 [of 102202910] (2.545247 %)

FP) num_correct_wrongly_fixed       = 1957 [of 102202910] (0.001915 %)

FN) num_wrong_stayed_wrong          = 12131 [of 102202910] (0.011870 %)

*) Recall (Sensitivity)             = 99.535824 %

*) Precision                        = 99.924825 %

*) FScore                           = 99.729945 %

*) Gain                             = 99.460942 %

-------------------------------------------------------------------------------

Evaluation Finished

 

 

 

scaffoldsにマッピングして観察した。上半分がエラー訂正前、下半分が訂正後。

f:id:kazumaxneo:20180221111236j:plain

クリップされた領域も表示されているので、ミスマッチが多いように見えるが、エラー訂正した後はそのミスマッチが多いリードが目に見えて減っている(面白いことに、このscaffoldsはエラー訂正前のリードをspadesでアセンブルしたものである。それでもエラーが減る。)

 

シミュレーションだと、もっと劇的な結果が出る。

まずwgsimで変異をわずかに発生させたリードをシミュレートした。それからエラー訂正有り無しでリファレンスにアライメントした。下図の上半分は変異入りリードを元のリファレンスにアライメントした結果で、下半分がkarectでエラー訂正してから元のリファレンスにアライメントした結果となる。wgsimで発生させたSNV、indelエラーが完全に消えている。

f:id:kazumaxneo:20180221173137j:plain

一方でホモの挿入変異は、当然そのまま残っている。

 

 

 

補足1

このあとspadesなどに持ち込む場合、--only-assemblerのフラグを立てることでエラーコレクションを無しでアセンブルを実行できます。一方MaSuRCAでは外部ツールでのエラーコレクションは推奨されていません。

 

補足2

バクテリアのリアルデータ(x30)でテストした。エラー補正後、ntcardでk-merカウントした結果が下になる(ntcard紹介)。

補正前

F1 69409538

F0 6243310

f1 2384162

f2 8925

f3 1407

f4 2304

f5 5057

f6 13572

f7 25671

f8 43404

補正後

F1 69409635

F0 3995728

f1 146017

f2 255

f3 640

f4 1472

f5 3584

f6 8834

f7 19268

f8 33607

固有のk-merのトータルを表すF0が大きく減少した。頻度が1、2回のf1とf2が大きく減少しているので、これが固有なk-merの母数が減った主な理由と思われる。

 

引用

Karect: accurate correction of substitution, insertion and deletion errors for next-generation sequencing data.

Allam A, Kalnis P, Solovyev V.

Bioinformatics. 2015 Nov 1;31(21):3421-8. doi: 10.1093/bioinformatics/btv415. Epub 2015 Jul 14.

 

ref.1

Evaluation of the impact of Illumina error correction tools on de novo genome assembly | BMC Bioinformatics | Full Text