Elections and the economy: Hibbs regression

Source: ElectionsEconomy/hibbs.Rmd

This page fits the incumbent-party vote regression with lapylace, matching the original Bayesian linear-regression target.

Setup

from pathlib import Path
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import lapylace as lp

root = Path("../../ROS-Examples")

def coef_medians(fit):
    return pd.Series(
        [np.median(fit.alpha_draws()), *np.median(fit.beta_draws(), axis=0)],
        index=["Intercept", *fit.columns],
    )

def fit_glm(formula, data, family=None, seed=1, prior_scale=2.5, intercept_scale=5, aux_scale=10, **kwargs):
    return lp.stan_glm(
        formula,
        data=data,
        family=family or lp.gaussian(),
        prior=lp.normal(0, prior_scale),
        prior_intercept=lp.normal(0, intercept_scale),
        prior_aux=lp.exponential(aux_scale),
        chains=2,
        parallel_chains=2,
        iter_warmup=300,
        iter_sampling=500,
        seed=seed,
        refresh=100,
        **kwargs,
    )

hibbs = pd.read_table(root / "ElectionsEconomy/data/hibbs.dat", sep=r"\s+")
hibbs.head()

	year	growth	vote	inc_party_candidate	other_candidate
0	1952	2.40	44.60	Stevenson	Eisenhower
1	1956	2.89	57.76	Eisenhower	Stevenson
2	1960	0.85	49.91	Nixon	Kennedy
3	1964	4.21	61.34	Johnson	Goldwater
4	1968	3.02	49.60	Humphrey	Nixon

fit = fit_glm("vote ~ growth", hibbs, seed=17701, prior_scale=10, intercept_scale=100, aux_scale=10)
coef = coef_medians(fit)
fit.summary(["alpha", "beta", "sigma"])

	Mean	MCSE	StdDev	5%	50%	95%	N_Eff	N_Eff/s	R_hat
alpha	46.33050	0.089919	1.819120	43.32570	46.36640	49.25970	409.2820	22737.90000	1.00371
beta[1]	3.02169	0.037608	0.762204	1.77057	3.03201	4.27065	410.7470	22819.30000	1.00520
sigma	3.98695	0.037430	0.814850	2.78476	3.92088	5.44362	473.8723	26326.23866	1.00149

grid = pd.DataFrame({"growth": np.linspace(hibbs.growth.min(), hibbs.growth.max(), 100)})
pred = fit.posterior_epred(grid)
fig, ax = plt.subplots()
ax.scatter(hibbs.growth, hibbs.vote, color="black")
ax.plot(grid.growth, pred.mean(axis=0), color="black")
ax.fill_between(grid.growth, np.quantile(pred,.1,axis=0), np.quantile(pred,.9,axis=0), color="0.85")
ax.set_xlabel("Economic growth")
ax.set_ylabel("Incumbent-party vote")

Text(0, 0.5, 'Incumbent-party vote')