🍰 Postre logístico

Plato probabilístico para la segunda. Too close to call: el postre le da un 45% de probabilidad de ganar las elecciones a Gustavo Petro

Autor/a

Afiliación

Recetas Electorales

Análisis independiente

“There is no inference without assumption, but do not choose your assumptions for the sake of inference. Build complex models one piece at a time. Be critical. Be kind.”
– Richard McElreath

El Postre 2022 estima qué tan probable es que gane Gustavo Petro. A diferencia del Ajiaco, el Postre solo hace una transformación logística a un modelo sencillo que se puede interpretar como un pronóstico probabilístico para la segunda vuelta.

Veamos cómo se prepara.

---

Postre Logístico 2022

Este nuevo Postre recoge la probabilidad de que Gustavo Petro sea presidente, estimada por el modelo que se describe abajo.

Las probabilidades se obtienen de muchas simulaciones de cada uno de los parámetros del modelo, y luego calcular la mediana de la densidad que resulta de esas simulaciones.

En estas primeras estimaciones la densidad es bimodal, porque las encuestas se contradicen: las telefónicas dan como ganador a Hernandez, mientras que las presenciales dan como ganador a Petro.

ImportanteCarrera apretada, según el modelo basado en las encuestas

Segun la más reciente estimación, Gustavo Petro tiene más o menos 45% de probabilidad de ganar

library(rethinking)
library(RColorBrewer)
library(lubridate)
library(tidyverse)
library(rstan)

# ponderador LSV ####
lsv_semaforo <- tribble(~encuestadora,~semaforo,
                          "CNC",7.1,
                          "Invamer",10,
                          "GAD3",7,
                          "AtlasIntel",7.5,
                          "MassiveCaller",5.6,
                          "YanHaas",5.5,
                          "Guarumo",8.9,
                          "TYSE",5.5,
                          "CELAG",5.1)

# Encuestas ####
postre2022_encuestas <- readr::read_csv("https://raw.githubusercontent.com/nelsonamayad/Elecciones-presidenciales-2022/main/Encuestas%202022/encuestas_2022.csv") %>% 
  dplyr::select(n,fecha,encuestadora,muestra,gustavo_petro,rodolfo_hernandez,tipo) %>%
  dplyr::filter(fecha>=lubridate::as_date("2022-05-29")) %>%
  dplyr::mutate(gp=gustavo_petro*muestra/100,
                rh=rodolfo_hernandez*muestra/100,
                encuestadora =factor(encuestadora)) %>%
  dplyr::mutate(dd = lubridate::as_date(as.character(today()), format="%Y-%m-%d") - lubridate::as_date(as.character(fecha), format="%Y-%m-%d")) %>%
  dplyr::mutate(dd = as.numeric(dd)) %>% 
  dplyr::mutate(gp_win = ifelse(gp>rh,1,0),
                rh_win = ifelse(gp<rh,1,0)) %>% 
  dplyr::mutate(tipo_1=ifelse(tipo=="presencial",1,0),
                tipo_2=ifelse(tipo=="telefonico y presencial",1,0),
                tipo_3=ifelse(tipo=="digital",1,0))

# Numero de encuestas disponibles
postre2022_disponibles <- readr::read_csv("https://raw.githubusercontent.com/nelsonamayad/Elecciones-presidenciales-2022/main/Encuestas%202022/encuestas_2022.csv") %>%
  dplyr::filter(lubridate::as_date(fecha) >= lubridate::as_date("2022-05-30")) %>%
  dplyr::tally() %>%
  as.numeric()

postre2022_encuestadoras <- readr::read_csv("https://raw.githubusercontent.com/nelsonamayad/Elecciones-presidenciales-2022/main/Encuestas%202022/encuestas_2022.csv") %>% 
  dplyr::filter(lubridate::as_date(fecha) >= lubridate::as_date("2022-05-30")) %>%
  dplyr::select(encuestadora) %>%
  dplyr::n_distinct() %>%
  as.numeric()

