# Immagini singole malattie: 
# 1- Clusters più accesi
# 2- Vlns con singolo colore e diverse sfumature
# 3- miRNAh ----> low solo riga continua e blocco sulla malattia.
# 4. Unico pannello
# 

library(Seurat)
library(ggplot2)
library(gridExtra)
library(Seurat)
library(ggplot2)
library(scales)
library(dplyr)

# Fix immagini singole - colori più accesi

objs <- list()
data2plot <- list()
for(i in c("L1","L2","L3","L5","L7")){
  objs[[i]] <- readRDS(file = paste0("/beegfs/scratch/ric.gentner/ric.gentner/scRNA_BALL/resultswithpediatric/",i,"/01-seurat/",i,"_final.rds"))
  sigs <- readRDS(file = paste0("/beegfs/scratch/ric.gentner/ric.gentner/scRNA_BALL/resultswithpediatric/",i,"/00-Signatures/",i,"_ModScores.rds")) 
  objs[[i]] <- AddMetaData(objs[[i]], sigs)
  data2plot[[i]] <- FetchData(object = objs[[i]], vars = c("UMAP_1","UMAP_2","tSNE_1","tSNE_2",colnames(objs[[i]]@meta.data)))
}

# Vlnplots

# Verde scuro: #64921E
# Verde chiaro: #EFF4E8

mypalette <- c("#FFE599", "#cc0000")
scaleFUN <- function(x) sprintf("%.1f", x)

colfunc <- list()
dimplots <- list()
diseases <- list(L1 = "PZ1227",
                 L2 = "PZ13000037",
                 L3 = "PZ114966",
                 L5 = "PZ170160",
                 L7 = "PZ5114")
vlnplots <- list()

