tried ED logo, but needs work

This commit is contained in:
Tanushree Tunstall 2022-01-18 16:53:16 +00:00
parent 00094f036a
commit 8f8a9db92c
3 changed files with 119 additions and 108 deletions

View file

@ -166,38 +166,42 @@ LogoPlotMSA <- function(msaSeq_mut
, method = msa_method , method = msa_method
, col_scheme = my_logo_col , col_scheme = my_logo_col
, seq_type = 'aa') + , seq_type = 'aa') +
scale_x_discrete(x_lab theme(legend.position = leg_pos
, breaks = i_extract , legend.direction = leg_dir
, labels = i_extract #, legend.title = element_blank()
#, limits = min(i_extract): max(i_extract)) , legend.title = element_text(size = leg_tts
, limits = factor(i_extract)) , colour = ytt_col)
, legend.text = element_text(size = leg_ts)
# further customisation
msa_mut_logo_P = p0 + theme(legend.position = leg_pos
, legend.direction = leg_dir
#, legend.title = element_blank()
, legend.title = element_text(size = leg_tts
, colour = ytt_col)
, legend.text = element_text(size = leg_ts)
, axis.text.x = element_text(size = x_ats , axis.text.x = element_text(size = x_ats
, angle = x_tangle , angle = x_tangle
, hjust = 1 , hjust = 1
, vjust = 0.4 , vjust = 0.4
, colour = xfont_bgc) , colour = xfont_bgc)
#, axis.text.y = element_blank() #, axis.text.y = element_blank()
, axis.text.y = element_text(size = y_ats , axis.text.y = element_text(size = y_ats
, angle = y_tangle , angle = y_tangle
, hjust = 1 , hjust = 1
, vjust = -1.0 , vjust = -1.0
, colour = yfont_bgc) , colour = yfont_bgc)
, axis.title.x = element_text(size = x_tts , axis.title.x = element_text(size = x_tts
, colour = xtt_col) , colour = xtt_col)
, axis.title.y = element_text(size = y_tts , axis.title.y = element_text(size = y_tts
, colour = ytt_col) , colour = ytt_col)
, plot.background = element_rect(fill = theme_bgc))+
, plot.background = element_rect(fill = theme_bgc)) xlab(x_lab)
if (missing(plot_positions)){
msa_mut_logo_P = p0
}else{
msa_mut_logo_P = p0 +
scale_y_continuous(expand = c(0,0.09)) +
scale_x_discrete(breaks = i_extract
, expand = c(0.09,0)
, labels = i_extract
, limits = factor(i_extract))
}
cat('\nDone: msa_mut_logo_P') cat('\nDone: msa_mut_logo_P')
#return(msa_mut_logoP) #return(msa_mut_logoP)
@ -210,15 +214,7 @@ LogoPlotMSA <- function(msaSeq_mut
, facet = "grid" , facet = "grid"
, method = msa_method , method = msa_method
, col_scheme = my_logo_col , col_scheme = my_logo_col
, seq_type = 'aa')+ , seq_type = 'aa') +
scale_x_discrete(x_lab
, breaks = i_extract
, labels = i_extract
#, limits = min(i_extract): max(i_extract))
, limits = factor(i_extract))
# further customisation
msa_wt_logo_P = p1 +
theme(legend.position = "none" theme(legend.position = "none"
, legend.direction = leg_dir , legend.direction = leg_dir
@ -240,8 +236,20 @@ LogoPlotMSA <- function(msaSeq_mut
, colour = ytt_col) , colour = ytt_col)
, plot.background = element_rect(fill = theme_bgc)) + , plot.background = element_rect(fill = theme_bgc)) +
ylab("") ylab("") + xlab("Wild-type position")
if (missing(plot_positions)){
msa_wt_logo_P = p1
}else{
msa_wt_logo_P = p1 +
scale_y_continuous(expand = c(0,0.09)) +
scale_x_discrete(breaks = i2_extract
, expand = c(0.09,0)
, labels = i2_extract
, limits = factor(i2_extract))
}
cat('\nDone: msa_wt_logo_P') cat('\nDone: msa_wt_logo_P')
#return(msa_wt_logoP) #return(msa_wt_logoP)
LogoPlotMSAL[['msa_wt_logoP']] <- msa_wt_logo_P LogoPlotMSAL[['msa_wt_logoP']] <- msa_wt_logo_P

View file

