Question: Calling SNPs with samtools: how to specify how many supporting observations are required?
1
gravatar for Lina F
2.2 years ago by
Lina F150
Boston, MA
Lina F150 wrote:

Hi all,

I am interested in calling SNPs for some alignment data that I have. I am looking at solutions with samtools mpileup and freebayes.

I noticed that freebayes has an option to specify how many observations have to support a variant before it is called (the -C option).

Is there something similar for the samtools mpileup option?

I've been using samtools as follows:

samtools mpileup -uf reference.fasta Sample_sorted.bam -d 8000 | bcftools call -mv -Oz -o variants.vcf

However, I don't know how many observations supported the resulting variants.

An example SNP call from the above command is as follows:

#CHROM  POS ID  REF ALT QUAL    FILTER  INFO    FORMAT  Sample_sorted.bam
Ecoli_genome    5159    .   A   G   228 .   DP=303;VDB=0.658475;SGB=-0.693147;RPB=1;MQB=1;MQSB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=0,1,120,171;MQ=60   GT:PL   1/1:255,255,0
Ecoli_genome    5162    .   A   T   228 .   DP=310;VDB=0.633206;SGB=-0.693147;MQSB=1;MQ0F=0;AC=2;AN=2;DP4=0,0,120,172;MQ=60 GT:PL   1/1:255,255,0
...

Any insight is appreciated!

ADD COMMENTlink modified 2.2 years ago by Brian Bushnell16k • written 2.2 years ago by Lina F150
4
gravatar for Brian Bushnell
2.2 years ago by
Walnut Creek, USA
Brian Bushnell16k wrote:

"DP4=0,1,120,171" tells you the number of reads at that position, in this format:

DP4=refPlus,refMinus,varPlus,varMinus

So, in that case (for the first variant) there was 1 read supporting the reference, and 120+171=291 supporting the variant.

Incidentally, if you're trying out variant callers, the BBMap package also has one callvariants.sh) that is extremely fast and easy to use. It has a flag "minreads" that can be set to specify the number of observations required to call a variant, as well as various other filters.

ADD COMMENTlink written 2.2 years ago by Brian Bushnell16k

Thanks for pointing me towards the DP4 field!

ADD REPLYlink written 2.2 years ago by Lina F150
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1743 users visited in the last hour