Replication and investigation of "What Are Trump's Chances Of Winning The GOP Primary?" <https://fivethirtyeight.com/features/trump-chances-to-win-republican-primary/>

primary_model.R

library(tidyverse)
library(janitor)
library(rvest)
library(rstanarm)
library(posterior)

# load data, add 2024 candidates

page <- read_html("https://fivethirtyeight.com/features/trump-chances-to-win-republican-primary/")
d <- html_table(page) |>
  pluck(1) |>
  clean_names() |>
  mutate(win = if_else(cycle < 2024L, str_length(winner), NA_integer_),
         aug_poll = parse_number(polling_average_in_august_beforeelection_year)/100,
         aug_poll = if_else(aug_poll == 0, 0.003, aug_poll),
         pop_vote = parse_number(popular_vote, na="TBD") / 100) |>
  select(cycle, party, candidate, aug_poll, pop_vote, win) |>
  bind_rows(
    data.frame(cycle=2024L, party="R", candidate="Ron DeSantis", aug_poll=0.146),
    data.frame(cycle=2024L, party="R", candidate="Vivek Ramaswamy", aug_poll=0.099),
    data.frame(cycle=2024L, party="R", candidate="Nikki Haley", aug_poll=0.053),
    data.frame(cycle=2024L, party="R", candidate="Mike Pence", aug_poll=0.045),
    data.frame(cycle=2024L, party="R", candidate="Chris Christie", aug_poll=0.037),
    data.frame(cycle=2024L, party="R", candidate="Tim Scott", aug_poll=0.032),
  ) |>
  add_count(cycle, party) |>
  arrange(desc(cycle), desc(aug_poll))
# count(d, cycle, party, wt=win) # sanity check

d_fit = filter(d, aug_poll >= 0.01)
d_new = select(d, cycle, party, n, aug_poll)


# fit the original model

m0_b = stan_glm(win ~ aug_poll, data=d_fit, family=binomial(),
                prior = normal(0, 1, autoscale=TRUE),
                iter=2000, warmup=1000, chains=1, refresh=0)
d |>
  mutate(.fitted = rvar(posterior_epred(m0_b, newdata=d_new)),
         q05 = quantile2(.fitted, 0.05),
         q95 = quantile2(.fitted, 0.95)) |>
  filter(cycle == 2024)

# fit the original linpred specification with the poisson trick

m0p_b = stan_glm(win ~ aug_poll, data=d_fit, family=poisson(),
                 prior = normal(0, 2, autoscale=TRUE),
                 iter=2000, warmup=1000, chains=1, refresh=0)
d |>
  mutate(.fitted = rvar(posterior_epred(m0p_b, newdata=d_new))) |>
  mutate(.fitted = .fitted / rvar_sum(.fitted),
         q05 = quantile2(.fitted, 0.05),
         q95 = quantile2(.fitted, 0.95),
         .by=c(cycle, party)) |>
  filter(cycle == 2024)

# alternative model

m2p_b = stan_glm(win ~ (log(n) + party) * aug_poll, data=d_fit, family=poisson(),
                 prior = normal(0, 1.5, autoscale=TRUE),
                 iter=2000, warmup=1000, chains=1, refresh=0)
d |>
  mutate(.fitted = rvar(posterior_epred(m2p_b, newdata=d_new))) |>
  mutate(.fitted = .fitted / rvar_sum(.fitted),
         q05 = quantile2(.fitted, 0.05),
         q95 = quantile2(.fitted, 0.95),
         .by=c(cycle, party)) |>
  filter(cycle == 2024)