Visualizing many hierarchial model parameters

General
1

I am fitting a model where unobserved parameters in a 2-simplex are transformed from a multivariate normal, then they generate the actual data. The relevant lines are not unlike

for (i in 1:N) {
  z[i] ~ multi_normal_cholesky(mu, L); // 2-dimensional mv normal
  a[i] = inv_logit(z[i][1]);           // transform to [0,1]
  b[i] = (1-a[i])*inv_logit(z[i][2]);  // transform to make a simplex
  // ... a[i] and b[i] unobserved, they are parameters to a multinomial
  // ... that generates observations, omitted here
}

The multivariate normal is just a modeling device, I care about a and b. N is around 1000, so the questions is how to visualize this. I thought of using HPD regions graphically. I could

  1. use the posterior mean of mu and L, and visualize the conditional a and b based on that, or
  2. I could use all the a and bs, and plot HPD regions for that.

I am wondering which one makes more sense conceptually. (1) does not capture posterior uncertainty, while I think that (2) would overstate it.

2

That’s the “logistic normal” if you use softmax. I think I mention it in the manual.

It’s explicitly the relative log odds.

Usually you want to present all the uncertainty, including the uncertainty from estimating the parameters and the residual sampling uncertainty.