I am trying to calculate p-values and ultimately FDR values for RI introns for 2 conditions of experiment I have(control and knockdown)

For every condition I have 3 biological replicates each: 3 control replicates and 3 knockdown replicates.

This is the header line of my input file.

chrom   Istart  Iend    strand  avg-L_R_counts-1    Sk_counts-1 avg-L_R_counts-2    Sk_counts-2

-1 and -2 are 2 sample conditions control and knockout respectively

avg-L_R_counts = for a particular intron, this is the average of reads 5bp up stream and downstream of the 5' end and at the 3'end of the intron

Sk_counts = This is the number of reads spanning exon-exon junction for that intron.

I have approximately 350 introns(this is yeast). I can get the Percent Spliced In values but I am interested in getting a p-value and then I can calculate the FDR values. Can anyone point me as to what approach I should take to calculate a p-value, or what test should I use?

  • 1
    $\begingroup$ What is your null hypothesis? What do you want to test? That one condition has more introns than the other? $\endgroup$
    – llrs
    Commented Apr 30, 2019 at 7:20
  • $\begingroup$ My null hypothesis is that there is no intron retention in either conditions $\endgroup$ Commented Apr 30, 2019 at 16:18

3 Answers 3


Here is how we do this:

First break your gene model down into non-overlapping chunks, such that any possible gene model can be built from a combination of those chunks. E.g.

transcript 1 |>>>>>>>>>|----------|>>>>>>>>|---|>>>>|
transcript 2 |>>>>>>|-------------|>>>>>>>>|---|>>>>|
transcript 3 |>>>>>>|-------------|>>>>>>>>>>>>>>>>>|
chunks       |+++1++|+2|+++++3++++|+++4++++|+5+|+6++|

this can be achieved using cgat gtf2gtf --method=genes-to-unique-chunks from the cgat-apps toolkit if you don't want to code it yourself.

I call
chunks 1, 4 and 6 constitutive exon
chunk 3 constitutive intron chunk 2 alternate intron
chunk 5 annotated retained intron.

Now use featureCounts to count the number of reads in each chunk.

Test the chunk counts for usage using DEXSeq.

Filter the results to just include chunks of type 2, 3 and 5 and recompute the FDRs.

  • $\begingroup$ If you have your transcripts in BED format, you can use bedops --partition to make chunks. $\endgroup$ Commented Oct 1, 2019 at 19:00

Did you look into literature? I am not an expert in this field, but a small search on google leads me to IRFinder (ref.), a tool that can assess intron retention. My advice is to see if this tool can help you, before reinventing the wheel...


You might take a look at the MISO algorithm for inspiration, or as a method for calculating confidence intervals of alternatively spliced isoforms: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3037023/


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.