2_combine_res_novax.R

######################################
####Combine res novax for FOISWEEP####
######################################

library(dplyr)
library(deSolve)
library(ggplot2)
library(tidyr)
#####Read in mod_scen for serologically triggered scenarios and no vax scenario
hpc_folder <- "/projects/blopman/vger/cliu/00_sens_varyfoi"
setwd(hpc_folder)

#setwd("C:/Users/cliu369/OneDrive - Emory University/Documents/DIssertation/Aim 3/COVIDVaxSero_sens/8_varyfoi")
#Read in lists of model results
mod_scen1 <- readRDS("mod_scen_novax1.RDS")
mod_scen2 <- readRDS("mod_scen_novax2.RDS")
mod_scen3 <- readRDS("mod_scen_novax3.RDS")
mod_scen4 <- readRDS("mod_scen_novax4.RDS")[1:32]
mod_scen5 <- readRDS("mod_scen_novax_foisp0.9.RDS")
#Combine the pieces into one
mod_scen <- c(mod_scen1,mod_scen2,mod_scen3,mod_scen4,mod_scen5)

##Working directories for scripts in Carol's home drive
setwd("/home/cliu369/COVIDVaxSero_sens/8_varyfoi")
##Read in initial sweep scenarios
#sweep <- readRDS("1_main/0_sweep_imm4_varyfoi_novax.RDS")
##Get the sweep unique to correctly identify from the r0 draws
sweep_unique <- readRDS("0_postprocess/scraps/00_sweep_unique.RDS")


#Number of model runs needing processing
n   <- length(mod_scen)

##Inputs: n, mod_scen, sweep associated with mod_scen

##Outputs:vax dose and deaths needed for nnt
source("0_postprocess/scraps/9_nnt_sero.R")
nnt <- nnt%>%
  select(-sweep_unique)%>%
  left_join(sweep_unique%>%select(sweep_unique:r03), by=c("r02","r03"))

##Output: Time series case, sero, and death
source("0_postprocess/scraps/9_case_sero_death_timeseries.R")
inc    <- inc%>%   left_join(sweep_unique%>%select(sweep_unique:r03), by =c("r02","r03"))
deaths <- deaths%>%left_join(sweep_unique%>%select(sweep_unique:r03), by =c("r02","r03"))
sero   <- sero%>%  left_join(sweep_unique%>%select(sweep_unique:r03), by =c("r02","r03"))

##Output: Time series immunity tiers and median value
#list1 <- which(sweep$sero_thresh %in% c(0))
list1 <- seq(1:n)
source("0_postprocess/scraps/9_imm_timeseries.R")
imm   <- imm%>%  left_join(sweep_unique%>%select(sweep_unique:r03), by =c("r02","r03"))

setwd(hpc_folder)
saveRDS(nnt, "0_comb_res/novax_nnt.RDS")
saveRDS(inc,    "0_comb_res/novax_inc.RDS")
saveRDS(sero,   "0_comb_res/novax_sero.RDS")
saveRDS(imm,    "0_comb_res/novax_imm.RDS")
saveRDS(imm_med,"0_comb_res/novax_imm_med.RDS")
saveRDS(deaths, "0_comb_res/novax_death.RDS")