Problem to use Bedtools after filtering uniquely mapped reads with samtools
 

Hi there,

I've aligned my ChIP-seq data with BWA and filtered for uniquely mapped reads with this command:

samtools view file.bam | grep "XT:A:U" > file.unique.bam

Then I tried to convert the file.unique.bam to .bed file with Bedtools but failed. The command I used is:

bamToBed -i file.unique.bam > file.bed

And I got error message like:

BgzfStream ERROR: read block failed - invalid block header
BamHeader ERROR: could not read magic number
BamReader ERROR: Could not load header data for file.unique.bam

Does anyone have this problem before? How can I solve it?

Thanks!

There are two reasons why it doesn't work:

1. "file.unique.bam" is not in the BAM format (look at "samtools view" for correct usage).
2. You stripped the SAM headers with your grep command.

You might need find a better way of filtering for uniquely mapped reads which preserves the SAM headers, or you add the SAM headers back to the header-stripped SAM file you created with your first command afterwards (look at "samtools view").

Then, convert it to BAM with SAMtools and feed it into BEDTools.

grab header
stick on the unique reads
convert to bam
also check out samtools reheader.