In this third component, we calibrate unknown model parameters by matching model outputs to specified calibration targets. Specifically, we calibrate the Sick-Sicker model to match survival, prevalence and the proportion who are Sicker, among all those afflicted (Sick+Sicker). We used a Bayesian calibration approach using the incremental mixture importance sampling (IMIS) algorithm (Steele, Raftery, and Emond 2006), which has been used to calibrate health policy models (Raftery and Bao 2010, Menzies et al. (2017), Rutter et al. (2018)). Bayesian methods allow us to quantify the uncertainty in the calibrated parameters even in the presence of non-identifiability (Alarid-Escudero et al. 2018). This analysis is coded in the 03_calibration.R file in the analysis
folder. The target data is stored in the 03_calibration_targets.RData file. Similar to component 02 ??, in the section 03.1 Load packages, we start by loading inputs and functions. In addition, we load the calibration targets data into the R workspace. In the next section, 03.2 Visualize targets, we plot each of the calibration targets with their confidence intervals.
In section 03.3 Run calibration algorithms, we set the parameters we need to calibrate to fixed values and test if the function calibration_out
that produces model outputs corresponding to the calibration targets works. This function takes a vector of parameters that need to be calibrated and a list with all parameters of decision model and computes model outputs to be used for calibration routines.
print.function(calibration_out) # print the functions
#> function(v_params_calib, l_params_all){ # User defined
#> # Substitute values of calibrated parameters in base-case with
#> # calibrated values
#> l_params_all <- update_param_list(l_params_all = l_params_all, params_updated = v_params_calib)
#>
#> # Run model with updated calibrated parameters
#> l_out_stm <- decision_model(l_params_all = l_params_all)
#>
#> ####### Epidemiological Output ###########################################
#> #### Overall Survival (OS) ####
#> v_os <- 1 - l_out_stm$m_M[, "D"]
#>
#> #### Disease prevalence #####
#> v_prev <- rowSums(l_out_stm$m_M[, c("S1", "S2")])/v_os
#>
#> #### Proportion of sick in S1 state #####
#> v_prop_S2 <- l_out_stm$m_M[, "S2"] / rowSums(l_out_stm$m_M[, c("S1", "S2")])
#>
#> ####### Return Output ###########################################
#> l_out <- list(Surv = v_os[c(11, 21, 31)],
#> Prev = v_prev[c(11, 21, 31)],
#> PropSicker = v_prop_S2[c(11, 21, 31)])
#> return(l_out)
#> }
#> <bytecode: 0x7fc4df626ef8>
#> <environment: namespace:darthpack>
This function is informed by two argument v_params_calib
and l_params_all
. The vector v_params_calib
contains the values of the three parameters of interest. The list l_params_all
contains all parameters of the decision model. The placeholder values are replaced by v_params_calib
and with these values the model is evaluated. Model evaluation takes place by running the decision_model
function, described in component 02. The result in a new list with output of the model corresponding to the parameter values in the v_params_calib
. With this new decision model output, the overall survival, disease prevalence and the proportion of Sicker in the Sick and Sicker states are calculated. The estimated values for these epidemiological outcomes at different timepoints are combined in a list called l_out
produced but the calibration_out
.
Once we make sure this code works, we specify the calibration parameters in section 03.3.1 Specify calibration parameters. These include setting the seed for the random number generation, specifying the number of random samples to obtain from the calibrated posterior distribution, the name of the input parameters and the range of these parameters that will inform the prior distributions of the calibrated parameters, and the name of the calibration targets: Surv
, Prev
, PropSick
.
In the next section, 03.3.2 Run IMIS algorithm, we calibrate the Sick-Sicker model with the IMIS algorithm. For this case-study, we assume a normal likelihood and uniform priors. For a more detailed description of IMIS for Bayesian calibration, different likelihood functions and prior distributions, we refer the reader to the tutorial for Bayesian calibration by Menzies et al. (Menzies et al. 2017). We use the IMIS
function from the IMIS
package that calls the functions likelihood
, sample.prior
and prior
, to draw samples from the posterior distribution (A. Raftery and Le Bao 2012). The functions are specified in the 03_calibration_functions.R file in the R
folder. For the IMIS
function, we specify the incremental sample size at each iteration of IMIS, the desired posterior sample size at the resample stage, the maximum number of iterations in IMIS and the number of optimizers which could be 0. The function returns a list, which we call l_fit_imis
, with the posterior samples, the diagnostic statistics at each IMIS iteration and the centers of Gaussian components (A. Raftery and Le Bao 2012). We store the posterior samples in the matrix m_calib_post
.
We then explore these posterior distributions in section 03.4 Exploring posterior distribution. We start by estimating the posterior mean, median and 95% credible interval, the mode and the maximum-a-posteriori (MAP). All for these summary statistics are combined in a dataframe called df_posterior_summ
. Table 1.1 shows the summary statistics of the posterior distribution.
Mean | 2.5% | 50% | 97.5% | MAP | |
---|---|---|---|---|---|
p_S1S2 | 0.1076682 | 0.0956305 | 0.107546 | 0.1201152 | 0.1073784 |
hr_S1 | 2.7150701 | 1.1119686 | 2.724407 | 4.3357750 | 2.4223530 |
hr_S2 | 9.6366087 | 7.2323662 | 9.660728 | 12.1053631 | 9.7262516 |
In section 03.4.2 Visualization of posterior distribution, we generate a pairwise scatter plot of the calibrated parameters (Figure 1.2) and a 3D scatter plot of the joint posterior distribution (Figure 1.1). These figures are saved in the figs directory.
Finally, the posterior distribution and MAP estimate from the IMIS calibration are stored in the file 03_imis_output.RData. Storing this data as an .RData file allows to import the data in following sections without needing to re-run the calibration component.
Alarid-Escudero, F, RF MacLehose, Y Peralta, KM Kuntz, and Enns EA. 2018. “Nonidentifiability in Model Calibration and Implications for Medical Decision Making.” Medical Decision Making 38 (7): 810–21. doi:10.1177/0272989X18792283.
Menzies, Nicolas A., Djøra I. Soeteman, Ankur Pandya, and Jane J. Kim. 2017. “Bayesian Methods for Calibrating Health Policy Models: A Tutorial.” PharmacoEconomics 35 (6). Springer International Publishing: 613–24. doi:10.1007/s40273-017-0494-4.
Raftery, A, and L Bao. 2010. “Estimating and Projecting Trends in HIV/AIDS Generalized Epidemics Using Incremental Mixture Importance Sampling.” Biometrics 66 (4): 1162–73.
Raftery, Adrian, and Le Bao. 2012. IMIS: Increamental Mixture Importance Sampling. https://CRAN.R-project.org/package=IMIS.
Rutter, C, J Ozik, M DeYoreo, and Collier N. 2018. “Microsimulation Model Calibration using Incremental Mixture Approximate Bayesian Computation.” arXiv, no. april: 1–20. https://arxiv.org/abs/1804.02090v3.
Steele, R, A Raftery, and M Emond. 2006. “Computing Normalizing Constants for Finite Mixture Models via Incremental Mixture Importance Sampling (IMIS).” Journal ofComputational and Graphical Statistics 15 (3): 712–34. doi:10.1198/106186006X132358.