I am trying to get WAIC from models run in JagsUI. I understand that I can use the loo package to extract WAIC using the log likelihood (like the code chunk below), but this can only be done with a stan model object of class stanfit.
LL_2 <- extract_log_lik(output, parameter_name = "log_lik", merge_chains = TRUE)
loo_2 <- loo(LL_2, save_psis = TRUE, cores = 4)
waic_2 <- waic(LL_2, save_psis = TRUE, cores = 4)
I am inexperienced in stan and it is not feasible for me to switch my models over to use stan; however, I am wondering if there is a a method for me to create an “empty” stanfit model object and then add the mcmc.list from my jagsUI model to this stanfit object. Then I imagine I might be able to use this object to extract the log likelihood and resulting WAIC?
Is there any way to create an empty stanfit object, or directly use the mcmc.list object from my jagsUI model to extract WAIC using the loo package?
Thank you for your help!
The loo function is flexible in that it takes a log-likelihood array, matrix, or function as input, so it should be generally applicable to Bayesian model fits, not just Stan. So if you’re able to extract the log-likelihood directly from your JagsUI object, that might be easier than trying to transform things into a Stan object.
The loo function documentation has details on what exactly the log likelihood should look like.
3 Likes
Thank you for your explanation and help!
I think the easiest way is to use CRAN: Package posterior. If you have the log_lik values in a mcmc.list then the following works
loo(as_draws_array(your_mcmc_list["log_lik"]))
I recommend learning to use posterior package, as it makes so many other posterior draw manipulations easier and has better convergence diagnostics than, e.g. coda package.
I also recommend to use LOO-CV instead of WAIC. See CV-FAQ for explanation.
Thank you for your suggestion – my main issue now is that I actually am unable to extract the log likelihood from my jagsUI model easily. I am running an integrated survival model which has a cormack-jolly-serber component and a known-fate component and I am struggling to determine how to extract the log-likelihoods, as they are not calculated automatically in JagsUI. Once I do this, I believe the package you suggested would work.
Just a note that most commonly for these sorts of models, the log-likelihoods of interest for most downstream applications will be the log-likelihoods conditional only on the parameters in the marginalized form of the model, and not conditional on any latent discrete parameters that may be present in your implementation. We discuss this a bit more in the context of occupancy models here:
See here https://www.biorxiv.org/content/10.1101/2023.10.26.564080v1.full and here flocker/R/log_lik_occupancy.R at main · jsocolar/flocker · GitHub for implementations in the context of various occupancy models.
2 Likes
I’m afraid that unless you can extract the log likelihoods, you’re not going to be able to evaluate things like WAIC.
I would suggest asking the developers of jagsUI how to do that.
It is well worth moving these CJS models over to Stan. HMC is way faster than Gibbs for these models and much more robust. I’d suggest checking out this paper (from fisheries and wildlife ecologists who work on ADMB/TMB, not the Stan team) for a comparison:
Faster estimation of Bayesian models in ecology using Hamiltonian Monte Carlo. Monnahan, Thorson, Branch. 2016. Methods in Ecology and Evolution.
I do a simple example of CJS in the User’s Guide and there’s a very early case study where I fit the Dorazio/Royle occupancy model.