I discovered tonight that I’ve spent years misunderstanding the example in the
SUG 1.13 (“Multivariate priors for hierarchical models”), and now have a couple questions:
The absence of a prior on
sigma(observation-level noise scale) is surely an error, yes? I can’t imagine that the implied
uniform(0,+Inf)could be expected to perform appropriately.
I feel like there are serious identification issues in the model at present due to the combination of three structures:
a) The data
Yconsist of single observations from each of
Nindividuals, modelled as
muis itself a vector of length
N, each entry of which encodes the dot-product between the corresponding row in the N \times K contrast matrix
xand the column in the K \times J coefficient matrix
betacorresponding to the group of that individual.
c) The K \times J coefficient matrix
betais modelled as a multivariate normal with the square of the scale vector
tauforming the diagonal.
Thus, I don’t think that
tau are actually identifiable with this structure. These quantities would be properly identified if we had multiple observations from each of the individuals (or even just one of the individuals), but that’s explicitly not how the data declarations are set up.
I actually think I long ago worked out (with help from other Stan gurus as I recall) code for a version of the model that doesn’t have this identification issue (old video walkthrough here, but the critical bit is the use of the
get_coefs() function in the transformed data section to transform the data to the coefficients that are then modelled as multivariate normal in the model).
Or am I missing something (again)?
(Tagging @andrewgelman but happy to receive input/correction from anyone)