library(splitstackshape)
library(dplyr)
library(ggplot2)
library(VennDiagram)
library(RColorBrewer)
library(gridExtra)
library(cowplot)
library(knitr)
library(kableExtra)
library(tidyverse)
library(Hmisc)
library(ggridges)
library(venn)
library(ggthemes)
library(stringr)
library(splitstackshape)
library(chron)
library(magicfor)#install
library(ggtext)#install
library(reshape2)
library(plyr)
library(formattable)
library(data.table)



#Elzar data analysis:
#https://man7.org/linux/man-pages/man2/getrusage.2.html
setwd("/home/asiimwe/projects/run_env/alpha_star_wasp_benchmarking/")
#elzar_runs <- read.table("benchmarking_stats_elzar.txt",comment.char = "", header = FALSE,fill=TRUE)
elzar_runs <- read.table("dataExtractions/benchmark_results_elzar_d10.txt",comment.char = "", header = FALSE,fill=TRUE)

#Data cleaning 
elzar_runs$V2 <- paste0(elzar_runs$V2, " ", elzar_runs$V3)
elzar_runs$V3 <- NULL

#elzar_runs <- elzar_runs %>% filter(V1 != "start_time" & V1 != "end_time")
#elzar_runs$V3 <- as.character(elzar_runs$V3)
#unique(elzar_runs$V3)
#elzar_runs$V3 <- NULL 

#Transposing dataset
#Creating sep variable"
elzar_runs_transposed <- as.data.frame(t(elzar_runs))
colnames(elzar_runs_transposed) <- unlist(elzar_runs_transposed[c(1),])
elzar_runs_transposed <- elzar_runs_transposed[-c(1),]

k <- 8
nr <- nrow(elzar_runs_transposed)
nc <- ncol(elzar_runs_transposed)
unames <- unique(names(elzar_runs_transposed))

a <- array(as.matrix(elzar_runs_transposed), c(nr, k, nc/k))
m <- matrix(aperm(a, c(1, 3, 2)),, k, dimnames = list(NULL, unames))
elzar_runs_transposed <- as.data.frame(m, stringsAsFactors = FALSE)

duplicated(elzar_runs_transposed[,1]) #need to filter given some runs were repeated on the same day - filter by time
head(elzar_runs_transposed$start_time)
elzar_runs_transposed <- cSplit(elzar_runs_transposed, "start_time", " ", direction = "wide", type.convert = "character")
elzar_runs_transposed <- cSplit(elzar_runs_transposed, "end_time", " ", direction = "wide", type.convert = "character")


#elzar_runs_transposed <- elzar_runs_transposed %>% filter(!duplicated(elzar_runs_transposed$jobname))


# #Capture sample, threads, wall clock and memory
# #Wall clock is reported in seconds                         
# elzar_runs_wallclock <- elzar_runs %>% select(V2) %>%  filter("ru_wallclock" )
# 
# 
# time <- "22:13:59.760"
# time2 <- "22:23:11.209"
# 
# library(chron)
# library(lubridate)
# x <- chron(times=time)
# y <- chron(times=time2)
# 
# y-x
# period_to_seconds(hms(y-x)) #551 exact reported wall clock is 551.510

elzar_runs_transposed$cwd_base <- elzar_runs_transposed$cwd
elzar_runs_transposed <- cSplit(elzar_runs_transposed, "cwd", "/", direction = "wide", type.convert="as.character")

#filtering to remove prior runs by name
unique(elzar_runs_transposed$cwd_07)
runs <- c("STAR_Run1",  "STAR_Run2",  "STAR_Run3",  "STAR_WASP_Run1",  "STAR_WASP_Run2",  "STAR_WASP_Run3",  "WASP_Run1",  "WASP_Run2",  "WASP_Run3")
elzar_runs_transposed <- elzar_runs_transposed %>% filter(cwd_07 %in% runs)
unique(elzar_runs_transposed$cwd_07)

elzar_runs_transposed[,c(11:16, 18)] <- NULL
elzar_runs_transposed <- na.omit(elzar_runs_transposed)
colnames(elzar_runs_transposed)[6:13] <- c("start_date", "start_time", "end_date", "end_time", "Path", "Run", "Sample", "Threads")
elzar_runs_transposed$wallclock <- as.numeric(as.character(elzar_runs_transposed$wallclock))
elzar_runs_transposed$ru_wallclock <- as.numeric(as.character(elzar_runs_transposed$ru_wallclock))

elzar_runs_transposed$maxrss <- gsub("G", "", elzar_runs_transposed$maxrss)
elzar_runs_transposed$maxrss <- as.numeric(as.character(elzar_runs_transposed$maxrss))
elzar_runs_transposed$ru_maxrss <- as.numeric(as.character(elzar_runs_transposed$ru_maxrss))



