Samtools mplileup detects wrong INDEL positions
1
0
Entering edit mode
2.3 years ago
vit.filippo ▴ 50

Dear all I've a big problem with samtools mpileup since I am trying to retriving base counts and indel but fails to locate indels. More in particular, it seems that INDELs are shifted by one position (i.e. if the deletion is located at position 2, samtools gives me that is at position 1).

I've tryied bam-readcount to do the same task and it calls INDEL at the right position (after check it with IGV).

What I've tryied so far is:

  • Using different samtools versions: nothing changes
  • Using different machines: nothing changes

I would like to ask you if it is a well known bug that others found or it's just me.

I know that I could pharse the output with awk shifting indels, but I prefer not to "touch" data in that sense. Thank you in advance for your help

samtools mpileup variant calling count nuceotides • 653 views
ADD COMMENT
0
Entering edit mode

Please add the command line, an example of the VCF and the evidence that made you think the indel is shifted (=reproducible example).

ADD REPLY
0
Entering edit mode
samtools mpileup -ABQ0 -f $REF -r {} -l $BED $input

I didn't produce a .vcf, but a .tab file with counts.

RIGHT output of bam-readcount

chr1     27105518    C   4069    3   4065    1   0                       
chr1     **27105519**    A   4122    4111    6   0   0         -A    4  
chr1     27105520    G   4183    1   1   4171    10

samtools WRONG parsed

chr1            27105517    C   2033    2   0   0   2009    1   0   1   NA   NA

chr1    **27105518**    C   2046    0   0   0   2022    0   0   1   NA    2:a|2:A

chr1            27105519    A   2068    0   3   0   2043    0   3   0   NA  NA

The point is that also the mpileup file reports the deletions at position chr1:27105518 (so the code is parsing correctly compared to readcount) but when I open the .bam file on IGV it's clear that the deletion is in position 27105519 as reported by readcount

ADD REPLY
0
Entering edit mode

Hello,

the output you are showing doesn't look like the mpileup format. Please show us all commands that lead to this output.

Which version of samtools are you using?

Also please note, that samtools mpileup is deprecated. Instead use bcftools mpileup.

fin swimmer

ADD REPLY
0
Entering edit mode

