geneL_normal  = c("pnca")
geneL_na      = c("gid", "rpob")
geneL_ppi2    = c("alr", "embb", "katg", "rpob")


if (tolower(gene)%in%geneL_na){
  lf_mcsm_na
  }

if (tolower(gene)%in%geneL_ppi2){
  lf_mcsm_ppi2
}

colnames(lf_duet)
table(lf_duet$param_type)

static_colsP = c("Lig Dist(Å)","ASA", "RSA","RD","KD","Log10(MAF)")

stability_suffix <- paste0(delta_symbol, delta_symbol, "G")
 
lf_commonP   = lf_duet[!lf_duet$param_type%in%c("DUET ΔΔG"),]
lf_commonP$param_type = levels(droplevels(lf_commonP$param_type))
table(lf_commonP$param_type); colnames(lf_commonP)
lf_commonP$outcome = lf_commonP$duet_outcome
lf_commonP$duet_outcome = NULL


lf_duet$outcome = lf_duet$duet_outcome
lf_duet$duet_outcome = NULL
lf_duetP                = lf_duet[!lf_duet$param_type%in%c(static_colsP, "outcome"),]
lf_duetP$param_type     = levels(droplevels(lf_duetP$param_type))
table(lf_duetP$param_type); colnames(lf_duetP)
colnames(lf_duetP)

lf_foldx$outcome = lf_foldx$foldx_outcome
lf_foldx$foldx_outcome = NULL
lf_foldxP               = lf_foldx[!lf_foldx$param_type%in%c(static_colsP,"outcome"),]
lf_foldxP$param_type    = levels(droplevels(lf_foldxP$param_type))
table(lf_foldxP$param_type) 
colnames(lf_foldxP)


lf_deepddg$outcome = lf_deepddg$deepddg_outcome
lf_deepddg$deepddg_outcome = NULL
lf_deepddgP             = lf_deepddg[!lf_deepddg$param_type%in%c(static_colsP, "outcome"),]
lf_deepddgP$param_type  = levels(droplevels(lf_deepddgP$param_type))
table(lf_deepddgP$param_type)
colnames(lf_deepddgP)


lf_dynamut2$outcome = lf_dynamut2$ddg_dynamut2_outcome
lf_dynamut2$ddg_dynamut2_outcome = NULL
lf_dynamut2P            = lf_dynamut2[!lf_dynamut2$param_type%in%c(static_colsP, "outcome"),]
lf_dynamut2P$param_type = levels(droplevels(lf_dynamut2P$param_type))
table(lf_dynamut2P$param_type)
colnames(lf_dynamut2P)


lf_consurf$outcome = lf_consurf$consurf_outcome
lf_consurf$consurf_outcome = NULL
lf_consurfP             = lf_consurf[!lf_consurf$param_type%in%c(static_colsP),]
lf_consurfP$param_type  = levels(droplevels(lf_consurfP$param_type))
table(lf_consurfP$param_type)
colnames(lf_consurfP)


lf_snap2$outcome = lf_snap2$snap2_outcome
lf_snap2$snap2_outcome = NULL
lf_snap2P               = lf_snap2[!lf_snap2$param_type%in%c(static_colsP),]
lf_snap2P$param_type    = levels(droplevels(lf_snap2P$param_type))
table(lf_snap2P$param_type)
colnames(lf_snap2P)


lf_provean$outcome = lf_provean$provean_outcome
lf_provean$provean_outcome = NULL
lf_proveanP             = lf_provean[!lf_provean$param_type%in%c(static_colsP),]
lf_proveanP$param_type  = levels(droplevels(lf_proveanP$param_type))
table(lf_proveanP$param_type)
colnames(lf_proveanP)

bar = rbind(colnames(lf_duetP)
               , colnames(lf_foldxP)
               , colnames(lf_deepddgP)
               , colnames(lf_dynamut2P)
               , colnames(lf_consurfP)
               , colnames(lf_snap2P)
               , colnames(lf_proveanP)
)
bar

lf_df_stabP = rbind((lf_duetP)
                      , (lf_foldxP)
                      , (lf_deepddgP)
                      , (lf_dynamut2P))

lf_df_consP = rbind((lf_consurfP)
                      , (lf_snap2P)
                      , (lf_proveanP))

table(lf_df_stabP$param_type)
# VERY USEFUL for seeing numbers for param types
table(lf_df_stabP$param_type,lf_df_stabP$outcome)
table(lf_df_consP$param_type,lf_df_consP$outcome)

#==============
# Plot:BP
#==============
stability_suffix <- paste0(delta_symbol, delta_symbol, "G")

# lf_bp(lf_df_stabP, p_title = paste0("Stability",stability_suffix)
#       , violin_quantiles = c(0.5))
# lf_bp(lf_df_consP, p_title = "Evolutionary Conservation"
# , violin_quantiles = c(0.5))

lf_bp2(lf_df_stabP, p_title = paste0("Stability",stability_suffix)
       , violin_quantiles = c(0.5), monochrome = F)

