heatmap.R

# HALLMARK GSEA
library('dplyr')
library('ggplot2')
library('reshape2')
library('RColorBrewer')
plot.dir="/DATA/31/molteni/vexas_mouse_4/results/heatmap/"
wdir="/DATA/31/molteni/vexas_mouse_4/results"
# 
# #heatmap custom cluster
# 
# 
plot.dir="/DATA/31/molteni/vexas_mouse_4/results/heatmap/"
wdir="/DATA/31/molteni/vexas_mouse_4/results"

#custom
Full_final <- readRDS("/DATA/31/molteni/vexas_mouse_4/results/Full/01-seurat/Full_final_annot.rds")
hallmark_ds <- c("A_B")
clusters <- unique(Full_final$celltypes_annot)
hallmark.full <- data.frame()
for (ds in hallmark_ds) {
  for (i in clusters) {
    s <- paste(ds,"markers_cluster", i,sep="_") #qui cambiare 1 con i e inserire ciclo for
    ff <- paste(wdir, "GSEA", "AvsB",s, "tables", paste(s, "NA-gsea-h.all.v7.2.sym_mouse.txt", sep = "-"), sep = "/")
    ds.t <- read.table(ff, header = T, sep = "\t")
    ds.t.filt <- ds.t[,c("ID", "setSize", "enrichmentScore", "NES", "pvalue", "p.adjust", "qvalues"), drop = FALSE]
    ds.t.filt$Dataset <- paste(ds,i,sep="_")
    if (nrow(hallmark.full) == 0) {
      hallmark.full <- ds.t.filt
    } else {
      hallmark.full <- rbind(hallmark.full, ds.t.filt)
    }
  }
}

hallmark.full.filt <- hallmark.full[, c("ID", "NES", "Dataset","p.adjust")]

##
all_hallmark <- hallmark.full.filt

all_hallmark$star <- cut(all_hallmark$p.adjust, breaks=c(-Inf, 0.001, 0.01, 0.05, Inf), label=c("***", "**", "*", ""))  # Create column of significance labels
all_hallmark$ID <- gsub(x = all_hallmark$ID, "HALLMARK_", "")
hallmark.order <- all_hallmark %>% group_by(ID) %>% summarise(Pos = sum(NES))
hallmark.order.terms <- hallmark.order[order(hallmark.order$Pos, decreasing = FALSE), "ID", drop = FALSE]
all_hallmark.nn <- melt(reshape2::dcast(all_hallmark, ID ~ Dataset, value.var="NES"), id.vars = c("ID"))
colnames(all_hallmark.nn) <- c("ID", "Dataset", "NES")
all_hallmark.pp <- melt(reshape2::dcast(all_hallmark, ID ~ Dataset, value.var="star"), id.vars = c("ID"))
colnames(all_hallmark.pp) <- c("ID", "Dataset", "STARS")
all_hallmark.tt <- merge(all_hallmark.nn,all_hallmark.pp)
all_hallmark.tt$ID <- factor(all_hallmark.tt$ID, levels = hallmark.order.terms$ID)
level_order <- c("A_B_HSC_MPP", "A_B_CMP", "A_B_GMP", "A_B_GMP_MEP", "A_B_Myelocytes", "A_B_Neutrophils", 
                 "A_B_MEP_Baso_Mast_cells", "A_B_Plasmacytoid_DC_prog", "A_B_Conventional_DC_prog", "A_B_Proliferating_myeloid_cells",
                 "A_B_Pro_monocyte", "A_B_Monocytes", "A_B_Dendritic_cells", "A_B_Macrophages", "A_B_Erythroid_cells",
                 "A_B_Undefined_1", "A_B_Undefined_2")

all_hallmark.tt$Dataset <- factor(all_hallmark.tt$Dataset, levels =level_order )
p.hallmark <- ggplot(all_hallmark.tt, aes(x = Dataset, y = ID)) +
  theme_bw(base_size = 12) +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  theme(axis.text.y=element_text(size=6))+
  geom_tile(aes(fill = NES), colour = "white") +
  geom_text(aes(label = STARS), color="black", size=3) +
  scale_fill_gradientn(colours = colorRampPalette(rev(brewer.pal(11,"RdBu")))(100),
                       limits = c(-4.5, 4.5),
                       na.value = "grey"
                       #guide = guide_colourbar(reverse = TRUE)
  ) +
  ylab("") + 
  xlab("") +
  coord_fixed(ratio = 0.4)

ggsave(filename = paste(plot.dir, paste("hallmark_AvsB", "heatmap.pdf", sep = "-"), sep = "/"), 
       plot = p.hallmark,
       width = 10, height = 10, dpi = 600)