for(i in c("L1","L2","L3","L5","L7")){
  if(i == "L2"){
    objs[[i]] <- SetIdent(object = objs[[i]], value = "RNA_snn_res.1.2")
    colfunc[[i]] <- rev(colorRampPalette(mypalette)(length(unique(Idents(objs[[i]]))) - 1))
    
    vlnplots[[i]] <- VlnPlot(object = objs[[i]], features = "dIntegrated_up", sort = T, cols = colfunc[[i]], 
                             pt.size = 0, idents = setdiff(levels(objs[[i]]$RNA_snn_res.1.2), "9")) + NoLegend() + 
      geom_boxplot(width=0.15, color="black",outlier.size = 0) + 
      ggtitle(label = paste0(diseases[[i]], " - signature score")) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), #plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
      ) 
    
    names(colfunc[[i]]) <- levels(vlnplots[[i]]$data$ident)
    colfunc[[i]][10] <- "lightgrey"
    names(colfunc[[i]])[10] <- "9"
    dimplots[[i]] <- DimPlot(object = objs[[i]], pt.size = 0.3, #group.by = "RNA_snn_res.1.2", 
                             label = T, cols = colfunc[[i]], order = T) + NoLegend() + 
      ggtitle(paste0(diseases[[i]], " - res 1.2")) +
      scale_y_continuous(labels = scaleFUN) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), #plot.title = element_text(hjust = 0.5, size = 12)
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
            ) 

    
  }
  if(i == "L7"){
    objs[[i]] <- SetIdent(object = objs[[i]], value = "RNA_snn_res.0.6")
    colfunc[[i]] <- rev(colorRampPalette(mypalette)(length(unique(Idents(objs[[i]]))) - 3))
    
    vlnplots[[i]] <- VlnPlot(object = objs[[i]], features = "dIntegrated_up", sort = T, cols = colfunc[[i]], 
                             pt.size = 0, idents = setdiff(levels(objs[[i]]$RNA_snn_res.0.6), c("8","9","10"))) + NoLegend() +
      geom_boxplot(width=0.15, color="black",outlier.size = 0) + 
      ggtitle(label = paste0(diseases[[i]], " - signature score")) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), #plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
      ) 
    
    names(colfunc[[i]]) <- levels(vlnplots[[i]]$data$ident)
    colfunc[[i]][9:11] <- "lightgrey"
    names(colfunc[[i]])[9:11] <- c("8","9","10")
    
    dimplots[[i]] <- DimPlot(object = objs[[i]], pt.size = 0.3, group.by = "RNA_snn_res.0.6", 
                             label = T, cols = colfunc[[i]], order = T) + NoLegend() + 
      ggtitle(paste0(diseases[[i]], " - res 0.6")) +
      scale_y_continuous(labels = scaleFUN) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), # plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
      ) 
  }
  if(i == "L1"){
    objs[[i]] <- SetIdent(object = objs[[i]], value = "RNA_snn_res.0.6")
    colfunc[[i]] <- c("#CC0000", "#DA412B", "#E16241", "#E98257","#F7C483" , "#FFE599")
    
    vlnplots[[i]] <- VlnPlot(object = objs[[i]], features = "dIntegrated_up", sort = T, cols = colfunc[[i]], 
                             pt.size = 0, idents = setdiff(levels(objs[[i]]$RNA_snn_res.0.6), c("6","7"))) + NoLegend() +
      geom_boxplot(width=0.15, color="black",outlier.size = 0) + 
      ggtitle(label = paste0(diseases[[i]], " - signature score")) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), #plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
      ) 
    
    names(colfunc[[i]]) <- levels(vlnplots[[i]]$data$ident)
    colfunc[[i]][7:8] <- "lightgrey"
    names(colfunc[[i]])[7:8] <- c("6","7")
    
    dimplots[[i]] <- DimPlot(object = objs[[i]], pt.size = 0.3, group.by = "RNA_snn_res.0.6", 
                             label = T, cols = colfunc[[i]], order = T) + NoLegend() + 
      ggtitle(paste0(diseases[[i]], " - res 0.6")) +
      scale_y_continuous(labels = scaleFUN) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), # plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
      ) 
  }
  if(i %in% c("L3","L5")){
    if(i == "L3"){
      colfunc[[i]] <- c("#CC0000","#DA412B","#E98257","#F7C483","#FFE599")
    }else{
      colfunc[[i]] <- rev(colorRampPalette(mypalette)(length(unique(Idents(objs[[i]])))))
    }
    objs[[i]] <- SetIdent(object = objs[[i]], value = "RNA_snn_res.0.6")
    vlnplots[[i]] <- VlnPlot(object = objs[[i]], features = "dIntegrated_up", sort = T, cols = colfunc[[i]], pt.size = 0) + NoLegend() +
      geom_boxplot(width=0.15, color="black",outlier.size = 0) + 
      ggtitle(label = paste0(diseases[[i]], " - signature score")) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), #plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
      )
    
    names(colfunc[[i]]) <- levels(vlnplots[[i]]$data$ident)

    
    dimplots[[i]] <- DimPlot(object = objs[[i]], pt.size = 0.3, group.by = "RNA_snn_res.0.6", 
                             label = T, cols = colfunc[[i]], order = T) + NoLegend() + 
      ggtitle(paste0(diseases[[i]], " - res 0.6")) +
      scale_y_continuous(labels = scaleFUN) +
      theme(axis.title = element_blank(), 
            axis.line = element_line(size = 0.2, color = "black"), 
            axis.text = element_text(size = 8), # plot.title = element_text(hjust = 0.5, size = 12),
            plot.title = element_blank(),
            axis.ticks = element_line(size = 0.2, color = "black")
            ) 

  }
}

png(filename = "Clusters_signature_distribution_res300.png", width = 6, height = 10, units = "in", res = 300)
grid.arrange(dimplots[["L1"]], vlnplots[["L1"]],
             dimplots[["L2"]], vlnplots[["L2"]],
             dimplots[["L3"]], vlnplots[["L3"]],
             dimplots[["L5"]], vlnplots[["L5"]],
             dimplots[["L7"]], vlnplots[["L7"]], nrow=5, ncol=2)
dev.off()

png(filename = "Clusters_signature_distribution_res200.png", width = 6, height = 10, units = "in", res = 200)
grid.arrange(dimplots[["L1"]], vlnplots[["L1"]],
             dimplots[["L2"]], vlnplots[["L2"]],
             dimplots[["L3"]], vlnplots[["L3"]],
             dimplots[["L5"]], vlnplots[["L5"]],
             dimplots[["L7"]], vlnplots[["L7"]], nrow=5, ncol=2)
dev.off()

save.image(file = "PaperPanels.RData")

##### UMAP - IMAGES #########

a <- sample(rainbow(50),5)
dgns <- readRDS("/beegfs/scratch/ric.gentner/ric.gentner/scRNA_BALL/resultswithpediatric/Diagnosis/01-seurat/Diagnosis_final.rds")

dgns <- SetIdent(object = dgns, value = "RNA_DonorID")

