Question: Does BBMap not include headers in BAM files?
0
gravatar for Sinji
2.3 years ago by
Sinji2.8k
UT Southwestern Medical Center
Sinji2.8k wrote:

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 • 871 views
ADD COMMENTlink modified 2.3 years ago • written 2.3 years ago by Sinji2.8k
1

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

ADD REPLYlink modified 2.3 years ago • written 2.3 years ago by genomax70k

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 REPLYlink modified 2.3 years ago • written 2.3 years ago by Sinji2.8k
1
gravatar for Sinji
2.3 years ago by
Sinji2.8k
UT Southwestern Medical Center
Sinji2.8k wrote:

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

ADD COMMENTlink written 2.3 years ago by Sinji2.8k
1
gravatar for Devon Ryan
2.3 years ago by
Devon Ryan91k
Freiburg, Germany
Devon Ryan91k wrote:

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 COMMENTlink modified 2.3 years ago • written 2.3 years ago by Devon Ryan91k

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 REPLYlink modified 2.3 years ago • written 2.3 years ago by Sinji2.8k

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

ADD REPLYlink written 2.3 years ago by Devon Ryan91k

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

ADD REPLYlink written 2.3 years ago by Sinji2.8k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 549 users visited in the last hour