Commit c0e03a31 authored by domingue's avatar domingue
Browse files

Limited number of plots when too many pathways are enriched

parent c93f3790
......@@ -75,7 +75,7 @@ dir.create("figure", showWarnings = FALSE)
# de_file <- "../de_results.txt"
# de_file <- file.path(Sys.getenv("PRJ_DATA"),"/data/alexaki/de_results.txt")
# category <- "None"
# category <- "H"
# category <- "C3"
# extra_sets <- "extra_sets.txt"
# species <- "Mus musculus"
# qcutoff <- 0.05
......@@ -226,10 +226,24 @@ fgseaRes <- fgsea(m_list, ranks, minSize = 15, maxSize = 500, nperm = 10000)
enriched_pathways <- fgseaRes[padj < qcutoff][["pathway"]]
clean_ranks <- ranks[!is.na(ranks)] # needed for plotting, NA cause error
#' We have have an overview of the pathway enrichment by plotting the normalized enrichment score (NES) and p-adjusted values for each of the enriched pathways. A positive NES indicates gene set enrichment at the top of the ranked list; a negative NES indicates gene set enrichment at the bottom of the ranked list.
## if too many pathways are enriched, only top 20 should be plotted in the report,otherwise it's a mess.
if(length(enriched_pathways) > 20){
msg <- paste0(length(enriched_pathways), " were found significantly enriched for a p-adjust value of ", qcutoff, ".In this report only the top 20 will be show for visualization reasons. Enrichment plots for all enriched pathways can be found in the folder figures/")
top_pathways <- fgseaRes[padj < qcutoff] %>%
arrange(padj) %>%
head(20) %>%
setDT()
} else {
top_pathways <- fgseaRes[padj < qcutoff]
}
#' `r msg`
#'
#' We can have an overview of the pathway enrichment by plotting the normalized enrichment score (NES) and p-adjusted values for each of the top enriched pathways. A positive NES indicates gene set enrichment at the top of the ranked list; a negative NES indicates gene set enrichment at the bottom of the ranked list.
#' Reference: https://www.gsea-msigdb.org/gsea/doc/GSEAUserGuideFrame.html
fgseaRes[padj < qcutoff] %>%
top_pathways %>%
mutate(pathway = fct_reorder(pathway, NES)) %>%
arrange(-padj) %>%
ggplot(aes(x = pathway, y = NES, fill = padj)) +
......@@ -241,7 +255,7 @@ fgseaRes[padj < qcutoff] %>%
coord_flip() +
theme_light()
for (pathway in enriched_pathways) {
for (pathway in top_pathways[["pathway"]]) {
p1 <- plotEnrichment(m_list[[pathway]], clean_ranks) +
labs(title = pathway)
print(p1)
......@@ -252,11 +266,11 @@ for (pathway in enriched_pathways) {
#
# tables of results
# --------------------------------------------------------------------------
#' Since the genes were ranked according to their differential expression significance and fold change, with the most significantly down-regulated genes at the top and up-regulated genes at the bottom of the list, the enriched gene sets provides us with some idea of the function of these genes. We can vislualize that in the table:
#' Since the genes were ranked according to their differential expression significance and fold change, with the most significantly down-regulated genes at the top and up-regulated genes at the bottom of the list, the enriched gene sets provides us with some idea of the function of these genes. We can visualize that in the table:
#'
dev.off()
topPathwaysUp <- fgseaRes[ES > 0 & padj < qcutoff][order(pval), pathway]
topPathwaysDown <- fgseaRes[ES < 0 & padj < qcutoff][order(pval), pathway]
topPathwaysUp <- top_pathways[ES > 0 & padj < qcutoff][order(pval), pathway]
topPathwaysDown <- top_pathways[ES < 0 & padj < qcutoff][order(pval), pathway]
topPathways <- c(topPathwaysUp, rev(topPathwaysDown))
p_table <- plotGseaTable(
......@@ -266,7 +280,7 @@ p_table <- plotGseaTable(
gseaParam = 0.5,
render = FALSE
)
#+ fig.width=7, fig.height=8
#+ fig.width=9, fig.height=8
grid.draw(p_table)
#' Often it is useful to extract the core members of high scoring gene sets that contribute to the enrichment score. As you may have noticed from the previous section, one of the results column was named “leadingEdge”; this contains the genes that contributed to the enrichment score.
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment