1988 Election: Multilevel Modeling

About the Data

We model individual vote choice in the 1988 US Presidential Election (Bush vs. Dukakis) using simulated survey data (following Gelman & Hill). The data has a natural hierarchical structure:

• Level 1 (Individuals): 6,003 survey respondents with demographic characteristics
• Level 2 (States): 50 states with different political climates (~120 respondents per state)
• Outcome: Binary vote choice (1 = Bush, 0 = Dukakis)
• Bush vote share: 53.5%

Individual-Level Predictors

Download the data: election_1988.csv

Step 5: Plot Odds Ratios

Exponentiate the coefficients and their credible intervals to plot on the odds ratio scale. The red dashed line at OR = 1 represents no effect.

# Exponentiate coefficients to get odds ratios
or_data = pd.DataFrame({
    'OR': np.exp(coef_data['mean']),
    'OR_low': np.exp(coef_data['0.025quant']),
    'OR_high': np.exp(coef_data['0.975quant'])
}, index=coef_data.index)
or_data = or_data.sort_values('OR', ascending=True)

# Green if CI excludes 1 (significant on OR scale)
y_pos = np.arange(len(or_data))
colors = ['#059669' if not ((row['OR_low'] < 1) and (row['OR_high'] > 1))
          else '#94a3b8' for _, row in or_data.iterrows()]

fig, ax = plt.subplots(figsize=(10, 6))
ax.hlines(y_pos, or_data['OR_low'], or_data['OR_high'],
          colors=colors, linewidth=2, alpha=0.7)
ax.scatter(or_data['OR'], y_pos, c=colors, s=80, zorder=5)
ax.axvline(x=1, color='#e11d48', linestyle='--', linewidth=1.5, alpha=0.7)

ax.set_yticks(y_pos)
ax.set_yticklabels(or_data.index)
ax.set_xlabel('Odds Ratio (vs reference category)')
ax.set_title('1988 Election: Odds Ratios for Bush Vote')
ax.grid(axis='x', alpha=0.3)
plt.tight_layout()
plt.savefig('election_odds_ratios.png', dpi=150, bbox_inches='tight')

Step 6: State Random Effects (Caterpillar Plot)

The caterpillar plot shows the estimated random intercept for each state, ordered from lowest to highest:

• States above zero: More pro-Bush than demographics alone predict
• States below zero: More pro-Dukakis than demographics alone predict
• Interval width: Wider intervals indicate more uncertainty (often smaller sample sizes)
• Shrinkage: Notice most estimates are close to zero, partial pooling in action

# Extract and sort state random effects
state_effects = py_result.summary_random['state'].copy()
state_effects = state_effects.sort_values('mean')

n_states = len(state_effects)
y_pos = np.arange(n_states)

# Blue if CI excludes zero (state is distinguishable from average)
colors = ['#3b82f6' if not ((row['0.025quant'] < 0) and (row['0.975quant'] > 0))
          else '#94a3b8' for _, row in state_effects.iterrows()]

fig, ax = plt.subplots(figsize=(12, 8))
ax.hlines(y_pos, state_effects['0.025quant'], state_effects['0.975quant'],
          colors=colors, linewidth=1.5, alpha=0.6)
ax.scatter(state_effects['mean'], y_pos, c=colors, s=30, zorder=5)
ax.axvline(x=0, color='#e11d48', linestyle='--', linewidth=1.5, alpha=0.7)

ax.set_ylabel('State (ordered by effect)')
ax.set_xlabel('State Random Effect (log-odds)')
ax.set_title('State-Level Variation in Bush Support (Caterpillar Plot)')
ax.set_yticks([])
plt.tight_layout()
plt.savefig('election_state_effects.png', dpi=150, bbox_inches='tight')

Step 7: State-Level Variance

How much do states vary after accounting for demographics? Convert precision (what INLA estimates) to variance and SD, then plot the distribution of state effects.

from scipy import stats

# Convert precision to variance and SD
prec_mean = py_result.summary_hyperpar['mean'].values[0]
state_var = 1 / prec_mean         # variance = 1 / precision
state_sd = np.sqrt(state_var)     # SD = sqrt(variance)

print(f"State effect precision: {prec_mean:.3f}")
print(f"State effect variance:  {state_var:.3f}")
print(f"State effect SD:        {state_sd:.3f}")

# Histogram of state effects with fitted normal overlay
state_means = py_result.summary_random['state']['mean'].values

fig, ax = plt.subplots(figsize=(8, 5))
ax.hist(state_means, bins=15, density=True, alpha=0.7,
        color='#3b82f6', edgecolor='white')
ax.axvline(0, color='#e11d48', linestyle='--', linewidth=2,
           label='Zero (no state effect)')

# Overlay the fitted N(0, state_sd) distribution
x = np.linspace(state_means.min() - 0.2, state_means.max() + 0.2, 100)
ax.plot(x, stats.norm.pdf(x, 0, state_sd), 'k-', linewidth=2,
        label=f'Fitted N(0, {state_sd:.2f})')

ax.set_xlabel('State Random Effect')
ax.set_ylabel('Density')
ax.set_title('Histogram of posterior mean state effects.\n'
             'The fitted normal distribution shows the estimated between-state variation.',
             fontsize=10, ha='center')
ax.legend()
plt.tight_layout()
plt.savefig('election_state_variance.png', dpi=150, bbox_inches='tight')

State effect precision: 5.992
State effect variance:  0.167
State effect SD:        0.409

True variance used in simulation: 0.2
True SD used in simulation:       0.447

The estimated state precision is 5.99, corresponding to a variance of 0.167 and a standard deviation of 0.409 on the log-odds scale. The true simulation variance was 0.2 (SD = 0.447), so the model recovers this well. In practical terms, a one-SD shift in the state effect translates to roughly a 50% relative change in odds between typical pro-Bush and pro-Dukakis states.

Key Takeaways

1. Multilevel models handle clustered data: voters within states are not independent
2. Random intercepts capture group-level variation: some states are more pro-Bush than demographics predict
3. Shrinkage improves estimates: states with little data borrow strength from the overall pattern
4. Race was the strongest predictor: Black voters had 75% lower odds of voting Bush (OR = 0.25)
5. Gender and age matter: women and 30–44-year-olds had significantly lower odds of supporting Bush
6. pyINLA syntax is simple: add 'random': [{'id': 'state', 'model': 'iid'}]

References

• Gelman, A., and Hill, J. (2006). Data Analysis Using Regression and Multilevel/Hierarchical Models. Cambridge University Press. (Chapters 11-14)
• Rue, H., Martino, S. and Chopin, N. (2009). "Approximate Bayesian inference for latent Gaussian models using integrated nested Laplace approximations." Journal of the Royal Statistical Society: Series B, 71(2), 319-392.
• ANES (American National Election Studies) for historical election survey data.