# Postre ulam Petro ####
postre2022_ulam_gp <- ulam(
  alist(
    #modelo
    gp ~ dbinom(muestra,p),
    logit(p) <- a[encuestadora] + b1*dd + b2*tipo_1,
    #priors
    a[encuestadora] ~ dnorm(a,s),
    a ~ dnorm(45,5),
    b1 ~ dnorm(0,5),
    b2 ~ dnorm(0,5),
    s ~ dcauchy(0,5)
  ),
  data=list(
    N = postre2022_disponibles,
    encuestadora = postre2022_encuestas$encuestadora,
    muestra = postre2022_encuestas$muestra,
    dd = postre2022_encuestas$dd,
    tipo_1 = postre2022_encuestas$tipo_1,
    gp= postre2022_encuestas$gp
    ),
  control=list(adapt_delta=0.96),
  iter=5000, 
  warmup=1000, 
  chains=4, 
  cores=4,
  log_lik=TRUE)

# Postre ulam Hernandez ####
postre2022_ulam_rh <- ulam(
  alist(
    #modelo
    rh ~ dbinom(muestra,p),
    logit(p) <- a[encuestadora] + b1*dd + b2*tipo_1,
    #priors
    a[encuestadora] ~ dnorm(a,s),
    a ~ dnorm(47,5),
    b1 ~ dnorm(0,5),
    b2 ~ dnorm(0,5),
    s ~ dcauchy(0,5)
  ),
  data=list(
    N = postre2022_disponibles,
    encuestadora = postre2022_encuestas$encuestadora,
    muestra = postre2022_encuestas$muestra,
    dd = postre2022_encuestas$dd,
    tipo_1 = postre2022_encuestas$tipo_1,
    rh= postre2022_encuestas$rh
    ),
  control=list(adapt_delta=0.96),
  iter=5000, 
  warmup=1000, 
  chains=4, 
  cores=4,
  log_lik=TRUE)

library(rethinking)
library(tidyverse)
library(RColorBrewer)

# Estimar prob 50% para incluir en grafico 
prob_win_gp <- postre2022_ulam_gp %>% 
  rethinking::link(n=20000) %>% 
  tibble::as_tibble() %>% 
  tidyr::pivot_longer(cols = everything(), names_to = "poll_n", values_to = "estimate") %>% 
  dplyr::summarize(m = median(estimate),
                   p10 = quantile(estimate,prob=0.1),
                   p90 = quantile(estimate,prob=0.9)) %>%
  dplyr::select(m) %>%
  as.numeric()

prob_win_rh <- postre2022_ulam_rh %>% 
  rethinking::link(n=20000) %>% 
  tibble::as_tibble() %>% 
  tidyr::pivot_longer(cols = everything(), names_to = "poll_n", values_to = "estimate") %>% 
  dplyr::summarize(m = median(estimate),
                   p10 = quantile(estimate,prob=0.1),
                   p90 = quantile(estimate,prob=0.9)) %>%
  dplyr::select(m) %>%
  as.numeric()


#postre logis
postre2022_ulam_gp %>% 
  rethinking::link(n=20000) %>% 
  tibble::as_tibble() %>% 
  tidyr::pivot_longer(cols = everything(), names_to = "poll_n", values_to = "estimate") %>% 
  ggplot()+
  geom_density(aes(x=estimate),fill="#800080",color="grey60", alpha=0.8)+
  geom_vline(xintercept = prob_win_gp, color="gold", linetype="dashed", size=1.2)+
  annotate("text",x=0.38,y=15,
           color="#800080",
           label=paste0("Probabilidad de ganar: ",round(100*prob_win_gp, digits = 0),"%"))+
  labs(x="Densidad de probabilidad de sacar mas de 50% de la votación",
       y="Frecuencia simulacion distribucion posterior",
       title=paste0("Postre Logístico para Gustavo Petro a 2022-06-10"), 
       subtitle=paste0("20.000 simulaciones del modelo logístico para la 2da vuelta","\nBasado en ",postre2022_disponibles," encuestas"), 
       caption = "Fuente: www.recetas-electorales.com")+
  theme(legend.position = "none",
        panel.background = element_blank())+
  xlim(c(0.35,0.55))

Postre logístico 2022 para Gustavo Petro

Ingredientes

Los únicos ingredientes de esta receta son las encuestas que han salido desde la primera vuelta. Como priors se toman los promedios y desviaciones estándar de cada candidato.

library(tidyverse)
library(lubridate)
library(kableExtra)

