#!/usr/bin/env Rscript  
#########################################################
# TASK: producing boxplots for dr and other muts

#########################################################
#=======================================================================
# working dir and loading libraries
getwd()
setwd("~/git/LSHTM_analysis/scripts/plotting")
getwd()

source("Header_TT.R")
#library(ggplot2)
#library(data.table)
#library(dplyr)
source("other_plots_data.R")

#=======
# output
#=======
dr_other_plots_combined = "dr_other_muts.svg"
plot_dr_other_plots_combined  =  paste0(plotdir,"/", dr_other_plots_combined)

########################################################################
#               end of data extraction and cleaning for plots          #
########################################################################
#my_comparisons <- list( c(dr_muts_col, other_muts_col) )
my_comparisons <- list( c("DM", "OM") )

#===========
# Plot1: PS
#===========

my_stat_ps = compare_means(param_value~mutation_info, group.by = "param_type"
                        , data = df_lf_ps,  paired = FALSE, p.adjust.method = "BH")

y_value = "param_value"

p1 = ggplot(df_lf_ps, aes(x = mutation_info
                    , y = eval(parse(text=y_value)) ))  + 
  facet_wrap(~ param_type
             , nrow = 1
             , scales = "free_y") + 
  geom_boxplot(fill = "white", outlier.colour = NA
               # position = position_dodge(width = 0.9)
               , width = 0.5) +
  geom_point(position = position_jitterdodge(dodge.width=0.01)
             , show.legend = FALSE
             , aes(colour = factor(duet_outcome))) +
  theme(axis.text.x = element_text(size = 15)
        , axis.text.y = element_text(size = 15
                                     , angle = 0
                                     , hjust = 1
                                     , vjust = 0)
        , axis.title.x = element_text(size = 15)
        , axis.title.y = element_text(size = 15)
        , plot.title = element_text(size = 20, hjust = 0.5)
        , strip.text.x = element_text(size = 15, colour = "black")
        , legend.title = element_text(color = "black", size = 20)
        , legend.text = element_text(size = 15)
        , legend.direction = "vertical") +
  labs(title = "DUET"
       , x = ""
       , y = "")+ 
  stat_compare_means(comparisons = my_comparisons
                     , method = "wilcox.test"
                     , paired = FALSE
                     #, label = "p.format")
                     , label = "p.signif")
p1

#===========
# Plot 2: Foldx
#===========
my_stat_foldx = compare_means(param_value~mutation_info, group.by = "param_type"
                              , data = df_lf_foldx,  paired = FALSE, p.adjust.method = "BH")

y_value = "param_value"

p2 = ggplot(df_lf_foldx, aes(x = mutation_info
                             , y = eval(parse(text=y_value)) ))  + 
  facet_wrap(~ param_type
             , nrow = 1
             , scales = "free_y") + 
  geom_boxplot(fill = "white", outlier.colour = NA
               # position = position_dodge(width = 0.9)
               , width = 0.5) +
  geom_point(position = position_jitterdodge(dodge.width=0.01)
             , show.legend = FALSE
             , aes(colour = factor(foldx_outcome))) +
  theme(axis.text.x = element_text(size = 15)
        , axis.text.y = element_text(size = 15
                                     , angle = 0
                                     , hjust = 1
                                     , vjust = 0)
        , axis.title.x = element_text(size = 15)
        , axis.title.y = element_text(size = 15)
        , plot.title = element_text(size = 20, hjust = 0.5)
        , strip.text.x = element_text(size = 15, colour = "black")
        , legend.title = element_text(color = "black", size = 20)
        , legend.text = element_text(size = 15)
        , legend.direction = "vertical") +
  labs(title = "Foldx"
       , x = ""
       , y = "")+ 
  stat_compare_means(comparisons = my_comparisons
                     , method = "wilcox.test"
                     , paired = FALSE
                     #, label = "p.format")
                     , label = "p.signif")
p2


#===========
# Plot 3: LIG
#===========

my_stat_lig = compare_means(param_value~mutation_info, group.by = "param_type"
                            , data = df_lf_lig,  paired = FALSE, p.adjust.method = "BH")

y_value = "param_value"

p3 = ggplot(df_lf_lig, aes(x = mutation_info
                           , y = eval(parse(text=y_value)) ))  + 
  facet_wrap(~ param_type
             , nrow = 1
             , scales = "free_y") + 
  geom_boxplot(fill = "white", outlier.colour = NA
               # position = position_dodge(width = 0.9)
               , width = 0.5) +
  geom_point(position = position_jitterdodge(dodge.width=0.01)
             , show.legend = FALSE
             , aes(colour = factor(ligand_outcome))) +
  theme(axis.text.x = element_text(size = 15)
        , axis.text.y = element_text(size = 15
                                     , angle = 0
                                     , hjust = 1
                                     , vjust = 0)
        , axis.title.x = element_text(size = 15)
        , axis.title.y = element_text(size = 15)
        , plot.title = element_text(size = 20, hjust = 0.5)
        , strip.text.x = element_text(size = 15, colour = "black")
        , legend.title = element_text(color = "black", size = 20)
        , legend.text = element_text(size = 15)
        , legend.direction = "vertical") +
  labs(title = "Ligand Affinity"
       , x = ""
       , y = "")+ 
  stat_compare_means(comparisons = my_comparisons
                     , method = "wilcox.test"
                     , paired = FALSE
                     #, label = "p.format")
                     , label = "p.signif")

p3

#===========================
# combine
#===========================
#svg(plot_or_combined, width = 32, height = 12)
svg("test.svg", width = 25, height = 12)
#theme_set(theme_gray()) # to preserve default theme

printFile = cowplot::plot_grid(plot_grid(p1, p2, p3
                                         , nrow = 3
                                         , align = 'h'
                                         , labels = c("", "", "")
                                         , label_size = 20))
print(printFile)
dev.off()