#Plotting Run1
geom_xspline <- function(mapping = NULL, data = NULL, stat = "xspline",
                         position = "identity", show.legend = NA,
                         inherit.aes = TRUE, na.rm = TRUE,
                         spline_shape=-0.25, open=TRUE, rep_ends=TRUE, ...) {
  layer(
    geom = GeomXspline,
    mapping = mapping,
    data = data,
    stat = stat,
    position = position,
    show.legend = show.legend,
    inherit.aes = inherit.aes,
    params = list(spline_shape=spline_shape,
                  open=open,
                  rep_ends=rep_ends,
                  ...)
  )
}

GeomXspline <- ggproto("GeomXspline", GeomLine,
                       required_aes = c("x", "y"),
                       default_aes = aes(colour = "black", size = 0.5, linetype = 1, alpha = NA)
)

stat_xspline <- function(mapping = NULL, data = NULL, geom = "line",
                         position = "identity", show.legend = NA, inherit.aes = TRUE,
                         spline_shape=-0.25, open=TRUE, rep_ends=TRUE, ...) {
  layer(
    stat = StatXspline,
    data = data,
    mapping = mapping,
    geom = geom,
    position = position,
    show.legend = show.legend,
    inherit.aes = inherit.aes,
    params = list(spline_shape=spline_shape,
                  open=open,
                  rep_ends=rep_ends,
                  ...
    )
  )
}

StatXspline <- ggproto("StatXspline", Stat,
                       
                       required_aes = c("x", "y"),
                       
                       compute_group = function(self, data, scales, params,
                                                spline_shape=-0.25, open=TRUE, rep_ends=TRUE) {
                         tf <- tempfile(fileext=".png")
                         png(tf)
                         plot.new()
                         tmp <- xspline(data$x, data$y, spline_shape, open, rep_ends, draw=FALSE, NA, NA)
                         invisible(dev.off())
                         unlink(tf)
                         
                         data.frame(x=tmp$x, y=tmp$y)
                       }
)


elzar_runs_transposed$Threads <- ordered(elzar_runs_transposed$Threads , levels = c("8threads", "16threads", "32threads"))
elzar_runs_transposed <- elzar_runs_transposed %>% filter(Sample != "HG00514" & Sample != "NA12878_Small" & Sample != "NA12878_RAMPAGE" & Sample != "NA12878_RAMPAGE_Rep")
unique(elzar_runs_transposed$Sample)

unique(elzar_runs_transposed$Run)
dim(STAR_Run1 <- elzar_runs_transposed %>% filter(Run == "STAR_Run1")) #96 we expect 16 * 3 = 48
dim(recnt_run_STAR_Run1 <- tail(STAR_Run1, n = 48))
all(recnt_run_STAR_Run1$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_STAR_Run1$Run)

dim(STAR_WASP_Run1 <- elzar_runs_transposed %>% filter(Run == "STAR_WASP_Run1")) #96
(recnt_run_STAR_WASP_Run1 <- tail(STAR_WASP_Run1, n = 48))
all(recnt_run_STAR_WASP_Run1$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_STAR_WASP_Run1$Run)

dim(WASP_Run1 <- elzar_runs_transposed %>% filter(Run == "WASP_Run1")) #84
dim(recnt_run_WASP_Run1 <- tail(WASP_Run1, n = 48))
all(recnt_run_WASP_Run1$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_WASP_Run1$Run)

dim(STAR_Run2 <- elzar_runs_transposed %>% filter(Run == "STAR_Run2"))#80
recnt_run_STAR_Run2 <- tail(STAR_Run2, n = 48)
all(recnt_run_STAR_Run2$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_STAR_Run2$Run)

dim(STAR_WASP_Run2 <- elzar_runs_transposed %>% filter(Run == "STAR_WASP_Run2")) #67
dim(recnt_run_STAR_WASP_Run2 <- tail(STAR_WASP_Run2, n = 48))
all(recnt_run_STAR_WASP_Run2$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_STAR_WASP_Run2$Run)

dim(WASP_Run2 <- elzar_runs_transposed %>% filter(Run == "WASP_Run2"))#60
dim(recnt_run_WASP_Run2 <- tail(WASP_Run2, n = 48))
all(recnt_run_WASP_Run2$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_WASP_Run2$Run)

dim(STAR_Run3 <- elzar_runs_transposed %>% filter(Run == "STAR_Run3"))#48
dim(STAR_WASP_Run3 <- elzar_runs_transposed %>% filter(Run == "STAR_WASP_Run3"))#48

