Brms: multivariate response, different distributions


I am interested in modelling survival (0/1) and reproduction (# offspring 0…3) as a multivariate response, combining Bernoulli with (truncated) Poisson distributions. Is it possible to do this with brms?


A model could look a follows:

bform <- bf(survival ~ ..., family = bernoulli()) +
  bf(offspring | trunc(<value>) ~ ..., family = poisson())

fit <- brm(bform, ...)

However, without correlated multilevel structure on the right-hand side of the model formulas, this model will just be equivalent to two univariate ones.


Thanks Paul. Is it possible to add a correlated multilevel structure on the rhs? Specifically, I want to model correlated year-effects to both survival and offspring. So, if I would add +(1|p|year) to both formulas, would a bivariate normal be estimated?




Great, thanks! Is it also possible to model the correlated year-effects as multivariate t-distributed or multivariate skew-normal rather than “ordinary” multivariate normal?


Student-t is possible albeit not documented and not yet officially supported until we better understand its implications. Go for (1 |p| gr(year, dist = "student")).


Thanks again. I will give the multivariate t a try because I have the impression that the year-effect empirical distribution has “fat tails”.

One more thing: when using the predict function, is it possible to marginalize/integrate over the “random” year effects, in order to obtain expected survival and offspring number?

BTW: thanks for creating bmrs. It looks like an enormously useful package that will save me a lot of Stan coding.


You can ignore random effects via re_formula = NA but marginalizing is a much more difficult thing (read: not yet possible).