After running a DRIMSeq pipeline and obtaining the genes that are differently used between the null and full models, how do I select the genes that are differently used in the different conditions?

According to the vignette, I would have to use the plotProportions function and say in the group_variable the variable of my interest. What if I want a list of the genes that are differently used specifically in the variables of my interest? Checking every gene manually sounds exhausting. Help would be much appreciated.


1 Answer 1


After running the test for differential transcript usage you get a list with gene- and transcript-level p-values that you can use to filter which differences are significant. The gene-level p-value tells you if there is differential usage of transcripts within a gene and the transcript-level p-value tells you wich transcripts actually are changing. Therefor you will want to select significant genes for further investigation.

The plots are mainly used to visualize top hits as detailed in the vignette.

res <- results(d)
# The data is ordered by p-value
res <- res[order(res$pvalue, decreasing = FALSE), ]
# The top hit is selected
top_gene_id <- res$gene_id[1]
plotProportions(d, gene_id = top_gene_id, group_variable = "group")

However you should consider running the two-stage test that is described on page 16 of the vignette you are using to reduce the number of false positive hits. From this you will get another list of adjusted p-values that can be used as above.

  • $\begingroup$ I'm very gratefull for your input. But what if I have more than one group_variable how do I know why p-values are related to that group_variable? $\endgroup$
    – BioStudent
    Commented Aug 19, 2019 at 15:39
  • $\begingroup$ When you run the test (d <- dmTest(d, coef = "groupKD", verbose = 1)) you are specifying the contrast you are investigating (in this case "groupKD"). If you for example have two variables that you want to look at independently you will have to make two DRIMSeq data sets and run dMTest two times. It would be best if you can specify your model.matrix in the question so I can give you more specific answers. $\endgroup$
    – PPK
    Commented Aug 19, 2019 at 16:18
  • $\begingroup$ Let's say I have a design matrix which is ~ treatment + celltype and my null model is ~1 . After running the pipeline with this design, how do I get a list of the genes whose transcripts are differentially used for the treatment or for the celltype ? I don't want to use two different models to get the DTU for each variable as the above additive model gave better gene precision estimates than models with individual variable. $\endgroup$
    – BioStudent
    Commented Aug 27, 2019 at 14:37
  • $\begingroup$ If you have two columns in your design matrix you can use the coef approach. for example dmTest(d, coef = "treatmentB, verbose = 1) will give you the p values for comparing treatment B vs A (in this example A is the reference level, if you want to change that you can refactor the columns in the column data). But please note that in this case the null model becomes ~ cellltype which is what you need to isolate the treatment from the celltype effect. $\endgroup$
    – PPK
    Commented Sep 10, 2019 at 11:24

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.