@ -7,83 +7,84 @@ source("~/git/LSHTM_analysis/config/gid.R")
# source("~/git/LSHTM_analysis/config/katg.R") # source("~/git/LSHTM_analysis/config/katg.R")
# source("~/git/LSHTM_analysis/config/rpob.R") # source("~/git/LSHTM_analysis/config/rpob.R")
#--------------------------------------------------- #---------------------------------------------------
source("~/git/LSHTM_analysis/scripts/plotting/get_plotting_dfs.R") #source("~/git/LSHTM_analysis/scripts/plotting/get_plotting_dfs.R")
#
################################ # ################################
# Logo plot with custom Y axis # # Logo plot with custom Y axis
# mainly OR # # mainly OR
# script: logoP.R # # script: logoP.R
################################ # ################################
LogoPlotCustomH (plot_df = merged_df3 # LogoPlotCustomH (plot_df = merged_df3
, x_axis_colname = "position" # , x_axis_colname = "position"
, y_axis_colname = "or_mychisq" # , y_axis_colname = "or_mychisq"
, symbol_colname = "mutant_type" # , symbol_colname = "mutant_type"
, y_axis_log = F # , y_axis_log = F
, log_value = log10 # , log_value = log10
, y_axis_increment = 5 # , y_axis_increment = 5
, rm_empty_y = F # , rm_empty_y = F
, my_logo_col = 'chemistry' # , my_logo_col = 'chemistry'
, x_lab = "Wild-type position" # , x_lab = "Wild-type position"
, y_lab = "Odds Ratio" # , y_lab = "Odds Ratio"
, x_ats = 12 # text size # , x_ats = 12 # text size
, x_tangle = 90 # text angle # , x_tangle = 90 # text angle
, y_ats = 22 # , y_ats = 22
, y_tangle = 0 # , y_tangle = 0
, x_tts = 19 # title size # , x_tts = 19 # title size
, y_tts = 22 # , y_tts = 22
#, leg_pos = c(0.05,-0.12) # #, leg_pos = c(0.05,-0.12)
, leg_pos = "top" # , leg_pos = "top"
, leg_dir = "horizontal" # , leg_dir = "horizontal"
, leg_ts = 15 # leg text size # , leg_ts = 15 # leg text size
, leg_tts = 16 # leg title size # , leg_tts = 16 # leg title size
) # )
#
#
######################################## # # ########################################
# Logo plot showing nsSNPs by positions # # # Logo plot showing nsSNPs by positions
# wild-type and mutant aa # # # wild-type and mutant aa
# script: logoP_snp.R # # # script: logoP_snp.R
######################################## # # ########################################
LogoPlotSnps(plot_df = merged_df3 # LogoPlotSnps(plot_df = merged_df3
, x_axis_colname = "position" # , x_axis_colname = "position"
, symbol_mut_colname = "mutant_type" # , symbol_mut_colname = "mutant_type"
, symbol_wt_colname = "wild_type" # , symbol_wt_colname = "wild_type"
, omit_snp_count = c(0)# can be 1, 2, etc. # , omit_snp_count = c(0)# can be 1, 2, etc.
, my_logo_col = "chemistry" # , my_logo_col = "chemistry"
, x_lab = "Wild-type position" # , x_lab = "Wild-type position"
, y_lab = "" # , y_lab = ""
, x_ats = 13 # text size # , x_ats = 13 # text size
, x_tangle = 90 # text angle # , x_tangle = 90 # text angle
, y_ats = 20 # , y_ats = 20
, y_tangle = 0 # , y_tangle = 0
, x_tts = 18 # title size # , x_tts = 18 # title size
, y_tts = 18 # , y_tts = 18
, leg_pos = "top" # can be top, left, right and bottom or c(0.8, 0.9) # , leg_pos = "top" # can be top, left, right and bottom or c(0.8, 0.9)
, leg_dir = "horizontal" # can be vertical or horizontal # , leg_dir = "horizontal" # can be vertical or horizontal
, leg_ts = 14 # leg text size # , leg_ts = 14 # leg text size
, leg_tts = 16 # leg title size # , leg_tts = 16 # leg title size
) # )
#
#
######################################## ########################################
# Logo plot MSA # Logo plot MSA
# Mutant and wild-type # Mutant and wild-type
# wild-type and mutant aa # wild-type and mutant aa
# Can select active site residues
# specify {plot_positions}
# To plot entire MSA, simply don't specify {plot_positions}
# script: logoP_msa.R # script: logoP_msa.R
######################################## # TODO perhaps: ED logo from Logolas
# BOTH WORK
#LogoPlotMSA(msa_seq, wt_seq)
########################################
LogoPlotMSA(msaSeq_mut = msa_seq LogoPlotMSA(msaSeq_mut = msa_seq
, msaSeq_wt = wt_seq , msaSeq_wt = wt_seq
, msa_method = 'bits' # or probability , msa_method = 'bits' # or probability
, my_logo_col = "chemistry" , my_logo_col = "taylor"
#, plot_positions = active_aa_pos #, plot_positions = active_aa_pos
, plot_positions , x_lab = "nsSNP position"
, x_lab = "Wild-type position"
, y_lab = "" , y_lab = ""
, x_ats = 13 # text size , x_ats = 10 # text size
, x_tangle = 90 # text angle , x_tangle = 90 # text angle
, y_ats = 13 , y_ats = 13
, y_tangle = 0 , y_tangle = 0
@ -93,4 +94,4 @@ LogoPlotMSA(msaSeq_mut = msa_seq
, leg_dir = "horizontal" #can be vertical or horizontal , leg_dir = "horizontal" #can be vertical or horizontal
, leg_ts = 16 # leg text size , leg_ts = 16 # leg text size
, leg_tts = 16 # leg title size , leg_tts = 16 # leg title size
) )

View file

@ -38,7 +38,9 @@ head(msa_seq)
# fasta file: wt # fasta file: wt
#================ #================
in_filename_fasta = paste0(tolower(gene), ".1fasta") #in_filename_fasta = paste0(tolower(gene), ".1fasta")
in_filename_fasta = paste0(tolower(gene), "2.2fasta") # gid
infile_fasta = paste0(indir, "/", in_filename_fasta) infile_fasta = paste0(indir, "/", in_filename_fasta)
cat("\nInput fasta file for WT: ", infile_fasta, "\n") cat("\nInput fasta file for WT: ", infile_fasta, "\n")