# Kable
readr::read_csv("https://raw.githubusercontent.com/nelsonamayad/Elecciones-presidenciales-2022/main/Encuestas%202022/encuestas_2022.csv") %>% 
  dplyr::select(n,fecha,encuestadora,muestra,gustavo_petro,rodolfo_hernandez) %>%
  dplyr::filter(fecha>=lubridate::as_date("2022-05-29")) %>% 
  tidyr::pivot_longer(cols=contains("_"),names_to = "candidato", values_to = "int_voto") %>% 
  dplyr::mutate(nombres = case_when(candidato=="gustavo_petro" ~ "Gustavo Petro",
                                    candidato=="rodolfo_hernandez" ~ "Rodolfo Hernandez")) %>%
  dplyr::group_by(nombres) %>% 
  dplyr::summarize(#promedio_prior = mean(int_voto),
                   promedio_ponderado_t = weighted.mean(int_voto,w=as.numeric(fecha)),
                   sd_prior = sd(int_voto)) %>%  
  kable("html", digits=1,caption = "Priors por candidato",
        col.names = c("Candidato",
                      "$\\mu_{t}$",
                      "$\\sigma_{t}$")) %>% 
  kable_styling(full_width = F) %>% 
  row_spec(0,bold=TRUE, background = "#FF4900", color = "white") %>%
  row_spec(1,bold=TRUE, background = "#800080", color = "white") %>%
  row_spec(2,bold=TRUE, background = "#F7BE0A", color = "black") %>%
  footnote(number = c(paste0("Fecha: 2022-05-21"),
                      paste0("Encuestas disponibles: ",postre2022_disponibles)))

Priors por candidato

Candidato	\(\mu_{t}\)	\(\sigma_{t}\)
Gustavo Petro	45.0	2.5
Rodolfo Hernandez	46.6	5.9
1 Fecha: 2022-05-21
2 Encuestas disponibles: 13

Postre 2022: Priors

Receta

El Postre modela directamente las encuestas como una distribución binomial donde cada encuesta es un ensayo:

\[ N^{Petro}_i \sim \textrm{Binomial}(Muestra_i, \pi^{Petro}_i) \] donde \(N_i\) es la muestra de cada encuesta y \(\pi^{Petro}_i\) es la proporción de la intención de voto para Petro. La proporción de intención de voto \(\pi^{Petro}_i\) se determina a través de una función logística (log-link function) y, como es costumbre en este recetario, efectos aleatorios por encuestadora, la duración entre la encuesta y la estimación y una dummy para el tipo de encuestas (i.e. telefónica, presencial) así:

\[\textrm{logit} (\pi^{Petro}_i) = \alpha_{encuestadora} + \beta_1*duracion + \beta_2*tipo\]

Preparación

Este es el postre completo, con los priors del caso:
\[ N^{Petro}_i \sim \textrm{Binomial}(N_i, \pi^{Petro}_i) \] \[\textrm{logit} (\pi^{Petro}_i) = \alpha_{encuestadora} + \beta_1*duracion + \beta_2*tipo\] \[\small\alpha_{encuestadora} \sim Normal(\alpha,\sigma) \] \[\small\alpha \sim Normal(44, 5) \] \[\small\beta \sim Normal(0, 5) \] \[\small\sigma \sim HalfCauchy(0,5) \]

Estimación con ulam

El modelo estimado con ulam de McElreath a continuación:

library(rethinking)

# Modelo logístico en map2stan:
postre2022_ulam <- ulam(
  alist(
    #modelo
    gp ~ dbinom(muestra,p),
    logit(p) <- a[encuestadora] + b1*dd + b2*tipo_1,
    #priors
    a[encuestadora] ~ dnorm(a,s),
    a ~ dnorm(45,5),
    b1 ~ dnorm(0,5),
    b2 ~ dnorm(0,5),
    s ~ dcauchy(0,5)
  ),
  data=list(
    N = postre2022_disponibles,
    encuestadora = postre2022_encuestas$encuestadora,
    muestra = postre2022_encuestas$muestra,
    dd = postre2022_encuestas$dd,
    tipo_1 = postre2022_encuestas$tipo_1,
    gp= postre2022_encuestas$gp
    ),
  control=list(adapt_delta=0.96),
  iter=1e4, 
  warmup=1000, 
  chains=4, 
  cores=4,
  log_lik=TRUE)