Does BBMap not include headers in BAM files?
2
0
Entering edit mode
4.3 years ago
Sinji ★ 3.1k

Been working with some ChIPseq data, and am trying to run it through ChIPQC from Bioconductor. I continually get this error:

[bam_header_read] EOF marker is absent. The input is probably truncated.
  [bam_header_read] invalid BAM binary header (this is not a BAM file).
  Error in value[[3L]](cond) : 
    failed to open BamFile: SAM/BAM header missing or empty
    file: '/usr/local/anaconda2/lib/R/bin/exec/R'
  Calls: ChIPQCsample ... tryCatch -> tryCatchList -> tryCatchOne -> <Anonymous>
  Execution halted

I then samtools view | head my bam file and sure enough I don't get any header information:

NS500579:38:HC3MLBGXX:3:22605:6021:20111    16  chr1    9980    2   49M *   0   0   GCCTACGCGTCCGCGTCCGCGTAACCCTAACCCTAACCCTAACCCTAAC   //////A//A///E/E/////6/A<E6EAEAEEEEEEEEEEEEEAAAAA   XT:A:R  NM:i:21 AM:i:2
NS500579:38:HC3MLBGXX:3:21402:17333:9302    16  chr1    9981    2   49M *   0   0   CATACGTGTGCTCTTCCGATTAACCCTAACCCTAACCCTAACCCTAACC   6//E/A//E///E/E///////EEEEEEEEEEAEEEEEEEEEEEAAAAA   XT:A:R  NM:i:20 AM:i:2
NS500579:38:HC3MLBGXX:1:11212:16057:4322    16  chr1    9981    2   49M *   0   0   GTTACTCATCAGATTACGAGTAACCCTAACCCTAACCCTAACCCTAACC   <//6//////////A/////E/////EAAE/AEEEEEEEEEEEEAA/AA   XT:A:R  NM:i:20 AM:i:2
NS500579:38:HC3MLBGXX:3:23606:5285:10927    16  chr1    10000   2   49M *   0   0   GTAACCGTTACCGTAACCGTAACCCTAACCCTAACCCTAACCCTAACCC   /////E/////E/////E/////E////EEA6EEEEEAEEEEEEAAAAA   XT:A:R  NM:i:5  AM:i:2
NS500579:38:HC3MLBGXX:3:13511:5485:4267 16  chr1    10003   2   49M *   0   0   ACGCTTCCGCTACCGCTCACCCTAACCCTAACCCTAACCCAAACCCTAA   /E/////A/AE//E//E//EEEEE/EAAEEEEEEAEEEEE/EEAAAAAA   XT:A:R  NM:i:8  AM:i:2
NS500579:38:HC3MLBGXX:1:23206:6162:17684    0   chr1    10016   3   49M *   0   0   CCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCAAACCCAAAC   AAAAAEEEEEEEEEEEEEEEEEEAEE<AE6EAE/E/E///E/A/</E//   XT:A:R  NM:i:2  AM:i:3
NS500579:38:HC3MLBGXX:1:22311:17702:3502    0   chr1    10023   3   49M *   0   0   CCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACC   AAAAAEEEEEEEEEEEEAEEEEE/E//AE/E//AE/<<///E////EE/   XT:A:R  NM:i:0  AM:i:3
NS500579:38:HC3MLBGXX:4:21408:9090:4962 16  chr1    10024   3   49M *   0   0   CTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCC   ////EEE////E/E///EEEE/E/E/E/EEEEAEAEEAEEAAEE//AAA   XT:A:R  NM:i:0  AM:i:3
NS500579:38:HC3MLBGXX:3:11512:25155:2180    0   chr1    10026   2   49M *   0   0   AACCCAAACCCTAACCCTAACCCTAACCCAACCCCTAACCCAAACCATA   AAAA//EEEAAEEAAAEAE/E/AEE//EE/////<A///E</A/EA/A/   XT:A:R  NM:i:5  AM:i:2
NS500579:38:HC3MLBGXX:4:11610:17453:12297   16  chr1    10050   3   49M *   0   0   AACCCTAACCCTAACCCAAACCCTAACCCTAACCCTAACCCTAACCCTA   //////A//////A/<//E//AE//6EAEAE//EE/EEEAAEEEAA/AA   XT:A:R  NM:i:1  AM:i:3

I figure it's something about the way i'm mapping and then converting the files, or maybe something specific with sambama, i'm not entirely sure, I can go ahead and test it out but I figure perhaps someone form Biostars can answer this quickly enough that I don't have to bother.

Here's my pipeline:

 bbmap.sh in=${read1} outm=${id}.mapped.bam outu=${id}.unmapped.bam keepnames=t trd sam=1.3 maxindel=1 ambig=random statsfile=${id}.alignmentReport.txt minid=${params.minid} usemodulo

 sambamba sort --tmpdir $baseDir -t ${params.threads} -o ${id}.sorted.mapped.bam ${id}.mapped.bam

 sambamba index -t ${params.threads} ${id}.sorted.mapped.bam

Thanks in advance!

EDIT: As Devon pointed out, I DO have headers.

bbmap • 1.4k views
ADD COMMENT
1
Entering edit mode

I doubt BBMap is the culprit here :) What version of samtools are you using?

ADD REPLY
0
Entering edit mode

Am not using samtools, instead using sambamba 0.6.5.

EDIT: Since BBMap does use samtools to convert to bam automatically, samtools version is 1.3.1.

ADD REPLY
1
Entering edit mode
4.3 years ago
Sinji ★ 3.1k

The culprit was my R script for ChIPQC doing some funny stuff with arguments. BBMap, Sambamba, and Samtools were not the problem.

ADD COMMENT
1
Entering edit mode
4.3 years ago

Your file has a header, otherwise you wouldn't see chromosome names with samtools view, where you forgot -h.

BTW, the error message suggests that you have a truncated file, since it didn't end with the magic number.

ADD COMMENT
0
Entering edit mode

Whoops. You're right. Any idea why ChIPQC isn't recognizing the header then?

EDIT: NVM, just saw your edit. Doesn't make any sense to me, the file maps just fine, and it works if I use Bowtie2 as my mapper. Odd.

ADD REPLY
0
Entering edit mode

My only guess is that sambamba produced a broken BAM file. Sort the BAM file with samtools and see what happens.

ADD REPLY
0
Entering edit mode

After sorting and indexing with samtools, it appears I still have the same problem.

ADD REPLY

Login before adding your answer.

Traffic: 2009 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