From the model graph you posted it looks like everything is implemented correctly. I would just suggest you make the district_to_block
and block_to_data
maps pm.Data
, in case you want to do out-of-sample prediction. You also don’t need the district_idx
at all – you can see that it’s not attached to anything in the graph.
For discussion on this type of “telescoping” hierarchy, see here and here. You will need to run prior predictive simulations and parameter recovery studies to make sure your model is identified. In general, I’m not sure these types of models are.
Wishart covariance priors are used in the literature because they have nice conjugate properties with Multivariate normal. PyMC does not take advantage of this. It is recommend to instead use an LKJCholeskyCov prior.