1. Goal of the Example
We want to show how a hierarchical linear model can:
- Give more realistic uncertainty for forecasts, and
- Fix problems that appear when we fit a simple (nonhierarchical) regression.
Application: forecasting U.S. presidential elections using data up to 1988 and then using the model to forecast the 1992 election.
2. Data: Units and Outcome
- Units of analysis: state–election pairs
- Each state in each of the 11 elections: 1948, 1952, …, 1988
- Exclusions:
- District of Columbia (always extremely Democratic)
- States where third-party candidates actually won
After cleaning: 511 units.
Outcome variable
For state $s$ in year $t$:
$y_{st} = \text{Democratic share of the two-party presidential vote in state } s \text{ at year } t.$
Two-party share means:
$y_{st} = \frac{\text{Dem votes}}{\text{Dem votes} + \text{Rep votes}}$
3. Exploratory Plots
Two key scatterplots of $y_{st}$ in two consecutive elections:
- (a) 1984 vs 1988: points lie close to a straight line → previous election is a strong predictor.
- (b) 1972 vs 1976: relation is much weaker and noisier, but:
- Big outlier: Georgia, Jimmy Carter’s home state (Dem in 1976).
- Other southern states are above the line when Carter runs.
This suggests:
- Previous vote is predictive,
- But we need extra variables:
- home state of candidates,
- home region,
- regional patterns (South vs others).
4. First Step: Nonhierarchical Regression
We start with a single linear regression on all 511 observations:
$y_{st} \mid \theta, X_{st} \sim N\!\big( X_{st} \beta,\ \sigma^{2} \big)$
where:
- $X_{st}$ = row vector of predictors for state $s$ in year $t$
- $\beta$ = regression coefficients for all predictors
- $\sigma^{2}$ = residual variance
Types of predictors (columns of $X$)
They capture three levels:
- Nationwide variables (same for all states within a year)
- Candidate support in September national polls
- Presidential approval in July
- Economic indicators (e.g., GNP growth), interacted with incumbent party
- These vary across years, not across states.
- Statewide variables
- State’s Democratic vote share in last election (deviation from national average)
- Democratic vote share two elections ago
- Home-state indicators for presidential and VP candidates ($-1, 0, 1$)
- State economic growth, state ideology, ideological compatibility
- Proportion Catholic (for 1960 effect)
- These vary across states and years.
- Regional/subregional variables
- Indicators for South, New England, North Central, West in particular years.
- These capture region-specific patterns in specific elections (e.g., South in 1964).