lf_bp2(lf_df_consP, p_title = "Evolutionary Conservation"
      , violin_quantiles = c(0.5), monochrome = F)

#HMMM: Bollocks!
lf_bp2(lf_commonP,  p_title = paste0("Residue level properties")
       , violin_quantiles = c(0.5)
      , monochrome = T) # doesn't plot stat bars

lf_bp(lf_commonP,  p_title = paste0("Residue level properties")
       , violin_quantiles = c(0.5)) #plots stat bars but incorrect result

lf_unpaired_stats(lf_duet)
wilcox.test(wf_duet$`Lig Dist(Å)`[wf_duet$mutation_info_labels=="R"]
            , wf_duet$`Lig Dist(Å)`[wf_duet$mutation_info_labels=="S"])

wilcox.test(wf_duet$ASA[wf_duet$mutation_info_labels=="R"]
            , wf_duet$ASA[wf_duet$mutation_info_labels=="S"])


# 1: variable
# 9: conserved

# CHECK THESE
foo = merged_df3[c("dst_mode", "mutation_info_labels", "consurf_colour_rev"
                   , "consurf_scaled"
                   , "consurf_score"
                   , "consurf_outcome"
                   , "snap2_score"
                   , "snap2_scaled"
                   , "snap2_outcome"
                   , "provean_score"
                   , "provean_scaled"
                   , "provean_outcome")]

################
# Affinity
################
# ligand
lf_mcsm_lig$outcome = lf_mcsm_lig$ligand_outcome
lf_mcsm_lig$ligand_outcome = NULL
colnames(lf_mcsm_lig)
table(lf_mcsm_lig$param_type)


lf_mcsm_lig$outcome = lf_mcsm_lig$ligand_outcome
lf_mcsm_ligP            = lf_mcsm_lig[!lf_mcsm_lig$param_type%in%c(static_colsP, "outcome"),]
#lf_mcsm_ligP$param_type = levels(droplevels(lf_mcsm_ligP))
table(lf_mcsm_ligP$param_type)
lf_mcsm_ligP$ligand_outcome = NULL
colnames(lf_mcsm_ligP)

if (tolower(gene)%in%geneL_na){
  lf_mcsm_na$outcome = lf_mcsm_na$mcsm_na_outcome
  #lf_mcsm_na$mcsm_na_outcome = NULL
  lf_mcsm_naP            = lf_mcsm_na[!lf_mcsm_na$param_type%in%c(static_colsP, "outcome"),]
  #lf_mcsm_naP$param_type = levels(droplevels(lf_mcsm_naP))
  table(lf_mcsm_naP$param_type)
  lf_mcsm_naP$mcsm_na_outcome = NULL
  colnames(lf_mcsm_naP)
}

if (tolower(gene)%in%geneL_ppi2){
  
  lf_mcsm_ppi2$outcome = lf_mcsm_ppi2$mcsm_ppi2_outcome
  colnames(lf_mcsm_ppi2)
  #lf_mcsm_ppi2$mcsm_ppi2_outcome = NULL
  lf_mcsm_ppi2P           = lf_mcsm_ppi2[!lf_mcsm_ppi2$param_type%in%c(static_colsP, "outcome"),]
  #lf_mcsm_ppi2P$param_type = levels(droplevels(lf_mcsm_ppi2P))
  table(lf_mcsm_ppi2P$param_type)
  lf_mcsm_ppi2P$mcsm_ppi2_outcome = NULL
  colnames(lf_mcsm_ppi2P)
  
}


bar = rbind(colnames(lf_mcsm_ligP)
            #, colnames(lf_mcsm_naP)
            , colnames(lf_mcsm_ppi2P))
bar

lf_df_affP = rbind((lf_mcsm_ligP)
                    , (lf_mcsm_ppi2P))

lf_bp(lf_df_affP,  p_title = paste0("Affinity changes")
       , violin_quantiles = c(0.5))
       #, monochrome = T) # doesn't plot stat bars

wilcox.test(wf_mcsm_lig$ASA[wf_mcsm_lig$mutation_info_labels=="R"]
            , wf_mcsm_lig$ASA[wf_mcsm_lig$mutation_info_labels=="S"])

wilcox.test(wf_mcsm_ppi2$ASA[wf_mcsm_ppi2$mutation_info_labels=="R"]
            , wf_mcsm_ppi2$ASA[wf_mcsm_ppi2$mutation_info_labels=="S"])

#===============================
p1 =   lf_bp2(lf_df_stabP, p_title = paste0("Stability",stability_suffix)
       , violin_quantiles = c(0.5), monochrome = F)

p2 = lf_bp2(lf_df_consP, p_title = "Evolutionary Conservation"
       , violin_quantiles = c(0.5), monochrome = F)

p3 = lf_bp2(lf_df_affP,  p_title = paste0("Affinity changes")
      , violin_quantiles = c(0.5), monochrome = F)

cowplot::plot_grid(p1,cowplot::plot_grid(p2, p3), nrow=2)