macでインフォマティクス

macでインフォマティクス

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

ユーザー定義の変異を再現可能なfastqのシミュレーター NEAT-genReads

 

NEAT-genReadsは2016年に発表されたfastqをシミュレートできるツール。変異のVCFファイルなどの情報も与えて現実に近いfastqを発生させることができる。fastq以外にポジコンとして使えるbamやVCFファイルも生成されるため、indel検出ツールの妥当性をポジコンと比較して考えることが可能になっている。

NEAT-genReadsにはすでに学習済みのプロファイルが用意されているが、ユーザーがbamやVCFファイルを指定することでユーザー定義のプロファイルでシーケンスbiasやエラー発生率を再現することができる。また領域をbedファイルで明示することで、特定の領域だけシミュレートすることができる(例えばexome解析のシミュレーション)。

 

8/29追記

 

インストール

Github

GitHub - zstephens/neat-genreads: NEAT read simulation tools

本体はpythonで記述されている。ダウンロードしてパスを通す。

 

依存

 

 

ラン

1、デフォルト条件でランすると平均カバレッジ10のfastqが作られる。

./genReads.py -r ref.fa -R 101 -o simulated_data
  • -R read length. Required.
  • -r ref.fa. Required.

ランが終わると"simulated_data_read1.fq" ができる。デフォルトのプロファイルに従い変異は導入されている。

 

 

 

2、インサートサイズが平均600で、カバレッジ x50、301bp x 2の ペアードエンドfastqを作る。アライメントのbamと変異が起きた部位を示したvcfファイルも出力させる。

./genReads.py --pe 600 50 -c 50 -R 301 -r ref.fa --bam --vcf -o simulated_data
  • --vcf Output golden VCF file
  • --bam Output golden BAM file
  • --pe Paired-end fragment length mean and standard deviation. To produce paired end data, one of --pe or --pe-model must be specified.
  • -c Average coverage across the entire dataset. Default: 10

"simulated_data_read1.fq" と "simulated_data4_read2.fq" ができる。 

 

bamとvcfをIGVに読み込ませて確認する。(igvがない人はこちらを参照)

samtools index simulated_data_golden.bam #index作成
igv -g ref.fa simulated_data_golden.bam,simulated_data_golden.vcf #IGV起動と読み込み

f:id:kazumaxneo:20170817191004j:plain

一番上が変異部位。ホモとヘテロの部位がある。defaultの条件では かなりの頻度でindel変異が入る。

 

 

 

3、ユーザーが指定したVCFと同じ変異を導入する。また、それ以外の変異の発生率をゼロにしてfastqを作る。

./genReads.py --pe 600 50 -c 50 -R 301 -v mutation_call.vcf -M 0 -r ref.fa -o simulated_data
  • -v Variants from this VCF will be inserted into the simulated sequence with 100% certainty.
  • -M The mutation rate model is rescaled to make this the average value. Must be between 0 and 0.3. These random mutations are inserted in addition to the once specified in the -v option.

 

入力のvcfと出力のvcfをIGVに読み込んで結果を確認。

samtools index simulated_data_golden.bam #index作成
igv -g ref.fa simulated_data_golden.bam,simulated_data_golden.vcf,mutation_call.vcf #IGV起動と読み込み

導入される変異がホモかヘテロかは入力のVCFファイルの GT:AD:DP:GQ:PLのGTで決定される。GTが1/1なら変異がホモ、0/1ならヘテロ、0/0ならリファレンスのホモになる。

 -参考-

ーtでbedを指定すれば、特定の領域だけシミュレートすることができます。exome解析やchip-seqのデータのシミュレートのほか、計算リソースを削減するために特定のクロモソームだけシミュレートしたい時に使えるオプションです。

 

vcfで情報を与えれば large indelも再現される。

f:id:kazumaxneo:20170829194446j:plain

large del(golden.bamをigvで開いて確認した)

 

 

 

4、既存データのGCバイアスに基づいてGC biasモデルを作り、それを使いfastqを作る。

bedtools genomecov -d -ibam normal.bam -g ref.fa > coverage
python2.7 computeGC.py -r ref.fa -i coverage -w 100 -o model
./genReads.py --gc-model GCcoverage -r ref.fa -R 101 -o simulated_data

はじめにbedtoolsでcoverageを計算させている(リンク)。

注意;リード長もモデルのbamと同じと仮定してシミュレートされるため、ラストのコマンドの-R の指定値は元のbamのリード長と揃える必要がある。

 

 

4と似た流れで変異、リードサイズ、シーケンスエラーのbiasを既存データからモデル化して利用することもできます。詳細はGitのHPで確認してください。

GitHub - zstephens/neat-genreads: NEAT read simulation tools

 

 

引用

Simulating Next-Generation Sequencing Datasets from Empirical Mutation and Sequencing Models

Zachary D. Stephens , Matthew E. Hudson, Liudmila S. Mainzer, Morgan Taschuk, Matthew R. Weber, Ravishankar K. Iyer

Published: November 28, 2016https://doi.org/10.1371/journal.pone.0167047