| estimate.mle {REffectivePred} | R Documentation |
Fit the Model
Description
Estimate the parameters of the model by maximizing the likelihood function (or, rather, by minimizing the negative log likelihood).
Usage
estimate.mle(
hessian = FALSE,
H.E = NULL,
H.W = NULL,
cases = NULL,
cfg = NULL,
ini_params = NULL,
params_limits = NULL,
restrictions = NULL,
restriction.starts = NULL,
ranges = NULL,
rt_func = 1,
silence.errors = FALSE,
fit.t.pred = NULL,
param_scale = NULL,
num.iter = NULL,
scenario = NULL,
adj.period = NULL,
population = NULL,
rho = NULL,
serial_mean = serial_mean,
serial_var = serial_var,
lt = NULL,
window_size = NULL,
verbose = FALSE
)
Arguments
hessian |
Logical. If TRUE, computes the variance-covariance matrix at the MLE. |
H.E |
Mobility metrics for category Retail & Entertainment. Currently unsupported. |
H.W |
Mobility metrics for category Workplaces. Currently unsupported. |
cases |
Vector of case counts. |
cfg |
The object that contains all variables from the configuration file.
This includes all function arguments except for |
ini_params |
Initial parameter values to be used in optimization. Includes the following sets of parameters in a vector, in this order:
|
params_limits |
Boundaries/limits of the ini_params. |
restrictions |
A numeric integer vector giving the severity of restrictions. Zero means no restriction, and higher numbers means greater severity/disruption. The ordered unique values should be consecutive integers starting from zero. Each number (other than 0) adds a new parameter to the fit. restrictions = NULL causes the function to use mobility data instead of the psi values (currently unsupported). |
restriction.starts |
A vector of same length as restrictions, of times when restrictions came into effect. Note: the first index time should be 1. |
ranges |
An vector of time ranges for the different waves. The waves ranges should be contiguous, with at least one unit of time between consecutive waves. |
rt_func |
The parametric form of function c(). Options are listed under function c_helper. |
silence.errors |
Logical. If TRUE, ignores certain errors to allow optimization to proceed. Not all errors can be ignored. |
fit.t.pred |
Time from which prediction is done. If use.actual.not.predicted is TRUE, values of |
param_scale |
Parameter scale. Passed as argument parscale to optim. |
num.iter |
Maximum number of iterations. Passed as argument maxit to optim. |
scenario |
A character string describing options to deal with restrictions. Currently unsupported. |
adj.period |
Delays in society adjusting. |
population |
total population size. |
rho |
A vector of under-reporting rates of the same length as cases. If a scalar is supplied, the vector will be constant with this value. |
serial_mean |
Mean of the serial interval on the log scale. |
serial_var |
Variance of the serial interval on the log scale. |
lt |
The length of cases. |
window_size |
The maximum value for the serial interval. |
verbose |
Logical. If TRUE, provides additional details while running the function. |
Value
A list of maximum likelihood estimates of the parameters. Includes:
a1
a2
a3
a4
nu
vt_params_est
Psi
betas
negative_log_lik
mle
hessian
SE
Examples
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)