次世代シークエンシング(NGS)に基づいて、単一ヌクレオチド変異(SNV)または挿入 - または欠失(indel)突然変異を呼び出すための様々なアプローチが開発されている。しかし、それらの大部分は、特定のタイプの突然変異に捧げられている。正常細胞における生殖系列SNV、癌/腫瘍細胞における体細胞SNV、またはindelsのみである。文献では、生殖系列および体細胞SNVs / indelsの両方に対する効率的かつ統合されたコーラーは、まだ広範に調査されていない。
SNVSnifferは、NGSデータからの生殖系列と体細胞SNVs / indelsの両方を同定する効率的かつ統合されたコーラーである。このアルゴリズムでは、ベイジアン確率モデルを使用してSNVを特定し、indelsを呼び出すため複数のungappedアライメントアプローチを調べることを提案する(一部略)。
多様なベンチマークを使用して、包括的なパフォーマンス評価を実施した。生殖系列変異コールのために、SNVSnifferは最先端のFaSD、GATKおよびSAMtoolsと比較して優れたスピードで競争力の高い正確性を示す。体細胞変異のコールのために、このアルゴリズムは、VarScan2、SomaticSniper、JointSNVMix2、MuTectよりも高速で、同等またはそれ以上の精度を達成する。
本ツールには、germlineとsomaticの変異を再現可能なシーケンスリードのシミュレータ機能も備わっています。ユーザー定義の条件で変異を発生させて、本ツールで解析し、評価する、という流れが本ツールのみで可能になっています。
公式サイト
http://snvsniffer.sourceforge.net/homepage.htm#latest
version2マニュアル
http://snvsniffer.sourceforge.net/version2.htm#latest
インストール
cent OSに導入した。
依存
- SAMtools
本体は上記公式から実行ファイルをダウンロードして使う。SNVSniffer 1.0はSNVのみがターゲットだが、SNVSniffer 2.0 はindelも検出可能になっている。
$ ./SNVSniffer
SNVSniffer (v2.0.4): identifying germline and somatic single nucleotide and indel mutations (released on March 07, 2016)
Command:
snp Identify single-sample SNVs/indels
somatic Identify germline/somatic SNVs/indels from tumor-normal pairs
gsim Simulate sample reads with germline SNVs/indels
ssim Simulate normal and tumor sample reads with somatic SNVs/indels
eval Evaluate predicted germline/somatic SNVs/indels against gold standard, only for evaluation purpose
ラン
5つのコマンドのうち、performance結果を評価するevalコマンド以外の4つを紹介する。テストデータがアップされているので、それを使ってテストランする。リファレンスは"Home / data / v1.0 / synthetic_tumors_error_0.01"にある(リンク)。
1、snp - Identify single-sample SNVs/indels
$ SNVSniffer snp
Usage: SNVSniffer snp [options] sam_header infile
* sam_header: a SAM header file associated with the BAM file
* infile: a mpileup/pileup/BAM file (can have no header)
Input:
-g <string> (reference genome file, [Required by BAM format])
-f <int> (input file format, default = 2)
0: mpileup format generated by SAMtools
1: pipeup format generated by MAQ
2: BAM file
Output:
-o <string> (output file name, default = STDOUT)
Base call and coverage:
-min_cov <int> (minimum coverage, default = 3)
-max_cov <int> (maximum coverage, default = 65536)
-min_bqual <int> (minimum Phred base quality score, default = 13)
-min_mapq <int> (minimum mapping quality score, default = 0)
-seq_err_rate <float> (sequencing error rate of the data, default = 0.01)
SNV calling:
-exe_mode <int> (execution model for SNV calling, default = 0)
0: very fast and lesst accurate
1: fast and accurate
2: slow and most accurate
-prior <int> (model used for genotype prior probabilities, default = 1)
0: equal probability
1: prior probabilities with no consideration of Ti/Tv ratio
2: prior probabilties by considering Ti/Tv ratio
-snp_rate <float> (SNP mutation rate for the species, default = 0.001)
-min_allele_freq <int> (minimum allele frequency, default = 0.2 [Important])
-stringency <int> (stringency level [0, 9] for low-confidence mutations, default = 6)
-homo_freq <float> (allelic frequency threshold for homozygous genotype, default = 0.75)
-relaxed <int> (using the relaxed confidence score calcuation (favours sensitivity), default = 1)
-min_locus_dist <int> (minimum distance interval between two neighboring SNP loci, default = 1 [0 disables it])
-use_strand_dist <bool> (use strand distribution information, default = 0)
-local_range <int> (set the local search range for strand distribution, default = 1000)
-call_snps <bool> (enable the calling of SNPs, default = 1)
-call_indel <bool> (enable the calling of insertions or deletions, default = 1)
-pvalue <float> (P value threshold for variant calling, default = 0.05)
ランにはSAMのヘッダー情報が必要。ヘッダーファイルを作成してからランする。
samtools view -H input.bam > header.sam
SNVSniffer snp -f 2 header.sam -g reference.fa input.bam > output.vcf
出力。
$ head -n 15 out.vcf
##fileformat=VCFv4.1
##source=SNVSniffer-v2.0.4
##INFO=<ID=DP,Number=1,Type=Integer,Description="Coverage">
##INFO=<ID=HC,Number=1,Type=Float,Description="High confidence">
##INFO=<ID=LC,Number=1,Type=Float,Description="Low confidence">
##INFO=<ID=VT,Number=1,Type=String,Description="Variant type: SNP, INS or DEL">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT SAMPLE
chr21 5010446 . T C . PASS DP=16;VT=SNP;HC=3.657 GT 0|1
chr21 5012778 . C G . PASS DP=26;VT=SNP;HC=4.350 GT 0|1
chr21 5013197 . C T . PASS DP=27;VT=SNP;HC=3.657 GT 0|1
chr21 5015438 . C A . PASS DP=45;VT=SNP;HC=4.350 GT 0|1
chr21 5016893 . C G . PASS DP=33;VT=SNP;HC=4.350 GT 0|1
chr21 5017812 . G A . PASS DP=31;VT=SNP;LC=3.657 GT 0|1
chr21 5018805 . G A . PASS DP=32;VT=SNP;HC=3.657 GT 0|1
2、somatic - Identify germline/somatic SNVs/indels from tumor-normal pairs
$ SNVSniffer somatic
Usage: SNVSniffer somatic [options] normal_sam_header tumor_sam_header normal tumor
1. normal_sam_header: SAM header file for normal
2. tumor_sam_header: SAM header file for tumor
3. normal: a mpileup/pileup/BAM file for normal
4. tumor: a mpileup/pileup/BAM file for tumor
Input:
-g <string> (reference genome file, required for BAM format input)
-f <int> (input file format, default = 2)
0: mpileup format generated by SAMtools
1: pipeup format generated by MAQ
2: BAM format
Output:
-o <string> (output file name, default = STDOUT)
-o_somatic <int> (output somatic mutations, default = 1)
-o_loh <int> (output loss of heterzygosity mutations, default = 1)
-o_germline <int> (output germline mutations, default = 0)
-o_unknown <int> (output unknown type mutations, default = 1)
-o_indel <int> (output indels, default = 1)
Base call and coverage:
-min_cov_normal <int> (minimum normal coverage, default = 3)
-min_cov_tumor <int> (minimum tumor coverage, default = 3)
-min_bqual <int> (minimum Phred base quality score, default = 13)
-min_mapq <int> (minimum mapping quality score, default = 0)
-seq_err_rate <float> (sequencing error rate of the data, default = 0.01)
SNV calling:
-exe_mode <int> (execution model for SNV calling, default = 0)
0: very fast
1: fast and more accurate
2: slow and most accurate
-prior <int> (model used for genotype prior probabilities, default = 1)
0: equal probability
1: prior probabilities with no consideration of Ti/Tv ratio
2: prior probabilties by considering Ti/Tv ratio
-somatic_rate <float> (somatic mutation rate, default = 0.01)
-tumor_purity <float> (estimated purity (tumor data) of tumor sample, default = 0 [0 means AUTO])
-min_allele_freq <float> (minimum allele frequency for the normal, default = 0.2)
-stringency <int> (stringency level [0, 9] for low-confidence mutations, default = 6)
-homo_freq <float> (allelic frequency threshold for homozygous genotype, default = 0.75)
-p_value <float> (P value threshold for variant calling per sample, default = 0.05 [0 means global])
-use_strand_dist <bool> (use strand distribution information, default = 0)
-local_range <int> (set the local search range for strand distribution, default = 1000)
最初にnormalとtumorのbamからheader.samを作る。
samtools view -H chr21_e15_c30_normal.realigned.bam> normal_header.sam
samtools view -H chr21_e15_c30_tumor.realigned.bam > tumor_header.sam
somaticコマンドを走らせる。tumorの純度を"-tumor_purity <float>"で指定しない時はSNVSnifferによって推定された値が使われる。
SNVSniffer -f 2 -g chr21.fa somatic normal_header.sam tumor_header.sam chr21_e15_c30_normal.realigned.bam chr21_e15_c30_tumor.realigned.bam -o out.vcf
3、gsim - Simulate sample reads with germline SNVs/indels
$ SNVSniffer gsim
SNVSniffer gsim: a diaploid germline mutation simulator derived from WGSIM in SAMtools
Usage: SVNSniffer gsim [options] reference.fa reads.base
reference.fa: reference file stored in FASTA/FASTQ format
reads.base: the base file name for reads and VCF-formated mutations
Options:
Output:
-z (compress all of the output using ZLIB)
Random numbers:
-p <int> (random number seed value for generating the simulated genome, default=11)
-q <int> (random number seed value for generating simulated reads from the simulated genome, default = 47)
Reads:
-c <int> (coverage of the genome, default = 0 [>0 distables option -n])
-n <int> (number of read pairs, default = 1000000)
-1 <int> (length of the first read, default = 100)
-2 <int> (length of the second read, default = 100)
-a <float> (disgard if the fraction of ambiguous bases higher than #FLOAT, default = 0.05)
Germline:
-e <float> (sequencing base error rate, default = 0.02)
-i <int> (average insert size, i.e. outer distance, between the two ends, default = 500)
-d <int> (standard deviation of insert size, default = 50)
-r <float> (rate of mutations, including substitutions and indels, default = 0.001)
-R <float> (fraction of indels, default = 0.15)
-X <float> (probability an indel is extended, default = 0.3)
-H <float> (homozygous variant ratio, default = 0.3333)
-h (haplotype mode (all reads are sequenced from a single sequence)
カバレッジ30で、SNVとindelが含まれたgermlineのfastqをシミュレートする。
SNVSniffer gsim reference.fa output -e 0.01 -c 30
ランが終わるとoutput_1.fqとoutput_2.fqができる。下のキャプチャ画像は bamに変換してIGVで表示してみたもの。
4、ssim - Simulate normal and tumor sample reads with somatic SNVs/indels
$ SNVSniffer ssim
SNPSniffer ssim: a diaploid somatic mutation simulator
Usage: SNVSniffer ssim [options] reference.fa reads.base
reference.fa: reference file stored in FASTA/FASTQ format
reads.base: the base file name for reads and VCF-formated mutations
Options:
Output:
-z (compress all of the output using ZLIB)
Random numbers:
u <int> (random number seed value for the simulated genome, default=11)
v <int> (random number seed value for generating somatic mutations on the simulated gneome, default=149)
w <int> (random number seed value for generating reads from the simulated tumor genomes, default=97)
Reads:
-c <int> (coverage of the genome, default = 0 [>0 distables option -n])
-n <int> (number of read pairs, default = 1000000)
-1 <int> (length of the first read, default = 100)
-2 <int> (length of the second read, default = 100)
-a <float> (disgard if the fraction of ambiguous bases higher than #FLOAT, default = 0.05)
-i <int> (average insert size, i.e. outer distance, between the two ends, default = 500)
Somatic:
-E <float> (sequencing base error rate from simulated tumor genome, default = 0.02)
-s <float> (somatic mutation rate, default = 0.01)
-t <float> (ratio of homozygous variants (the rest are heterozygous), default = 0.3333)
-p <float> (probabiltiy of observing tumor reads at any somatic site [bionomial distribution], default = 0.9)
-q <float> (fraction of somatic indels, default=0.15)
-x <float> (probability an somatic indel is extended, default = 0.3)
-h <float> (homozygous rate for somatic variants, default = 0.3333)
Germline:
-e <float> (sequencing base error rate from simulated normal genome, default = 0.02)
-d <int> (standard deviation of insert size, default = 50)
-r <float> (rate of mutations, including substitutions and indels, default = 0.001)
-R <float> (fraction of indels, default = 0.15)
-X <float> (probability an indel is extended, default = 0.3)
-H <float> (homozygous variant rate, default = 0.3333)
使い方はgsimと同様。
SNVSniffer ssim reference.fa reads_base -e 0.01 -c 30
参考
SNVSnifferの解説ページに、SNVSniffer以外のツールについても、normalとtumorを比較しtumor specificな変異を検出するコマンドが記載されている。覚書としてここに載せておきます。
準備
正常組織。
#解説ページではbwa memでなくbwa alnでマッピングしている。
bwa aln -t 4 reference.fa normal_1.fa.gz >tmp1.sai
bwa aln -t 4 reference.fa normal_2.fq.gz >tmp2.sai
#saiからpaired.samを作り、bamに変換。
bwa sampe -r '@RG id:gatech SM:illumina' -a 700 reference.fa tmp1.sai tmp2.sai normal_1.fq.gz normal_2.fq.gz | samtools view -bhS -@ 10 - >normal.bam
#sortとindex
samtools sort -@ 10 normal.bam > normal_sorted.bam
samtools index normal_sorted.bam
#pileupは時間がかかる
samtools mpileup -s -C 50 -f reference.fa normal_sorted.bam >normal.pileup
腫瘍サンプルのマッピングとpileupも同様の手順で進める。
tumorの変異コール。
Varscan
java -Xmx4g -jar VarScan.jar somatic normal.mpileup tumor.mpileup out.vcf -output-vcf 1 --strand-filter --tumor-purity 0.9 --min-coverage-normal 3 --min-coverage-tumor 3
JointSNVMix ( paper)
java -Xmx4g -jar VarScan.jar somatic normal.mpileup tumor.mpileup out.vcf -output-vcf 1 --strand-filter --tumor-purity 0.9 --min-coverage-normal 3 --min-coverage-tumor 3
SomaticSniper (HP)
bam-somaticsniper -F vcf -f reference.fa tumor.sorted.bam normal.sorted.bam out.vcf
Mutect
java -Xmx16g -jar muTect.jar --analysis_type MuTect --reference_sequence reference.fa --input_file:normal normal.sorted.bam --input_file:tumor tumor.sorted.bam --vcf out.vcf
SNVSniffer
SNVSniffer somatic normal_header.sam tumor_header.sam normal.sorted.pileup tumor.sorted.pileup -o out.vcf
引用
SNVSniffer: an integrated caller for germline and somatic single-nucleotide and indel mutations
Liu Y, Loewer M, Aluru S, Schmidt B
BMC Syst Biol. 2016 Aug 1;10 Suppl 2:47