mydata <- FetchData(object = dgns, vars = c("UMAPh_1","UMAPh_2", "orig.ident"))
plotdata <- mydata
patients <- c("DIAGNOSIS","PZ1227","PZ13000037","PZ114966","PZ170160","PZ5114")
pn <- 1
plotdata$plotnum <- patients[[(pn)]]
for (id in unique(mydata$orig.ident)) {
  pn <- pn + 1
  pdata <- subset(mydata, orig.ident == id)
  pdata$plotnum <- patients[[(pn)]]
  plotdata <- rbind(plotdata, pdata)
}

png("Panel_DGNs_singleDiseases_distribution.png", width = 9,height = 6, units = "in", res = 300)
ggplot(plotdata, aes(x = UMAPh_1, y = UMAPh_2, col = orig.ident)) +
  #theme_bw() +
  theme(legend.position = "none",
        #panel.grid.major = element_blank(),
        #panel.grid.minor = element_blank(),
        #panel.background = element_rect(fill = "#FAFAFA"),
        panel.background = element_rect(fill = "white"),
        strip.background = element_rect(fill = "#FAFAFA"),
        axis.line = element_line(colour = "black", size = 0.2),
        axis.ticks = element_line(size = 0.2),
        axis.text = element_text(size = 8)
        #strip.text = element_blank()
        ) +
  geom_point(size = .1, alpha = .5) +
  xlab("") +
  ylab("") +
  facet_wrap(.~plotnum, ncol = 3)
dev.off()

#### PLOTTING SINGLE-R SIGNATURES DISTRIBUTION ####

mydata.full <- FetchData(object = dgns, vars = c("UMAPh_1","UMAPh_2", colnames(dgns@meta.data)))
singler.sig <- grep("Database", x = grep("SingleR", x = colnames(dgns@meta.data), value = T), invert = T, value = T)

u <- as.data.frame(prop.table(table(dgns$RNA_DonorID, dgns$SingleR_BlueprintEncodeData_labels), 1))
u$dataset <- "SingleR_BlueprintEncodeData_labels"

for (id in singler.sig[2:length(singler.sig)]) {
  
  g <- as.data.frame(prop.table(table(dgns@meta.data[,"RNA_DonorID"], 
                                      dgns@meta.data[,id]), 1))
  g$dataset <- id
  u <- rbind(u,g)
}

u$dataset <- gsub(u$dataset, pattern = "SingleR_", replacement = "")
u$dataset <- gsub(u$dataset, pattern = "_labels", replacement = "")
u$dataset <- gsub(u$dataset, pattern = "BlueprintEncodeData", replacement = "BPE")
u$dataset <- gsub(u$dataset, pattern = "HumanPrimaryCellAtlasData", replacement = "HPCA")
u$dataset <- gsub(u$dataset, pattern = "SingleRrefined_", replacement = "r")
u$dataset <- gsub(u$dataset, pattern = "MonacoImmuneData", replacement = "Monaco")
u$dataset <- gsub(u$dataset, pattern = "NovershternHematopoieticData", replacement = "Nover")


# re-classify labels with % < 2% as "<2%"
u_filt <- u
u_filt$Var2 <- as.character(u_filt$Var2)
u_filt$Var2[u_filt$Freq < 0.02] <- "<2%"
u_filt$Var2 <- as.factor(u_filt$Var2)
disease_distr_list <- split(u_filt, u$dataset)

dname <- list(BPE = "Blueprint Encode",
              HPCA = "Human primary cell atlas",
              Monaco = "Monaco Immune Data",
              Nover = "Novershtern Hematopoietic Data")


