S4 class to store data from differentially expression analysis. It should be compatible with different package and stores the information in a way the methods will work with all of them.

DEGSet(resList, default)

DEGSet(resList, default)

DEGSetFromEdgeR(object, ...)

DEGSetFromDESeq2(object, ...)

# S4 method for TopTags
DEGSetFromEdgeR(object, default = "shrunken",
  extras = NULL)

# S4 method for DESeqResults
DEGSetFromDESeq2(object, default = "shrunken",
  extras = NULL)

Arguments

resList

List with results as elements containing log2FoldChange, pvalues and padj as column. Rownames should be feature names. Elements should have names.

default

The name of the element to use by default.

object

Different objects to be transformed to DEGSet.

...

Optional parameters of the generic.

extras

List of extra tables related to the same comparison.

Details

For now supporting only DESeq2::results() output. Use constructor degComps() to create the object.

The list will contain one element for each comparison done. Each element has the following structure:

  • DEG table

  • Optional table with shrunk Fold Change when it has been done.

To access the raw table use deg(dgs, "raw"), to access the shrunken table use deg(dgs, "shrunken") or just deg(dgs).

Examples

library(DESeq2) dds <- makeExampleDESeqDataSet(betaSD = 1) colData(dds)[["treatment"]] <- sample(colData(dds)[["condition"]], 12) design(dds) <- ~ condition + treatment dds <- DESeq(dds)
#> estimating size factors
#> estimating dispersions
#> gene-wise dispersion estimates
#> mean-dispersion relationship
#> final dispersion estimates
#> fitting model and testing
res <- degComps(dds, combs = c("condition"))
#> Doing 1 element(s).
#> Doing results() for each element.
#> Doing lcfSrink() for each element.
deg(res[[1]])
#> log2 fold change (MAP): condition B vs A #> Wald test p-value: condition B vs A #> DataFrame with 1000 rows and 6 columns #> baseMean log2FoldChange lfcSE stat pvalue padj #> <numeric> <numeric> <numeric> <numeric> <numeric> <numeric> #> gene928 132.31 3.468 0.3472 9.489 2.324e-21 2.138e-18 #> gene952 70.21 3.625 0.3834 8.875 6.976e-19 3.209e-16 #> gene644 455.54 2.266 0.2618 8.634 5.940e-18 1.822e-15 #> gene515 372.73 2.804 0.3344 8.348 6.967e-17 1.602e-14 #> gene368 37.82 -2.763 0.3557 -7.532 4.992e-14 9.185e-12 #> ... ... ... ... ... ... ... #> gene906 1.3820 -0.35406 0.5346 -1.09121 0.2752 NA #> gene943 0.0000 NA NA NA NA NA #> gene956 0.6292 -0.10818 0.3952 -0.15208 0.8791 NA #> gene963 0.4705 0.09311 0.3518 -0.01693 0.9865 NA #> gene995 0.1578 -0.11622 0.3514 -0.41438 0.6786 NA
deg(res[[1]], tidy = "tibble")
#> # A tibble: 1,000 x 7 #> gene baseMean log2FoldChange lfcSE stat pvalue padj #> * <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> #> 1 gene928 132 3.47 0.347 9.49 2.32e⁻²¹ 2.14e⁻¹⁸ #> 2 gene952 70.2 3.63 0.383 8.88 6.98e⁻¹⁹ 3.21e⁻¹⁶ #> 3 gene644 456 2.27 0.262 8.63 5.94e⁻¹⁸ 1.82e⁻¹⁵ #> 4 gene515 373 2.80 0.334 8.35 6.97e⁻¹⁷ 1.60e⁻¹⁴ #> 5 gene368 37.8 -2.76 0.356 -7.53 4.99e⁻¹⁴ 9.19e⁻¹² #> 6 gene409 68.1 -2.05 0.283 -7.07 1.56e⁻¹² 2.39e⁻¹⁰ #> 7 gene130 404 1.81 0.258 6.95 3.76e⁻¹² 4.94e⁻¹⁰ #> 8 gene948 898 2.15 0.300 6.89 5.61e⁻¹² 6.45e⁻¹⁰ #> 9 gene60 115 -1.97 0.280 -6.83 8.74e⁻¹² 8.93e⁻¹⁰ #> 10 gene605 133 1.80 0.285 6.13 8.97e⁻¹⁰ 8.25e⁻ ⁸ #> # ... with 990 more rows