I thank you in advance for your precious hepl! I handle the .mpileup file with an already tested perl script (https://github.com/riverlee/pileup2base/blob/master/pileup2baseindel.pl ). Before parsing the .mpileup I made sure that the script worked well and everything was good exept for deletions. I need to you samtools instead of bcftools since i need to count % of nucleoides for each position and it seems that this is the fastest way. I've tryied samtools 0.1.19, 0.1.18, 1.3 and 1.9 and ALL gave me same bug with ALL the deletions were shifted. Count and insertion positions are still ok.

ADD REPLY
0
Entering edit mode

To rule out that this perl script is the problem, can you please show the exact command line for mpileup and its raw output for the position. It is very unlikely that samtools is reporting incorrect indel positions given how established it is and the number of users running it to scan for SNVs for almost a decade now.

ADD REPLY
0
Entering edit mode

The command is the following:

samtools mpileup -ABQ0 -f $REF -r {} -l $BED $input

bam-readcount reports that the deletion is at position 27105519 (according to IGV), but samtools says that it is at position 27105518 The output of the .mpileup is the following:

chr1    27105518        C       3409    ..,........,.,..,,,..,.,,,,.,,..,,.,.,.,.,.....,,,,..,,...,.,.,,.....,,,,,...,,...,,,..,,.......,,...,,,.,..,,,,,.,...,,,.,,.,,,....,,,,...,,,,..,,......,,,,,.,....,,,,..,,..,,.,......,,,,,,..,,,..,....,,,..,,.,....
.,,,,.....,,,,,.,...,,,,.....,,,,...,...,,,,..,,......,,,,,,...,,,,..,,,,...,,,...,,,.,,.,,,....,,,,....,,,,,.......,,,,,,,....,,,...,,.......,,,,,,,.......,,,,,,.....,,,,,,,....,,,,,.....,,,,.....,,.........,,,,,,,,.,.,.........,,,,,,,,,,..........,,,,,,
,,,,....,,,...,,,..........,,,,,,,,,,...........,,,,,,,,,,.....,,,,,,...,,,.....,,,,,.....,,,,,.......,,,,,,,,,,............,,,,,,,,,,,,...................,,,,,,,,,,,,,,,,,...............,,,,,,,,,,,,,,,,,...............,,,,,,,,,,,,.............,,,,,,,,,,,
,......,,,,,,......,,,,,,,....,,...,,,...........-1A.,,,,,,,,,,,,,-1a,.....,,,,,..................,,,,,,,,,,,,,,,,,,............-1A.....,,,,,,,,,-1a,,,,,,....................,,,,,,,,,,,,,,,,,,........,,,,,,,.............,,,,,,,,,,.........................
,,,,,,,,,,,,,,,,,,,,,,,,,,.........,,,,,,,,,,..........,,,,,,,,,.......................,,,,,,,,,,,,,,,,,,,..................,,,,,,,,,,,,,,,,..........,,,,,,,,,........,,,,,,,,...................,,,,,,,,,,,,,,,,.......,,,,,..................,,,,,,,,,,,,,,,
,.........,,,,,,,,,,..........,,,,,,,,,.............,,,,,,,,,,,,,................,,,,,,,,,,,,,,..................,,,,,,,,,,,,,,,........................,,,,,,,,,,,,,,,,,,,,,,,............,,,,,,,,,,,,,,...............,,,,,,,,,,,,,,,,...,,,,...........,,,,,
,,,,,...................,,,,,,,,,,,,,,,,,,,...............,,,,,,,,,,,,,,..............,,,,,,,,,,,,,............,,,,,,,,,,,,............,,,,,,,,,,,,........,,,,,,,,,..................,,,,,,,,,,,,,,,,,,.....................................,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,...................,,,,,,,,,,,,,,,,,,..........................,,,,,,,,,,,,,,,,,,,,,,,,,.............,,,,,,,,,,,,,.............,,,,,,,,,,,,.............,,,,,,,,,,,,,,.............................,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,........................................,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,.....,,,,,....................,,,,,,,,,,,,,,,,,,,,........,,,,,,,,,.......................................................,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,..............,,,,,,,,,,,,,...................,,,,,,,,,,,,,,,,,,,..........................,,,,,,,,,,,,,,,,,,,,,,,,........................,,,,,,,,,,,,,,,,,,,,,........................,,,,,,,,,,,,,,,,,,,,,,...................,,,,,,,,,,,,,,,,,,,,,.......................,,,,,,,,,,,,,,,,,,,,,,......................,,,,,,,,,,,,,,,,,,,,,,,,,,..........................,,,,,,,,,,,,,,,,,,,,,,,,,.........................,,,,,,,,,,,,,,,,,,,,,,.....................,,,,,,,,,,,,,,,,,,.....................,,,,,,,,,,,,,,,,,,,,,.....................,,,,,,,,,,,,,,,,,,,,,..........................,,,,,,,,,,,,,,,,,,,,,,,,,,,.........................,,,,,,,,,,,,,,,,,,,,,,,,,..............................,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,...............,,,,,,,,,,,,,,.............................,,,,,,,,,,,,,,,,,,,,,,,,,,,................,,,,,,,,,,,,,,,,.................................,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,...................,,,,,,,,,,,,,,,,,,,,,,,.....................,,,,,,,,,,,,,,,,,,.........................,,,,,,,,,,,,,,,,,,,,,,,,,,.......................,,,,,,,,,,,,,,,,,,,,,,,,,,..................,,,,,,,,,,,,,,,,,,,....................,,,,,,,,,,,,,,,,^S.^M.^S.^S.^S.^S.^S.^S.^M.^M.^S,^M,^S,^S,^S,^S,^S,^S,^M,^M,

I cut the output ...

ADD REPLY
4
Entering edit mode
2.3 years ago

Hello,

than the perl script misinterpret the mpileup output. Deletions (-) and insertion (+) are always described from the next position on.

See:

fin swimmer

ADD COMMENT
0
Entering edit mode

Thank you very much for your help!

ADD REPLY
0
Entering edit mode

enter image description here

ADD REPLY

Login before adding your answer.

Traffic: 2473 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6