Skip to content

Instantly share code, notes, and snippets.

@padpadpadpad
Last active September 5, 2023 13:44
Show Gist options
  • Save padpadpadpad/77aa9c661653aed421bd497343cf6767 to your computer and use it in GitHub Desktop.
Save padpadpadpad/77aa9c661653aed421bd497343cf6767 to your computer and use it in GitHub Desktop.
Very simple example of using nimble, ggdist, and tidybayes.
# nimble and tidybayes example
# load packages
librarian::shelf(tidybayes, nimble, tidyverse, bayesplot, ggdist, nlist)
# run example from runMCMC
code <- nimbleCode({
mu ~ dnorm(0, sd = 1000)
sigma ~ dunif(0, 1000)
for(i in 1:10) {
x[i] ~ dnorm(mu, sd = sigma)
}
})
Rmodel <- nimbleModel(code)
Rmodel$setData(list(x = c(2, 5, 3, 4, 1, 0, 1, 3, 5, 3)))
Rmcmc <- buildMCMC(Rmodel)
Cmodel <- compileNimble(Rmodel)
Cmcmc <- compileNimble(Rmcmc, project = Rmodel)
inits <- function() list(mu = rnorm(1,0,1), sigma = runif(1,0,10))
samplesList <- runMCMC(Cmcmc, niter = 10000, nchains = 3, inits = inits, samplesAsCodaMCMC = TRUE)
samplesList
# check model parameters
nlist::pars(samplesList)
# check caterpillars
mcmc_trace(samplesList)
# grab draws - uses tidybayes
d <- gather_draws(samplesList, mu, sigma)
head(d)
# make plot of distribution of parameter
ggplot(d, aes(.value, group = .chain)) +
geom_density() +
facet_wrap(~.variable, scales = 'free')
# use ggdist to make nice plot
ggplot(d, aes(y = .variable, x = .value)) +
stat_halfeye()
ggplot(d, aes(y = .variable, x = .value)) +
stat_pointinterval()
# get summary of data
d_summary <- group_by(d, .variable) %>%
median_qi()
d_summary
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment