BBB 20 - Prevendo a final

Tentando prever a final do BBB 20 com modelo

Olá a todos!

Depois de um longo período conturbado, retorno (ainda não com a frequência que gostaria) o blog. Aproveitando essa quarentena por conta do COVID-19 (que posteriormente analisarei aqui), raspei e analisei alguns dados do BBB 20.

O intuito do post de hoje é utilizar um modelo logístico para determinar a final do reality show, ou seja, os três participantes que sairão por último, independente da ordem. Um modelo de regressão logística, grosso modo, difere-se da regressão linear pois o seu output, ou seja, a variável independente, ficará entre 0 e 1, sendo utilizado então para variáveis categóricas, em geral binárias. Sem nos ater diretamente à parte estatística, deixando isso para outro post, vamos aos dados.

O banco foi raspado da Wikipedia, de todas as páginas de cada edição do programa. Lá, é possível encontrar o nome, nascimento, origem, profissão, quando saiu, qual a colocação no reality, e o número de votos. O banco completo e atualizado a cada semana pode ser baixado aqui.

pacman::p_load(tidyverse, rio, hrbrthemes)
bbb <- import("bbb.csv") %>% 
    mutate(final = as.factor(final),
         vencedor = as.factor(vencedor),
         ano = as.factor(ano))

Pronto, os pacotes tidyverse e rio foram carregados, e o banco também. Para abrir o banco, utilizei a função import() do rio. Também aproveitei e converti a variável binária final e vencedor para fator, assim como o ano. Removi a variável profissão e localidade, pois não serão utilizadas nessa análise.

summary(bbb)
##      nome            nascimento                                 profissao  
##  Length:316         Length:316         São Paulo, São Paulo          : 20  
##  Class :character   Class :character   Modelo                        : 15  
##  Mode  :character   Mode  :character   Rio de Janeiro, Rio de Janeiro: 12  
##                                        Belo Horizonte, Minas Gerais  : 11  
##                                        Estudante                     :  9  
##                                        Goiânia, Goiás                :  7  
##                                        (Other)                       :242  
##   localidade         resultado         n_votos_total     final     vencedor  
##  Length:316         Length:316         Min.   : 0.000   0   :247   0   :285  
##  Class :character   Class :character   1st Qu.: 3.000   1   : 57   1   : 19  
##  Mode  :character   Mode  :character   Median : 6.000   NA's: 12   NA's: 12  
##                                        Mean   : 7.595                        
##                                        3rd Qu.:11.000                        
##                                        Max.   :31.000                        
##                                                                              
##       ano     
##  2002   : 24  
##  2014   : 20  
##  2020   : 20  
##  2011   : 19  
##  2009   : 18  
##  2010   : 17  
##  (Other):198

Quem foram os 5 mais votados de todas as edições?

bbb %>% 
  top_n(5, n_votos_total) %>% 
  ggplot(aes(y = n_votos_total, x = reorder(nome, n_votos_total)))+
  geom_col(aes(fill = nome), show.legend = F)+
  coord_flip()+
  theme_ipsum_tw()+
  labs(title = "Mais votados do BBB",
       subtitle = "Atualizado em 24/03/2020",
       x = "Participante",
       y = "Nº de votos")+
  paletteer::scale_fill_paletteer_d("vapoRwave::floralShoppe")

Quantos votos os vencedores tiveram no total, dentro da casa?

bbb %>% 
  filter(vencedor == 1) %>% 
  ggplot(aes(x = reorder(nome, as.integer(ano)), y = n_votos_total, fill = nome))+
  geom_col(show.legend = FALSE)+
  theme_ipsum_tw()+
  coord_flip()+
  labs(title = "Vencedores x Votação",
       subtitle = "Atualizado em 24/03/2020",
       x = "Vencedores",
       y = "Nº de votos")+
    paletteer::scale_fill_paletteer_d("dichromat::SteppedSequential.5")

Podemos perceber que muitos vencedores receberam uma quantidade razoável de votos. Vamos olhar a distribuição.

bbb %>% 
  filter(vencedor == 1) %>% 
  ggplot(aes(x = n_votos_total))+
  geom_density(fill = "blue4", alpha = 0.5, size = 0.3)+
  theme_ipsum_tw()+
  labs(x = "Nº de votos", y = "Densidade",
       title = "Densidade da distribuição dos votos dos vencedores",
       subtitle = "Atualizado em 24/03/2020")

bbb %>% 
  filter(final == 1) %>% 
  ggplot(aes(x = n_votos_total))+
  geom_density(fill = "blue4", alpha = 0.5, size = 0.3)+
  theme_ipsum_tw()+
  labs(x = "Nº de votos", y = "Densidade",
       title = "Densidade da distribuição dos votos dos finalistas",
       subtitle = "Atualizado em 24/03/2020")

Ou seja, observamos que há normalidade na curva de distribuição dos votos dos vencedores, mas a distribuição dos votos dos finalistas é bi-modal: há muitos finalistas com poucos votos, queridos desde sempre, e alguns com muitos votos.

Vamos pensar agora nos modelos. Façamos dois modelos: um que observe a probabilidade de ser finalista, com base no número de votos recebidos na edição, e outro que analise a probabilidade de ser vencedor. Excluo de ambos os participantes que tiveram 0 votos: tais participantes ou desistiram ou foram para algum paredão no início por conta alguma prova especial, sem ser indicado por ninguém.

O modelo se dará em função de ir a final ~ número de votos recebidos ao longo da edição e vencer ~ número de votos recebidos ao longo da edição. Isso se justifica pelo fato de haver o pensamento de que um participante muito votado (e que voltou dos paredões) tem mais chance de chegar à final e vencer o programa. Um grande exemplo disso foi o participante Marcelo Dourado.

Para fazer o modelo, vamos utilizar a função glm(), que representa modelos lineares generalizados, o que inclui o logit, o probit, e toda a família.

mod_finalistas <- glm(final ~ n_votos_total,
            family = "binomial",
            data = bbb %>%
              filter(n_votos_total > 0))

mod_vencedor <- glm(vencedor ~ n_votos_total,
            family = "binomial",
            data = bbb %>%
              filter(n_votos_total > 0))
summary(mod_finalistas)
## 
## Call:
## glm(formula = final ~ n_votos_total, family = "binomial", data = bbb %>% 
##     filter(n_votos_total > 0))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.3191  -0.6569  -0.5221  -0.4119   2.2402  
## 
## Coefficients:
##               Estimate Std. Error z value             Pr(>|z|)    
## (Intercept)   -2.54947    0.29750  -8.570 < 0.0000000000000002 ***
## n_votos_total  0.12507    0.02576   4.856            0.0000012 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 290.46  on 296  degrees of freedom
## Residual deviance: 265.80  on 295  degrees of freedom
##   (9 observations deleted due to missingness)
## AIC: 269.8
## 
## Number of Fisher Scoring iterations: 4
summary(mod_vencedor)
## 
## Call:
## glm(formula = vencedor ~ n_votos_total, family = "binomial", 
##     data = bbb %>% filter(n_votos_total > 0))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.0625  -0.3386  -0.2304  -0.1721   2.8378  
## 
## Coefficients:
##               Estimate Std. Error z value            Pr(>|z|)    
## (Intercept)   -4.79427    0.60681  -7.901 0.00000000000000277 ***
## n_votos_total  0.19643    0.04042   4.860 0.00000117302755526 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 141.23  on 296  degrees of freedom
## Residual deviance: 114.78  on 295  degrees of freedom
##   (9 observations deleted due to missingness)
## AIC: 118.78
## 
## Number of Fisher Scoring iterations: 6

O que podemos observar em ambos é que há alta significância estatística, ou seja, a princípio, faz algum sentido o que queremos verificar… Todavia, ao analisarmos os coeficientes, não podemos fazê-lo diretamente, como em uma regressão linear. Isso se dá por conta da função de ligação do logit, que é diferente. Logo, ao observarmos o coeficiente de n_votos_total no modelo de vencedor, cujo valor é 0,19, não podemos interpretar como “o aumento de 1 voto no participante significa um aumento de 19% de chances dele vencer o jogo” pois a translação da função de ligação não é direta. Uma solução para que possamos interpretar diretamente é passar para razão de chances (odds-ratio) e então para probabilidade. Isso pode ser feito da seguinte maneira:

Se a função de ligação do logit utiliza log, ao fazermos o inverso, ou seja, exponenciar, obtemos a razão de chances. Logo, para isso, podemos fazer:

# Odds-Ratio de mod_finalistas
exp(coef(mod_finalistas))
##   (Intercept) n_votos_total 
##    0.07812338    1.13322368
# Odds-ratio de mod_vencedor
exp(coef(mod_vencedor))
##   (Intercept) n_votos_total 
##   0.008277014   1.217046189

Todavia, razão de chances não é algo muito utilizado no Brasil, sendo mais verificado na estatística inglesa e americana. Podemos transformar de odds-ratio para probabilidade da seguinte forma:

\[Probabilidade = \frac{\text{Razão de Chances}}{1 + \text{Razão de Chances}}\]

E para isso, podemos construir uma função no R chamada ConvLogit(), que pode ser redigida da seguinte forma:

