library(REffectivePred)
## Read in the data
path_to_data <- system.file("extdata/NY_OCT_4_2022.csv", package = "REffectivePred")
data <- read.csv(path_to_data)
head(data)
cases <- diff(c(0, data$cases)) # Convert cumulative cases into daily cases
lt <- length(cases) # Length of cases
Time <- as.Date(data$date, tryFormats = c("%d-%m-%Y", "%d/%m/%Y"))
navigate_to_config() # Open the config file, make any necessary changes here.
path_to_config <- system.file("config.yml", package = "REffectivePred") # Read config file
cfg <- load_config() # Build the cfg object
##### Option 1: populate the global environment with args to pass to function.
population <- cfg$population # Population size
window_size <- cfg$window.size
adj.period <- cfg$adj.period
fit.t.pred <- cfg$fit.t.pred # Time of prediction
not.predict <- cfg$not.predict
rt.func.num <- cfg$rt.func.num # choose which Rt function you want to use
num.iter <- cfg$num.iter
silence.errors <- cfg$silence.errors
predict.beyond <- cfg$predict.beyond
curve_params <- as.double(unlist(cfg$curve_params))
vt_params <- as.double(unlist(cfg$vt_params)) # The vt initial values, starting at wave 2
restriction_levels <- as.double(unlist(cfg$restriction_levels)) # Psi, u, and v parameters
betas <- as.double(unlist(cfg$betas)) # betas
ini_params <- c(curve_params, vt_params, restriction_levels, betas)
restrictions_params <- cfg$restrictions_params
restriction_st_params <- cfg$restriction_st_params
param_scale <- abs(ini_params) / 10
waves_list <- ranges_to_waves(cfg$waves_list)
params_limits <- cfg$params_limits
num_waves <- cfg$num_waves
waves <- waves_1d_list(num_waves, waves_list)
rho <- eval(parse(text = cfg$rho))
serial_mean <- cfg$serial_mean
serial_var <- cfg$serial_var
est <- estimate.mle(
ini_params = ini_params,
params_limits = params_limits,
restrictions = restrictions_params,
restriction.starts = restriction_st_params,
ranges = waves,
rt_func = rt.func.num,
silence.errors = silence.errors,
fit.t.pred = fit.t.pred,
param_scale = param_scale,
num.iter = num.iter,
cases = cases,
scenario = NULL,
H.E = NULL,
H.W = NULL,
adj.period = adj.period,
population = population,
rho = rho,
serial_mean = serial_mean,
serial_var = serial_var,
lt = lt,
window_size = window_size,
hessian = FALSE
)
print(est)
##### Option 2: pass the cfg object instead.
est <- estimate.mle(
cases = cases,
cfg = cfg,
hessian = FALSE
)
print(est)
Run the code above in your browser using DataLab