Question: Incorrect Reference Base In 'Samtools Mpileup' Output
gravatar for Sergei Lebedev
8.1 years ago by
Saint-Petersburg, Russia
Sergei Lebedev30 wrote:

Hello, I'm trying to call variants using samtools and bcftools and (for unknown reasons) I get incorrect base in the REF column of the resulting VCF file; example:

#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  ../algn/algn_400K.sorted.bam
...    226300  .       T       C       119     .       DP=465;VDB=0.0418;AF1=0.5;AC1=1;DP4=135,80,131,105;MQ=20;FQ=105;PV4=0.13,0.27,1,1       GT:PL:GQ    0/1:149,0,132:99

But when I check the above position in the reference, the actual base seems to be G:

In [18]: from Bio import SeqIO
In [19]: [ref] = SeqIO.parse(open("MG1655-K12.fasta"), "fasta")
In [20]: ref[226300] 
Out[20]: 'G'

Here's a summary of what I did to obtain the VCF:

$ samtools mpileup -uf MG1655-K12.fasta algn_400K.sorted.bam > 400K.raw.bcf
$ bcftools view -bvcg 400K.raw.bcf > 400K.vcf

Any hints on why this might happen?

vcf samtools bcftools • 3.4k views
ADD COMMENTlink written 8.1 years ago by Sergei Lebedev30
gravatar for Pierre Lindenbaum
8.1 years ago by
France/Nantes/Institut du Thorax - INSERM UMR1087
Pierre Lindenbaum129k wrote:

The bases in the VCF are numbered starting with 1 as the first base, NOT '0' and I suppose that SeqIO uses '0' for the first base, NOT 1.

ADD COMMENTlink written 8.1 years ago by Pierre Lindenbaum129k

Thanks, that's exactly the case!

In [23]: ref[226299]
Out[23]: 'T'
ADD REPLYlink written 8.1 years ago by Sergei Lebedev30
gravatar for Stephwen
8.1 years ago by
Stephwen140 wrote:

Perhaps we could further complete the information present in this table : and save it on some wiki somewhere (seqanswers?), as this kind of question seems to arise periodically.

ADD COMMENTlink modified 8.1 years ago • written 8.1 years ago by Stephwen140

Perhaps a better solution would be to unify the tools and make the table smaller? ;) Gregor

ADD REPLYlink written 8.1 years ago by Gregor Rot440
Please log in to add an answer.


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