R/Prehos.R

In cyhsuTN/COVID19: Prediction of Cases Volumne and Hospitalization by SIR Model

#' Prediction of hospitalized, ICU and ventilated cases.
#' @description Prediction of hospitalized, ICU and ventilated cases based on SIR model.
#' @import stats
#' @importFrom utils tail head

#' @param obj Input. Object from function \emph{fitSIR}.
#' @param inihos Input. Initial number of hospitalized cases.
#' @param iniicu Input. Initial number of ICU cases.
#' @param iniven Input. Initial number of ventilated cases.
#' @param hosrate Input. Hospitalization rate of infected people (percentage between 0 to 100).
#' @param icurate Input. ICU rate of infected people (percentage between 0 to 100).
#' @param venrate Input. Ventilated rate of infected people (percentage between 0 to 100).
#' @param outhosdays Input. Hospital Length of Stay (days).
#' @param outicudays Input. ICU Length of Stay (days).
#' @param outvendays Input. Vent Length of Stay (days).
#' @param hms Input. Hospital market share (percentage between 0 to 100).


#' @examples ## To predicte 100 days from today (dayFT=100).
#' @examples casevolumne <- fitSIR(susceptible=4119405, Infected=3733, inihos=14,
#' @examples      hosrate=2.5, hms=15, inidbt=4, mrt=14, sdis=30, dayFT=100)
#' @examples hospitalization <- Prehos(casevolumne, inihos=14, iniicu=0, iniven=0,
#' @examples      hosrate=2.5, icurate=0.75, venrate=0.5, outhosdays=7, outicudays=9,
#' @examples      outvendays=10, hms=15)
#' @examples head(hospitalization, 21) ## show the first 20 days
#'

#' @export


Prehos <- function(obj, inihos=14, iniicu=0, iniven=0, hosrate=2.5, icurate=0.75, venrate=0.5,
                   outhosdays=7, outicudays=9, outvendays=10, hms=15) {

  hrate <- (hms/100) * (hosrate/100)
  irate <- (hms/100) * (icurate/100)
  vrate <- (hms/100) * (venrate/100)

  ts <- length(obj$result$susceptible)-1

  allnewinflected <- obj$result$susceptible[1] - tail(obj$result$susceptible,-1)
  reps <- cbind(1:ts, (allnewinflected * hrate) + inihos,
                (allnewinflected * irate) + iniicu,
                (allnewinflected * vrate) + iniven)

  reps <- rbind(c(0, inihos, iniicu, iniven), reps)

  set.seed(20200326)
  rdrop.hos <- cumsum(c(rmultinom(n=1, size=inihos, prob=rep(1/outhosdays,outhosdays))))
  rdrop.icu <- cumsum(c(rmultinom(n=1, size=iniicu, prob=rep(1/outicudays,outicudays))))
  rdrop.ven <- cumsum(c(rmultinom(n=1, size=iniven, prob=rep(1/outvendays,outvendays))))

  minus.h <- head(c(0, rdrop.hos, reps[-1,2]), -outhosdays)
  minus.i <- head(c(0, rdrop.icu, reps[-1,3]), -outicudays)
  minus.v <- head(c(0, rdrop.ven, reps[-1,4]), -outvendays)

  reps[,2] <- reps[,2] - minus.h
  reps[,3] <- reps[,3] - minus.i
  reps[,4] <- reps[,4] - minus.v

  reps[reps<0] <- 0

  colnames(reps) <- c("days.from.today","hosp","icu","vent")
  ceiling(data.frame(reps))

}