dim(WASP_Run3 <- elzar_runs_transposed %>% filter(Run == "WASP_Run3"))#64
dim(recnt_run_WASP_Run3 <- tail(WASP_Run3, n = 48))
all(recnt_run_WASP_Run3$Sample %in% unique(elzar_runs_transposed$Sample))
unique(recnt_run_WASP_Run3$Run)


#merging dfs
dim(elzar_runs_transposed_filt <- rbind(recnt_run_STAR_Run1,recnt_run_STAR_WASP_Run1, recnt_run_WASP_Run1, recnt_run_STAR_Run2, recnt_run_STAR_WASP_Run2,
                                    recnt_run_WASP_Run2,  STAR_Run3,STAR_WASP_Run3, recnt_run_WASP_Run3)) #should be 9 * 48 = 432

write.csv(elzar_runs_transposed_filt, file="Downstream_Analysis/ref_bias/elzar_runs_filtered")

run1_subset <- elzar_runs_transposed_filt %>% filter(Run == "STAR_Run1" | Run ==  "STAR_WASP_Run1" | Run =="WASP_Run1")
unique(run1_memory$Run)

run2_subset <- elzar_runs_transposed_filt %>% filter(Run == "STAR_Run2" | Run ==  "STAR_WASP_Run2" | Run =="WASP_Run2")
unique(run2_subset$Run)

run3_subset <- elzar_runs_transposed_filt %>% filter(Run == "STAR_Run3" | Run ==  "STAR_WASP_Run3" | Run =="WASP_Run3")
unique(run3_subset$Run)


# Wall Clock
elzar_runs_transposed %>% 
  ggplot(aes(x = reorder(Sample, wallclock), y = wallclock, group=Run,color=factor(Run))) +
  geom_point(color="white") +
  geom_smooth(se=FALSE, linetype="dashed", size=0.5) +
  scale_color_manual(values = c("gray50", "black", "gray80")) +
  facet_wrap(~Threads) + 
  labs(y = "Wall Clock", x="") + scale_y_continuous() +
  theme_bw() + theme(legend.title = element_blank()) + 
  theme(strip.background =element_rect(fill="white", colour = "white"))+
  theme(strip.text = element_text(colour = 'black',size=10), strip.text.x = element_markdown(hjust = 0.5)) +
  theme(axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0.5, size=10)) 

# run_wall clock
elzar_runs_transposed %>% 
  ggplot(aes(x = reorder(Sample, ru_wallclock), y = ru_wallclock, group=Run,color=factor(Run))) +
  geom_point(color="white") +
  geom_smooth(se=FALSE, linetype="dashed", size=0.5) +
  scale_color_manual(values = c("gray50", "black", "gray80")) +
  facet_wrap(~Threads) + 
  labs(y = "Wall Clock", x="") + scale_y_continuous() +
  theme_bw() + theme(legend.title = element_blank()) + 
  theme(strip.background =element_rect(fill="white", colour = "white"))+
  theme(strip.text = element_text(colour = 'black',size=10), strip.text.x = element_markdown(hjust = 0.5)) +
  theme(axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0.5, size=10)) 

# Memory
elzar_runs_transposed %>% 
  ggplot(aes(x = reorder(Sample, maxrss), y = maxrss, group=Run,color=factor(Run))) +
  geom_point(color="white") +
  geom_smooth(se=FALSE, linetype="dashed", size=0.5) +
  scale_color_manual(values = c("firebrick3", "dodgerblue4", "orange")) +
  facet_wrap(~Threads) + 
  labs(y = "Memory (GB)", x="") + scale_y_continuous() +
  theme_bw() + theme(legend.title = element_blank()) + 
  theme(strip.background =element_rect(fill="white", colour = "white"))+
  theme(strip.text = element_text(colour = 'black', size = 10), strip.text.x = element_markdown(hjust = 0.5)) +
  theme(axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0.5, size=10)) 

# ru_maxrss
elzar_runs_transposed %>% 
  ggplot(aes(x = reorder(Sample, ru_maxrss), y = ru_maxrss, group=Run,color=factor(Run))) +
  geom_point(color="white") +
  geom_smooth(se=FALSE, linetype="dashed", size=0.5) +
  scale_color_manual(values = c("firebrick3", "dodgerblue4", "orange")) +
  facet_wrap(~Threads) + 
  labs(y = "Memory (GB)", x="") + scale_y_continuous() +
  theme_bw() + theme(legend.title = element_blank()) + 
  theme(strip.background =element_rect(fill="white", colour = "white"))+
  theme(strip.text = element_text(colour = 'black', size = 10), strip.text.x = element_markdown(hjust = 0.5)) +
  theme(axis.text.x = element_text(angle = 90,hjust = 1,vjust = 0.5, size=10))