Beta-binomial distributions using BRMS in R - How to use the custom_family argument

Interfaces brms
1

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://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
}

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
)

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’ve attached an example of the data I use. 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.
Example_data.csv (4.9 KB)

Best regards,

2

Hi Ida,

You’re missing the call to stanvars:

added_fun <- stanvar(scode = stan_beta_binomial2, block = "functions")

fit1ab = brm(Disease ~ (Ants*Treatment*Plant_sort*Pests) + (1|Block/Plant_ID),
             data=mydata, family = beta_binomial2, stanvars=added_fun)
3

Unfortunately not, I get a similar error when including/excluding stanvars. The error I get when including the stanvars call is:

“SYNTAX ERROR, MESSAGE(S) FROM PARSER:
No matches for:
beta_binomial2_lpmf(int, int, real, real)
Function beta_binomial2_lpmf not found.
error in ‘model1d78509b7147_file1d78cfe3b16’ at line 5, column 63
3:
4: real beta_binomial2_lpmf(int y, real mu, real phi, int T) {
5: return beta_binomial2_lpmf(y | T, mu * phi, (1 - mu) * phi);
^
6: }
Error in stanc(model_code = paste(program, collapse = “\n”), model_name = model_cppname, :
failed to parse Stan model ‘file1d78cfe3b16’ due to the above error.
SYNTAX ERROR, MESSAGE(S) FROM PARSER:
No matches for:
beta_binomial2_lpmf(int, int, real, real)
Function beta_binomial2_lpmf not found.
error in ‘model1d786c397d18_file1d78cfe3b16’ at line 5, column 63
3:
4: real beta_binomial2_lpmf(int y, real mu, real phi, int T) {
5: return beta_binomial2_lpmf(y | T, mu * phi, (1 - mu) * phi);
^
6: }
Error in stanc(model_code = paste(program, collapse = “\n”), model_name = model_cppname, :
failed to parse Stan model ‘file1d78cfe3b16’ due to the above error.”

4

That error is because of how you’ve specified the distribution:

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);
}

In the vignette the model specification is:

real beta_binomial2_lpmf(int y, real mu, real phi, int T) {
return beta_binomial_lpmf(y | T, mu * phi, (1 - mu) * phi);
}
int beta_binomial2_rng(real mu, real phi, int T) {
return beta_binomial_rng(T, mu * phi, (1 - mu) * phi);
}

In other words, the beta_binomial2_lpmf is supposed to call the beta_binomial_lpmf function, not the beta_binomial2_lpmf

1 Like
5

Note also that your brms formula is not consistent with the vignette either, as you need to use the vint function to specify the total number of trials:

fit2 <- brm(
  incidence | vint(size) ~ period + (1|herd), data = cbpp, 
  family = beta_binomial2, stanvars = stanvars
)
2 Likes
6

Thank you! I’ll look into how I best implement the vint() argument. Thank you very much for taking your time to help me, I’ll send plenty of good karma in your direction! :-)

2 Likes