2024-2: 5b.4 GO (Gene ontology) term enrichment
- Vicki Thomson
4a. GO term enrichment
GO terms are functional groupings of genes. We can look at enrichment of GO terms in a similar way as we did with KEGG pathways.
NOTE: GO terms are hierarchical, necessitating that we initially provide one of three top level terms in ont <- "..": these three terms are Biological Processes (ont <- "BP") Molecular Functions (ont <- "MF") and Cellular Components (ont <- "CC"). We should run the two following code blocks (enrichment and plotting) three times, providing each ont <- "..".
Copy and paste this section into R.
#### 11. GO term enrichment ####
# Fist select which of the top level terms you want to examine. You can run the analysis for each top-level term by changing this to another term then re-running this code cell and then the following GO code cells. The three terms are "BP", "MF", and "CC".
ont <- "BP"
# Then run the enrichGO function
gg <- enrichGO(gene = as.character(gene_list$ENTREZID), OrgDb = org.Mm.eg.db, ont = ont, pvalueCutoff = 0.05, qvalueCutoff = 0.2)
# Now you can save this as a table of enriched GO terms
write.csv(as.data.frame(gg), paste0("./results/Enriched_GO_terms_", ont, ".csv"))
4b. Plotting enriched GO terms
As with the KEGG pathways, we can plot a bar plot, dot plot and network plot.
In addition, we’re adding a fourth plot type here: a directed acyclic graph.
As previously mentioned, GO terms are hierarchical, with individual terms being linked to parent or child levels. The previous plots combine all levels under each of three top-level terms. A directed acyclic graph allows a visualisation of how the enriched GO terms fit into the GO hierarchy.
Copy and paste this section into R.
#### 12. Plotting GO results ####
# Bar plot
p <- barplot(gg, showCategory = 14, font.size = 14)
p
# Output as tiff
tiff_exp <- paste0("./results/", ont, "_GO_bar.tiff")
ggsave(file = tiff_exp, dpi = 300, compression = "lzw", device = "tiff", plot = p, width = 20, height = 20, units = "cm")
# Output as pdf
pdf_exp <- paste0("./results/", ont, "_GO_bar.pdf")
ggsave(file = pdf_exp, device = "pdf", plot = p, width = 20, height = 20, units = "cm")
# Dot plot
p <- dotplot(gg, showCategory = 8, font.size = 14)
p
# Output as tiff
tiff_exp <- paste0("./results/", ont, "_GO_dot.tiff")
ggsave(file = tiff_exp, dpi = 300, compression = "lzw", device = "tiff", plot = p, width = 20, height = 20, units = "cm")
# Output as pdf
pdf_exp <- paste0("./results/", ont, "_GO_dot.pdf")
ggsave(file = pdf_exp, device = "pdf", plot = p, width = 20, height = 20, units = "cm")
# Network plot
p <- cnetplot(gg, showCategory = 8, colorEdge = TRUE, node_label = "category")
p
# Output as tiff
tiff_exp <- paste0("./results/", ont, "_GO_network.tiff")
ggsave(file = tiff_exp, dpi = 300, compression = "lzw", device = "tiff", plot = p, width = 20, height = 20, units = "cm")
# Output as pdf
pdf_exp <- paste0("./results/", ont, "_GO_network.pdf")
ggsave(file = pdf_exp, device = "pdf", plot = p, width = 20, height = 20, units = "cm")
# Directed acyclic graph
p <- goplot(gg, showCategory = 30)
p
tiff_exp <- paste0("./results/", ont, "_GO_DAG.tiff")
ggsave(file = tiff_exp, dpi = 300, compression = "lzw", device = "tiff", plot = p, width = 20, height = 20, units = "cm")
# Output as pdf
pdf_exp <- paste0("./results/", ont, "_GO_DAG.pdf")
ggsave(file = pdf_exp, device = "pdf", plot = p, width = 20, height = 20, units = "cm")
