[Samtools-help] Rare variants calling

[Jody Phelan]

Hi all,

I am trying to detect rare variants on a haploid organism using samtools
mpileup + bcftools call (htslib v1.4).
I have a few variants which I am expecting to see however they are not
called as variant.

When I use the command like this:

samtools mpileup -ugf ref.fa 1.bam | bcftools call --ploidy 1 -Am

I get an output like this:
Chromosome 1473246 . A G 24.8295 .
DP=80;VDB=0.44194;SGB=-0.693021;RPB=0.139427;MQB=1;MQSB=1;BQB=0.993938;MQ0F=0;AC=0;AN=1;DP4=28,23,12,15;MQ=60
GT:PL 0:255,255

From "DP4" I can see that there is quite a good 0.65:0.35 ratio for ref and
alt alleles.

Now when I remove the -A flag and add the -v flag form bcftools call like
so:
samtools mpileup -ugf ref.fa 1.bam | bcftools call --ploidy 1 -vm

I do not get this variant. So my question: Is samtools not able to call
mixed variants?
I thought the -m flag in bcftools was used for this purpose?
If I use --ploidy ==1 does it always pick the major allele?

Any help would be appreciated,
Thanks,
Jody

[i***@gmail.com]

I would suggest you use GATK rather then samtools

Post by Jody Phelan
Hi all,
I am trying to detect rare variants on a haploid organism using samtools mpileup + bcftools call (htslib v1.4).
I have a few variants which I am expecting to see however they are not called as variant.
samtools mpileup -ugf ref.fa 1.bam | bcftools call --ploidy 1 -Am
Chromosome 1473246 . A G 24.8295 . DP=80;VDB=0.44194;SGB=-0.693021;RPB=0.139427;MQB=1;MQSB=1;BQB=0.993938;MQ0F=0;AC=0;AN=1;DP4=28,23,12,15;MQ=60 GT:PL 0:255,255
From "DP4" I can see that there is quite a good 0.65:0.35 ratio for ref and alt alleles.
samtools mpileup -ugf ref.fa 1.bam | bcftools call --ploidy 1 -vm
I do not get this variant. So my question: Is samtools not able to call mixed variants?
I thought the -m flag in bcftools was used for this purpose?
If I use --ploidy ==1 does it always pick the major allele?
Any help would be appreciated,
Thanks,
Jody
------------------------------------------------------------------------------
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
_______________________________________________
Samtools-help mailing list
https://lists.sourceforge.net/lists/listinfo/samtools-help

[Peter Johansson]

Hi Jody,

Post by Jody Phelan
Hi all,
I am trying to detect rare variants on a haploid organism using
samtools mpileup + bcftools call (htslib v1.4).
I have a few variants which I am expecting to see however they are not
called as variant.
samtools mpileup -ugf ref.fa 1.bam | bcftools call --ploidy 1 -Am
Chromosome 1473246 . A G 24.8295 .
DP=80;VDB=0.44194;SGB=-0.693021;RPB=0.139427;MQB=1;MQSB=1;BQB=0.993938;MQ0F=0;AC=0;AN=1;DP4=28,23,12,15;MQ=60
GT:PL 0:255,255
From "DP4" I can see that there is quite a good 0.65:0.35 ratio for
ref and alt alleles.
Now when I remove the -A flag and add the -v flag form bcftools call
samtools mpileup -ugf ref.fa 1.bam | bcftools call --ploidy 1 -vm
I do not get this variant. So my question: Is samtools not able to
call mixed variants?
I thought the -m flag in bcftools was used for this purpose?
If I use --ploidy ==1 does it always pick the major allele?

I'm not sure what you expect here. With '--ploidy 1' the caller can
choose between A or G, no heterozygous, so it calls the most likely,
which typically is the one with most reads.

If you have a mixture of cell populations, it might work with --ploidy 2
to pick up homozygous mutations in one of the population appearing as hets.

Cheers,
Peter