From e60b4c54923a8622ff5daab9f239c5f51a1453c2 Mon Sep 17 00:00:00 2001 From: Tanushree Tunstall Date: Tue, 6 Oct 2020 17:47:24 +0100 Subject: [PATCH] output corr plots with coloured dots --- scripts/plotting/corr_PS_LIG.R | 165 +++++---------------------------- 1 file changed, 24 insertions(+), 141 deletions(-) diff --git a/scripts/plotting/corr_PS_LIG.R b/scripts/plotting/corr_PS_LIG.R index b6ea1d7..30c3033 100644 --- a/scripts/plotting/corr_PS_LIG.R +++ b/scripts/plotting/corr_PS_LIG.R @@ -14,7 +14,7 @@ getwd() source("Header_TT.R") require(cowplot) source("combining_dfs_plotting.R") -source("my_pairs_panel.R") # with lower panel turned off + # should return the following dfs, directories and variables # PS combined: @@ -51,84 +51,29 @@ cat(paste0("Variables imported:" #======= # output #======= -# can't combine by cowplot because not ggplots -#corr_plot_combined = "corr_combined.svg" -#plot_corr_plot_combined = paste0(plotdir,"/", corr_plot_combined) # PS corr_ps = "corr_PS.svg" plot_corr_ps = paste0(plotdir,"/", corr_ps) -corr_ps_duet_col = "corr_PS_duet_coloured.svg" -plot_corr_ps_duet_col = paste0(plotdir,"/", corr_ps_duet_col) - -corr_upper_ps = "corr_upper_PS.svg" -plot_corr_upper_ps = paste0(plotdir,"/", corr_upper_ps) - # LIG corr_lig = "corr_LIG.svg" plot_corr_lig = paste0(plotdir,"/", corr_lig) -corr_upper_lig = "corr_upper_LIG.svg" -plot_corr_upper_lig = paste0(plotdir,"/", corr_upper_lig) - #################################################################### # end of loading libraries and functions # ######################################################################## #%%%%%%%%%%%%%%%%%%%%%%%%% -df_ps = merged_df3_comp -df_lig = merged_df3_comp_lig +df_ps = merged_df3 +df_lig = merged_df3_lig #%%%%%%%%%%%%%%%%%%%%%%%%% - rm( merged_df2, merged_df2_comp, merged_df2_lig, merged_df2_comp_lig, my_df_u, my_df_u_lig) - -######################################################################## -# end of data extraction and cleaning for plots # ######################################################################## -#============================ -# adding foldx scaled values -# scale data b/w -1 and 1 -#============================ -n = which(colnames(df_ps) == "ddg"); n - -my_min = min(df_ps[,n]); my_min -my_max = max(df_ps[,n]); my_max - -df_ps$foldx_scaled = ifelse(df_ps[,n] < 0 - , df_ps[,n]/abs(my_min) - , df_ps[,n]/my_max) -# sanity check -my_min = min(df_ps$foldx_scaled); my_min -my_max = max(df_ps$foldx_scaled); my_max - -if (my_min == -1 && my_max == 1){ - cat("PASS: foldx ddg successfully scaled b/w -1 and 1" - , "\nProceeding with assigning foldx outcome category") -}else{ - cat("FAIL: could not scale foldx ddg values" - , "Aborting!") -} - - -#================================ -# adding foldx outcome category -# ddg<0 = "Stabilising" (-ve) -#================================= - -c1 = table(df_ps$ddg < 0) -df_ps$foldx_outcome = ifelse(df_ps$ddg < 0, "Stabilising", "Destabilising") -c2 = table(df_ps$ddg < 0) - -if ( all(c1 == c2) ){ - cat("PASS: foldx outcome successfully created") -}else{ - cat("FAIL: foldx outcome could not be created. Aborting!") - exit() -} - -table(df_ps$foldx_outcome) +################################ +# Data for Correlation plots: PS +################################# #====================== # adding log cols @@ -160,9 +105,6 @@ corr_data_ps = df_ps[cols_to_select] dim(corr_data_ps) -#p_italic = substitute(paste("-Log(", italic('P'), ")"));p_italic -#p_adjusted_italic = substitute(paste("-Log(", italic('P adjusted'), ")"));p_adjusted_italic - # assign nice colnames (for display) my_corr_colnames = c("DUET" @@ -171,7 +113,7 @@ my_corr_colnames = c("DUET" , "Log (OR)" , "-Log (P)" - , "AF" + , "MAF" , "duet_outcome" , drug) @@ -194,15 +136,12 @@ head(my_corr_ps) # deep red: #ae301e #--------------------------------------- -# generate corr PS plot 1: both panels +# generate corr PS plot: both panels #--------------------------------------- cat("Corr plot PS DUET with coloured dots:", plot_corr_ps) svg(plot_corr_ps, width = 15, height = 15) -#cat("Corr plot PS DUET with coloured dots:",plot_corr_ps_duet_col) -#svg(plot_corr_ps_duet_col, width = 15, height = 15) - -OutPlot1 = pairs.panels(my_corr_ps[1:(length(my_corr_ps)-1)] +pairs.panels(my_corr_ps[1:(length(my_corr_ps)-1)] , method = "spearman" # correlation method , hist.col = "grey" ##00AFBB , density = TRUE # show density plots @@ -211,12 +150,10 @@ OutPlot1 = pairs.panels(my_corr_ps[1:(length(my_corr_ps)-1)] , rug = F , breaks = "Sturges" , show.points = T - #, bg = c("#f8766d", "#00bfc4")[unclass(factor(my_corr_ps$duet_outcome))] # can't use colour as duet and foldx are opposite - #, pch = 21 # for bg - , pch = 19 + , bg = c("#f8766d", "#00bfc4")[unclass(factor(my_corr_ps$duet_outcome))] # foldx colours are reveresed + , pch = 21 # for bg , jitter = T , alpha = 1 - #, points(pch = 19, col = c("#f8766d", "#00bfc4")) , cex = 1.8 , cex.axis = 2 , cex.labels = 4 @@ -225,40 +162,15 @@ OutPlot1 = pairs.panels(my_corr_ps[1:(length(my_corr_ps)-1)] ) -print(OutPlot1) dev.off() - -#---------------------------------------------- -# generate corr PS plot 2: upper_panel only -#---------------------------------------------- -cat("Corr plot upper PS:", plot_corr_upper_ps) -svg(plot_corr_upper_ps, width = 15, height = 15) - -OutPlot1_upper = my_pp(my_corr_ps[1:(length(my_corr_ps)-1)] # no lower panel - , method = "spearman" # correlation method - , hist.col = "grey" ##00AFBB - , density = TRUE # show density plots - , ellipses = F # show correlation ellipses - , stars = T - , rug = F - , breaks = "Sturges" - , show.points = F - , cex = 3 - , cex.axis = 1.6 - , cex.labels = 4 - , cex.cor = 1 - , smooth = F - -) - -print(OutPlot1_upper) -dev.off() +corr_ps_rho = corr.test(my_corr_ps[1:5], method = "spearman")$r +corr_ps_p = corr.test(my_corr_ps[1:5], method = "spearman")$p ################################################################################################ -#=========================== +################################### # Data for Correlation plots: LIG -#=========================== +################################## table(df_lig$ligand_outcome) df_lig$log10_or_mychisq = log10(df_lig$or_mychisq) @@ -289,7 +201,7 @@ my_corr_colnames = c("Ligand Affinity" , "Log (OR)" , "-Log (P)" - , "AF" + , "MAF" , "ligand_outcome" , drug) @@ -308,13 +220,13 @@ my_corr_lig = corr_data_lig[start:(end-offset)] head(my_corr_lig) #--------------------------------------- -# generate corr LIG plot 1: both panels +# generate corr LIG plot: both panels #--------------------------------------- cat("Corr LIG plot:", plot_corr_lig) svg(plot_corr_lig, width = 15, height = 15) # uncomment as necessary -OutPlot2 = pairs.panels(my_corr_lig[1:(length(my_corr_lig)-1)] +pairs.panels(my_corr_lig[1:(length(my_corr_lig)-1)] , method = "spearman" # correlation method , hist.col = "grey" ##00AFBB , density = TRUE # show density plots @@ -323,12 +235,9 @@ OutPlot2 = pairs.panels(my_corr_lig[1:(length(my_corr_lig)-1)] , rug = F , breaks = "Sturges" , show.points = T - #, bg = c("#f8766d", "#00bfc4")[unclass(factor(my_corr_lig$ligand_outcome))] # can't use colour as duet and foldx are opposite - #, pch = 21 # for bg - , pch = 19 + , bg = c("#f8766d", "#00bfc4")[unclass(factor(my_corr_lig$ligand_outcome))] + , pch = 21 # for bg , jitter = T - #, alpha = .05 - #, points(pch = 19, col = c("#f8766d", "#00bfc4")) , cex = 2 , cex.axis = 2 , cex.labels = 4 @@ -336,35 +245,9 @@ OutPlot2 = pairs.panels(my_corr_lig[1:(length(my_corr_lig)-1)] , smooth = F ) -print(OutPlot2) -dev.off() - - -#--------------------------------------- -# generate corr LIG plot 2: upper panels -#--------------------------------------- -cat("Corr LIG plot:", plot_corr_upper_lig) -svg(plot_corr_upper_lig, width = 15, height = 15) - -# uncomment as necessary -OutPlot2_upper = my_pp(my_corr_lig[1:(length(my_corr_lig)-1)] - , method = "spearman" # correlation method - , hist.col = "grey" ##00AFBB - , density = TRUE # show density plots - , ellipses = F # show correlation ellipses - , stars = T - , rug = F - , breaks = "Sturges" - , show.points = F - #, alpha = .05 - #, points(pch = 19, col = c("#f8766d", "#00bfc4")) - , cex = 3 - , cex.axis = 1.6 - , cex.labels = 4 - , cex.cor = 1 - , smooth = F -) - -print(OutPlot2_upper) dev.off() +corr_lig_rho = corr.test(my_corr_lig[1:4], method = "spearman")$r +corr_lig_p = corr.test(my_corr_lig[1:4], method = "spearman")$p ####################################################### + +