ConvLogit <- function(modelo){
  oddsratio <- exp(modelo)
  probabilidade <- oddsratio / (1 + oddsratio)
  return(probabilidade)
}

E então podemos utilizar essa função para testar, podendo seus resultados serem interpretados diretamente:

ConvLogit(coef(mod_finalistas))
##   (Intercept) n_votos_total 
##    0.07246237    0.53122591
ConvLogit(coef(mod_vencedor))
##   (Intercept) n_votos_total 
##   0.008209068   0.548949406

Agora, para testar, vamos analisar a probabilidade tanto de vencerem quanto de chegarem à final, dos participantes do atual BBB 20. Para isso, vamos filtrar pelo ano, e alterar as variáveis VENCEDOR e FINAL de NA para 1,

bbb20 <- bbb %>% 
  filter(ano == 2020) %>% 
  mutate(final = ifelse(is.na(final), 1, 0),
         vencedor = ifelse(is.na(vencedor), 1, 0))

bbb20 %>% 
  select(nome, n_votos_total, final, vencedor) %>% 
  mutate(prob_final = ConvLogit(predict(mod_finalistas, bbb20)),
         prob_vitoria = ConvLogit(predict(mod_vencedor, bbb20))) %>% 
  filter(final == 1) %>% 
  select(-final, -vencedor) -> bbb20_probs

bbb20_probs %>% 
  kable(col.names = c("Nome", "Nº total de votos", "Probabilidade de Final", "Probabilidade de Vitória")) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = F, position = "center") 
Nome Nº total de votos Probabilidade de Final Probabilidade de Vitória
Alexandre da Silva Santana (Babu) 26 0.6686887 0.5775742
Daniel Xavier Lenhardt 2 0.0911782 0.0121114
Felipe Antoniazzi Prior 31 0.7904421 0.7849850
Flayslane Raiane Pereira da Silva 14 0.3103387 0.1146298
Gabriela Piedade Martins (Gabi) 0 0.0724624 0.0082091
Gizelly Bicalho Abreu 5 0.1274013 0.0216230
Ivy Moraes Barbosa 3 0.1020853 0.0147015
Manoela Latini Gavassi Francisco (Manu) 3 0.1020853 0.0147015
Marcela Olmedo McGowan 0 0.0724624 0.0082091
Mariana Decânio Gonzalez (Mari) 8 0.1752402 0.0383145
Rafaella Freitas Ferreira de Castro (Rafa Kalimann) 0 0.0724624 0.0082091
Thelma Regina Maria dos Santos Assis 3 0.1020853 0.0147015

Lembrando que a probabilidade nesse caso não é cumulativa, logo, não poderão somar 100%, pois se referem à chances individuais.

Analisando os 3 participantes com mais chance de chegarem à final do BBB 20:

bbb20_probs %>% 
  select(-n_votos_total) %>% 
  top_n(3, prob_final) %>% 
  arrange(-prob_final) %>% 
  kable(col.names = c("Nome", "Probabilidade de Final", "Probabilidade de Vitória")) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), 
                full_width = F, position = "center") 
Nome Probabilidade de Final Probabilidade de Vitória
Felipe Antoniazzi Prior 0.7904421 0.7849850
Alexandre da Silva Santana (Babu) 0.6686887 0.5775742
Flayslane Raiane Pereira da Silva 0.3103387 0.1146298

O que percebemos é que Prior, Babu e Flayslane são os mais cotados para os finalistas. Todavia, algo pode acontecer e os três irem juntos para o paredão, o que já alteraria todo o conjunto, removendo algum deles da análise.

Continuarei atualizando o banco com o passar das semanas, recalculando as probabilidades para, ao final, verificarmos a capacidade do modelo. Colocarei abaixo o palpite alterado a cada semana:

24/03/2020 - Final: Prior, Babu, Flayslane

Atualização: Como visto, errei na previsão. Todavia, Babu chegou muito perto da final, o que me faria considerar um erro muito, muito pequeno. O que observamos é que em modelos, levamos em consideração variáveis que são correlacionadas ao fato que queremos prever, mas tais modelos só ajudam a diminuir a incerteza: tudo pode acontecer no meio do caminho, e o modelo, com sua margem de erro, se adequa a isso. Optei por não ir atualizando semanalmente, tanto pelo trabalho quanto pela vontade de manter a mesma previsão até o final sem incorporar viés. De qualquer maneira, no âmbito pessoal, apesar da saída de Babu ter sido triste para mim, a vitória da Thelma foi excelente. =)

Avatar
Mateus C. Pestana
Doutorando e Mestre em Ciência Política

Interessado em ciência de dados, ciência política, política russa e métodos de impressão 3d.

Relacionados