SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
RNA-Seq: MITIE: Simultaneous RNA-Seq-based Transcript Identification and Quantificati Newsbot! Literature Watch 0 08-28-2013 02:00 AM
RNA-Seq: rQuant.web: a tool for RNA-Seq-based transcript quantitation. Newsbot! Literature Watch 0 06-17-2010 02:00 AM

Reply
 
Thread Tools
Old 01-05-2014, 01:40 PM   #1
lnzzz
Junior Member
 
Location: France

Join Date: Jan 2014
Posts: 4
Question RNA-seq : 5 to 3 transcript coverage

Hello everybody,

I would like to get the read distribution along transcripts in my RNA-seq libraries. Basically, I would like to see if I have a 5' or 3' shortening of transcripts.


I already tried without success to different tools :
1) RNA-SeQC. I have some issues to use it on my tophat bam files
2) RSeQC (geneBody_coverage.py). For this one, I get this error :

ImportError: dlopen(/Library/Python/2.7/site-packages/RSeQC-2.3.7-py2.7-macosx-10.8-intel.egg/csamtools.so, 2): Symbol not found: ___ks_insertsort_heap
Referenced from: /Library/Python/2.7/site-packages/RSeQC-2.3.7-py2.7-macosx-10.8-intel.egg/csamtools.so
Expected in: flat namespace
in /Library/Python/2.7/site-packages/RSeQC-2.3.7-py2.7-macosx-10.8-intel.egg/csamtools.so


It seems to be an installation problem but my informatic skills are too weak to know how to deal with this problem .


I'm blocked. I would very much appreciate your help. Could you help me to use RNA-SeQC or RSeQC? Perhabs do you know a best strategy to get the 5 to 3 transcript coverage?

Thank you!!!

ln


PS: I wish you all the best for the new year
lnzzz is offline   Reply With Quote
Old 01-05-2014, 04:17 PM   #2
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,978
Default

You can use the coverageBed from the bedtools program to get the coverage information: http://bedtools.readthedocs.org/en/l.../coverage.html
GenoMax is offline   Reply With Quote
Old 01-07-2014, 01:25 AM   #3
lnzzz
Junior Member
 
Location: France

Join Date: Jan 2014
Posts: 4
Default

Thank you very much.
I tried coverageBed. I got a file with the read coverage for each position in the genome. I succeeded in extracting the distribution along one particular transcript. However, what I really need, it is an estimation of the global read distribution along all transcripts. I don't really how I can get it with coverageBed .
lnzzz is offline   Reply With Quote
Old 01-07-2014, 03:36 AM   #4
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,978
Default

Can you post the command you used for coverageBed?

BTW: Did you try the test dataset for RSeQC with your local install? Is that generating an error? Looking at your file paths it appears that you are using OS X (10.8?)

Last edited by GenoMax; 01-07-2014 at 03:44 AM.
GenoMax is offline   Reply With Quote
Old 01-07-2014, 04:38 AM   #5
lnzzz
Junior Member
 
Location: France

Join Date: Jan 2014
Posts: 4
Default

Yes I used a bed file with the coordinate of my genomes. Here is the command I use for coverage:
/Users/bedtools2/bin/coverageBed -d -s -split -abam 621_hits.bam -b TAIR10.bed >coverage.txt

This command gave me a very big file (20Go)... As this file obtained with coveragebed was huge (20Go), I also use the genomeCoverageBed and obtained the depth at each position:

/Users/bedtools2/bin/genomeCoverageBed -bg -split -trackline -ibam 621_hits.bam -g TAIR10.bed > coverage2.txt

What I would like is the global read distribution along transcript. That is to say, what is the read percentage in the first 10% bases of the transcripts, in the next 10%… (Enclosed an example of what I could get for one transcript and I would to obtain for all the transcripts.)
Attached Images
File Type: png readdensity.png (80.1 KB, 28 views)
lnzzz is offline   Reply With Quote
Old 01-07-2014, 06:49 AM   #6
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,978
Default

The -d option reports coverage at each position. Can you try the coverageBed without the -d and -split? Your transcriptome BED file is start/stop positions of exons? If you are set on the 10% intervals then you may need to create a custom BED file.

Last edited by GenoMax; 01-07-2014 at 06:51 AM.
GenoMax is offline   Reply With Quote
Old 01-07-2014, 07:03 AM   #7
lnzzz
Junior Member
 
Location: France

Join Date: Jan 2014
Posts: 4
Default

Yes, I think I need to create custom bed file.

Thank you for your help
lnzzz is offline   Reply With Quote
Old 01-07-2014, 09:25 AM   #8
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default

Quick and dirty way
This works best if you align to a list of transcripts, instead of genome. Sure, it's not quite as accurate as aligning to genome with TopHat, but you don't need exact figures, just a ballpark.

1) get a list of your transcripts, and how long each one is. (If you align to a list of transcripts, samtools idxstats will tell you this)
2) take line of your sam file, and associate it with a transcript (If you align to a list of transcripts, each line will already have that info)
3) Go through each line of the .sam, and change the alignment position to a corrected integer position that is the position / total length of the transcript
4) Bin up all your new positions.

Last edited by swbarnes2; 01-07-2014 at 09:28 AM.
swbarnes2 is offline   Reply With Quote
Reply

Tags
rna-seq, transcript coverage

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 09:59 PM.


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