How to get mu and shape for Randomized Quantile Residuals

Interfaces brms
1

I came across this paper about Randomized Quantile Residuals and have found the idea pretty useful:

I’m using this idea on a mixed effects model with a negative binomial error distribution. On page 9 of the paper linked above, you’ll see the formula I want to calculate in R after the Stan run has been completed. For each observation, I want the associated mu and shape parameters for the negative binomial distribution.

Since brms calculates mu in the Model block, I don’t think Stan will print it (page 101 of the Stan 2.17.0 manual). When I coded models myself, I’d get around this by calculating mu in the Transformed Parameters block. For that matter, if I were to have the shape parameter vary by a class variable, the shape parameter also ends up calculated in the Model block.

So my question: is there a simple way in R to take a brms object and have it spit out the associated mu and shape parameter for each observation?

Thanks,
Paul

  • Operating System: windows
  • brms Version: 2.2
2

For mu, you can call posterior_linpred with transform = TRUE.

3

posterior_linpred (or equivalently fitted) will always return the mean by default, which is indead mu in your case, but for some other families the mean may actually be not the same as mu (lognormal is a good example).

For your purpose, go for

fitted(<your model>, dpar = "mu")
fitted(<your model>, dpar = "shape")
4

Exactly what I was looking for. Thanks for the quick reply.

5

A follow-up, and please tell me if I’m correct:

If I have a simple model with a negative binomial error distribution:

model_Formula <- bf(y ~ x1)

the fitted function was giving me an error for shape:

Error: Distributional parameter 'shape' was not predicted.

This only seemed to be a problem when I had not added a specific formula for shape in the formula (i.e. shape was not calculated in the Model block). In this case, I believe posterior_summary and/or posterior_samples gives me the equivalent results. I think this makes sense because, in this first example, there is only one shape parameter for the entire model.

But if I had a formula for shape, something like:

model_Formula <- bf(y ~ x1,
                    shape ~ (1|x2))

then fitted worked as described in earlier replies. Do I have this right?

6

you are exactly right.

7

In the dev version of brms from github, argument dpar in fitted will now work even if the parameter has no associated formula and is thus constant across observations.