Hierachical models

Combine variation (and covariates) from multiple levels (analagous to mixed models). Need to preserve exchangability, usually in terms of conditional indepdence (given parameters). Call model distributions of model parameters hyperpriors. Need to be very careful when using non-proper hyperpriors.

Often two posterior predictive distributions of interest:

• distribution of future observations from this experiment
• distribution of future observations from all possible experimenters (ie. hyperparmeters varied as well)

Computation is harder than before because we have more parameters. Easiest when population distribution is conjugate to likelihood. Usual steps:

• write `p(phi, theta | y) propto p(phi) p(theta | phi) p(y | theta)`
• determine `p(theta | phi, y)` analytically – easy for conjugate models
• derive marginal `p(phi | y)` by integrating out `theta`, or `p(phi, theta|y) / p(theta|y, phi)` (but need to be careful because normalising “constant” might not be constant)

Steps for computation are usually:

• draw `phi` from `p(phi | y)`
• draw `theta` `p(theta | phi, y)` – usually factors into product of independents
• draw `tilde y` from appropriate posterior distribution

Beta-binomial

`y_j ~ "Bin"(n_j, theta_j)`, `theta_j ~ "Beta"(alpha, beta)`. `p(alpha, beta)` non informative (will need to check integrability of posterior). `p(alpha, beta) propto 1` doesn’t work, `p(alpha/(alpha + beta), (alpha + beta)^(-1/2)) propto 1` does.

Normal hierachical model

(equivalent to one-way random effects). Set up: `J` indepdendent experiments, want to estimate mean `theta_j` from each, sampling model: `y_ij | theta_j ~ N(theta_j, sigma^2)` (`sigma^2` known), and `bar y_j | theta_j ~ N(theta_j, sigma_j^2)`. This sampling model is quite general because of CLTs. What type of posterior estimates for `theta_j`’s might be reasonable? `hat theta_j = bar theta_j` if `n_n` large, `hat theta_j = bar y` reasonable if we believe all means are the same, or more general estimator `hat theta_j = lambda_j bar y_j + (1 - lambda_j) bar y`. `(1 - lambda_j)` called a shrinkage factor.

• `bar_j | theta ~ N(theta_j, sigma_j^2`
• `theta_1, ..., theta_j | mu, tau ~ product_j N(mu, tau)`
• so `p(theta_1, ..., theta_j) = int product_j N(theta_j; mu, tau^2) p(mu, tau^2) dmu dtau^2`

Joint prior can be written as `p(mu, tau) = p(mu | tau) p(tau)`. Will consider a conditionally flat prior so that `p(mu, tau) propr p(tau)`

For `p(mu | tau) propto 1`, `p(theta_j | mu, tau) = N(hat theta_j, V_j)` where `hat theta_j = (sigma_j^2 mu + tau^2 bar y_j)/(sigma_j^2 + tau^2)`, `1/(V_j) = 1/(sigma_j^2) + 1/(tau^2)`

`p(mu, tau | y) propto p(tau) product_j N(bar y_j; mu, sigma_j^2 + tau^2)`. Now `p(tau | y) propto p(tau) V_mu^(-1/2) product_j N(bar y_j; hat mu, sigma_j^2 + tau^2)`. What should we used for `p(tau)`? A conjugate is always a safe choice so, `p(tau) ~ Inv-chi^2(nu_0, tau_0^2)`