Let:
- $n = 511$ observations
- $p$ = number of predictors in Table (≈ 20)
The preliminary model:
$y \mid \beta, \sigma^{2}, X \sim N(X\beta,\ \sigma^{2} I_{n})$
We can fit this with a noninformative prior:
$p(\beta, \sigma^{2}) \propto \frac{1}{\sigma^{2}}$
so posterior means/intervals match classical regression.
5. Forecasting with the Simple Model
Given a posterior for $(\beta, \sigma^{2})$, we could:
- For each MCMC draw $(\beta^{(m)}, \sigma^{2(m)})$,
- Use 1992 explanatory variables $\tilde{X}$ (national, regional, state predictors available before the election),
- Generate a predictive draw:
$\tilde{y}^{(m)} \sim N\!\big( \tilde{X}\beta^{(m)},\ \sigma^{2(m)} I_{50} \big)$
for the 50 states in 1992.
From these draws we can estimate:
- Probability Democrats win each state
- Probability Democrats win the Electoral College
- Expected number of Democratic states, etc.
But before trusting predictions, we must check the model.
6. Problem: Model Ignores Year-Level Correlation
The simple regression treats all 511 observations as independent, conditional on $X$:
- In reality, observations in the same year are correlated:
- If the Democrat exceeds expectations in Massachusetts in 1980,
they tend to exceed expectations in other states in that same year.
- If the Democrat exceeds expectations in Massachusetts in 1980,
So the assumption “conditional independence given $X$” is not realistic.
Consequence:
- We might underestimate the uncertainty in the national total (national average vote share) for each year.
7. Posterior Predictive Check for Nationwide Variation
We build a test statistic that measures how much national vote outcomes vary across elections.
For each year $t = 1,\dots,11$:
- Let $S_t$ = set of states in year $t$, $|S_t| \approx 50$.
- Residuals:
$e_{st} = y_{st} – X_{st} \beta$
- Nationwide average residual in year $t$:
$\bar{e}_t = \frac{1}{|S_t|} \sum_{s \in S_t} e_{st}$
Define the test statistic $T(y,\beta)$ as:
$T(y,\beta) = \sqrt{ \frac{1}{11} \sum_{t=1}^{11} \bar{e}_t^{2} }$
This measures year-to-year variability in national outcomes not explained by the model.
Posterior predictive check
- Draw $(\beta^{(m)}, \sigma^{2(m)})$ from posterior.
- Compute $T(y,\beta^{(m)})$ (realized test statistic).
- Simulate replicated data:
$y_{\text{rep}}^{(m)} \sim N\!\big( X\beta^{(m)},\ \sigma^{2(m)} I_{n} \big)$
and compute $T(y_{\text{rep}}^{(m)},\beta^{(m)})$.
Plot pairs:
$\big(T(y,\beta^{(m)}),\ T(y_{\text{rep}}^{(m)},\beta^{(m)})\big)$
for many draws (e.g. 200) in a scatterplot.
Result (nonhierarchical model)
- All points lie far below the 45° line.
- That means:
$T(y,\beta) \gg T(y_{\text{rep}},\beta) \ \text{typically}$
Interpretation:
- The actual election results show more year-to-year national variability than the model expects.
- The model is overconfident: it underestimates uncertainty in national outcomes.
So we need a better model.
8. Extending to a Hierarchical (Varying-Coefficients) Model
We introduce extra random effects for:
- Each election year (national swings)
- Each region × year (regional swings within each year)
Let:
- States indexed by $s$
- Regions indexed by $r(s) \in \{1,2,3,4\}$ (e.g., Northeast, Midwest, West, South)
- Years indexed by $t$
Model:
$y_{st} \mid \beta, \gamma, \delta, \sigma^{2} \sim N\!\big( X_{st}\beta + \gamma_{r(s)t} + \delta_t,\ \sigma^{2} \big)$
where:
- $X_{st}\beta$ = main regression part (from Table variables)
- $\gamma_{r(s)t}$ = region × year effect
- $\delta_t$ = year effect (national swing not captured by other covariates)
Hierarchical priors for varying coefficients
Regions:
suppose the 4 regions are indexed$r = 1,2,3,4$, and region 4 = South. Then:
$\gamma_{rt} \mid \tau_{\gamma1}^{2}, \tau_{\gamma2}^{2} \sim \begin{cases} N(0,\ \tau_{\gamma1}^{2}), & r=1,2,3 \ (\text{non-South}), \\ N(0,\ \tau_{\gamma2}^{2}), & r=4 \ (\text{South}). \end{cases}$
Year effects:
$\delta_t \mid \tau_{\delta}^{2} \sim N(0,\ \tau_{\delta}^{2})$
We also have residual variance:
$\sigma^{2}$
Hyperparameters:
$\tau_{\gamma1},\ \tau_{\gamma2},\ \tau_{\delta},\ \sigma > 0$
Hyperprior (weak):
$p(\beta, \sigma, \tau_{\gamma1}, \tau_{\gamma2}, \tau_{\delta}) \propto 1$
or something similar (essentially uniform on these scales; equivalently, noninformative).
Block notation
We can stack all coefficients into one big vector:
$\theta = \begin{pmatrix} \beta \\ \gamma \\ \delta \end{pmatrix}$
and define a prior:
$\theta \mid \tau_{\gamma1},\tau_{\gamma2},\tau_{\delta} \sim N\!\left( 0,\, \Sigma_{\theta} \right)$
where:
- Components of $\beta$ have variance $\infty$ (noninformative).
- Components corresponding to non-South $\gamma_{rt}$ have variance $\tau_{\gamma1}^{2}$.
- South $\gamma_{rt}$ have variance $\tau_{\gamma2}^{2}$.
- $\delta_t$ have variance $\tau_{\delta}^{2}$.
This is exactly a mixed-effects hierarchical model:
- Fixed effects: main β’s (Table)
- Random effects: year and region×year coefficients
9. Forecasting 1992 with the Hierarchical Model
For 1992, we have:
- New states × predictors: $\tilde{X}_{st}$ for 50 states
- New year $t = 1992$ (call it $t = 12$ conceptually)
- New region×year effects for 1992:
- $\gamma_{r,1992}$ for each region $r$
- $\delta_{1992}$ (year effect)
These new random effects are not in the original β-vector fitted from 1948–1988 data.
So for each posterior draw:
- Draw hyperparameters:
$(\beta^{(m)}, \sigma^{2(m)}, \tau_{\gamma1}^{(m)}, \tau_{\gamma2}^{(m)}, \tau_{\delta}^{(m)})$
- Draw new year effect:
$\delta_{1992}^{(m)} \sim N\!\left(0,\ \big(\tau_{\delta}^{(m)}\big)^{2}\right)$
- Draw new region×year effects:
$\gamma_{r,1992}^{(m)} \sim \begin{cases} N\!\left(0,\ \big(\tau_{\gamma1}^{(m)}\big)^{2}\right), & r \neq 4, \\ N\!\left(0,\ \big(\tau_{\gamma2}^{(m)}\big)^{2}\right), & r = 4. \end{cases}$
- For each state $s$ (in region $r(s)$), define mean:
$\mu_{s,1992}^{(m)} = \tilde{X}_{s,1992}\beta^{(m)} + \gamma_{r(s),1992}^{(m)} + \delta_{1992}^{(m)}$
- Generate predictive outcome:
$\tilde{y}_{s,1992}^{(m)} \sim N\!\big( \mu_{s,1992}^{(m)},\ \sigma^{2(m)} \big)$
From all draws $\{\tilde{y}^{(m)}\}$, we can:
- Estimate predictive mean and interval for each state
- Compute probability the Democrat wins each state
- Convert to Electoral College probabilities (e.g., Dem win probability ≈ 85%)
They report:
- Predictive standard errors for state vote shares: 5–6 percentage points
→ more realistic than the overconfident nonhierarchical model.
10. Posterior Inference for the Hierarchical Model
They fit the model using:
- EM algorithm to find a sensible starting point (posterior mode or approximation)
- Vector Gibbs sampler to simulate from the joint posterior:
$p(\beta, \gamma, \delta, \sigma, \tau_{\gamma1}, \tau_{\gamma2}, \tau_{\delta} \mid y)$
They run:
- 10 parallel chains
- 500 iterations each
- Potential scale reduction factors:
$\hat{R} < 1.1$
for all parameters → convergence seems good.
Key results:
- Coefficients for Table variables are similar to nonhierarchical regression.
- Standard deviations $\tau_{\gamma1}, \tau_{\gamma2}, \tau_{\delta}$ are not very precisely determined → but that uncertainty is propagated into predictions (a strength of Bayesian modeling).
11. Model Checking for the Hierarchical Model
Use the same test statistic as before:
$T(y,\beta) = \sqrt{ \frac{1}{11} \sum_{t=1}^{11} \bar{e}_t^{2} }$
Now:
- Simulate $y_{\text{rep}}$ from the hierarchical model, including year and region×year effects.
- Compute pairs:
$\big(T(y,\beta^{(m)}),\ T(y_{\text{rep}}^{(m)},\beta^{(m)})\big)$
Result:
- Points are scattered around the 45° line, not all below it.
- The observed value $T(y,\beta)$ is typical under the model.
Conclusion:
- The hierarchical model captures the nationwide year-to-year variability much better.
- Uncertainty for national outcomes is now more realistic.
12. Why the Hierarchical Model Is Better
The hierarchical model improves the election forecast in three key ways:
- Captures correlation within years and regions
- Year effects $\delta_t$ model nationwide shifts.
- Region×year effects $\gamma_{rt}$ model regional swings.
- Observations in same year/region are no longer treated as independent.
- Generalizes between two extremes
- If $\tau_{\delta}, \tau_{\gamma1}, \tau_{\gamma2} \to 0$: all year and region×year effects are 0 → pure nonhierarchical model.
- If $\tau_{\delta}, \tau_{\gamma1}, \tau_{\gamma2} \to \infty$: extremely large variability in these effects.
- The hierarchical model learns intermediate values from data → partial pooling.
- More realistic predictive uncertainty
- Predictions include extra uncertainty from:
- residual noise $\sigma^{2}$
- uncertain future year effect $\delta_{1992}$
- uncertain future region×year effects $\gamma_{r,1992}$
- uncertainty in hyperparameters $\tau_{\cdot}$
- Result: wider but honest predictive intervals for state and national outcomes.
- Predictions include extra uncertainty from:
