(brms) Loo pit overlay differences when mu and sigma are in the transformed parameters block and loglik/yrep are in the generated quantities

Interfaces brms
1

Hello all, I have been having a problem with my loo-pit diagnostics.

I have a stan model and a brms model that are supposed to be equal. The point is that it’s an exercise in writing stan models and confirming that they’re written correctly by comparing the results to those of brms.
My default assumption is that when model results are different then brms is correct.

Both my stan model and brms models are supposed to be approximately identical except for the fact that my stan model creates mu and sigma terms in the transformed parameters block which I then use in the model block for the likelihood and the generated quantities block to create log_lik and yrep.

The two models produce essentially equal parameter estimates and the pp_check(brms.fit)/ppc_dens_overlay(stan.fit) produce plots that look identical to me.
The loo() results are approximately the same, but not exactly.

The pp_check(brms.fit, type=“loo_pit_overlay”)/ppc_loo_pit_overlay() plots produce wildly different results: the brms loo_pit_overlay result indicates that there are no glaring model problems whereas the stan model results in a classic frowny face.

I assume that this is an error on my part and that there are 2 possibilities:

  1. I am calculating yrep incorrectly.
  2. I am calculating log_lik incorrectly.

I cannot share my data, however, I have attached code to reproduce a similar effect using the Oats data from nlme in the attached file.
This code fits two models: a brms model and a stan model that uses slightly modified brms code (mu and sigma are moved to the tparms block and log_lik and yrep are created in the generated quantities block).

This example shows that the resulting loo pit overlay plots are different and that it does not seem to be a difference in pp_check vs. ppc_loo_pit_overlay.

This is the brms loo-pit output:

image
image1735×1061 46.3 KB

Here is the stan loo-pit output where I calculate log_lik and yrep in the generated quantities:

image
image1723×1063 43.9 KB

For the life of me I cannot figure out why these plots are so different when they should be about the same (if not identical).

Thanks!

Software

  • Windows 10 64bit
  • R 4.2.1
  • rstan version 2.35.0.9000 (Stan version 2.35.0)
  • brms version 2.21.10
  • bayesplot version 1.11.1.9000

pit_ex.R (5.1 KB)

2

This is due to a very unfortunate design choice in RStan to make extract() by default to return the draws in randomly permutated order. Thus, your yrep and lw don’t have the same order and loo-pit output is nonsense. You can use extract(, permute=FALSE), but the another very unfortunate design choice makes the output of the extract to have different structure. If I use

yrep.stan <- as_draws_matrix(rstan::extract(fit2, pars="yrep", permute=FALSE))

I get the same plot as from brms.

3 Likes
3

I would have never guessed that in a million years.

Thanks for the clear explanation!