<p>All introns will have some level of retention and I guess no intron will be 100% retained. In the absence of a comparator, you will need to work out how much intron retention an intron needs to be called retained. This is a human judgement, there is no automated way of deciding this. </p>

<p>What you need to do this is some sort of PSI calculation (percent-spliced-in). MISO will provide this, but it doesn't do replicates, so I guess you'll need to run it on each sample and then average the PSIs and determine if it surpasses some sort of threshold you have set. </p>

CC BY-SA 3.0
2018-01-22T08:46:14.243Z
I second this, but would suggest using [SUPPA](https://github.com/comprna/SUPPA#psi-calculation-for-each-event) which handles replicates out of the box and is actively developed. As for deciding which introns are retained, I expect that most biological conditions most of them will have a low PSI (which a nice, easily interpretable metric), so I have in the past plotted the distribution and used some upper bound of the distribution to be on the safe side. That said, some condition to use as a reference would be better.
fridaymeetssunday
CC BY-SA 3.0
2018-01-09T15:30:13.830Z
2018-01-09T15:30:13.830Z
Ian Sudbery
<p>IRFinder can be used to detect the IR events and can easily reach your goal.</p>

<p>Here is an example. IRFinder can calculate the IR ratio to detect IR event. mimirna.centenary.org.au/irfinder/example1.html. (You just need run IRFinder for one condition and each replicate will output a folder). Two conditions can be used to detect the differential IR events. If you have only one condition, you can select the IR events from the common events of the three replicates or select IR events exist at least two of three replicates based on the main output files( named IRFinder-IR-dir.txt). The IR events can be selected by filtering the IR ratio, Intron coverage etc</p>

CC BY-SA 3.0
2018-01-22T09:04:04.960Z
Hi Renhua, this answer doesn't seem to really address the issue. Could you please elaborate how IRFinder helps when there is one condition with replicates? How would that be done with IRFinder? Many thanks for your time
llrs
CC BY-SA 3.0
2018-01-22T06:07:18.690Z
2018-03-04T08:02:22.913Z
2018-03-04T08:02:22.913Z
Renhua Song
<p>The method you linked to is interesting, though I suspect that Ian's answer is both easier and generally sufficient to get what you want.</p>

<p>Having said that, the general idea behind the method you linked to is to use DESeq2 to test for differentially expressed introns. To generate the control samples, they're effectively taking the number of reads across all introns in a gene and distributing them over the introns. They weight the distribution such that longer introns receive more reads (they use an effective length for this, but I'm skeptical that this makes that much of a difference). DESeq2 is then giving the p-values you're after. Unless you can get the authors to send you there code, you'll have to write a bit of code to do all of this...but it shouldn't be too bad.</p>

