Set-up DESeq2 analysis - Answer Key

Author

Will Gammerdinger, Noor Sohail

Published

September 5, 2025

Exercise 1

Using the code above as a template, run the over-representation analysis on the significantly down-regulated genes from the pseudobulk analysis.

# Extract significant results for down-regulated
sigDown <- dplyr::filter(dge_deseq2_noNAs, padj < 0.05, log2FoldChange < 0)
sigDown_genes <- as.character(sigDown$gene)

# Run GO enrichment analysis 
egoDown <- enrichGO(gene = sigDown_genes, 
                    universe = all_genes,
                    keyType = "SYMBOL",
                    OrgDb = org.Mm.eg.db, 
                    ont = "BP", 
                    pAdjustMethod = "BH", 
                    qvalueCutoff = 0.05, 
                    readable = TRUE)

# Output results from GO analysis to a table
cluster_summaryDown <- data.frame(egoDown)

# Save results
write.csv(cluster_summaryDown, "results/clusterProfiler_VSM_TNvsCold7_downregulated.csv")

How many significant terms do you find?

# Number of significantlt down-regulated terms
nrow(cluster_summaryDown)

[1] 108

What are some of the prominent biological processes that are observed?

head(cluster_summaryDown$Description)

[1] "cytoplasmic translation"                     
[2] "translation at presynapse"                   
[3] "translation at synapse"                      
[4] "translation at postsynapse"                  
[5] "regulation of mRNA splicing, via spliceosome"
[6] "regulation of RNA splicing"

Exercise 2

Now that we have run through functional analysis with the results from Pseudobulk DE, let’s see what results we derive from the DGE lists from our FindMarkers DE analysis.

Create a significant DE genes data frame from the FindMarkers results with an added fold change criteria to reduce the gene list size. You can do this by running the code below:

# Create significant DE genes from FindMarkers()
sig_fc_dge <- dge_vsm %>% dplyr::filter(p_val_adj < 0.05, abs(avg_log2FC) > 1)

Use this gene list to run an over-representation analysis. Be sure to separate genes into up- and down-regulated first. Also keep in mind that the background gene dataset is different than for the DESeq2 analysis.

# Create background dataset for hypergeometric testing using all tested genes for significance in the results
all_genes_fm <- dge_vsm$gene

# Extract significant results for up- and down-regulated
sigUp_fm <- dplyr::filter(sig_fc_dge, avg_log2FC > 0)
sigUp_fm_genes <- sigUp_fm$gene
sigDown_fm <- dplyr::filter(sig_fc_dge, avg_log2FC < 0)
sigDown_fm_genes <- sigDown_fm$gene

# Run GO enrichment analysis 
egoUp_fm <- enrichGO(gene = sigUp_fm_genes, 
                     universe = all_genes_fm,
                     keyType = "SYMBOL",
                     OrgDb = org.Mm.eg.db, 
                     ont = "BP", 
                     pAdjustMethod = "BH", 
                     qvalueCutoff = 0.05, 
                     readable = TRUE)

egoDown_fm <- enrichGO(gene = sigDown_fm_genes, 
                       universe = all_genes_fm,
                       keyType = "SYMBOL",
                       OrgDb = org.Mm.eg.db, 
                       ont = "BP", 
                       pAdjustMethod = "BH", 
                       qvalueCutoff = 0.05, 
                       readable = TRUE)

# Output results from GO analysis to a table
cluster_summaryUp_fm <- data.frame(egoUp_fm)
cluster_summaryDown_fm <- data.frame(egoDown_fm)

What are the top terms enriched among up-regulated genes?

# Print top up-regulated terms from FindMarkers()
head(cluster_summaryUp_fm$Description)

[1] "extracellular matrix organization"            
[2] "extracellular structure organization"         
[3] "external encapsulating structure organization"
[4] "regulation of cell adhesion"                  
[5] "blood vessel development"                     
[6] "cell-substrate adhesion"

What are the top terms enriched among down-regulated genes?

# Print top down-regulated terms from FindMarkers()
head(cluster_summaryDown_fm$Description)

[1] "muscle tissue development"          "fat cell differentiation"          
[3] "brown fat cell differentiation"     "striated muscle tissue development"

How do these results compare with what we observed from the Pseudobulk DE functional analysis?

In the pseudobulk results, we also saw terms related to the extracellular matrix and cell adhesion upregulated. However, the downregulated results are different: for FindMarkers we see muscle and fat development, whereas for pseudobulk we saw RNA splicing and translation downregulated.

