Using loo for brms non-linear model incorporating Stan functions

Interfaces brms
1

I’m working with some models of binary outcome data where I want the probability parameter of a Bernoulli distribution to depend on several Gamma survival functions assessed at a specified threshold value. To achieve this, I’m using brms’ non-linear functionality, so that the model formula looks something like this (simplified):

bf(response ~ 1-gamma_cdf(exp(threshold), exp(alpha), exp(beta))),
alpha ~ condition,
beta ~ condition,
threshold ~ 1,
nl = TRUE,
family = bernoulli(link="identity"))

The gamma_cdf function comes from Stan, since as far as I’m aware Stan has no pgamma function that corresponds to the one in R.

The models are fitting the data more or less effectively, and there are clear differences between the data I’ve (manually) simulated from the posteriors of most of them that allow me to rule out some models via eye test. However, there are a couple of models that look to be fitting the data roughly equally well, so I’d like to use loo to get a quantitative estimate of their relative performance—if only to confirm my judgement that the performance is fairly similar. Unfortunately when I try to do this (and also when I try to use posterior_predict and conditional_effects) I get an error message telling me:

Error in gamma_cdf(exp(threshold), exp(alpha), exp(beta)) : 
  could not find function "gamma_cdf"
 Most likely this is because you used a Stan function in the non-linear model formula that is not defined in R. If this is a user-defined function, please run 'expose_functions(., vectorize = TRUE)' on your fitted model and try again.

I’m just wondering if anyone who’s more experienced with this sort of thing might have a suggestion regarding how I can circumvent this problem so I can perform the comparison I’m interested in. If so, I’d love to hear it!

1 Like
2

To answer my own question nearly a year later, the solution seems as simple as creating the relevant function in R. So, to be able to use loo with my model containing gamma_cdf, I just had to write an R function called gamma_cdf (which, in this case, simply called R’s pgamma function).

Not sure why I didn’t think of this in the first place, but I hope this helps anyone else who happens to get stuck on the same problem!

2 Likes