#########################################################
# 1: Installing and loading required packages
#########################################################

#source("../Header_TT.R")
install.packages("qqman")
library(qqman)

source("combining_dfs_plotting.R")
#mcsm_data: raw file, 225, 15
#merged_df2 = 2201, 35
#merged_df3 =  205, 35 ("Can't trust non-numerical params')

#===============================================
# PLOTS: DUET vs GWAS: non-numerical
# lineage, country_code, etc
# merged_df2: 1592, 35
#===============================================

#########################
# Data for plot
#########################
df = merged_df2
#df = merged_df2_comp 


#========================
# Plot 1a: Lineage barplot
# x = lineage y = No of samples
# col = Lineage
# fill = lineage
#========================
table(df$lineage)

# subset only lineages1-4
sel_lineages = c("lineage1"
                 , "lineage2"
                 , "lineage3"
                 , "lineage4"
                 #, "lineage5"
                 #, "lineage6"
                 #, "lineage7"
)

# uncomment as necessary
df_lin = subset(df, subset = lineage %in% sel_lineages )
table(df_lin$lineage)

#%%%%%%%%%%%%%%%%%%%%%%%%%
# REASSIGNMENT
df <- df_lin
#%%%%%%%%%%%%%%%%%%%%%%%%%

#%%%%%%%%%%%%%%%%%%%%%%%%
# REASSIGNMENT
df2 = df
#%%%%%%%%%%%%%%%%%%%%%%%%
df2 = df2%>%
  add_count(country_code)

str(df2$country_code); str(df2$n)

n = which(colnames(df2) == "n")
colnames(df2)[n] = "count_country"

table(df2$count_country>100 & df$country_code!= "")
df3 = subset(df2, df2$count_country>100 & df2$country_code != "")


#%%%%%%%%%%%%%%%%%%%%%%%%
# REASSIGNMENT
df = df3
#%%%%%%%%%%%%%%%%%%%%%%%%

sample = sum(table(unique(df$id))); sample
table(df$country_code)
tab = sum(table(df$country_code)); tab


View(table(df$country_code))
View(t1)

############## begin plot
g = ggplot(df, aes(x = lineage))
g + geom_bar(aes(fill = lineage)) +
  theme( axis.text.x = element_text(size = 13
                                    , angle = 90
                                    , hjust = 1
                                    , vjust = 0.4)
         , 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) ) +
  labs(title = "Lineage"
       , x = "Lineage"
       , y = "No of samples")


#========================
# Plot 2: DUET, lineage, country_code and or_mychisq
# x = lineage y = DUET
# col = Lineage
# fill = country_code
#========================
### begin plot
g = ggplot(df, aes(x = country_code
                   , y = duet_scaled))
g + geom_point(aes(col = lineage
                   , size = or_mychisq)) +
  theme(axis.text.x = element_text(size = 13
                                    , angle = 90
                                    , hjust = 1
                                    , vjust = 0.4)
         , 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) ) +
  labs(title = "DUET, country_code, lineage, or_mychisq"
       , x = "Lineage"
       , y = "DUET (PS)")


#############
#========================
# Plot 3: DUET, lineage, or_mychisq
# x = lineage y = DUET
# col = Lineage
# fill = country_code
#========================

### begin plot
table(df$lineage)

g = ggplot(df_lin, aes(x = lineage
                   , y = duet_scaled))
g + geom_point(aes(col = lineage
                   , size = or_mychisq)) +
  theme(axis.text.x = element_text(size = 13
                                   , angle = 90
                                   , hjust = 1
                                   , vjust = 0.4)
        , 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) ) +
  labs(title = "DUET, lineage, or_mychisq"
       , x = "Lineage"
       , y = "DUET (PS)")

#========================
# Plot 4-5: Distributions
# ggrdiges
#========================


#==================================================
my_ats = 15 # axis text size
my_als = 20 # axis label size

my_labels = c('Lineage 1', 'Lineage 2', 'Lineage 3', 'Lineage 4'
              #, 'Lineage 5', 'Lineage 6', 'Lineage 7'
              )
names(my_labels) = c('lineage1', 'lineage2', 'lineage3', 'lineage4'
                     # , 'lineage5', 'lineage6', 'lineage7'
                     )
      

#========================
# Plot 4: Distribution
# x = duet_scaled
# y = country
# fill = country_code
# facet = lineage
#========================
# works neatly!

p1 = ggplot(df, aes(x = duet_scaled
                    , y = country_code))+
  
  #printFile=geom_density_ridges_gradient(
  geom_density_ridges_gradient(aes(fill = country_code)
                               , jittered_points = TRUE
                               , scale = 3
                               , size = 0.3 ) +
  facet_wrap( ~lineage
              , scales = "free"
              , switch = 'x'
              , labeller = labeller(lineage = my_labels)
              ) +
  coord_cartesian( xlim = c(-1, 1)) +
  #scale_fill_gradientn(colours = c("#f8766d", "white", "#00bfc4")
  #                     , name = "DUET" ) + 
  theme(axis.text.x = element_text(size = my_ats
                                   , angle = 90
                                   , hjust = 1
                                   , vjust = 0.4)
        
        #, axis.text.y = element_blank()
        , axis.title.x = element_blank()
        , axis.title.y = element_blank()
        , axis.ticks.y = element_blank()
        , plot.title = element_blank()
        , strip.text = element_text(size = my_als)
        , legend.text = element_text(size = my_als-5)
        , legend.title = element_text(size = my_als)
  ) 

p1


#========================
# Plot 5: Distribution
# x = duet_scaled
# y = country_code
# fill = lineage
# facet = NONE
#========================
# no facet wrap

p2 = ggplot(df, aes(x = duet_scaled
                    , y = country_code))+
  
  geom_density_ridges_gradient(aes(fill = factor(lineage))
                               , scale = 3
                               , size = 0.3 ) +
  coord_cartesian( xlim = c(-1, 1)) +
  #scale_fill_gradientn(colours = c("#f8766d", "white", "#00bfc4")
  #                     , name = "DUET" ) + 
  #scale_fill_continuous(colours = c("darkgreen", "pink", "orange", "brown")
  #                     , name = "lineage" ) + 
  theme(axis.text.x = element_text(size = my_ats
                                   , angle = 90
                                   , hjust = 1
                                   , vjust = 0.4)
        
        #, axis.text.y = element_blank()
        , axis.title.x = element_blank()
        , axis.title.y = element_blank()
        , axis.ticks.y = element_blank()
        , plot.title = element_blank()
        , strip.text = element_text(size = my_als)
        , legend.text = element_text(size = my_als-5)
        , legend.title = element_text(size = my_als)
  ) 

p2


#===============
# lineage only
#================
#svg(plot_lineage_duet)
p3 = ggplot(df, aes(x = duet_scaled
                    , y = duet_outcome))+
  geom_density_ridges_gradient(aes(fill = ..x..)
                               , jittered_points = TRUE
                               , scale = 3
                               , size = 0.3 ) +
  facet_wrap( ~lineage
              , scales = "free"
              #, switch = 'x'
              , labeller = labeller(lineage = my_labels) ) +
  coord_cartesian( xlim = c(-1, 1)) +
  scale_fill_gradientn(colours = c("#f8766d", "white", "#00bfc4")
                       , name = "DUET" ) + 
  theme(axis.text.x = element_text(size = my_ats
                                   , angle = 90
                                   , hjust = 1
                                   , vjust = 0.4)
        
        , axis.text.y = element_blank()
        , axis.title.x = element_blank()
        , axis.title.y = element_blank()
        , axis.ticks.y = element_blank()
        , plot.title = element_blank()
        , strip.text = element_text(size = my_als)
        , legend.text = element_text(size = my_als-5)
        , legend.title = element_text(size = my_als)
  ) 

print(p3)