I am working on alternative splicing (AS) events on 4 different tomato species. I am trying to find
interesting AS events.
Lets say for example, the AS event is
exon-skipping (ES). After mapping RNA-Seq reads to tomato genome, I look for
intron coordinates) where they are
spliced normally -
normal junctions (NJ) and where the 3prime end Exon is skipped. So, for every junction, I have a count of reads that map to the junction normally (exactly where the intron is and supposed splicing should occur) and count of reads where at the same junction an ES event had occurred (the 3prime exon is skipped). At the end I have a table like this for each junction (I already remove where there is no ES event in ALL 4 species).
S1 S2 S3 S4 ES 10 0 27 0 NJ 95 20 50 380
Then I do a
fisher-test on this 2*4 table and correct for multiple-testing using
Benjamini-Hochberg method (from R
multtest package) to obtain those events that are significantly different across species.
Now of course the question is, what if 1) the gene where this junction (or intron) belongs is over- (or under-) expressed between these species. For ex: S2 has only a total of 20 reads mapped. 2) How about the number of reads for these species that was sequenced? 3) what about gene length? ( as the transcript abundance is also found to be positively correlated with gene length).
So, I have to somehow normalize this data. So far, with the exception of RPKM (which I am not convinced as an appropriate measure), all other methods were about finding differential expression of genes (and demand 2 or more samples), for ex: quantile normalization, TMM, the edgeR package etc. However, I would like to normalize gene expression in each of these samples.
Does anyone have an idea how to go about it?
Thank you! Arun.