count reads on exons and introns for each transcripts

[krapulaxdoctor]

Hi,

I would like to ask for an advice:

I have an RNA-seq bam-file and i have to count the exonic reads and intronic reads for each annotated refseq transcript.

So, in the end I would like to have a table of something like this:

name - - - reands on exons - - - reads on introns

nm_00001 - - - 3214 - - - 1025

Could somebody please tell me a simple way to count reads on summa exons and summa introns for each transcripts? ( for each refseq gene ?)

I was thinking about this method:

First splitting the bam.file into only exonic reads and intronic reads using:

bedtools intersect -split -bed -a bam.file -b exonscoordinates.bed > exonicreads.bed
bedtools intersect -split -bed -a bam.file -b intronscoordinates.bed > intronicreads.bed

Then count reads on each gene from each intron/exon.bed file separately.
(however i am afraid that the exon-exon junction reads will be counted twice by this way... right?)

Thank you in advance for the help.
If i happened to miss an existing thread about this/similar problem, please let me know.

[Mchicken]

I would recommend HTSeq count (http://www-huber.embl.de/users/ander...doc/count.html). With -t you can specify the feature type on which you count (e.g. intron or exon). By making two runs, one for exons and one for introns, you should get the results you want.

[Michael.Ante]

Denote, transcripts of the same gene share at least one common exon or parts of the same. The same holds true for introns / overlapping parts of two introns.
Therefore, you cannot assign a read of a shared feature unambiguously to a certain transcript just by counting. That is one reason why programs like cufflinks or eXpress are assigning abundances to transcript with statistical routines.

If you still want to continue the read-count per exon approach, I would recommend to use DEXSeq (which has a python script for counting exon reads).

[krapulaxdoctor]

Thank you for the quick responses.
I will look after the mentioned programs.

Is there any software available that does count reads on introns the same way as it is mentioned in the case of cufflinks/eXpress?

[dpryan]

As in, is there something that will do some form of expectation maximization to get at an estimated count/TPM/FPKM? Not directly. Are you trying to probe intron retention? If so, I just talked a colleague who's doing that into just using Salmon and mapping against the transcriptome with/without introns. That seems to work well.