plot_dat_dist <- list()
for(myset in c("BPE","HPCA","Monaco","Nover")){
  if(myset %in% c("BPE","HPCA")){
    plot_dat_dist[[myset]] <- ggplot(data = disease_distr_list[[myset]], 
                                     mapping = aes(Var1, Freq*100, fill = Var2)) + 
      scale_fill_brewer(palette = "Paired") + 
      ggtitle(label = dname[[myset]]) + 
      geom_bar(stat = "identity", width = 0.5) + 
      guides(fill = guide_legend(nrow = 2)) +
      #guides(colour = guide_legend(title = "", override.aes = list(size=20))) +
      theme(legend.position = "bottom",
            plot.title = element_text(size = 7, hjust = 0.5, face = "plain"),
            panel.background = element_rect(fill = "white"),
            axis.line = element_line(colour = "black", size = 0.2),
            axis.title = element_blank(), 
            plot.margin = margin(0.2, 0.2,0.2,0.2, unit = "cm"),
            axis.text = element_text(size = 5),
            legend.title = element_blank(),
            legend.text = element_text(size = 4),
            legend.spacing.y = unit(0.1, 'cm'),
            legend.spacing.x = unit(0.1, 'cm'),
            legend.key.size = unit(0.2, 'cm'), 
            legend.key = element_rect(fill = "white"),
            legend.background = element_blank(),
            #legend.key.height = unit(0, 'cm'),
            legend.key.width = unit(0.2, 'cm'),
            legend.box.margin=margin(-10,0,0,0)
      ) #+ facet_wrap(.~Var1, ncol = 5)
  }
  else{
    plot_dat_dist[[myset]] <- ggplot(data = disease_distr_list[[myset]], 
                                     mapping = aes(Var1, Freq*100, fill = Var2)) + 
      scale_fill_brewer(palette = "Paired") + 
      ggtitle(label = dname[[myset]]) + 
      geom_bar(stat = "identity", width = 0.5) + 
      guides(fill = guide_legend(nrow = 2)) +
      #guides(colour = guide_legend(title = "", override.aes = list(size=20))) +
      theme(legend.position = "bottom",
            plot.title = element_text(size = 7, hjust = 0.5, face = "plain"),
            panel.background = element_rect(fill = "white"),
            axis.line = element_line(colour = "black", size = 0.2),
            axis.title = element_blank(), 
            plot.margin = margin(0, 0.2,0.2,0.2, unit = "cm"),
            axis.text = element_text(size = 5),
            legend.title = element_blank(),
            legend.key = element_rect(fill = "white"),
            legend.text = element_text(size = 4),
            legend.spacing.y = unit(0.1, 'cm'),
            legend.spacing.x = unit(0.1, 'cm'),
            legend.key.size = unit(0.2, 'cm'), 
            legend.background = element_blank(),
            #legend.key.height = unit(0, 'cm'),
            legend.key.width = unit(0.2, 'cm'),
            legend.box.margin=margin(-10,0,0,0)
      ) #+ facet_wrap(.~Var1, ncol = 5)
    
  }
}

png(filename = "Panel_cell_type_distribution_across_diseases_res_200.png",
    width = 6, height = 4.5, units = "in", res = 200)
grid.arrange(plot_dat_dist[["BPE"]], plot_dat_dist[["HPCA"]], plot_dat_dist[["Monaco"]], plot_dat_dist[["Nover"]])
dev.off()

png(filename = "Panel_cell_type_distribution_across_diseases_res_400.png",
    width = 6, height = 4.5, units = "in", res = 400)
grid.arrange(plot_dat_dist[["BPE"]], plot_dat_dist[["HPCA"]], plot_dat_dist[["Monaco"]], plot_dat_dist[["Nover"]])
dev.off()

############## SHINYGO - PLOTS #############

hallmarks.sgo.high <- read.csv(file = "shinyGOinputs/high_hallmarks.csv")
hallmarks.sgo.high$condition <- "High"
hallmarks.sgo.high$Pathway <- gsub(x = hallmarks.sgo.high$Pathway, pattern = "HALLMARK ", replacement = "")
hallmarks.sgo.low <- read.csv(file = "shinyGOinputs/low_hallmark.csv")
hallmarks.sgo.low$condition <- "Low"
hallmarks.sgo.low$Pathway <- gsub(x = hallmarks.sgo.low$Pathway, pattern = "HALLMARK ", replacement = "")

order_h_l <- hallmarks.sgo.high[order(x = hallmarks.sgo.high$Fold.Enrichment, decreasing = T),]
order_l_h <- hallmarks.sgo.low[order(x = hallmarks.sgo.low$Fold.Enrichment, decreasing = F),]

hallmarks.sgo <- rbind(hallmarks.sgo.high, hallmarks.sgo.low)


df <- transform(hallmarks.sgo, Score=ifelse(as.character(condition) %in% c("Low"), -Fold.Enrichment, Fold.Enrichment))
df <- droplevels.data.frame(df)
df$Pathway <- factor(x = df$Pathway, levels = rev(c(order_h_l$Pathway, c(setdiff(y = order_h_l$Pathway, x =order_l_h$Pathway)))))
colnames(df)[which(colnames(df) == "Score")] <- "Fold Enrichment"

