SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
GATK "MAPQ should be 0 for unmapped read" complaint for mapped read efoss Bioinformatics 12 10-18-2012 01:19 AM
In SAM, how can I know a read is mapped ambiguously? sorrychen Bioinformatics 18 08-29-2012 01:30 PM
Tophat suppress the output of read pairs which map to different chromosome? ruping RNA Sequencing 1 10-14-2011 11:56 AM
target mapped read percentage eren Illumina/Solexa 0 08-13-2011 10:06 PM
bowtie: mapped read not in paired alignment bowtietester Bioinformatics 3 04-19-2010 04:17 AM

Reply
 
Thread Tools
Old 08-23-2011, 07:55 AM   #1
davidehs
Junior Member
 
Location: Maryland

Join Date: Jun 2011
Posts: 4
Default How to get the output for each mapped read from TopHat?

For the Tophat output file "junctions.bed", every line represents a single junction when mapping short sequencing reads to the genome. Each junction consists of two connected BED blocks, where each block is as long as the maximal overhang of any read spanning the junction (from TopHat manual).

For the following line,

chr20 257975 259040 JUNC00000002 6 - 257975 259040 255,0,0 2 49,70 0,995

it means 6 reads are mapped to this junction site. The limits of the two blockes, 257975 to 258024 (Block1) and 258970 to 259040 (Block2), were defined by any of the reads that spanning farthest at this junction site. Here is the problem. I want to know the start and end position for each of the 6 reads. So, the one line above may need to be expanded to 6 lines, with each line specifying the mapping position for each read.

The present Tophat can not generate the results I need. It seems that tophat produces some temporary files when running. Is there any way I can generate the results I need from these temporary file? Does anybody has similar experience, or knows any other program that can give the results I want?

Any suggestions will be greatly appreciated!
davidehs is offline   Reply With Quote
Old 08-23-2011, 09:42 AM   #2
ShaunMahony
Member
 
Location: University Park, PA

Join Date: Apr 2008
Posts: 27
Default

Parse the BAM output file for spliced reads that overlap your blocks.
ShaunMahony is offline   Reply With Quote
Old 08-23-2011, 02:35 PM   #3
davidehs
Junior Member
 
Location: Maryland

Join Date: Jun 2011
Posts: 4
Default

Hi ShaunMahony,

Thanks very much for your reply.

I am not familiar with BAM. Could you please give a bit more detailed explanation?

I appreciate it.


David
davidehs is offline   Reply With Quote
Old 08-23-2011, 03:35 PM   #4
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default

I don't know TopHat, but BEDTools can take a bam file, and intersect it with a .bed file, and output the .bam entries that overlap the .bed file intervals.

Something like:

BEDTools/bin/intersectBed -abam TopHat.bam -b junctions.bed > junctions.bam
swbarnes2 is offline   Reply With Quote
Old 08-23-2011, 05:11 PM   #5
ShaunMahony
Member
 
Location: University Park, PA

Join Date: Apr 2008
Posts: 27
Default

swbarnes2 -- yes that is a good idea, but the bed files produced by tophat describe splice junctions. So doing a simple overlap with a BAM could return all reads that are between the junction coordinates (including over introns and skipped exons), not just the reads that overlap the exon-exon boundary. I haven't used BEDTools, so I don't know how it will handle junction bed files.

davidehs -- I don't know of any other software tool to do exactly what you want, but it would be easy to write a simple script to parse the BAM file for junction reads that overlap block coordinates of interest (converting the BAM to SAM using samtools would allow you to text parse in perl or python, etc).

On the other hand, if you are just interested in a handful of junctions, why don't you just view the BAM file in a genome browser (UCSC?) and that should show you where all the reads start.
ShaunMahony is offline   Reply With Quote
Old 08-24-2011, 04:07 AM   #6
davidehs
Junior Member
 
Location: Maryland

Join Date: Jun 2011
Posts: 4
Default

It seems that I have to spend some time to get familiar with SAM, BAM, and BEDTools. They look pretty tough. Hope I can get through them soon.

Thanks for you guys help.

David
davidehs 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 10:55 PM.


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