macでインフォマティクス

macでインフォマティクス

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

SNVやindel変異を再現できるfastqのシミュレーターwgsim

#標準出力修正

 

wgsimはfastqをシミュレートできるツールである。Wgsimは、SNPと挿入/欠失多型と二倍体ゲノムをシミュレートできる(シーケンスエラーを再現したり、diploidゲノムの多型を想定して、一定の確率で変異を入れることができる)。もともとSAMtoolsのパッケージでリリースされたが、その後スタンドアローンになった(ref.1)。

wgsimはARTなどのツールでは不可能な300bp以上の配列を発生させることでも可能である。ただしシーケンスプロファイル等はないので、シーケンスのbias等を再現するには向いていない点には注意が必要。論文にはまだなっていない(リンク)。

 

インストール

Github

https://github.com/lh3/wgsim

git clone https://github.com/lh3/wgsim.git
cd wgsim/
gcc -g -O2 -Wall -o wgsim wgsim.c -lz -lm

 

ラン

ref.faの配列を使い、デフォルトの条件でfastqを発生させる。

wgsim ref.fasta R1.fastq R2.fastq > mutation
  • -e base error rate [0.020]
  • -d outer distance between the two ends [500]
  • -s standard deviation [50]
  • -N number of read pairs [1000000]
  • -1 length of the first read [70]
  • -2 length of the second read [70]
  • -r rate of mutations [0.0010]
  • -R fraction of indels [0.15]
  • -X probability an indel is extended [0.30]
  • -S seed for random generator [-1]
  • -A disgard if the fraction of ambiguous bases higher than FLOAT [0.05]
  • -h haplotype mode

カレントディレクトリにR1.fastqとR2.fastqができる。デフォルトだと70-bpでインサートサイズが500、 SDが50、リード数が1000000などとなっている。また、上記のデフォルトに設定された確率でindelやSNV、シーケンスエラーが導入されている。カバレッジはリード数xリード長をゲノムサイズで割って計算できる(ペアエンドはx2)。

f:id:kazumaxneo:20170816201633j:plain

カバレッジにはある程度ばらつきができる。 

 

リード長150-bpのfastqを30000リード数発生させる。

wgsim -1 150 -2 150 -N30000 ref.fasta R1.fastq R2.fastq > mutation

1-kbp以上のfastqを発生させることも可能であるが、シミュレーションには相当な時間がかかる。ロングリードを作るなら、pacbioかnanoporeのシミュレーターを使うことをお勧めします。

 

 indelやSNVについては以下のパラメータで調整する。

  • -e0 => シーケンスエラーがゼロになる。indelやSNVは起きる。
  • -r0 => indelがゼロになる。シーケンスエラーやSNVは起きる。

 よってindelもSNVもないfastqをシミュレートするならこのようになる。

wgsim -e0 -r0 ref.fasta R1.fastq R2.fastq > mutation

 

 

-hをつけると2倍体ゲノムのシミュレートになる。

f:id:kazumaxneo:20171230134738j:plain

 

変異の頻度ですが、ヒトゲノムだとヘテロ接合の部位の平均存在率は0.001という数値がアフリカ系人種から出ているそうです(-rのデフォルトは0.001)。

https://www.biostars.org/p/6177/

 

 

引用

https://www.biostars.org/p/137165/

 

https://popmodels.cancercontrol.cancer.gov/gsr/packages/wgsim/

 

ref.1

https://bmcresnotes.biomedcentral.com/articles/10.1186/1756-0500-7-40