bedtools multicov

[daniel....@gmail.com]

Hi,

I have a strange issue where multicov works fine on two bams independently, but fails when I try both together. e.g.:

bedtools multicov -bams 1.bam -bed test.bed

gives:

"chr1 1 10200 test 12"

bedtools multicov -bams 2.bam -bed test.bed

gives:

"chr1 1 10200 test 3"

but 

bedtools multicov -bams 1.bam 2.bam -bed test.bed

gives:

"Could not open input BAM files."

The bam files are different read lengths - could this be an issue? Files to replicate the problem are in the attached tar archive

Daniel

[Aaron Quinlan]

Hi Daniel,

Sorry to take so long to respond; I have been traveling and this took me a bit of time to track down.

Basically, this is failing because bamtools, the underlying BAM library upon which bedtools depends, enforces that each bam has the same set and order of chromosomes in the BAM headers.  The headers in your two files are almost identical, but one file has chrM before chrX and the other has chrX before chrM.

The following works for me.

samtools view -H 1.bam > header

samtools reheader header 2.bam > 2.fixed.bam

samtools index 2.fixed.bam

bedtools multicov -bams 1.bam 2.fixed.bam -bed test.bed

Best,

- Aaron

quinlanlab.org

--
You received this message because you are subscribed to the Google Groups "bedtools-discuss" group.
To unsubscribe from this group and stop receiving emails from it, send an email to bedtools-discu...@googlegroups.com.
For more options, visit https://groups.google.com/groups/opt_out.
<multicovTest.tgz>

[daniel....@gmail.com]

Hi Aaron,

Awesome thanks for the response. Tried reheader'ing and it also works for me, so problem solved. Just out of interest, is tracking this down and dealing with it difficult partly because bamtools tends not to give informative error messages when this type of thing happens? I ask because because there is some work locally to try and improve error reporting in samtools + bamtools, so I could pass this on if it were helpful. In fact, if you have a list of your favourite errors that could be reported better from bamtools let me know and I can pass it on. Although perhaps it should be up to the end user to make sure their bams are sorted correctly ... ;)

Best

Dan