Pick a different cell type and run FindMarkers(). Name the population you chose and print the top 6 genes in the results.
# Vector of unique cell typescell_types <- seurat_rctd$first_type %>%unique() %>%sort()# Store DGE resultslist_dge <-list()for (ct in cell_types) {# Subset to cell type seurat_sub <-subset(seurat_rctd, first_type == ct)Idents(seurat_sub) <-"orig.ident"# Calculate DGEs with wilcox test dge <-FindMarkers(seurat_sub,ident.1 ="P5CRC",ident.2 ="P5NAT") dge$gene <-rownames(dge)# Store full DGE results list_dge[[ct]] <- dge# Filter significant genes dge_sig <- dge %>%subset(p_val_adj <0.05)# Get the gene names and get the first 6 values genes <- dge_sig %>%pull(gene) %>%head(6)# Print top 6 genesprint(ct)print(genes)}
Run GSEA on the cell type that you calculated DE results for in the first exercise. What is the top pathway that appears?
# Use a specific collection; C5 GO signaturesm_t2g <-msigdbr(species ="Homo sapiens",collection ="C5") %>% dplyr::select(gs_name, gene_symbol)for (ct innames(list_dge)) { dge <- list_dge[[ct]]# Filter significant genes dge_sig <-subset(dge, p_val_adj <0.05)# Extract fold changes foldchanges <- dge_sig$avg_log2FC# Name each fold change with the corresponding gene symbolnames(foldchanges) <- dge_sig$gene# Sort fold changes in decreasing order foldchanges <-sort(foldchanges, decreasing =TRUE)# Run GSEA msig_GSEA <-GSEA(foldchanges,TERM2GENE = m_t2g,verbose =FALSE)# Extract the GSEA results msigGSEA_results <- msig_GSEA@result# Order by NES and pull first row's ID top_id <- msigGSEA_results %>%arrange(desc(abs(NES))) %>%head(1) %>%pull(ID)print(paste(ct, ":", top_id))}
---title: "Differential Gene Expression - Answer Key"author: - Noor Sohaildate: "2026-04-29"license: "CC-BY-4.0"editor_options: markdown: wrap: 72---```{r}#| label: load_libraries_data#| echo: false# Load libraries and data# Load libraries and datalibrary(Seurat)library(tidyverse)library(EnhancedVolcano)# ORA and GSEA librarieslibrary(tidyverse)library(clusterProfiler)library(org.Hs.eg.db)library(msigdbr)set.seed(12345)seurat_rctd <- qs2::qs_read("intermediate/12_seurat_RCTD.qs")```# Exercise 11. Pick a different cell type and run `FindMarkers()`. Name the population you chose and print the top 6 genes in the results.```{r}#| label: dge_cell_types# Vector of unique cell typescell_types <- seurat_rctd$first_type %>%unique() %>%sort()# Store DGE resultslist_dge <-list()for (ct in cell_types) {# Subset to cell type seurat_sub <-subset(seurat_rctd, first_type == ct)Idents(seurat_sub) <-"orig.ident"# Calculate DGEs with wilcox test dge <-FindMarkers(seurat_sub,ident.1 ="P5CRC",ident.2 ="P5NAT") dge$gene <-rownames(dge)# Store full DGE results list_dge[[ct]] <- dge# Filter significant genes dge_sig <- dge %>%subset(p_val_adj <0.05)# Get the gene names and get the first 6 values genes <- dge_sig %>%pull(gene) %>%head(6)# Print top 6 genesprint(ct)print(genes)}```# Exercise 22. Run an ORA analysis on the cell type that you calculated DE results for in the previous exercise. What is the top pathway that appears?```{r}#| label: ora_cell_typesfor (ct innames(list_dge)) { dge <- list_dge[[ct]]# Background genes all_genes <-as.character(dge$gene)# Significant up-regulated genes sigUp <- dplyr::filter(dge, p_val_adj <0.05, avg_log2FC >0) sigUp_genes <-as.character(sigUp$gene)# GO enrichment egoUp <-enrichGO(gene = sigUp_genes,universe = all_genes,keyType ="SYMBOL",OrgDb = org.Hs.eg.db,ont ="BP",pAdjustMethod ="BH",qvalueCutoff =0.05,readable =TRUE) cluster_summaryUp <-data.frame(egoUp)print(paste(ct, ":", egoUp$Description[1]))}```# Exercise 33. Run GSEA on the cell type that you calculated DE results for in the first exercise. What is the top pathway that appears?```{r}#| label: gsea_cell_types# Use a specific collection; C5 GO signaturesm_t2g <-msigdbr(species ="Homo sapiens",collection ="C5") %>% dplyr::select(gs_name, gene_symbol)for (ct innames(list_dge)) { dge <- list_dge[[ct]]# Filter significant genes dge_sig <-subset(dge, p_val_adj <0.05)# Extract fold changes foldchanges <- dge_sig$avg_log2FC# Name each fold change with the corresponding gene symbolnames(foldchanges) <- dge_sig$gene# Sort fold changes in decreasing order foldchanges <-sort(foldchanges, decreasing =TRUE)# Run GSEA msig_GSEA <-GSEA(foldchanges,TERM2GENE = m_t2g,verbose =FALSE)# Extract the GSEA results msigGSEA_results <- msig_GSEA@result# Order by NES and pull first row's ID top_id <- msigGSEA_results %>%arrange(desc(abs(NES))) %>%head(1) %>%pull(ID)print(paste(ct, ":", top_id))}```***[Back to Lesson >>](13_dge_pathway.qmd)[Back to Schedule](../schedule/schedule.qmd)