Skip to content

Instantly share code, notes, and snippets.

@JimGrange
Created November 8, 2021 15:32
Show Gist options
  • Save JimGrange/3aaf870d388128f3b6e492ed2fcf4ac4 to your computer and use it in GitHub Desktop.
Save JimGrange/3aaf870d388128f3b6e492ed2fcf4ac4 to your computer and use it in GitHub Desktop.
mod_1 <- brm(value ~ questionnaire * condition,
data = idealised_data,
seed = 42,
cores = 4)
mod_2 <- brm(value ~ questionnaire + condition,
data = idealised_data,
seed = 42,
cores = 4)
mod_3 <- brm(value ~ questionnaire,
data = idealised_data,
seed = 42,
cores = 4)
mod_4 <- brm(value ~ condition,
data = idealised_data,
seed = 42,
cores = 4)
all_waic <- waic(mod_1,
mod_2,
mod_3,
mod_4)
weight <- function(d_waic, all_d_waic){
num <- exp(-0.5 * d_waic)
den <- sum(exp(-0.5 * all_d_waic[1]),
exp(-0.5 * all_d_waic[2]),
exp(-0.5 * all_d_waic[3]),
exp(-0.5 * all_d_waic[4]))
return(num / den)
}
all_waic <- c(749.7, 755.4, 754.5, 753.6)
d_waic <- all_waic - all_waic[which.min(all_waic)]
# model 1 is superior
weight(d_waic[1], d_waic)
weight(d_waic[2], d_waic)
weight(d_waic[3], d_waic)
weight(d_waic[4], d_waic)
# compare model 1 to model 4 in terms of Kullback–Leibler discrepancy
weight(d_waic[1], d_waic) / weight(d_waic[4], d_waic)
# evidence ratio as the normalized probability that Model 1 is to be
# preferred over Model 4
weight(d_waic[1], d_waic) / (weight(d_waic[1], d_waic) +
weight(d_waic[4], d_waic))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment