testDTU {satuRn} | R Documentation |
Function to test for differential transcript usage (DTU)
testDTU(object, contrasts, plot = FALSE, sort = FALSE)
object |
A 'SummarizedExperiment' instance containing a list of objects of the 'StatModel' class as obtained by the 'fitDTU' function of the 'satuRn' package. |
contrasts |
'numeric' matrix specifying one or more contrasts of the linear model coefficients to be tested. The rownames of the matrix should be equal to the names of parameters of the model that are involved in the contrast. The column names of the matrix will be used to construct names to store the results in the rowData of the SummarizedExperiment. |
plot |
'boolean(1)' Logical, defaults to FALSE. If set to TRUE, a plot of the histogram of the empirical z-scores and the standard normal distribution will be displayed. |
sort |
'boolean(1)' Logical, defaults to FALSE. If set to TRUE, the output of the topTable test function is sorted according to the empirical p-values. |
An updated 'SummarizedExperiment' that contains the 'Dataframes' that display the significance of DTU for each transcript in each contrast of interest.
Jeroen Gilis
data(sumExp_example, package = "satuRn") library(SummarizedExperiment) sumExp <- fitDTU( object = sumExp_example, formula = ~ 0 + group, parallel = FALSE, BPPARAM = BiocParallel::bpparam(), verbose = TRUE ) group <- as.factor(colData(sumExp)$group) design <- model.matrix(~ 0 + group) colnames(design) <- levels(group) L <- matrix(0, ncol = 2, nrow = ncol(design)) rownames(L) <- colnames(design) colnames(L) <- c("Contrast1", "Contrast2") L[c("VISp.L5_IT_VISp_Hsd11b1_Endou", "ALM.L5_IT_ALM_Tnc"), 1] <- c(1, -1) L[c("VISp.L5_IT_VISp_Hsd11b1_Endou", "ALM.L5_IT_ALM_Tmem163_Dmrtb1"), 2] <- c(1, -1) sumExp <- testDTU(object = sumExp, contrasts = L, plot = FALSE, sort = FALSE)