Tweet
Hi everyone,
I am currently using bcftools to call variants and saw some strange behavior. It would be perfect if someone could explain this to me, because I did not find any
hint from reading the paper and the bcftools documentation. The question is, looking at the mpileup I can see that roughly one third of the reads has a
C instead of a T. With a coverage of several hundred, this should not be a problem to call. However, if I do not use the option -C50 this variant is not called by bcftools. On the other hand I do not want to use -C50, because then I miss variants that are close to other variants in the genome.
Any idea?
The result is the following (I show just some part of the pileup columns, the coverage is high):
I am currently using bcftools to call variants and saw some strange behavior. It would be perfect if someone could explain this to me, because I did not find any
hint from reading the paper and the bcftools documentation. The question is, looking at the mpileup I can see that roughly one third of the reads has a
C instead of a T. With a coverage of several hundred, this should not be a problem to call. However, if I do not use the option -C50 this variant is not called by bcftools. On the other hand I do not want to use -C50, because then I miss variants that are close to other variants in the genome.
Any idea?
Code:
echo mpileup $samtools mpileup $fn_bam -f $fn_genome -E -r $chr:$((pos-$off))-$((pos+$off)) echo "mpileup | bcftools" $samtools mpileup $fn_bam -f $fn_genome -E -r $chr:$((pos-$off))-$((pos+$off)) -u | $bcftools view -cg '-' | grep -v "^#" echo "mpileup -C50 | bcftools" $samtools mpileup $fn_bam -C50 -f $fn_genome -E -r $chr:$((pos-$off))-$((pos+$off)) -u | $bcftools view -cg '-' | grep -v "^#"
Code:
mpileup 1 88212463 T 293 ,$,,,.,..,.,.,..,,,,,,,....,,...,......,.,..,...,,....,,.......,,.,....,..,,,,..,...,,,......,,..,...,....,.,,,.,,,,.,.,.......,,.....,,..,,,..,.......,.,,,,..,...,......,,,...,,.,.......,,.,,..,,... 1 88212464 T 293 ,$,$,.,..,.,.,..,,,,,,,....,,...,......,.,..,...,,....,,.......,,.,....,...,,,..,...,,,.....,,..,...,....,.,,,.,,,,.,.,.......,,.....,,..,,,..,.......,.,,,,..,...,......,,,...,,.,.......,,.,,..,,..., 1 88212465 T 296 ,$.,$..,.,C,.C,,,,,,,....,cC..cCC..C.,.,..,..,c...C,,.....C.,,.c....,..Ccc,,.C,..C,,c...C..,cC.,C..,....c.,,,.,,,,C,.c...C..C,,.....,cC.,c,..cCC...C.,.,,,c..,...,.CCC..,,c...c,C,.C.C.Cc,.,,.C,,..C,,. 1 88212466 C 299 .$.$.$,.,.,..,,,,,,,....,,...,......,.,..,...,,....,,.......,,.,....,...,,,..,...,,,.....,,..,...,....,.,,,.,,,,,.,.......,,.....,,..,,,..,.......,.,,,,..,...,......,,,...,,.,.......,,.,,..,,...,,.., 1 88212467 C 298 ,$,.,..,,,,,,,....,,...,......,.,..,...,,....,,.......,,.,....,...,,,..,...,,,......,,..,...,....,.,,,.,,,,.,.,.......,,.....,,..,,,..,.......,.,,,,.,...,......,,,...,,.,......,,.,,..,,...,,..,,...., mpileup | bcftools 1 88212463 . T . 283 . DP=302;AF1=0;AC1=0;DP4=184,109,0,0;MQ=20;FQ=-282 PL 0 1 88212464 . T . 283 . DP=306;AF1=0;AC1=0;DP4=186,107,0,0;MQ=20;FQ=-282 PL 0 1 88212465 . T . 144 . DP=313;VDB=6.428646e-02;RPB=7.216996e-01;AF1=0;AC1=0;DP4=146,86,42,22;MQ=20;FQ=-141;PV4=0.77,0.43,1,1 PL 0 1 88212466 . C . 283 . DP=314;AF1=0;AC1=0;DP4=193,106,0,0;MQ=20;FQ=-282 PL 0 1 88212467 . C . 283 . DP=316;AF1=0;AC1=0;DP4=191,107,0,0;MQ=20;FQ=-282 PL 0 mpileup -C50 | bcftools 1 88212463 . T . 283 . DP=292;AF1=0;AC1=0;DP4=178,105,0,0;MQ=43;FQ=-282 PL 0 1 88212464 . T . 283 . DP=295;AF1=0;AC1=0;DP4=179,103,0,0;MQ=43;FQ=-282 PL 0 1 88212465 . T C 60 . DP=299;VDB=6.578340e-03;RPB=-5.410734e-01;AF1=0.5;AC1=1;DP4=141,83,38,20;MQ=43;FQ=63;PV4=0.76,0.34,0.0003,1 GT:PL:GQ 0/1:90,0,255:93 1 88212466 . C . 283 . DP=299;AF1=0;AC1=0;DP4=183,101,0,0;MQ=43;FQ=-282 PL 0 1 88212467 . C . 283 . DP=299;AF1=0;AC1=0;DP4=179,102,0,0;MQ=43;FQ=-282 PL 0