--- title: "Running a subnational level multi-country model" output: rmarkdown::html_vignette: toc: true vignette: > %\VignetteEngine{knitr::rmarkdown} %\VignetteIndexEntry{Running a subnational level multi-country model} %\usepackage[UTF-8]{inputenc} --- ```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE) ``` ## Load your library ```{r, include=TRUE, message=FALSE, eval=FALSE} library(mcmsupply) set.seed(1209) ``` ## Get your data ```{r, include=TRUE, message=FALSE, eval=FALSE} cleaned_subnatdata <- get_data(national=FALSE) ``` ## Get the JAGS model inputs from the cleaned data ```{r, include=TRUE, message=FALSE, eval=FALSE} pkg_data <- get_modelinputs(startyear=1990, endyear=2025.5, nsegments=12, raw_data = cleaned_subnatdata) ``` ## Run JAGS model and get posterior point estimates with uncertainty. ```{r message=FALSE, eval=FALSE} mod <- run_jags_model(jagsdata = pkg_data, jagsparams = NULL, n_iter = 80000, n_burnin = 10000, n_thin = 35) ``` # Check the model diagnostics We use this to evaluate the convergence of the model parameters. We should expect to see R-hat values of approximately 1.05. The plot function will give you a visual summary for each parameter monitored. ```{r, include=TRUE, message=FALSE, eval=FALSE} plot(mod$JAGS) print(mod$JAGS) ``` Using the ggplot2 and tidybayes R packages, we will check the trace plots to assess the convergence of individual parameters. We expect to see a 'caterpillar' like appearance of the chains over the iterations. ```{r, include=TRUE, message=FALSE, eval=FALSE} sample_draws <- tidybayes::tidy_draws(mod$JAGS$BUGSoutput$sims.matrix) var <- sample_draws %>% dplyr::select(.chain, .iteration, .draw,`P[1,2,1,1]`) %>% dplyr::mutate(chain = rep(1:2, each=mod$JAGS$BUGSoutput$n.keep)) %>% dplyr::mutate(iteration = rep(1:mod$JAGS$BUGSoutput$n.keep, 2)) ggplot2::ggplot(data=var) + ggplot2::geom_line(ggplot2::aes(x=iteration, y=`P[1,2,1,1]`, color=as.factor(chain))) ``` ## Plot posterior point estimates with uncertainty ```{r, include=TRUE, message=FALSE, eval=FALSE} plots <- plot_estimates(jagsdata = pkg_data, model_output = mod) ``` ## Pull out estimates that you are particularly interested in ```{r, include=TRUE, message=FALSE, eval=FALSE} estimates_2018 <- pull_estimates(model_output = mod, country = 'Nepal', year=2018) head(estimates_2018) ```