Jonah Piovia-Scott writes,

I am Ph.D. student studying ecology and evolution at the University of California, Davis. A group of us are trying to work through your latest book (co-authored with Jennifer Hill) to get a better understanding of bayesian hierarchical modeling. During the course of this week’s discussion, a question came up concerning the calculation of standard errors associated with multilevel (partial pooling) parameter estimates.

Specifically, figure 12.3b shows estimated intercepts from a multilevel model plotted against within-group sample size (data are from the radon example, so the groups are counties). Equation 12.1 gives us an approximation of how those multilevel intercept parameters are calculated, but we couldn’t figure out how the standard errors that appear on figure 12.3b are calculated.

It looks like you use one of the built-in functions in lmer, but we were curious about what actually goes into that calculation. Any clarification you could provide would be greatly appreciated.

My reply: Yes, we used the standard errors from lmer (in particular, se.coef() or se.ranef()). The general formula for the s.e. of these coefficients is the same matrix formula for classical regression s.e.’s, but with the data matrices augmented to account for the prior distribution (which, in this case, is itself estimated from the data). We discuss in detail in chapter 15 of Bayesian Data Analysis.

If you want to get the quicker take on this, see Section 18.3 of Gelman and Hill–in particular, p.394 gives the basic idea for the model with no individual-level predictors, and pp.396-397 sets up the matrices for the more general regression model.