<p>BTW, you can use the <code>GenomicFeatures</code> package from bioconductor, specifically the <code>exonsBy(..., by="gene)</code> follows by <code>gaps()</code>.</p>

CC BY-SA 3.0
2018-01-09T20:59:22.203Z
Hi Devon, 
Thanks for your help.The paper compares insilco replicate counts with actual counts and that is what is input for DESeq2. From your explanation, for a gene, they took all the reads which map introns and then lets say we get 1000 reads for all the possible introns for this gene and then for each intron, they give a weightage like 1000/(len of intron). I am looking at GenomicFeatures and I have already created the intron file coordinates which you mentioned at seqanswers. If possible and if you have time can you guide me in writing the above code. Thanks for your help always
Regards
user3138373
CC BY-SA 3.0
2018-01-09T21:10:13.987Z
It's not 1000/(len of intron), but rather 1000/(len of intron / sum of all intron lengths in gene).
Devon Ryan
CC BY-SA 3.0
2018-01-09T21:13:50.130Z
I only have time to provide vague guidance on the code. You can use `findOverlaps()` (from GenomicRanges in bioconductor) for the counts. Aside from that it's mostly `split()` and `mendoapply()` that you'll need.
Devon Ryan
CC BY-SA 3.0
2018-01-09T21:15:37.253Z
so that would give us insilco replicate1 count for a particular intron and similarly insilco replicate2 count for 2nd replicate and then we compare these 2 with actual counts we are getting in rep1 and rep2 via DESeq2
user3138373
CC BY-SA 3.0
2018-01-09T21:17:10.597Z
Yup, exactly. Once you have the counts its just down to using DESeq2 as normal.
Devon Ryan
CC BY-SA 3.0
2018-01-11T11:14:24.577Z
D'oh! I've actaully done this before (that'll teach me not to look at the links in the question). I used feature counts to count the reads in introns, then a python script to calculated the effective length from mappability tracks downloaded from UCSC, and then used an R script to find the "detained introns".  See these gists: https://gist.github.com/IanSudbery/df92ebf4f84d76027394750a6516db1b , https://gist.github.com/IanSudbery/eb36ba5d59b4455f9cd7b2085447eb3d. I offer no guarantees with these, they were pulled from my archives.
Ian Sudbery
CC BY-SA 3.0
2018-01-11T11:15:09.063Z
It doesn't work too bad. Better on nuclear RNAseq than cytoplasmic, but it is very sensitive. Finds introns with only a small amount of retention.
Ian Sudbery
CC BY-SA 3.0
2018-01-09T20:43:04.013Z
2018-01-09T20:43:04.013Z
Devon Ryan
<p>I have 3 biological replicates of RNA seq data for a particular condition. I want to find out intron retention events from those biological replicates for a given condition. There is no comparison I am doing here. Just from the three replicates, I am trying to find out what introns are retained in my sample. There is only one condition here. Tools like rMATS require 2 conditions. Is there any tool I can use to find out intron retention events given that I only have one condition (with replicates).</p>

<p>Ideally I would like to use a method similar to that used <a href="http://genesdev.cshlp.org/content/suppl/2014/12/29/29.1.63.DC1/SuppMaterial.pdf" rel="nofollow noreferrer">in this paper</a>, though I can't fully grasp their methodology. This would allow me to get p-values for enrichment of specific introns.</p>

<rna-seq><rna-splicing>
Intron retention events for a condition sample with 3 replicates
# Intron retention events for a condition sample with 3 replicates

(Asked by: username_0 on 2018-01-08)

I have 3 biological replicates of RNA seq data for a particular condition. I want to find out intron retention events from those biological replicates for a given condition. There is no comparison I am doing here. Just from the three replicates, I am trying to find out what introns are retained in my sample. There is only one condition here. Tools like rMATS require 2 conditions. Is there any tool I can use to find out intron retention events given that I only have one condition (with replicates).

Ideally I would like to use a method similar to that used [in this paper](http://genesdev.cshlp.org/content/suppl/2014/12/29/29.1.63.DC1/SuppMaterial.pdf), though I can't fully grasp their methodology. This would allow me to get p-values for enrichment of specific introns.

**Tags:** rna-seq, rna-splicing

---

## 3 answers

### Answer by username_1 on 2018-01-09. Score: 1

All introns will have some level of retention and I guess no intron will be 100% retained. In the absence of a comparator, you will need to work out how much intron retention an intron needs to be called retained. This is a human judgement, there is no automated way of deciding this.

What you need to do this is some sort of PSI calculation (percent-spliced-in). MISO will provide this, but it doesn't do replicates, so I guess you'll need to run it on each sample and then average the PSIs and determine if it surpasses some sort of threshold you have set.

> **Comments:**
> **username_2**: I second this, but would suggest using [SUPPA](https://github.com/comprna/SUPPA#psi-calculation-for-each-event) which handles replicates out of the box and is actively developed. As for deciding which introns are retained, I expect that most biological conditions most of them will have a low PSI (which a nice, easily interpretable metric), so I have in the past plotted the distribution and used some upper bound of the distribution to be on the safe side. That said, some condition to use as a reference would be better.  

---

### Answer by username_3 on 2018-01-22. Score: 0

IRFinder can be used to detect the IR events and can easily reach your goal.

Here is an example. IRFinder can calculate the IR ratio to detect IR event. mimirna.centenary.org.au/irfinder/example1.html. (You just need run IRFinder for one condition and each replicate will output a folder). Two conditions can be used to detect the differential IR events. If you have only one condition, you can select the IR events from the common events of the three replicates or select IR events exist at least two of three replicates based on the main output files( named IRFinder-IR-dir.txt). The IR events can be selected by filtering the IR ratio, Intron coverage etc

> **Comments:**
> **username_4**: Hi Renhua, this answer doesn't seem to really address the issue. Could you please elaborate how IRFinder helps when there is one condition with replicates? How would that be done with IRFinder? Many thanks for your time  

---

### Answer by username_5 on 2018-01-09. Score: 0

The method you linked to is interesting, though I suspect that Ian's answer is both easier and generally sufficient to get what you want.

Having said that, the general idea behind the method you linked to is to use DESeq2 to test for differentially expressed introns. To generate the control samples, they're effectively taking the number of reads across all introns in a gene and distributing them over the introns. They weight the distribution such that longer introns receive more reads (they use an effective length for this, but I'm skeptical that this makes that much of a difference). DESeq2 is then giving the p-values you're after. Unless you can get the authors to send you there code, you'll have to write a bit of code to do all of this...but it shouldn't be too bad.

BTW, you can use the `GenomicFeatures` package from bioconductor, specifically the `exonsBy(..., by="gene)` follows by `gaps()`.

> **Comments:**
> **username_0**: Hi Devon,  Thanks for your help.The paper compares insilco replicate counts with actual counts and that is what is input for DESeq2. From your explanation, for a gene, they took all the reads which map introns and then lets say we get 1000 reads for all the possible introns for this gene and then for each intron, they give a weightage like 1000/(len of intron). I am looking at GenomicFeatures and I have already created the intron file coordinates which you mentioned at seqanswers. If possible and if you have time can you guide me in writing the above code. Thanks for your help always Regards  
> **username_5**: It's not 1000/(len of intron), but rather 1000/(len of intron / sum of all intron lengths in gene).  
> **username_5**: I only have time to provide vague guidance on the code. You can use `findOverlaps()` (from GenomicRanges in bioconductor) for the counts. Aside from that it's mostly `split()` and `mendoapply()` that you'll need.  
> **username_0**: so that would give us insilco replicate1 count for a particular intron and similarly insilco replicate2 count for 2nd replicate and then we compare these 2 with actual counts we are getting in rep1 and rep2 via DESeq2  
> **username_5**: Yup, exactly. Once you have the counts its just down to using DESeq2 as normal.  
> **username_1**: D'oh! I've actaully done this before (that'll teach me not to look at the links in the question). I used feature counts to count the reads in introns, then a python script to calculated the effective length from mappability tracks downloaded from UCSC, and then used an R script to find the "detained introns".  See these gists: https://gist.github.com/IanSudbery/df92ebf4f84d76027394750a6516db1b , https://gist.github.com/IanSudbery/eb36ba5d59b4455f9cd7b2085447eb3d. I offer no guarantees with these, they were pulled from my archives.  
> **username_1**: It doesn't work too bad. Better on nuclear RNAseq than cytoplasmic, but it is very sensitive. Finds introns with only a small amount of retention.  

---

