Samtools mplileup detects wrong INDEL positions

Question

[Gama313]

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

[ATpoint]

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

[Gama313]

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

[finswimmer]

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

[Gama313]

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.

[ATpoint]

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.

[Gama313]

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 ...

Answer 1

[finswimmer]

Hello,

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

See:

• http://samtools.sourceforge.net/pileup.shtml
• https://en.wikipedia.org/wiki/Pileup_format#Column_5:_The_bases_string

fin swimmer

[Gama313]

Thank you very much for your help!

[ATpoint]