hallmarks_plot <- ggplot() + 
  geom_bar(data=df,aes(x=Pathway, y=`Fold Enrichment`, fill=condition), stat="identity", width = 0.2) + 
  scale_fill_manual(values = c("#D4070F","#0066E7")) +
  geom_point(data=df,aes(x=Pathway, y=`Fold Enrichment`, colour=condition, size = nGenes), stat="identity") +
  scale_color_manual(values = alpha(c("#D4070F","#0066E7"), alpha = 1)) +
  geom_hline(yintercept=0) +
  coord_flip() +
  scale_y_continuous(labels=abs,limits=c(-7.5,7.5)) +
  theme(legend.position = "right",
        plot.title = element_text(size = 7, hjust = 0.5, face = "plain"),
        panel.background = element_rect(fill = "white"),
        axis.line = element_line(colour = "black", size = 0.3),
        axis.title = element_blank(), 
        axis.text.y = element_text(size = 10),
        axis.title.x = element_text(size = 10)#,
        #axis.text = element_text(size = 7)
)


#### KEGGS 

kegg.sgo.high <- read.csv(file = "shinyGOinputs/high_kegg.csv")
kegg.sgo.high$condition <- "High"
kegg.sgo.low <- read.csv(file = "shinyGOinputs/low_kegg.csv")
kegg.sgo.low$condition <- "Low"

order_h_l <- kegg.sgo.high[order(x = kegg.sgo.high$Fold.Enrichment, decreasing = T),]
order_l_h <- kegg.sgo.low[order(x = kegg.sgo.low$Fold.Enrichment, decreasing = F),]

kegg.sgo <- rbind(kegg.sgo.high, kegg.sgo.low)

df <- transform(kegg.sgo, Score=ifelse(as.character(condition) %in% c("Low"), -Fold.Enrichment, Fold.Enrichment))
df <- droplevels.data.frame(df)
df$Pathway <- factor(x = df$Pathway, levels = rev(c(order_h_l$Pathway, c(setdiff(y = order_h_l$Pathway, x =order_l_h$Pathway)))))
colnames(df)[which(colnames(df) == "Score")] <- "Fold Enrichment"

kegg_plot <- ggplot() + 
  geom_bar(data=df,aes(x=Pathway, y=`Fold Enrichment`, fill=condition), stat="identity", width = 0.2) + 
  scale_fill_manual(values = c("#D4070F","#0066E7")) +
  geom_point(data=df,aes(x=Pathway, y=`Fold Enrichment`, colour=condition, size = nGenes), stat="identity") +
  scale_color_manual(values = alpha(c("#D4070F","#0066E7"), alpha = 0.8)) +
  geom_hline(yintercept=0) +
  coord_flip() +
  scale_y_continuous(labels=abs,limits=c(-10,10)) +
  theme(legend.position = "right",
        plot.title = element_text(size = 7, hjust = 0.5, face = "plain"),
        panel.background = element_rect(fill = "white"),
        axis.line = element_line(colour = "black", size = 0.3),
        axis.title = element_blank(), 
        axis.text.y = element_text(size = 10),
        axis.title.x = element_text(size = 10)#,
        #axis.text = element_text(size = 7)
  )

### TF ####

TFfact.sgo.high <- read.csv(file = "shinyGOinputs/high_TFfact.csv")
TFfact.sgo.high$condition <- "High"
TFfact.sgo.low <- read.csv(file = "shinyGOinputs/low_TFfact.csv")
TFfact.sgo.low$condition <- "Low"

order_h_l <- TFfact.sgo.high[order(x = TFfact.sgo.high$Fold.Enrichment, decreasing = T),]
order_l_h <- TFfact.sgo.low[order(x = TFfact.sgo.low$Fold.Enrichment, decreasing = F),]

TFfact.sgo <- rbind(TFfact.sgo.high, TFfact.sgo.low)

df <- transform(TFfact.sgo, Score=ifelse(as.character(condition) %in% c("Low"), -Fold.Enrichment, Fold.Enrichment))
df <- droplevels.data.frame(df)
df$Pathway <- factor(x = df$Pathway, levels = rev(c(order_h_l$Pathway, c(setdiff(y = order_h_l$Pathway, x =order_l_h$Pathway)))))
colnames(df)[which(colnames(df) == "Score")] <- "Fold Enrichment"

