Skip to contents

Wrapper Function for DESeq2 Analysis

wrap_dds.Rd

This function performs differential expression analysis using DESeq2 based on the provided DESeqDataSet (dds) object. It calculates the dispersion values from the dds object and then performs inference on the log-fold change (LFC) values using the specified parameters.

Usage

wrap_dds(dds, lfcThreshold, altHypothesis, correction_method = "BH")

Arguments

dds

A DESeqDataSet object containing the count data and experimental design.

lfcThreshold

The threshold for minimum log-fold change (LFC) to consider differentially expressed.

altHypothesis

The alternative hypothesis for the analysis, indicating the direction of change. Options are "greater", "less", or "two.sided".

correction_method

The method for p-value correction. Default is "BH" (Benjamini-Hochberg).

Value

A list containing the dispersion values and the results of the differential expression analysis. The dispersion values are calculated from the dds object and named according to sample names. The inference results include adjusted p-values and log2 fold changes for each gene.

Examples

N_GENES = 100
MAX_REPLICATES = 5
MIN_REPLICATES = 5
BASAL_EXP = 10
## --init variable
input_var_list <- init_variable( name = "genotype", sd = 1, level = 3) %>%
                   init_variable(name = "environment", sd = 0.9 , level = 2) 
#> Variable name should not contain digits, spaces, or special characters.
#> If any of these are present, they will be removed from the variable name.
#> Variable name should not contain digits, spaces, or special characters.
#> If any of these are present, they will be removed from the variable name.
mock_data <- mock_rnaseq(input_var_list, N_GENES, MIN_REPLICATES, 
                         MAX_REPLICATES, basal_expression = BASAL_EXP)
#> Building mu_ij matrix
#> k_ij ~ Nbinom(mu_ij, dispersion)
#> Counts simulation: Done
dds <- DESeq2::DESeqDataSetFromMatrix(mock_data$counts , 
                   mock_data$metadata, ~ genotype + environment)
#> converting counts to integer mode
dds <- DESeq2::DESeq(dds, quiet = TRUE)
#> -- note: fitType='parametric', but the dispersion trend was not well captured by the
#>    function: y = a/x + b, and a local regression fit was automatically substituted.
#>    specify fitType='local' or 'mean' to avoid this message next time.
result <- wrap_dds(dds, lfcThreshold = 1, altHypothesis = "greater")
#> INFO: The dispersion values from DESeq2 were reparametrized to their reciprocals (1/dispersion).
#> INFO: The <coeff_threshold> is adjusted by multiplying it by log(2) to convert it to log2 scale for DESeq2.
#> INFO: The log2-fold change estimates from DESeq2 were converted to the natural logarithm scale.