SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Tophat uniquely mapped reads mrfox Bioinformatics 2 05-23-2013 05:58 AM
How to get uniquely mapped reads from Tophat subeet Bioinformatics 10 11-28-2012 05:56 AM
not uniquely mapped reads unidodo RNA Sequencing 2 04-22-2011 02:07 PM
DEGseq and uniquely mapped reads PFS Bioinformatics 1 07-14-2010 09:37 AM
cufflinks and non-uniquely mapped reads clariet Bioinformatics 1 05-08-2010 11:13 AM

Reply
 
Thread Tools
Old 12-14-2011, 04:11 PM   #1
eilosei
Member
 
Location: New York

Join Date: Nov 2011
Posts: 19
Unhappy 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!
eilosei is offline   Reply With Quote
Old 12-15-2011, 03:57 AM   #2
arvid
Senior Member
 
Location: Berlin

Join Date: Jul 2011
Posts: 156
Default

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.
arvid is offline   Reply With Quote
Old 12-21-2011, 03:51 PM   #3
jbrwn
Member
 
Location: Denver, CO

Join Date: Mar 2011
Posts: 37
Default

Quote:
Originally Posted by eilosei View Post
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!
grab header
Code:
samtools view -H .bam > new.sam
stick on the unique reads
Code:
samtools view .bam | grep ... >> new.sam
convert to bam
Code:
samtools view -Sb -o unique.bam new.sam
also check out samtools reheader.
jbrwn is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 12:35 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO