samtools mpileup using -l or -r

[vplagnol]

Dear all,

I am having some strange issues with samtools mpileup. I am running exactly the same command using either:

-r 17:45249328-45249329
or
-l chr17c.bed

where chr17c.bed is a one line BED file:
17 45249328 45249329

The latter, which uses a bed file, is incredibly slow. We also seem to observe a bug at some rare locations of the genome but we haven't gotten to the bottom of this.

Is there a real difference between using -l or -r? Why do I observe a different behavior? Can anyone give me dome details on the way samtools works in this case?

Thank you in advance, and happy to provide more details if it helps,

Vincent


PS: the full command line looks like
samtools mpileup -ugf /ugi/home/shared/vincent/reference_genome/fasta/human_g1k_v37.fasta -D -l chr17c.bed -b my_220_bamfiles.tab

[dakl]

Hi,

I see the same thing. Looping over a set of ≈40 positions (in hg19) using -r for each position is super fast. Giving a BED or list file with the full set of 40 positions takes ≈ 20 min to process. -r finishes in seconds.

I would guess that there is a bug somewhere in the code controlling the -l parameter not using the extremely fast random access code that -r is using... just a non-qualified guess hope for a fix in a future version.

I see this regularly using samtools v0.1.12 and v0.1.18

BR
Daniel

[chrisamiller]

This thread explains the -r vs -l issue. The -l is not using the bam index:

http://sourceforge.net/mailarchive/forum.php?thread_name=443287A31897D446961707538481382A01765A%40ICTS-S-MBX7.luna.kuleuven.be&forum_name=samtools-help

[dakl]

To quote Heng in that thread: " Ideally, samtools should be more clever about that."

+1 to that. It would be awesome if -l also used the index.