
rm(list = ls())

## plot
library(tidyverse)
library(corrplot)
library(ggcor)
library(ggnewscale)
library(egg)

##
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]

data = dplyr::full_join(hrs, tme, by = "ID")
colnames(data)

##  hypoxia survival plot
library(survival)
library(survminer)
library(ggthemes)

surv = dplyr::filter(data, Cohort %in% c("CGGA_325", "CGGA_693"))
surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 

surv$group = paste(surv$GRS, surv$TME, sep = "/")
table(surv$group)

surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))

fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p1 <- ggsurvplot(fit,
                 title = "Meta_CGGA",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p1



## time-ROC
set.seed(2022)

surv = dplyr::filter(data, Cohort %in% c("CGGA_325", "CGGA_693"))
cox <- coxph(Surv(OS.time, OS) ~ TMEscore + GRSscore, data = surv)
surv$GRSTMEscore=predict(cox, type="risk", newdata=surv)

library(timeROC)
result <-with(surv, timeROC(T=OS.time,
                            delta=OS,
                            marker=GRSTMEscore,
                            cause=1,
                            times=c(12,36,60),
                            iid = TRUE))
plot_dat = data.frame(fpr = as.numeric(result$FP),
                      tpr = as.numeric(result$TP),
                      time = rep(as.factor(c(12,36,60)), each = nrow(result$TP)))

library(ggplot2)
p2 <- ggplot() + 
  geom_line(data = plot_dat, aes(x = fpr, y = tpr, color = time), size = 1) + 
  scale_color_manual(name = NULL,values = c("#00BFC4", "#7CAE00", "#F8766D"),
                     labels = paste0("AUC of ",c(1,3,5)," year: ",
                                     format(round(result$AUC,3),nsmall = 3)))+
  geom_line(aes(x=c(0,1), y=c(0,1)), color = "black")+
  theme_article()+
  theme(axis.text = element_text(color = "black"),
        legend.position = c(0.78,0.125))+
  scale_x_continuous(expand = c(0.005,0.005))+
  scale_y_continuous(expand = c(0.005,0.005))+
  labs(x = "1 - Specificity",
       y = "Sensitivity")
p2

## CGGA_325
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::full_join(hrs, tme, by = "ID")
surv = dplyr::filter(data, Cohort %in% c("CGGA_325"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))


fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p3 <- ggsurvplot(fit,
                 title = "CGGA_325",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p3

## CGGA_693
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::full_join(hrs, tme, by = "ID")
surv = dplyr::filter(data, Cohort %in% c("CGGA_693"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))


fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p4 <- ggsurvplot(fit,
                 title = "CGGA_693",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p4

## IDH WT
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::inner_join(hrs, tme, by = "ID")
cli = read.csv("CGGA.cli.csv")
dat1 = dplyr::inner_join(cli, data, by="ID")
colnames(dat1)
table(dat1$IDH_mutation_status)
surv = dplyr::filter(dat1, IDH_mutation_status %in% c("Wildtype"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))

fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p5 <- ggsurvplot(fit,
                 title = "IDH Wildtype",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p5

## IDH MT
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::inner_join(hrs, tme, by = "ID")
cli = read.csv("CGGA.cli.csv")
dat1 = dplyr::inner_join(cli, data, by="ID")
colnames(dat1)
table(dat1$IDH_mutation_status)
surv = dplyr::filter(dat1, IDH_mutation_status %in% c("Mutant"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))

fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p6 <- ggsurvplot(fit,
                 title = "IDH Mutant",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p6

## WHO II
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::inner_join(hrs, tme, by = "ID")
cli = read.csv("CGGA.cli.csv")
cli = read.csv("CGGA.cli.csv")
dat1 = dplyr::inner_join(cli, data, by="ID")
colnames(dat1)
table(dat1$Grade)
surv = dplyr::filter(dat1, Grade %in% c("WHO II"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))

fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p7 <- ggsurvplot(fit,
                 title = "WHO II",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p7


## WHO III
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::inner_join(hrs, tme, by = "ID")
cli = read.csv("CGGA.cli.csv")
cli = read.csv("CGGA.cli.csv")
dat1 = dplyr::inner_join(cli, data, by="ID")
colnames(dat1)
table(dat1$Grade)
surv = dplyr::filter(dat1, Grade %in% c("WHO III"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))

fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p8 <- ggsurvplot(fit,
                 title = "WHO III",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p8

## WHO IV
hrs = read.csv("0105 Glycolysis.lasson.csv")
tme = read.csv("00202 TME.lasson.csv")
tme = tme[,-c(2,3,4)]
data = dplyr::inner_join(hrs, tme, by = "ID")
cli = read.csv("CGGA.cli.csv")
cli = read.csv("CGGA.cli.csv")
dat1 = dplyr::inner_join(cli, data, by="ID")
colnames(dat1)
table(dat1$Grade)
surv = dplyr::filter(dat1, Grade %in% c("WHO IV"))

surv$GRS <- ifelse(surv$GRSscore > median(surv$GRSscore),'GRSHigh', 'GRSLow')
surv$TME <- ifelse(surv$TMEscore > median(surv$TMEscore),'TMELow', 'TMEHigh') 
surv$group = paste(surv$GRS, surv$TME, sep = "/")
surv$GRSTME = forcats::fct_collapse(surv$group, Mixed = c("GRSHigh/TMEHigh", "GRSLow/TMELow"))
surv$GRSTME = factor(surv$GRSTME, 
                     levels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"),
                     labels = c("GRSLow/TMEHigh", "Mixed", "GRSHigh/TMELow"))

fit <- survfit(Surv(OS.time, OS) ~ GRSTME, data = surv)
p9 <- ggsurvplot(fit,
                 title = "WHO IV",
                 pval = TRUE, 
                 #conf.int = TRUE,
                 risk.table = F, 
                 risk.table.col = "strata",
                 xlab = "Time (months)",
                 ylab = "Overall Survival",
                 #ylab = "Disease Free Survival (%)",
                 linetype = "strata",
                 #surv.median.line = "hv", 
                 ggtheme = theme_article(),
                 #legend.labs = c("Low TME", "High TME"),
                 palette = c("#00BFC4", "#7CAE00", "#F8766D"),
                 legend = c(0.45,0.22),
                 legend.title = "",
                 xlim = c(0,150),
                 #break.x.by = 24,
                 surv.scale = "percent")
p9


### Figure 
###################
library(patchwork)
library(ggplotify)

p = (p1$plot|p2|p3$plot)/(p4$plot|p5$plot|p6$plot)/(p7$plot|p8$plot|p9$plot) 




