Passing a variable as part of a prior in brms

#1

I am trying to pass a computed variable as part of a prior, e.g.

s <- ...
prs <- c(prior(normal(0, s), class='b', coef='age'), ...)

But I get an error message during Stan compilation. How should I do this? Note that I don’t want this to represent a hyperprior; it’s just a computed constant for the standard deviation of the prior.

Here’s the code and error message:

prs <- c(prior(student_t(3,0,10), class='Intercept'),
         prior(normal(0, s1), class='b', coef='sexmale'),
         prior(normal(0, s2), class='b', coef='age'))

fit_brms <- brm(response ~ sex + age, data=sex.age.response,
              family=bernoulli("logit"), prior=prs, iter=5000)

SYNTAX ERROR, MESSAGE(S) FROM PARSER:

variable "s1" does not exist.
  error in 'model1d075d3cb4f5_file1d0711f5317d' at line 29, column 37
  -------------------------------------------------
    27:   vector[N] mu = Xc * b + temp_Intercept; 
    28:   // priors including all constants 
    29:   target += normal_lpdf(b[1] | 0, s1); 
                                            ^
    30:   target += normal_lpdf(b[2] | 0, s2); 
  -------------------------------------------------

Error in stanc(model_code = paste(program, collapse = "\n"), model_name = model_cppname,  : 
  failed to parse Stan model 'file1d0711f5317d' due to the above error.
  • Operating System: xubuntu 17.04
  • brms Version: 2.2.0
#2

Hi @harrelfe, can you share the error message you get during compilation? That will point me in the right direction. This is also something that should we should look into catching before passing to Stan if possible, so that there would be an informative error message from R.

#3

Function stanvar() is what you are looking for.

@jonah I guess the error message just says that s cannot be found.

#4

Sorry I neglected to include that. I edited the original post to do so.

#5

This worked. Thank Paul! Here’s the resulting code.

# Solve for SD such that sex effect has only a 0.025 chance of
# being above 5 (or being below -5)

s1 <- 5 / qnorm(0.975)

# Solve for SD such that 10-year age effect has only 0.025 chance
# of being above 20

s2 <- (20 / qnorm(0.975)) / 10   # divide by 10 since ratio on 10*b scale
stanvars <- stanvar(s1, name='s1') + stanvar(s2, name='s2')

prs <- c(prior(student_t(3,0,10), class='Intercept'),
         prior(normal(0, s1), class='b', coef='sexmale'),
         prior(normal(0, s2), class='b', coef='age'))

fit_brms <- brm(response ~ sex + age, data=sex.age.response,
              family=bernoulli("logit"), prior=prs, iter=5000,
              stanvars=stanvars)