Proper use of sum_to_zero_vector in nested multilevel models

[spinkney]

@aseyboldt

@spinkney If you add a normal_sum_to_zero(sigma), I’d consider defining it such that the marginal variance is \sigma^2 (1 - \tfrac{1}{N}). That is how pymc and numpyro define the ZeroSumNormal distribution, and having slightly different definitions around in the libraries could be annoying in the long term. I think that’s also usually what’s needed in applications of it. (At least it is in pretty much all cases where I’m using it, not sure how well that generalizes). It also feels quite natural to me if you look at the eigenvalues of the covariance matrix. One of those eigenvalues is zero, all others are \sigma in this parametrization.

This is something the Stan devs will need to discuss. I can see the benefit of what you’re saying although it seems harder to look at the code and know what’s going on. One benefit, and you mention this, is that we’d end up “backing-out” that sqrt(N * inv(N - 1)) to recover the population stuff, which is annoying. Initially the discussion was around how to say that each marginal is explicitly given the intended prior distribution. When a user writes z ~ sum_to_zero_normal(sigma) where z is a sum-to-zero vector it’s convenient to know that marginally you’re getting a standard deviation of sigma.