Modelizando resultados de partidos de fútbol: Regresión de Poisson

Este post está plenamente inspirado en este artículo.

En los últimos días he estado intentando leer información sobre aplicaciones de estadística y machine learning en el campo de los deportes, concretamente en el fútbol. De mis primeras impresiones del llamado soccer/football analytics, extraigo algunas conclusiones preliminares:

  1. Hay mucho camino recorrido, pero queda por recorrer: es relativamente fácil encontrar todo tipo de papers y artículos relacionados con todos los subámbitos del análisis en el fútbol: predicción de resultados, rendimiento de jugadores, scouting, la familia xG y modelos de predicción de éstos, etc, aunque no he encontrado manuales/libros de referencia que cubran todo lo anterior en detalle.

  2. Todo pasa por pagar: los datos abiertos en el mundo del fútbol no están a la orden del día. Hay cosas, sí, pero difícilmente tratables en un entorno online. Las principales plataformas son de pago. Desde luego, hay gente haciendo esfuerzos en democratizar este campo, pero en algún momento hay que pasar por caja (o buscar alternativas mediante web scrapping, lo cual no he indagado).

En este post lo que voy a hacer es utilizar uno de los modelos más básicos de predicción de resultados, la regresión de Poisson. Específicamente, lo que se modeliza es el número de goles que marcarán los equipos en un partido en concreto. La regresión de Poisson se encuadra dentro de los modelos lineales generalizados y se utiliza para modelizar datos de conteo. Cuando lo estudié en la carrera recuerdo que comentaban que una de sus primeras aplicaciones prácticas fue a finales del siglo XIX para estudiar el número de soldados prusianos muertos por coz de mulos. De esas cosas de las que uno no se olvida 😄.

El porqué de utilizar la regresión de Poisson viene dado porque la distribución de Poisson parece muy conveniente para entender cómo se distribuyen los goles en un partido. La distribución de Poisson es una distribución de probabilidad discreta que describe la probabilidad de ocurrencia de un número de eventos independientes en un periodo de tiempo. Resalto lo de independiente porque precisamente esta es una de las debilidades que le encuentro a este modelo, ya que se antoja difícil pensar que pueda existir independencia entre goles en un partido concreto. Este modelo solo cuenta con un parámetro, \(\lambda\), que corresponde a la media y varianza de la distribución. Y este es otro potencial problema del modelo de Poisson: supongamos que la media de la distribución es 1.5 goles. Este modelo asume que la varianza también es 1.5, lo cual es claramente limitante. Para corregir este hecho, una opción es optar por un modelo quasi-poisson que permite el ajuste de un parámetro de dispersión. Otra opción es ir a por otra distribución, como Binomial Negativa o Weibull.

Pasemos al código. En primer lugar, vamos a cargar las librerías necesarias:

pkgs <- c("engsoccerdata", "ggplot2", "purrr", "tidyr", "MASS", "magrittr", "AER", "dplyr")
invisible(lapply(pkgs, require, character.only = TRUE))

Ahora vamos a leer, a través de la librería engsoccerdata, los datos de los partidos de la temporada 2016-2017 de la Primera División en España:

df <- engsoccerdata::spain %>% 
  as_tibble() %>% 
  dplyr::filter(Season == 2016) %>% 
  select(home, visitor, hgoal, vgoal)

head(df)
## # A tibble: 6 x 4
##   home                visitor            hgoal vgoal
##   <chr>               <chr>              <int> <int>
## 1 Deportivo La Coruna SD Eibar               2     1
## 2 Malaga CF           CA Osasuna             1     1
## 3 FC Barcelona        Real Betis             6     2
## 4 Granada CF          Villarreal CF          1     1
## 5 Sevilla FC          Espanyol Barcelona     6     4
## 6 Atletico Madrid     CD Alaves              1     1

Una cosa que podemos hacer es ver el número de goles medio por localización (local o visitante). La ventaja de jugar de local que los aficionados al fútbol dan por sentada (y que los datos indican corroborar) ha sido bien estudiada (y se puede leer sobre ella aquí).

df %>% 
  summarise(`Goles local` = mean(hgoal), 
            `Goles visitante` = mean(vgoal)) %>% 
  pivot_longer(cols = c('Goles local', 'Goles visitante'), names_to = "tipo", values_to = "goles") %>% 
  ggplot(aes(x = tipo, y = goles)) + geom_bar(stat="identity", width=0.5) + 
  labs(title = "Media de goles por partido", y = "Goles", x = "Tipo") +
  theme_minimal()

Ahora vamos a construir el dataset que alimenta a la regresión de Poisson. Lo que haremos será pasar de un formato ancho a un largo, pasando a crear una variable dummy indicando si los goles de la columna goals hacen referencia al equipo local o visitante. Mejor si lo vemos. Por ejemplo, para la primera fila, vemos un partido en el que el Barcelona como visitante (home = 0) marcó 5 goles al Leganés. En la tercera fila vemos el caso en el que el Málaga no marcó ningún gol a la Real Sociedad en su partido como local, etc