--- title: "Set-up DESeq2 analysis - Answer Key" authors: "Will Gammerdinger, Noor Sohail" date: "Friday, September 5, 2025" editor_options: markdown: wrap: 72 --- ```{r} #| label: load_data_libraries #| echo: false # Load libraries library(tidyverse) library(clusterProfiler) library(org.Mm.eg.db) library(msigdbr) library(DT) dge_deseq2 <- read.csv("results/deseq2_VSM_cold7_vs_TN.csv") # Untested genes have padj = NA, so let's keep genes with padj != NA dge_deseq2_noNAs <- filter(dge_deseq2, padj != "NA" ) # Create background dataset for hypergeometric testing using all tested genes for significance in the results all_genes <- as.character(dge_deseq2_noNAs$gene) dge_vsm <- read.csv("results/findmarkers_vsm_cold7_vs_TN.csv") # Make sure there is a column called gene dge_vsm <- dplyr::rename(dge_vsm, c("gene"="X")) ``` # Exercise 1 Using the code above as a template, run the over-representation analysis on the significantly down-regulated genes from the pseudobulk analysis. ```{r} #| label: ORA_sigDown # Extract significant results for down-regulated sigDown <- dplyr::filter(dge_deseq2_noNAs, padj < 0.05, log2FoldChange < 0) sigDown_genes <- as.character(sigDown$gene) # Run GO enrichment analysis egoDown <- enrichGO(gene = sigDown_genes, universe = all_genes, keyType = "SYMBOL", OrgDb = org.Mm.eg.db, ont = "BP", pAdjustMethod = "BH", qvalueCutoff = 0.05, readable = TRUE) # Output results from GO analysis to a table cluster_summaryDown <- data.frame(egoDown) # Save results write.csv(cluster_summaryDown, "results/clusterProfiler_VSM_TNvsCold7_downregulated.csv") ``` * How many significant terms do you find? ```{r} #| label: sigDown_terms # Number of significantlt down-regulated terms nrow(cluster_summaryDown) ``` * What are some of the prominent biological processes that are observed? ```{r} #| label: view_prominent_BPs head(cluster_summaryDown$Description) ``` # Exercise 2 Now that we have run through functional analysis with the results from Pseudobulk DE, let's see what results we derive from the DGE lists from our [FindMarkers DE analysis](03_DEanalysis_using_FindMarkers.qmd). Create a significant DE genes data frame from the FindMarkers results with an added fold change criteria to reduce the gene list size. You can do this by running the code below: ```{r} #| label: create_sig_DE_findmarkers # Create significant DE genes from FindMarkers() sig_fc_dge <- dge_vsm %>% dplyr::filter(p_val_adj < 0.05, abs(avg_log2FC) > 1) ``` 1. Use this gene list to run an over-representation analysis. Be sure to separate genes into up- and down-regulated first. Also keep in mind that the background gene dataset is different than for the DESeq2 analysis. ```{r} #| label: run_ORA # Create background dataset for hypergeometric testing using all tested genes for significance in the results all_genes_fm <- dge_vsm$gene # Extract significant results for up- and down-regulated sigUp_fm <- dplyr::filter(sig_fc_dge, avg_log2FC > 0) sigUp_fm_genes <- sigUp_fm$gene sigDown_fm <- dplyr::filter(sig_fc_dge, avg_log2FC < 0) sigDown_fm_genes <- sigDown_fm$gene # Run GO enrichment analysis egoUp_fm <- enrichGO(gene = sigUp_fm_genes, universe = all_genes_fm, keyType = "SYMBOL", OrgDb = org.Mm.eg.db, ont = "BP", pAdjustMethod = "BH", qvalueCutoff = 0.05, readable = TRUE) egoDown_fm <- enrichGO(gene = sigDown_fm_genes, universe = all_genes_fm, keyType = "SYMBOL", OrgDb = org.Mm.eg.db, ont = "BP", pAdjustMethod = "BH", qvalueCutoff = 0.05, readable = TRUE) # Output results from GO analysis to a table cluster_summaryUp_fm <- data.frame(egoUp_fm) cluster_summaryDown_fm <- data.frame(egoDown_fm) ``` What are the top terms enriched among up-regulated genes? ```{r} #| label: up_terms # Print top up-regulated terms from FindMarkers() head(cluster_summaryUp_fm$Description) ``` What are the top terms enriched among down-regulated genes? ```{r} #| label: down_terms # Print top down-regulated terms from FindMarkers() head(cluster_summaryDown_fm$Description) ``` 2. How do these results compare with what we observed from the Pseudobulk DE functional analysis? In the pseudobulk results, we also saw terms related to the extracellular matrix and cell adhesion upregulated. However, the downregulated results are different: for FindMarkers we see muscle and fat development, whereas for pseudobulk we saw RNA splicing and translation downregulated.