Reference
Contents
Index
UCIWWEIHR.uciwweihr_model_params1
UCIWWEIHR.uciwweihr_model_params2
UCIWWEIHR.uciwweihr_sim_params
UCIWWEIHR.ChainsCustomIndex
UCIWWEIHR.NegativeBinomial2
UCIWWEIHR.calculate_quantiles
UCIWWEIHR.calculate_quantiles_without_chain
UCIWWEIHR.create_uciwweihr_model_params2
UCIWWEIHR.create_uciwweihr_sim_params
UCIWWEIHR.eihr_ode!
UCIWWEIHR.generate_colors
UCIWWEIHR.generate_logit_normal_random_walk
UCIWWEIHR.generate_random_walk
UCIWWEIHR.generate_simulation_data_agent
UCIWWEIHR.generate_simulation_data_uciwweihr
UCIWWEIHR.is_time_varying_above_n
UCIWWEIHR.mcmcdiags_vis
UCIWWEIHR.non_time_varying_param_vis
UCIWWEIHR.ode_solution_vis
UCIWWEIHR.optimize_many_MAP
UCIWWEIHR.optimize_many_MAP2
UCIWWEIHR.optimize_many_MAP2_wrapper
UCIWWEIHR.power
UCIWWEIHR.predictive_param_vis
UCIWWEIHR.repeat_last_n_elements
UCIWWEIHR.save_plots_to_docs
UCIWWEIHR.startswith_any
UCIWWEIHR.subset_desired_ode_from_gq
UCIWWEIHR.uciwweihr_fit
UCIWWEIHR.uciwweihr_gq_pp
UCIWWEIHR.uciwweihr_init_param
UCIWWEIHR.uciwweihr_likelihood_helpers
UCIWWEIHR.uciwweihr_model
UCIWWEIHR.uciwweihr_visualizer
UCIWWEIHR.uciwweihr_model_params1
— Typeuciwweihr_model_params1
Struct for holding parameters used in the UCIWWEIHR ODE compartmental model. Use create_uciwweihr_model_params
to create an instance of this struct. With hard coded sigmawastewater and sigmahosp.
Fields
E_init_sd::Float64=50.0
: Standard deviation for the initial number of exposed individuals.log_E_init_mean::Int64=200
: Mean for the initial number of exposed individuals, on log scale.I_init_sd::Float64=20.0
: Standard deviation for the initial number of infected individuals.log_I_init_mean::Int64=100
: Mean for the initial number of infected individuals, on log scale.H_init_sd::Float64=5.0
: Standard deviation for the initial number of hospitalized individuals.log_H_init_mean::Int64=20
: Mean for the initial number of hospitalized individuals, on log scale.gamma_sd::Float64=0.02
: Standard deviation for the rate of incubation.log_gamma_mean::Float64=log(1/4)
: Mean for the rate of incubation on log scale.nu_sd::Float64=0.02
: Standard deviation for the rate of leaving the infected compartment.log_nu_mean::Float64=log(1/7)
: Mean for the rate of leaving the infected compartment on the log scale.epsilon_sd::Float64=0.02
: Standard deviation for the rate of hospitalization recovery.log_epsilon_mean::Float64=log(1/5)
: Mean for the rate of hospitalization recovery on the log scale.rho_gene_sd::Float64=0.02
: Standard deviation for the rho prior.log_rho_gene_mean::Float64=log(0.011)
: Mean for the row prior on log scale.sigma_ww::Float64=log(0.1)
: Standard deviation for the normal distribution for wastewater data. Not infered.sigma_hosp::Float64=500.0
: Standard deviation for the negative binomial distribution for hospital data. Not infered.Rt_init_sd::Float64=0.3
: Standard deviation for the initial value of the time-varying reproduction number.Rt_init_mean::Float64=0.2
: Mean for the initial value of the time-varying reproduction number.sigma_Rt_sd::Float64=0.2
: Standard deviation for normal prior of log time-varying reproduction number standard deviation.sigma_Rt_mean::Float64=-3.0
: Mean for normal prior of log time-varying reproduction number standard deviation.w_init_sd::Float64=0.1
: Standard deviation for the initial value of the time-varying hospitalization rate.w_init_mean::Float64=log(0.35)
: Mean for the initial value of the time-varying hospitalization rate.sigma_w_sd::Float64=0.2
: Standard deviation for normal prior of log time-varying hospitalization rate standard deviation.sigma_w_mean::Float64=-3.5
: Mean for normal prior of time-varying hospitalization rate standard deviation.
UCIWWEIHR.uciwweihr_model_params2
— Typeuciwweihr_model_params2
Struct for holding parameters used in the UCIWWEIHR ODE compartmental model. Use create_uciwweihr_model_params
to create an instance of this struct. With prior on sigmawastewater and sigmahosp.
Fields
E_init_sd::Float64=50.0
: Standard deviation for the initial number of exposed individuals.log_E_init_mean::Int64=200
: Mean for the initial number of exposed individuals, on log scale.I_init_sd::Float64=20.0
: Standard deviation for the initial number of infected individuals.log_I_init_mean::Int64=100
: Mean for the initial number of infected individuals, on log scale.H_init_sd::Float64=5.0
: Standard deviation for the initial number of hospitalized individuals.log_H_init_mean::Int64=20
: Mean for the initial number of hospitalized individuals, on log scale.gamma_sd::Float64=0.02
: Standard deviation for the rate of incubation.log_gamma_mean::Float64=log(1/4)
: Mean for the rate of incubation on log scale.nu_sd::Float64=0.02
: Standard deviation for the rate of leaving the infected compartment.log_nu_mean::Float64=log(1/7)
: Mean for the rate of leaving the infected compartment on the log scale.epsilon_sd::Float64=0.02
: Standard deviation for the rate of hospitalization recovery.log_epsilon_mean::Float64=log(1/5)
: Mean for the rate of hospitalization recovery on the log scale.rho_gene_sd::Float64=0.02
: Standard deviation for the rho prior.log_rho_gene_mean::Float64=log(0.011)
: Mean for the row prior on log scale.sigma_ww_sd::Float64=nothing
: Standard deviation for the normal prior of the log standard deviation of the wastewater data. Ifnothing
, the sigma_ww is used.log_sigma_ww_mean::Float64=nothing
: Mean for the normal prior of the log standard deviation of the wastewater data. Ifnothing
, the sigma_ww is used.sigma_hosp_sd::Float64=nothing
: Standard deviation for the normal prior of the log standard deviation of the hospital data. Ifnothing
, the sigma_hosp is used.sigma_hosp_mean::Float64=nothing
: Mean for the normal prior of the log standard deviation of the hospital data. Ifnothing
, the sigma_hosp is used.Rt_init_sd::Float64=0.3
: Standard deviation for the initial value of the time-varying reproduction number.Rt_init_mean::Float64=0.2
: Mean for the initial value of the time-varying reproduction number.sigma_Rt_sd::Float64=0.2
: Standard deviation for normal prior of log time-varying reproduction number standard deviation.sigma_Rt_mean::Float64=-3.0
: Mean for normal prior of log time-varying reproduction number standard deviation.w_init_sd::Float64=0.1
: Standard deviation for the initial value of the time-varying hospitalization rate.w_init_mean::Float64=log(0.35)
: Mean for the initial value of the time-varying hospitalization rate.sigma_w_sd::Float64=0.2
: Standard deviation for normal prior of log time-varying hospitalization rate standard deviation.sigma_w_mean::Float64=-3.5
: Mean for normal prior of time-varying hospitalization rate standard deviation.
UCIWWEIHR.uciwweihr_sim_params
— Typeuciwweihr_sim_params
Struct for holding parameters used in the UCIWWEIHR ODE compartmental model simulation. Use create_uciwweihr_sim_params
to create an instance of this struct.
Fields
time_points::Int64
: Number of time points for the simulation.seed::Int64
: Seed for random number generation.E_init::Int64
: Initial number of exposed individuals.I_init::Int64
: Initial number of infected individuals.H_init::Int64
: Initial number of hospitalized individuals.gamma::Float64
: Rate of incubation.nu::Float64
: Rate of leaving the infected compartment.epsilon::Float64
: Rate of hospitalization recovery.rho_gene::Float64
: Contribution of infected individual's pathogen genome into wastewater.sigma_ww::Float64
: standard deviation of the log concentration of pathogen genome in wastewater.sigma_hosp::Float64
: Standard deviation for the negative binomial distribution for hospital data.Rt::Union{Float64, Vector{Float64}}
: Initial value or time series of the time-varying reproduction number.sigma_Rt::Float64
: Standard deviation for random walk of time-varying reproduction number.w::Union{Float64, Vector{Float64}}
: Initial value or time series of the time-varying hospitalization rate.sigma_w::Float64
: Standard deviation for random walk of time-varying hospitalization rate.rt_init::Float64
: Initial value of the time-varying reproduction number, NOT USER SPECIFIEDcreate_uciwweihr_params
TAKES CARE OF THIS.w_init::Float64
: Initial value of the time-varying hospitalization rate, NOT USER SPECIFIEDcreate_uciwweihr_params
TAKES CARE OF THIS.
UCIWWEIHR.ChainsCustomIndex
— MethodChainsCustomIndex(c::Chains, indices_to_keep::BitMatrix)
Reduce Chains object to only wanted indices.
Function created by Damon Bayer.
UCIWWEIHR.NegativeBinomial2
— MethodCreate a re-parametrized negative binomial distribution in terms of mean and overdispersion.
Arguments
μ
: Mean of the distribution.ϕ
: Overdispersion parameter.
Returns
A Distributions.NegativeBinomial
distribution object.
UCIWWEIHR.calculate_quantiles
— Methodcalculate_quantiles(df, chain, var_prefix, quantiles)
Calculate quantiles for a given chain and variable prefix. Quantiles can be any user desired quantile.
UCIWWEIHR.calculate_quantiles_without_chain
— Methodcalculate_quantiles_without_chain(df, var_prefix, quantiles)
Calculates quantiles for a given variable prefix without considering the chain.
UCIWWEIHR.create_uciwweihr_model_params2
— Methodcreate_uciwweihr_model_params1(; kwargs...) or create_uciwweihr_model_params2(; kwargs...)
Creates a uciwweihr_sim_params1
or uciwweihr_sim_params2
struct with the option to either hard code sigmawastewater and sigmahosp or provide a prior.
Arguments
kwargs...
: Named arguments corresponding to the fields inuciwweihr_sim_params1
oruciwweihr_sim_params2
.
Returns
params::uciwweihr_sim_params1
orparams::uciwweihr_sim_params2
: A struct with simulation parameters.
UCIWWEIHR.create_uciwweihr_sim_params
— Methodcreateuciwweihrsim_params(; kwargs...)
Creates a uciwweihr_sim_params
struct with the option to either use a predetermined Rt
and w
or generate them as random walks.
Arguments
kwargs...
: Named arguments corresponding to the fields inuciwweihr_sim_params
.
Returns
params::uciwweihr_sim_params
: A struct with simulation parameters.
UCIWWEIHR.eihr_ode!
— Methodeihr_ode!(du, u, p, t)
Calculate the ordinary differential equations (ODEs) for the EIHR model.
Parameters:
du
: Array{Float64,1} - The derivative of the state variables.u
: Array{Float64,1} - The current state variables.p
: Tuple{Float64,Float64,Float64,Float64,Float64} - The model parameters (alpha, gamma, nu, w, epsilon).t
: Float64 - The current time.
Returns:
du
: Array{Float64,1} - The derivative of the state variables.
UCIWWEIHR.generate_colors
— Methodgenerate_ribbon_colors(number_of_colors)
Generates a vector with colors for ribbons in plots.
UCIWWEIHR.generate_logit_normal_random_walk
— Methodgenerate_logit_normal_random_walk(time_points::Int64, sigma::Float64, init_val::Float64)
Generates a logit-normal random walk time series.
Arguments
time_points::Int64
: Number of time points.sigma::Float64
: Standard deviation of the random walk in logit space.init_val::Float64
: Initial value of the random walk on the probability scale.
Returns
walk::Vector{Float64}
: Generated random walk on the probability scale.
UCIWWEIHR.generate_random_walk
— Methodgenerate_random_walk(time_points::Int64, sigma::Float64, init_val::Float64)
Generates a random walk time series.
Arguments
time_points::Int64
: Number of time points.sigma::Float64
: Standard deviation of the random walk.init_val::Float64
: Initial value of the random walk.
Returns
walk::Vector{Float64}
: Generated random walk.
UCIWWEIHR.generate_simulation_data_agent
— FunctionGenerating Simulation Data for Agent Based Model
To generate simulation data using the agent based model, you can use the generate_simulation_data_agent
function defined in the UCIWWEIHR.jl
package. This function allows you to customize various parameters for the simulation. NOT FINISHED, STILL NEEDS WW AND RT
Function Signature
Arguments
- seed::Int64: Seed for random number generation. Default value is 1.
- pop_size::Int64: Size of the population. Default value is 1000.
- I_init::Int64: Initial number of infected individuals. Default value is 200.
- H_init::Int64: Initial number of hospitalized individuals. Default value is 20.
- beta::Float64: Transmission rate. Default value is 0.001.
- gamma::Float64: Rate of exposed individuals becoming infectious. Default value is 1/4.
- nu::Float64: Rate of infected individuals recovering or getting hospitalized. Default value is 1/7.
- epsilon::Float64: Rate of hospitalized individuals recovering. Default value is 1/5.
- w_init::Float64: Probability of an infected individual becoming hospitalized. Default value is 0.35.
Returns
- df::DataFrame: A DataFrame containing the simulation data with columns
Time
,S
,E
,I
,H
, andR
.
UCIWWEIHR.generate_simulation_data_uciwweihr
— Methodgenerate_simulation_data(params::UCIWWEIHRParams)
Generates simulation data for the UCIWWEIHR ODE compartmental model.
Arguments
params::uciwweihr_sim_params
: Struct containing parameters for the simulation.
Returns
df::DataFrame
: A DataFrame containing the simulation data with columnsobstimes
,log_ww_conc
,hosp
,rt
, andwt
.
UCIWWEIHR.is_time_varying_above_n
— Methodis_time_varying_above_n(name, n)
Checks if the time varying parameter is above a given time point.
UCIWWEIHR.mcmcdiags_vis
— Methodmcmcdiags_vis(...)
Default visualizer for results of the UCIWWEIHR model, includes posterior/priors of generated quantities and posterior predictive samples for forecasting. Forecasting plots will have the observed data alongside.
Arguments
gq_samples
: Generated quantities samples from the posterior/prior distribution, index 2 in uciwweihrgqpp output.desired_params
: A list of lists of parameters to visualize. Each list will be visualized in a separate plot. Default is [["Einit", "Iinit", "Hinit"], ["gamma", "nu", "epsilon"], ["rhogene", "tau"], ["sigma_hosp"]].actual_non_time_varying_vals::uciwweihr_sim_params
: A uciwweihrsimparams object of actual non-time varying parameter values if user has access to them. Default is nothing.save_plots::Bool=false
: A boolean to indicate if user wants to save the plots as pngs into a plots folder.plot_name_to_save
: A string to indicate the name of the plot to save. Default is "mcmcdiagnosisplots".
UCIWWEIHR.non_time_varying_param_vis
— Methodnon_time_varying_param_vis(...)
Used in the uciwweihr_visualizer
to create visuals for non-time varying parameters.
Arguments
build_params::uciwweihr_model_params
: A struct of model parameters used to buildgq_samples
, used only if user desired priors next to posteriors.data_hosp
: Hospitalization data, used only if user desired priors next to posteriors.data_wastewater
: Wastewater data, if model does not use this do not specify this, if user desires priors next to plot (do not specify if you do not want prior plots).obstimes_hosp
: An array of time points for hospital data, used only if user desired priors next to posteriors.obstimes_wastewater
: An array of time points for wastewater data, used only if user desired priors next to posteriors.param_change_times
: An array of time points where the parameters change, used only if user desired priors next to posteriors.seed
: An integer to set the seed for reproducibility, used only if user desired priors next to posteriors.forecast
: A boolean to indicate if user wants to forecast, used only if user desired priors next to posteriors.forecast_weeks
: An integer to indicate the number of weeks to forecast, used only if user desired priors next to posteriors.gq_samples
: Generated quantities samples from the posterior/prior distribution, index 2 in uciwweihrgqpp output.desired_params
: A list of lists of parameters to visualize. Each list will be visualized in a separate plot. Default is any parameter not in this list : ["alphat", "wt", "rtvals", "loggenes_mean", "H"]actual_non_time_varying_vals::uciwweihr_sim_params
: A uciwweihrsimparams object of actual non-time varying parameter values if user has access to them. Default is nothing.save_plots::Bool=false
: A boolean to indicate if user wants to save the plots as pngs into a plots folder.plot_name_to_save
: A string to indicate the name of the plot to save. Default is "mcmcnontimevaryingparameterplots".
UCIWWEIHR.ode_solution_vis
— Methodode_solution_vis(...)
Visualizer for ODE solutions produced by the UCIWWEIHR model. If actual ode solutions are provided, they will be included in the plots.
Arguments
gq_samples
: Generated quantities samples from the posterior/prior distribution, index 2 in uciwweihrgqpp output.actual_non_time_varying_vals::uciwweihr_sim_params
: A uciwweihrsimparams object of actual non-time varying parameter values if user has access to them. Default is nothing.save_plots::Bool=false
: A boolean to indicate if user wants to save the plots as pngs into a plots folder.plot_name_to_save
: A string to indicate the name of the plot to save. Default is "mcmcdiagnosisplots".
UCIWWEIHR.optimize_many_MAP
— Functionoptimize_many_MAP(model, n_reps = 100, top_n = 1, verbose = true)
Try n_reps different initializations to get MAP estimate.
Function by Damon Bayer
UCIWWEIHR.optimize_many_MAP2
— Functionoptimize_many_MAP2(model, n_reps = 100, top_n = 1, verbose = true)
Try n_reps different initializations to get MAP estimate.
Modified by Christian Bernal Zelaya
UCIWWEIHR.optimize_many_MAP2_wrapper
— Methodoptimize_many_MAP2_wrapper(...)
Wrapper function for optimizemanyMAP2 that uses model based on inputs to function.
Created by Christian Bernal Zelaya
UCIWWEIHR.power
— Methodpower(a,b)
Raise `a` to the `b` power
UCIWWEIHR.predictive_param_vis
— Methodpredictive_param_vis(...)
Used in the uciwweihr_visualizer
to create visuals for wastewater data and hospitalization data.
Arguments
pp_samples
: A DataFrame of posterior or prior predictive samples.data_wastewater
: An array of actual wastewater values if user has access to them assumed, using time scale of observed time points. Default is nothing.data_hosp
: An array of actual hospitalization values if user has access to them assumed, , using time scale of observed time points. Default is nothing.forecast_weeks
: An integer of the number of weeks forecasted. Default is 0.vars_to_pred
: A list of variables to predict. Default is ["datawastewater", "datahosp"].quantiles
: A list of quantiles to calculate for ploting uncertainty. Default is [0.5, 0.8, 0.95].bayes_dist_type
: A string to indicate if user is using Posterior or Prior distribution. Default is "Posterior".save_plots::Bool=false
: A boolean to indicate if user wants to save the plots as pngs into a plots folder.plot_name_to_save
: A string to indicate the name of the plot to save. Default is "mcmcpredparameter_plots".
UCIWWEIHR.repeat_last_n_elements
— Methodrepeat_last_n_elements(x::Vector{T}, n::Int, w::Int) where T
Modifies a given array so that the last n elements are repeated w times.
UCIWWEIHR.save_plots_to_docs
— Methodsave_plots_to_docs(plot, filename; format = "png")
Saves plots to docs/plots directory.
UCIWWEIHR.startswith_any
— Methodstartswith_any(name, patterns)
Checks if the name of time varying paramter starts with any of the patterns.
UCIWWEIHR.subset_desired_ode_from_gq
— Methodsubset_desired_ode_from_gq(gq, desired_var::Regex)
Subsets the generated quantities dataframe from the uciwweihr_gq_pp
output to only include the desired variable. desired_var should be put in the following form : r"I(\[|\])"
for example.
UCIWWEIHR.uciwweihr_fit
— Methoduciwweihr_fit(...)
This is the sampler for the bayesian semi-parametric model for the wastewater EIHR compartmental model. The defaults for this fuction will follow those of the default simulation in generatesimulationdatawweihr.jl function.
Arguments
data_hosp
: An array of hospital data.data_wastewater
: An array of pathogen genome concentration in localized wastewater data. If this is not avaliable, the model used will be one that only uses hospital data.obstimes
: An array of timepoints for observed hosp/wastewater.priors_only::Bool=false
: A boolean to indicate if only priors are to be sampled.n_samples::Int64=500
: Number of samples to be drawn.n_chains::Int64=1
: Number of chains to be run.n_discard_initial::Int64=0
: Number of samples to be discarded.seed::Int64=2024
: Seed for the random number generator.params::uciwweihr_model_params#
: A struct containing parameters for the model. # is either 1 or 2 for prior or hardcoded sigmaww and sigmahosp.init_params
: Initial parameters for the model. If supplied,uciwweihr_init_param
can supply these values.return_bool
: A boolean to indicate if the model is to use the return statement. Only set to false if only forecast is desired
Returns
- Samples from the posterior or prior distribution.
UCIWWEIHR.uciwweihr_gq_pp
— Methoduciwweihr_gq_pp(...)
This generates quantities and a posterior predictive distribution for the bayesian semi-parametric model for the wastewater EIHR compartmental model. The defaults for this fuction will follow those of the default simulation in generatesimulationdatawweihr.jl function.
Arguments
samples
: Samples from the posterior/prior distribution.data_hosp
: An array of hospital data.data_wastewater
: An array of pathogen genome concentration in localized wastewater data. If this is not avaliable, the model used will be one that only uses hospital data.obstimes
: An array of timepoints for observed hosp/wastewater.param_change_times
: An array of timepoints where the parameters change.seed::Int64=2024
: Seed for the random number generator.params::uciwweihr_model_params
: A struct containing parameters for the model.forecast::Bool=false
: A boolean to indicate if forecasting is to be done.forecast_weeks::Int64=4
: Number of weeks to forecast.return_bool::Bool=true
: A boolean to indicate if the model is to use the return statement. Only set to false if only forecast is desiredgq_bool::Bool=true
: A boolean to indicate if the model is to generate quantities. Only set to false if only forecast is desired
Returns
- Samples from the posterior or prior distribution.
UCIWWEIHR.uciwweihr_init_param
— Methoduciwweihr_init_param(...)
Gets initial parameters values for the UCIWWEIHR model. Only need to run once.
Arguments
data_hosp
: An array of hospital data.data_wastewater
: An array of pathogen genome concentration in localized wastewater data. If this is not avaliable, the model used will be one that only uses hospital data.obstimes_hosp
: An array of timepoints for observed hosp.obstimes_wastewater
: An array of timepoints for observed wastewater.param_change_times
: An array of timepoints where the parameters change.n_chains::Int64=1
: Number of chains to be run.seed::Int64=2024
: Seed for the random number generator.params::uciwweihr_model_params
: A struct containing parameters for the model.verbose_optimize::Bool=false
: A boolean to indicate if the optimization process is to be verbose.
Returns
- Samples from the posterior or prior distribution.
UCIWWEIHR.uciwweihr_likelihood_helpers
— Methoduciwweihr_likelihood_helpers(...)
asdf
Arguments
- asdf
UCIWWEIHR.uciwweihr_model
— Methoduciwweihr_model(...)
This is the bayesian semi-parametric model for the wastewater EIHR compartmental model. The defaults for this fuction will follow those of the default simulation in generatesimulationdatawweihr.jl function.
Arguments
data_hosp
: An array of hospital data.data_wastewater
: An array of pathogen genome concentration in localized wastewater data. If this is not avaliable, the model used will be one that only uses hospital data.obstimes_hosp
: An array of timepoints for observed hospital data.obstimes_wastewater
: An array of timepoints for observed wastewater data.param_change_times
: An array of timepoints where the parameters change.params::uciwweihr_model_params
: A struct containing parameters for the model.return_bool
: A boolean to indicate if the model is to use the return statement. Only set to false if only forecast is desired
UCIWWEIHR.uciwweihr_visualizer
— Methoduciwweihr_visualizer(...)
Default visualizer for results of the UCIWWEIHR model, includes posterior/priors of generated quantities and posterior predictive samples for forecasting. Forecasting plots will have the observed data alongside.
Arguments
build_params::uciwweihr_model_params
: A uciwweihrmodelparams object, if user desires priors next to plot (do not specify if you do not want prior plots).data_hosp
: An array of hospital data used for model fitting, if user desires priors next to plot (do not specify if you do not want prior plots).data_wastewater
: An array of wastewater data used for model fitting, if model does not use this do not specify this, if user desires priors next to plot (do not specify if you do not want prior plots).obstimes
: An array of timepoints for observed hosp/wastewater, if user desires priors next to plot (do not specify if you do not want prior plots).param_change_times
: An array of timepoints where the parameters change, if user desires priors next to plot (do not specify if you do not want prior plots).seed
: Seed for the random number generator, if user desires priors next to plot (do not specify if you do not want prior plots).forecast
: A boolean to indicate if forecasting is to be done, if user desires priors next to plot (do not specify if you do not want prior plots).forecast_weeks
: Number of weeks to forecast, if user desires priors next to plot (do not specify if you do not want prior plots).pp_sampeles
: Posterior predictive samples from the posterior/prior distribution, index 1 in uciwweihrgqpp output.gq_samples
: Generated quantities samples from the posterior/prior distribution, index 2 in uciwweihrgqpp output.obs_data_hosp
: An array of hospital data, data used for model fitting or extened timeseries for evaluation of forecast.obs_data_wastewater
: An array of wastewater data, data used for model fitting or extened timeseries for evaluation of forecast.actual_rt_vals
: An array of actual Rt values if user has access to them assumed to be on a daily scale. This typically will come from some simulation. Default is nothing.actual_w_t
: An array of actual w_t values if user has access to them assumed to be on a daily scale. This typically will come from some simulation. Default is nothing.actual_E_ode_sol
: An array of actual E ODE values if user has access to them assumed to be on a daily scale. This typically will come from some simulation. Default is nothing.actual_I_ode_sol
: An array of actual I ODE values if user has access to them assumed to be on a daily scale. This typically will come from some simulation. Default is nothing.actual_H_ode_sol
: An array of actual H ODE values if user has access to them assumed to be on a daily scale. This typically will come from some simulation. Default is nothing.actual_non_time_varying_vals::uciwweihr_sim_params
: A uciwweihrsimparams object of actual non-time varying parameter values if user has access to them. Default is nothing.desired_params
: A list of lists of parameters to visualize. Each list will be visualized in a separate plot. Default is [["Einit", "Iinit", "Hinit"], ["gamma", "nu", "epsilon"], ["rhogene", "tau", "df"], ["sigma_hosp"]].time_varying_params
: A list of time varying parameters to visualize. Default is ["rtvals", "wt"].var_to_pred
: A list of variables to predict. Default is ["datawastewater", "datahosp"].quantiles
: A list of quantiles to calculate for ploting uncertainty. Default is [0.5, 0.8, 0.95].bayes_dist_type
: A string to indicate if user is using Posterior or Prior distribution ("Posterior" / "Prior").mcmcdaigs::Bool=true
: A boolean to indicate if user wants to visualize mcmc diagnosis plots and Effective Sample Size(ESS).time_varying_plots::Bool=true
: A boolean to indicate if user wants to visualize time varying parameters.non_time_varying_plots::Bool=true
: A boolean to indicate if user wants to visualize non-time varying parameters.ode_sol_plots::Bool=true
: A boolean to indicate if user wants to visualize ODE solutions.pred_param_plots::Bool=true
: A boolean to indicate if user wants to visualize posterior (or prior) predictive parameter results.save_plots::Bool=false
: A boolean to indicate if user wants to save the plots as pngs into a plots folder.plot_name_to_save_mcmcdiag
: A string to indicate the name of the plot to save for MCMC diagnostics. Default is "mcmcdiagnosisplots".plot_name_to_save_time_varying
: A string to indicate the name of the plot to save for time varying parameters. Default is "mcmctimevaryingparameterplots".plot_name_to_save_non_time_varying
: A string to indicate the name of the plot to save for non-time varying parameters. Default is "mcmcnontimevaryingparameterplots".plot_name_to_save_ode_sol
: A string to indicate the name of the plot to save for ODE solutions. Default is "mcmcodesolution_plots".plot_name_to_save_pred_param
: A string to indicate the name of the plot to save for posterior (or prior) predictive parameter results. Default is "mcmcpredparameter_plots".