get_weekly_local.R

#' Extract information of local and weekly estimates from simulation 
#'
#'
#' @param sim An \code{icm} object returned by \link{simulate_seiqhrf}. 
#' @param market.share between 0 and 1, percentage of local hospital beds in 
#'        the simulated unit (e.g. state)
#' @param icu_percent between 0 and 1, percentage of patients that should go to 
#'        ICU among the ones that need hospitalization
#' @param start_date Epidemic start date. Default is 'na', if not provided will 
#'        plot week numbers, if provided will plot the first day (Sunday) of the
#'        week.
#' @param time_limit Number of days to include. Default = 90.
#' @param total_population True population size, needed only if simulation size 
#'        is smaller than the true population size due to computational cost 
#'        etc.
#' 
#' @return 
#' \itemize{
#' \item \code{plot:} A \code{ggplot} object, bar charts of count of patients 
#'              requiring hospitalization and ICU respectively
#' \item \code{result:} A dataframe
#'    \itemize{\item \code{week:}  week number from input \code{sim},
#'             \item \code{hosp:} the number of patients that require hospitalization locally,
#'             \item \code{icu:} the number of patients that require ICU locally. }
#
#' }
#' 
#' @importFrom tidyr pivot_wider
#' 
get_weekly_local <- function(sim, 
                             market.share = .04,
                             icu_percent = .1, 
                             start_date = 'na',
                             time_lim = 90,
                             total_population = NULL){
  
  hosp <- sim$df$h.num
  
  if(!is.null(total_population)){
    if(total_population < max(sim$df$s.num)) 
      stop("total Population should be larger than simulated size")
    cat("Scalling w.r.t total population")
    hosp <- hosp*total_population/max(sim$df$s.num)
  } 
  
  if(market.share < 0 || market.share > 1) stop("Market share has to be between 
                                                0 and 1")
  if(icu_percent < 0 || icu_percent > 1) stop("ICU percentage has to be between
                                              0 and 1")
  
  hosp[is.na(hosp)] <- 0
  hosp <- hosp[1: time_lim]
  
  hosp_week <- split(hosp, ceiling(seq_along(hosp)/7))
  hosp_sum_week <- unlist(lapply(hosp_week, sum))
  t_sz <- length(hosp_sum_week)
  
  hosp_wk_df <- data.frame(wk = rep(seq_along(hosp_sum_week), 2), 
                           group = rep(c("general", "icu"), 
                                       each = t_sz),
                           hosp_icu = c(hosp_sum_week - 
                                           (hosp_sum_week*icu_percent), 
                                         hosp_sum_week*icu_percent))
  
  if(class(start_date) == 'Date'){
    
    hosp_wk_df <- data.frame(append(hosp_wk_df,
                                    list(Date=start_date + 
                                           (7 * (hosp_wk_df$wk - 1))),
                                    after=match("wk", names(hosp_wk_df))))
    
    gg <- ggplot(data=hosp_wk_df, aes(x = Date, y = hosp_icu, fill = group)) +
      geom_bar(stat="identity") + theme_bw() +
      scale_x_date(date_breaks = "1 week", date_labels = "%m-%d") + 
      labs(y="Weekly Hospital Load (sum over week)", x = "Week")
    
  }else{

    gg <- ggplot(data=hosp_wk_df, aes(x = wk, y = hosp_icu, fill = group)) +
      geom_bar(stat="identity") + theme_bw() +
      labs(y="Weekly Hospital Load (sum over week)", x = "Week") + 
      scale_x_continuous(breaks = seq(0,t_sz,5), labels= seq(0,t_sz,5))
  }
  
  res <- hosp_wk_df %>% tidyr::pivot_wider(names_from = group, values_from = hosp_icu)
                    
  return(list("plot" = gg, "result" = res))
 
}