This function will extract the output of DESeq2::results() and DESeq2::lfcShrink() for multiple comparison using:

degComps(dds, combs = NULL, contrast = NULL, alpha = 0.05,
  skip = FALSE, type = "normal", pairs = FALSE, fdr = "default")

Arguments

dds

DESeq2::DESeqDataSet obcject.

combs

Optional vector indicating the coefficients or columns fom colData(dds) to create group comparisons.

contrast

Optional vector to specify contrast. See DESeq2::results().

alpha

Numeric value used in independent filtering in DESeq2::results().

skip

Boolean to indicate whether skip shrinkage. For instance when it comes from LRT method.

type

Type of shrinkage estimator. See DESeq2::lfcShrink().

pairs

Boolean to indicate whether create all comparisons or only use the coefficient already created from DESeq2::resultsNames().

fdr

type of fdr correction. default is FDR value, lfdr-stat is for local FDR using the statistics of the test, lfdr-pvalue is for local FDR using the p-value of the test

Value

DEGSet with unSrunken and Srunken results.

Details

  • coefficients

  • contrast

  • Multiple columns in colData that match coefficients

  • Multiple columns in colData to create all possible contrasts

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", 2), contrast = list("treatment_B_vs_A", c("condition", "A", "B")))
#> Doing 3 element(s).
#> Doing results() for each element.
#> Doing lcfSrink() for each element.
#> using 'normal' for LFC shrinkage, the Normal prior from Love et al (2014). #> #> Note that type='apeglm' and type='ashr' have shown to have less bias than type='normal'. #> See ?lfcShrink for more details on shrinkage type, and the DESeq2 vignette. #> Reference: https://doi.org/10.1093/bioinformatics/bty895
#> using 'normal' for LFC shrinkage, the Normal prior from Love et al (2014). #> #> Note that type='apeglm' and type='ashr' have shown to have less bias than type='normal'. #> See ?lfcShrink for more details on shrinkage type, and the DESeq2 vignette. #> Reference: https://doi.org/10.1093/bioinformatics/bty895
#> using 'normal' for LFC shrinkage, the Normal prior from Love et al (2014). #> #> Note that type='apeglm' and type='ashr' have shown to have less bias than type='normal'. #> See ?lfcShrink for more details on shrinkage type, and the DESeq2 vignette. #> Reference: https://doi.org/10.1093/bioinformatics/bty895
res <- degComps(dds,contrast = list("treatment_B_vs_A"), fdr="lfdr-stat")
#> Doing 1 element(s).
#> Doing results() for each element.
#> Step 1... determine cutoff point #> Step 2... estimate parameters of null distribution and eta0 #> Step 3... compute p-values and estimate empirical PDF/CDF #> Step 4... compute q-values and local fdr #>
#> Doing lcfSrink() for each element.
#> using 'normal' for LFC shrinkage, the Normal prior from Love et al (2014). #> #> Note that type='apeglm' and type='ashr' have shown to have less bias than type='normal'. #> See ?lfcShrink for more details on shrinkage type, and the DESeq2 vignette. #> Reference: https://doi.org/10.1093/bioinformatics/bty895