df <-
  bind_rows(
    tibble(goals = df$hgoal,
           team = df$home,
           opponent = df$visitor,
           home = 1),
    tibble(goals = df$vgoal,
           team = df$visitor,
           opponent = df$home,
           home = 0)
  )

set.seed(123)
sample_n(df, size = 7, replace = FALSE)
## # A tibble: 7 x 4
##   goals team            opponent        home
##   <int> <chr>           <chr>          <dbl>
## 1     5 FC Barcelona    CD Leganes         0
## 2     0 Sporting Gijon  Granada CF         0
## 3     0 Malaga CF       Real Sociedad      1
## 4     0 Athletic Bilbao Real Betis         0
## 5     2 Athletic Bilbao Sporting Gijon     1
## 6     1 Sporting Gijon  Real Sociedad      1
## 7     3 Valencia CF     Celta Vigo         1

Ahora corremos el modelo con, evidentente, el número de goles como variable dependiente y el resto como independientes:

poisson_model = glm(goals ~ home + team + opponent, family=poisson(link = "log"), data = df)
summary(poisson_model)
## 
## Call:
## glm(formula = goals ~ home + team + opponent, family = poisson(link = "log"), 
##     data = df)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -2.58114  -1.08917  -0.08988   0.52872   2.55228  
## 
## Coefficients:
##                              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                 -0.082140   0.210395  -0.390 0.696233    
## home                         0.262681   0.060331   4.354 1.34e-05 ***
## teamAtletico Madrid          0.263810   0.182247   1.448 0.147746    
## teamCA Osasuna              -0.233764   0.209740  -1.115 0.265046    
## teamCD Alaves               -0.257162   0.208162  -1.235 0.216686    
## teamCD Leganes              -0.376340   0.216171  -1.741 0.081695 .  
## teamCelta Vigo               0.024552   0.194513   0.126 0.899554    
## teamDeportivo La Coruna     -0.192708   0.205461  -0.938 0.348282    
## teamEspanyol Barcelona      -0.072102   0.198386  -0.363 0.716275    
## teamFC Barcelona             0.779852   0.166017   4.697 2.63e-06 ***
## teamGranada CF              -0.533738   0.228715  -2.334 0.019615 *  
## teamMalaga CF               -0.067417   0.198397  -0.340 0.734002    
## teamReal Betis              -0.237627   0.208210  -1.141 0.253752    
## teamReal Madrid              0.693264   0.168459   4.115 3.87e-05 ***
## teamReal Sociedad            0.116893   0.189477   0.617 0.537285    
## teamSD Eibar                 0.062679   0.191852   0.327 0.743893    
## teamSevilla FC               0.270120   0.182872   1.477 0.139649    
## teamSporting Gijon          -0.205913   0.206830  -0.996 0.319461    
## teamUD Las Palmas            0.029351   0.194526   0.151 0.880065    
## teamValencia CF              0.075968   0.191888   0.396 0.692179    
## teamVillarreal CF            0.045875   0.191808   0.239 0.810972    
## opponentAtletico Madrid     -0.450375   0.245761  -1.833 0.066866 .  
## opponentCA Osasuna           0.771875   0.184329   4.187 2.82e-05 ***
## opponentCD Alaves           -0.011112   0.215852  -0.051 0.958942    
## opponentCD Leganes           0.230796   0.203750   1.133 0.257324    
## opponentCelta Vigo           0.474132   0.194532   2.437 0.014798 *  
## opponentDeportivo La Coruna  0.341073   0.199331   1.711 0.087065 .  
## opponentEspanyol Barcelona   0.147401   0.208193   0.708 0.478945    
## opponentFC Barcelona        -0.090313   0.224583  -0.402 0.687583    
## opponentGranada CF           0.625343   0.188476   3.318 0.000907 ***
## opponentMalaga CF            0.242925   0.203781   1.192 0.233225    
## opponentReal Betis           0.387314   0.197390   1.962 0.049742 *  
## opponentReal Madrid          0.002919   0.218612   0.013 0.989347    
## opponentReal Sociedad        0.215221   0.205481   1.047 0.294914    
## opponentSD Eibar             0.173818   0.207266   0.839 0.401680    
## opponentSevilla FC           0.146021   0.209218   0.698 0.485217    
## opponentSporting Gijon       0.506337   0.192946   2.624 0.008684 ** 
## opponentUD Las Palmas        0.544334   0.192000   2.835 0.004582 ** 
## opponentValencia CF          0.417084   0.196818   2.119 0.034079 *  
## opponentVillarreal CF       -0.262377   0.231630  -1.133 0.257323    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 1006.47  on 759  degrees of freedom
## Residual deviance:  760.16  on 720  degrees of freedom
## AIC: 2241.6
## 
## Number of Fisher Scoring iterations: 5