TFfact_plot <- ggplot() + 
  geom_bar(data=df,aes(x=Pathway, y=`Fold Enrichment`, fill=condition), stat="identity", width = 0.2) + 
  scale_fill_manual(values = c("#D4070F","#0066E7")) +
  geom_point(data=df,aes(x=Pathway, y=`Fold Enrichment`, colour=condition, size = nGenes), stat="identity") +
  scale_color_manual(values = alpha(c("#D4070F","#0066E7"), alpha = 1)) +
  geom_hline(yintercept=0) +
  coord_flip() +
  scale_y_continuous(labels=abs,limits=c(-12.5,12.5)) +
  theme(legend.position = "right",
        plot.title = element_text(size = 7, hjust = 0.5, face = "plain"),
        panel.background = element_rect(fill = "white"),
        axis.line = element_line(colour = "black", size = 0.3),
        axis.title = element_blank(), 
        axis.text.y = element_text(size = 10),
        axis.title.x = element_text(size = 10)#,
        #axis.text = element_text(size = 7)
  )

# GOs

GObiopro.sgo.high <- read.csv(file = "shinyGOinputs/high_GObiopro.csv")
GObiopro.sgo.high$condition <- "High"
GObiopro.sgo.low <- read.csv(file = "shinyGOinputs/low_GObiopro.csv")
GObiopro.sgo.low$condition <- "Low"

order_h_l <- GObiopro.sgo.high[order(x = GObiopro.sgo.high$Fold.Enrichment, decreasing = T),]
order_l_h <- GObiopro.sgo.low[order(x = GObiopro.sgo.low$Fold.Enrichment, decreasing = F),]

GObiopro.sgo <- rbind(GObiopro.sgo.high, GObiopro.sgo.low)

df <- transform(GObiopro.sgo, Score=ifelse(as.character(condition) %in% c("Low"), -Fold.Enrichment, Fold.Enrichment))
df <- droplevels.data.frame(df)
df$Pathway <- factor(x = df$Pathway, levels = rev(c(order_h_l$Pathway, c(setdiff(y = order_h_l$Pathway, x =order_l_h$Pathway)))))
colnames(df)[which(colnames(df) == "Score")] <- "Fold Enrichment"

GObiopro_plot <- ggplot() + 
  geom_bar(data=df,aes(x=Pathway, y=`Fold Enrichment`, fill=condition), stat="identity", width = 0.2) + 
  scale_fill_manual(values = c("#D4070F","#0066E7")) +
  geom_point(data=df,aes(x=Pathway, y=`Fold Enrichment`, colour=condition, size = nGenes), stat="identity") +
  scale_color_manual(values = alpha(c("#D4070F","#0066E7"), alpha = 1)) +
  geom_hline(yintercept=0) +
  coord_flip() +
  scale_y_continuous(labels=abs,limits=c(-4,4)) +
  theme(legend.position = "right",
        plot.title = element_text(size = 7, hjust = 0.5, face = "plain"),
        panel.background = element_rect(fill = "white"),
        axis.line = element_line(colour = "black", size = 0.3),
        axis.title = element_blank(),
        axis.text.y = element_text(size = 8), 
        axis.title.x = element_text(size = 7)#,
        #axis.text = element_text(size = 7)
  )


# png(filename = "ShinyGO_enrichments_miRNA_high_low_singlecell_markers.png", width = 12, height = 8, units = "in", res = 300)
# grid.arrange(kegg_plot, hallmarks_plot , GObiopro_plot, TFfact_plot, ncol = 2, nrow = 2)
# dev.off()

# Single plots

png(filename = "GObiopro_enrichments_miRNA_high_low_singlecell_markers.png", width = 12, height = 8, units = "in", res = 300)
grid.arrange(GObiopro_plot, ncol = 1, nrow = 1)
dev.off()

png(filename = "hallmarks_enrichments_miRNA_high_low_singlecell_markers.png", width = 12, height = 8, units = "in", res = 300)
grid.arrange(hallmarks_plot, ncol = 1, nrow = 1)
dev.off()

png(filename = "kegg_enrichments_miRNA_high_low_singlecell_markers.png", width = 12, height = 8, units = "in", res = 300)
grid.arrange(kegg_plot, ncol = 1, nrow = 1)
dev.off()

png(filename = "TF_enrichments_miRNA_high_low_singlecell_markers.png", width = 12, height = 8, units = "in", res = 300)
grid.arrange(TFfact_plot, ncol = 1, nrow = 1)
dev.off()