I’m trying to compute a beta-binomial model using the R package brms. However, beta-binomial distributions are not native to this package, so I’m trying to use the custom_family argument, to create the family myself. I’ve found a vignette for the purpose: https://github.com/paul-buerkner/custom-brms-families/commit/32c9f2bfc92abade64c91fe601a3e6370070d30b AND https://cran.r-project.org/web/packages/brms/vignettes/brms_customfamilies.html
However, after trying the code in various ways, I get the following syntax error:
"SYNTAX ERROR, MESSAGE(S) FROM PARSER:
*No matches for:
beta_binomial2_lpmf(int, real, real, int)
Function beta_binomial2_lpmf not found.
error in ‘model445074614f85_file44502d6f3’ at line 71, column 65
-------------------------------------------------
69: if (!prior_only) {
70: for (n in 1:N) {
71: target += beta_binomial2_lpmf(Y[n] | mu[n], phi, trials[n]);
^
72: }
-------------------------------------------------
Error in stanc(model_code = paste(program, collapse = “\n”), model_name = model_cppname, :
failed to parse Stan model ‘file44502d6f3’ due to the above error.
SYNTAX ERROR, MESSAGE(S) FROM PARSER:
No matches for:
beta_binomial2_lpmf(int, real, real, int)
Function beta_binomial2_lpmf not found.
error in ‘model44506c6352ff_file44502d6f3’ at line 71, column 65
-------------------------------------------------
69: if (!prior_only) {
70: for (n in 1:N) {
71: target += beta_binomial2_lpmf(Y[n] | mu[n], phi, trials[n]);
^
72: }
-------------------------------------------------
Error in stanc(model_code = paste(program, collapse = “\n”), model_name = model_cppname, :
failed to parse Stan model ‘file44502d6f3’ due to the above error."
The code I’m running when I get this error is as follows:
"log_lik_beta_binomial2 <- function(i, draws) {
mu <- draws$dpars$mu[, i]
phi <- draws$dpars$phi
N <- draws$data$trials[i]
y <- draws$data$Y[i]
beta_binomial2_lpmf(y, mu, phi, N)
}
predict_beta_binomial2 <- function(i, draws, …) {
mu <- draws$dpars$mu[, i]
phi <- draws$dpars$phi
N <- draws$data$trials[i]
beta_binomial2_rng(mu, phi, N)
}
fitted_beta_binomial2 <- function(draws) {
mu <- draws$dpars$mu
trials <- draws$data$trials
trials <- matrix(trials, nrow = nrow(mu), ncol = ncol(mu), byrow = TRUE)
mu * trials
}
# definition of the custom family
beta_binomial2 <- custom_family(
“beta_binomial2”,
dpars = c(“mu”, “phi”),
links = c(“logit”, “log”),
lb = c(0, 0),
ub = c(1, NA),
type = “int”, vars = “trials[n]”,
log_lik = log_lik_beta_binomial2,
predict = predict_beta_binomial2,
fitted = fitted_beta_binomial2
)
# additionally required Stan code
stan_beta_binomial2 <- "
real beta_binomial2_lpmf(int y, real mu, real phi, int T) {
return beta_binomial2_lpmf(y | T, mu * phi, (1 - mu) * phi);
}
int beta_binomial2_rng(real mu, real phi, int T) {
return beta_binomial2_rng(T, mu * phi, (1 - mu) * phi);
}
"
fit1ab=brm(Disease~Ants+Treatment+Plant_sort+Pests+(AntsTreatmentPlant_sortPests)+(1|Block/Plant_ID),data=mydata,family = beta_binomial2)"*
I’m new to brms and not very experienced in problem solving in R, so if anyone can help as to what I should do differently it would be greatly appreciated.
Best regards,