9: Bulk RNA-Seq Analysis with DESeq2

Chapter 9: Bulk RNA-Seq Analysis: Differential Expression with DESeq2

DESeq2 is the gold standard package for differential gene expression analysis of RNA-seq data. It fits negative binomial generalized linear models (GLMs) with dispersion shrinkage to account for biological variability across replicates.

1. Experimental Design & DESeqDataSet Construction

library(DESeq2)

# Construct DESeqDataSet object
dds <- DESeqDataSetFromMatrix(
    countData = counts_matrix,
    colData   = col_data,
    design    = ~ Batch + Condition
)

2. Library Size Normalization (Median of Ratios)

DESeq2 calculates size factors by taking the median ratio of each sample's counts to the geometric mean across all samples:

dds <- estimateSizeFactors(dds)
print(sizeFactors(dds))

3. Running Differential Analysis & Results Extraction

# Run full DESeq pipeline
dds <- DESeq(dds)

# Extract differential results (Treated vs Control)
res <- results(dds, contrast = c("Condition", "Treated", "Control"), alpha = 0.05)
summary(res)

# Filter for significantly upregulated genes
sig_up <- subset(res, padj < 0.05 & log2FoldChange > 1.5)
head(sig_up[order(sig_up$padj), ])

4. Exploratory Data Analysis: PCA Projection

Visualizing sample clustering on variance-stabilized transformed (VST) counts:

vsd <- vst(dds, blind = FALSE)
plotPCA(vsd, intgroup = c("Condition", "Batch"))