En primer lugar, se puede ver que la estimación de home es de ~0.26, indicando que hay un factor \(\exp(0.26) \approx 1.3\) de marcar por jugar en casa, en términos medios. Sin embargo, mirando las estimaciones de los parámetros team* se puede ver que hay equipos mucho más goleadores que otros (destacan Madrid y Barcelona como goleadores y Granada y Leganés por “no goleadores”). Los parámetros opponent* vienen a indicar la dificultad de hacerle goles a un determinado equipo, en este caso valores negativos indican dificultad para marcarle goles a dicho oponente. Aquí destaca el Atlético de Madrid, y no es sorprendente, fue el equipo que menos goles encajó (27).

Antes de pasar a la segunda parte del análisis, antes se comentó que uno de los problemas de este modelo podría ser la sobredispersión, que no parece ser problema en este caso:

dispersiontest(poisson_model)
## 
##  Overdispersion test
## 
## data:  poisson_model
## z = -3.5546, p-value = 0.9998
## alternative hypothesis: true dispersion is greater than 1
## sample estimates:
## dispersion 
##  0.8516443

Ahora estamos en disposición de crear predicciones para un partido y un equipo concreto. Supongamos que queremos saber cuál es el número esperado de goles de un partido Real Sociedad - Málaga:

home_goals_avg <- predict(poisson_model,
                          data.frame(home=1, team="Real Sociedad", 
                                     opponent="Malaga CF"), type="response")

away_goals_avg <- predict(poisson_model, 
                          data.frame(home=0, team="Malaga CF", 
                                     opponent="Real Sociedad"), type="response")

# Goles esperados 
tibble("Real Sociedad" = home_goals_avg, "Málaga" = away_goals_avg)
## # A tibble: 1 x 2
##   `Real Sociedad` Málaga
##             <dbl>  <dbl>
## 1            1.72   1.07

Ahora veamos cual sería la probabilidad estimada de que los equipos marquen X goles:

tst <- tibble(prob = c(dpois(0:10, home_goals_avg), dpois(0:10, away_goals_avg)), 
              team = c(rep("Real Sociedad", 11), rep("Malaga CF", 11)),
              index = rep(c(0:10),2))

ggplot(data=tst, aes(x=as.factor(index), y=prob, color=team, group=team)) +
  geom_line() + 
  geom_point() +
  labs(title = "Probabilidad estimada de marcar", x = "Goles", color = "Equipo") +
  theme_minimal()

Con esto podemos calcular la probabilidad conjunta de dos eventos. Por ejemplo, la probabilidad de que la Real Sociedad marque un gol y el Málaga cero sería: \(P(G = 1|Real Sociedad) \times P(G = 0|Málaga) = 0.308 \times 0.343 = 0.1056 \approx 10.6\%\). Podemos calcular una serie de resultados con la misma lógica:

tst1 <- tibble(prob = c(dpois(0:5, home_goals_avg)), 
               team = c(rep("Real Sociedad", 6)),
               index = 0:5)

tst2 <- tibble(prob = c(dpois(0:5, away_goals_avg)), 
               team = c(rep("Málaga CF", 6)),
               index = 0:5)

tst3 <- merge(x = tst1, y = tst2, by = NULL) %>% 
  mutate(prob = prob.x*prob.y)


ggplot(tst3, aes(x = as.factor(index.x), y = as.factor(index.y), fill=prob)) +
  geom_tile() +
  theme_minimal() +
  geom_text(aes(label = scales::percent(prob)), color = 'black') +
  scale_fill_distiller(palette = "Blues", direction = +1, labels = scales::percent) +
  labs(title = "Posibles resultados",
       x = "Goles Real Sociedad",
       y = "Goles Málaga")

Para este partido, la mayor probabilidad estimada se da en 1-1, aunque los resultados “colindantes” son muy similares.

Para terminar, vemos cuál sería el resultado estimado en términos de victorias (local o visitante) o empate.

rs_ma <- dpois(0:10, home_goals_avg) %o% dpois(0:10, away_goals_avg)
home_win <- sum(rs_ma[lower.tri(rs_ma)])
away_win <- sum(rs_ma[upper.tri(rs_ma)])
draw <- sum(diag(rs_ma))
 
tibble(`Victoria local` = round(100*sum(rs_ma[lower.tri(rs_ma)]), 2),
       `Victoria visitante` = round(100*sum(rs_ma[upper.tri(rs_ma)]), 2),
       'Empate' = round(100*sum(diag(rs_ma)),2)) %>% 
  pivot_longer(cols = c('Victoria local', 'Victoria visitante', 'Empate'), names_to = "res", values_to = "prob") %>% 
  ggplot(aes(x = res, y = prob)) + geom_bar(stat="identity", width=0.5) + 
  labs(title = "Resultado", y = "Probabilidad estimada (%)", x = "Resultado") +
  theme_minimal()