gs_heatmap {GeneTonic} | R Documentation |
Plot a heatmap for the selected gene signature on the provided data, with the possibility to compactly display also DE only genes
gs_heatmap( se, res_de, res_enrich, annotation_obj = NULL, geneset_id = NULL, genelist = NULL, FDR = 0.05, de_only = FALSE, cluster_rows = TRUE, cluster_columns = FALSE, center_mean = TRUE, scale_row = FALSE, anno_col_info = NULL, plot_title = NULL )
se |
A |
res_de |
A |
res_enrich |
A |
annotation_obj |
A |
geneset_id |
Character specifying the gene set identifier to be plotted |
genelist |
A vector of character strings, specifying the identifiers
contained in the row names of the |
FDR |
Numeric value, specifying the significance level for thresholding adjusted p-values. Defaults to 0.05. |
de_only |
Logical, whether to include only differentially expressed genes in the plot |
cluster_rows |
Logical, determining if rows should be clustered, as
specified by |
cluster_columns |
Logical, determining if columns should be clustered, as
specified by |
center_mean |
Logical, whether to perform mean centering on the row-wise |
scale_row |
Logical, whether to standardize by row the expression values |
anno_col_info |
A character vector of names in |
plot_title |
Character string, to specify the title of the plot,
displayed over the heatmap. If left to |
A plot returned by the ComplexHeatmap::Heatmap()
function
library("macrophage") library("DESeq2") library("org.Hs.eg.db") library("AnnotationDbi") # dds object data("gse", package = "macrophage") dds_macrophage <- DESeqDataSet(gse, design = ~line + condition) rownames(dds_macrophage) <- substr(rownames(dds_macrophage), 1, 15) dds_macrophage <- estimateSizeFactors(dds_macrophage) vst_macrophage <- vst(dds_macrophage) # annotation object anno_df <- data.frame( gene_id = rownames(dds_macrophage), gene_name = mapIds(org.Hs.eg.db, keys = rownames(dds_macrophage), column = "SYMBOL", keytype = "ENSEMBL"), stringsAsFactors = FALSE, row.names = rownames(dds_macrophage) ) # res object data(res_de_macrophage, package = "GeneTonic") res_de <- res_macrophage_IFNg_vs_naive # res_enrich object data(res_enrich_macrophage, package = "GeneTonic") res_enrich <- shake_topGOtableResult(topgoDE_macrophage_IFNg_vs_naive) res_enrich <- get_aggrscores(res_enrich, res_de, anno_df) gs_heatmap(vst_macrophage, res_de, res_enrich, anno_df, geneset_id = res_enrich$gs_id[1], cluster_columns = TRUE, anno_col_info = "condition")