Carga de datos

Lo primero que tenemos que hacer para comenzar a analizar el dataset de datos 2015 de los equipos de futbol, es cargarlos.

# Carga de datos
datos_2015 <- read.csv2("datos_2015.csv", encoding = "latin1")
datos_2015$valor <- as.numeric(as.character(datos_2015$valor))

Análisis exploratorio de los datos 2015

Al cargar la base de datos vemos que tiene una dimensión de 30 observaciones o filas y de 49 variables o columnas. Las columas tiene los siguientes nombres:

names(datos_2015)
##  [1] "club"               "Pts_2012_13"        "Pts_2013_14"       
##  [4] "Pts_2014"           "inversion_abs"      "inversion_relativa"
##  [7] "valor"              "libertadores"       "sudamericana"      
## [10] "ascenso"            "Pos"                "Pts"               
## [13] "PJ"                 "PG"                 "PE"                
## [16] "PP"                 "GF"                 "GC"                
## [19] "Dif"                "Pos_Fecha_1"        "Pos_Fecha_2"       
## [22] "Pos_Fecha_3"        "Pos_Fecha_4"        "Pos_Fecha_5"       
## [25] "Pos_Fecha_6"        "Pos_Fecha_7"        "Pos_Fecha_8"       
## [28] "Pos_Fecha_9"        "Pos_Fecha_10"       "Pos_Fecha_11"      
## [31] "Pos_Fecha_12"       "Pos_Fecha_13"       "Pos_Fecha_14"      
## [34] "Pos_Fecha_15"       "Pos_Fecha_16"       "Pos_Fecha_17"      
## [37] "Pos_Fecha_18"       "Pos_Fecha_19"       "Pos_Fecha_20"      
## [40] "Pos_Fecha_21"       "Pos_Fecha_22"       "Pos_Fecha_23"      
## [43] "Pos_Fecha_24"       "Pos_Fecha_25"       "Pos_Fecha_26"      
## [46] "Pos_Fecha_27"       "Pos_Fecha_28"       "Pos_Fecha_29"      
## [49] "Pos_Fecha_30"

Así vemos que la base de datos tiene el nombre del club de futbol correspondiente en la variable club, y luego los puntos que obtuvieron esos clubes en los torneos 2012/13, 2013/14 y 2014 (cuando no participó de ellos se indica 0 puntos), luego las variables relacionadas con la inversión del club, absoluta y relativa, y el valor del equipo. Luego la variable libertadores y sudamericana sobre si jugaron esas copas durante 2015, la variable ascenso sobre si esos equipos ascendieron en el 2014 (y el del 2015 fue su primer torneo en la primera); luego la Posición final en el torneo 2015 (Pos), con la cantidad de Puntos (Pts); Partidos Jugados (PJ), ganados (PG), empatados (PE), perdidos (PP), goles a favor (GF), goles en contra (GC), y diferencia de goles (Dif).

Luego, 30 columnas con las posiciones de cada equipo en las distintas fechas.

Podemos ver quienes cómo fueron los resultados, en la tabla de posiciones, ordenada por posición final:

datos_ordenados <- datos_2015 %>% 
  arrange(Pos) %>% 
  select(club, Pts, Pos, PG, PE, PP, GF, GC, Dif, libertadores, sudamericana)

knitr::kable(datos_ordenados)  %>% 
  kable_styling(bootstrap_options = "striped", full_width = F)
club Pts Pos PG PE PP GF GC Dif libertadores sudamericana
Boca Juniors 64 1 20 4 6 49 26 23 1 0
San Lorenzo 61 2 18 7 5 44 20 24 1 0
Rosario Central 59 3 16 11 3 47 26 21 0 0
Racing 57 4 16 9 5 40 23 17 1 0
Independiente 54 5 14 12 4 44 22 22 0 1
Belgrano 51 6 14 9 7 33 23 10 0 1
Estudiantes 51 7 14 9 7 34 28 6 1 0
Banfield 50 8 14 8 8 38 32 6 0 0
River Plate 49 9 13 10 7 46 33 13 1 1
Tigre 46 10 12 10 8 32 25 7 0 1
Quilmes 45 11 13 6 11 38 37 1 0 0
Gimnasia y Esgrima 44 12 12 8 10 41 38 3 0 0
Lanús 42 13 10 12 8 33 29 4 0 1
Unión 41 14 9 14 7 38 37 1 0 0
Aldosivi 40 15 11 7 12 37 40 -3 0 0
Newell’s Old Boys 40 16 10 10 10 27 30 -3 0 0
San Martín (SJ) 37 17 8 13 9 32 34 -2 0 0
Olimpo 36 18 8 12 10 23 26 -3 0 0
Colón 34 19 7 13 10 26 31 -5 0 0
Argentinos Juniors 33 20 8 9 13 30 38 -8 0 0
Defensa y Justicia 32 21 8 8 14 27 31 -4 0 0
Godoy Cruz 32 22 8 8 14 32 40 -8 0 0
Huracán 30 23 6 12 12 29 37 -8 1 1
Sarmiento 30 24 7 9 14 24 34 -10 0 0
Temperley 30 25 6 12 12 19 29 -10 0 0
Nueva Chicago 29 26 7 8 15 29 38 -9 0 0
Vélez 29 27 7 8 15 27 37 -10 0 0
Arsenal 27 28 7 6 17 25 44 -19 0 1
Atlético de Rafaela 23 29 4 11 15 29 51 -22 0 0
Crucero del Norte 14 30 3 5 22 21 55 -34 0 0

De esta tabla se desprende que Boca Juniors quedó en el primer lugar, y Crucero del Norte en el último.

Podemos analizar las trayectorias de los primeros y de los últimos cinco a lo largo de las distintas fechas.

xfecha <- datos_2015 %>% 
  arrange(Pos) %>% 
  slice(1:5) %>% 
  select(club, Pos_Fecha_1:Pos_Fecha_30) %>% 
  gather(key = fecha, value = "Posicion", Pos_Fecha_1:Pos_Fecha_30)
xfecha <- xfecha %>% 
  mutate(fecha =  as.numeric(str_replace(fecha, pattern =  "Pos_Fecha_", replacement = "")))

ggplot(data = xfecha) +
  geom_path(aes(y = Posicion, x = fecha, group = club, color = club), size = 1.2)+
  labs(title = "Trayectoria de los primeros 5 equipos", 
       subtitle = "Torneo 2015", x = "Fecha", y = "Posición")+
  scale_y_reverse()+
  theme_bw()

Y la trayectoria de los cinco últimos

xfecha2 <- datos_2015 %>% 
  arrange(desc(Pos)) %>% 
  slice(1:5) %>% 
  select(club, Pos_Fecha_1:Pos_Fecha_30) %>% 
  gather(key = fecha, value = "Posicion", Pos_Fecha_1:Pos_Fecha_30)
xfecha2 <- xfecha2 %>% 
  mutate(fecha =  as.numeric(str_replace(fecha, pattern =  "Pos_Fecha_", replacement = "")))

ggplot(data = xfecha2) +
  geom_path(aes(y = Posicion, x = fecha, group = club, color = club), size = 1.2)+
  labs(title = "Trayectoria de los últimos 5 equipos", 
       subtitle = "Torneo 2015", x = "Fecha", y = "Posición")+
  scale_y_reverse()+
  theme_bw()

Con respecto a los puntos de los distintos equipos, podemos ver el siguiente análisis

summary(datos_2015$Pts)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   14.00   30.50   40.00   40.33   49.75   64.00
sd(datos_2015$Pts)
## [1] 12.14661

Así, el promedio de puntos en el torneo fue de 40.33 y el desvío estandar de 12.15. Y podemos hacer un gráfico. Parece una distribución normal, pero cuando ampliamos la cantidad de breaks no es tan acampanadada…

par(mfrow = c(1,2))
hist(datos_2015$Pts, breaks = 5, xlab = "Puntos", main = "Histograma breaks = 5",
     col = "lightblue")
hist(datos_2015$Pts, breaks = 10, xlab = "Puntos", main = "Histograma breaks = 10",
     col = "lightblue")

También podemos hacer un plot con la dispersión de puntos del promedio (marcamos también los dos desvíos estandar del promedio)

ggplot(data = datos_ordenados)+
  geom_point (aes(x = Pts, y = reorder(club, Pts)))+
  geom_vline(xintercept =  mean(datos_2015$Pts), color = "red")+
  geom_vline(xintercept =  c(mean(datos_2015$Pts)+ 2*sd(datos_2015$Pts),
                             mean(datos_2015$Pts)- 2*sd(datos_2015$Pts)),
             color = "red", linetype = 2)+
   geom_vline(xintercept =  c(mean(datos_2015$Pts)+ 1*sd(datos_2015$Pts),
                             mean(datos_2015$Pts)- 1*sd(datos_2015$Pts)),
             color = "red", linetype = 3)+
  labs(title = "Puntajes del torneo 2015",
       subtitle = "Promedio y 1 y 2 desvío estándares", x = "Puntos", y = "Club")+
  theme_bw()

Regresión lineal

Usar una reg. lineal para explicar la variable Pts (puntos) a partir de las variables Pts_2012_13, Pts_2013_14,Pts_2014, inversion_abs (inversión absoluta), inversion_relativa (inversión relativa), valor (valor en millones deeuros), libertadores (si participó en el 2015 en la Copa Libertadores), sudamericana (si participó en el 2015 en laCopa Sudamericana), ascenso (si ascendió).

Comenzaré por regresiones lineales simples. Primero voy a realizar tres regresiones para analizar si existe relación entre los puntos obtenidos en el torneo 2015, y los puntos obtenidos por el mismo equipo en los torneos 2012_13, 2013_14, y 2014. Para ello, deberé primero filtrar de la regresión los equipos que no jugaron esos torneos, por lo que figuran con 0 puntos (por ello los convierto a NA).

datos_regresiones <- datos_2015 %>% 
  select(club, Pts, Pts_2012_13, Pts_2013_14, Pts_2014) %>% 
  na_if(0)

regresion_2012 <- lm(Pts ~ Pts_2012_13, data = datos_regresiones)
regresion_2013 <- lm(Pts ~ Pts_2013_14, data = datos_regresiones)
regresion_2014 <- lm(Pts ~ Pts_2014, data = datos_regresiones)

par(mfrow = c(2,2))
plot(y = datos_regresiones$Pts, x = datos_regresiones$Pts_2012_13,
     ylab = "Torneo 2015", xlab = "Torneo 2012/2013")
  abline(regresion_2012, lwd = 3, col = "red")
plot(y = datos_regresiones$Pts, x = datos_regresiones$Pts_2013_14,
     ylab = "Torneo 2015", xlab = "Torneo 2013/2014")
  abline(regresion_2013, lwd = 3, col = "red")
plot(y = datos_regresiones$Pts, x = datos_regresiones$Pts_2014,
     ylab = "Torneo 2015", xlab = "Torneo 2014")
  abline(regresion_2014, lwd = 3, col = "red")

Luego de ver que estos gráficos señalen que probablemente no exista correlación, armo una regresión lineal múltiple, y para ver si el modelo es posible, grafico los residudos

regresion_multiple <- lm(Pts ~ Pts_2012_13 + Pts_2013_14 + Pts_2014, 
                         data = datos_regresiones)

tab_model(regresion_multiple, 
          show.se = TRUE,  show.df = TRUE,  show.stat = TRUE,  show.ci = FALSE,  show.fstat = TRUE)
  Pts
Predictors Estimates std. Error Statistic p
(Intercept) 27.20 31.46 0.86 0.408
Pts_2012_13 -0.11 0.38 -0.30 0.774
Pts_2013_14 0.07 0.47 0.15 0.881
Pts_2014 0.70 0.54 1.30 0.223
Observations 14
R2 / adjusted R2 0.147 / -0.109 
par(mfrow = c(2,2))
plot(regresion_multiple)

Lamentablemente, estos diagnosticos nos indican que este modelo lineal no es satisfactorio. Especialmente la distribución de los residuos en el plot 1, donde la distribución no es regular al rededor de la línea, indicándonos que probablemente no sea una regresióin lineal, o que este modelo no encaja. Lo mismo sucede con el estadístico \(R^2\), que al ser muy bajo (0.1468036) nos señala que el modelo tiene problemas.

Otras variables 1. Show me the money

Podemos pensar otra regresión con otras variables. Por ejemplo, ¿existe alguna relación en el valor del equipo y los puntos que obtiene? Es decir, ¿la inversión en jugadores de futbol repercute en cuántos puntos obtiene el club en el torneo? Analicemos eso.

regresion_valor <- lm(Pts ~ valor, data = datos_2015)
plot(y = datos_2015$Pts, x = datos_2015$valor, 
     ylab = "Puntos torneo 2015", xlab = "Valor del club")
abline(regresion_valor, lwd = 3, col = "red")

tab_model(regresion_valor, 
          show.se = TRUE,  show.df = TRUE,  show.stat = TRUE,  show.ci = FALSE,  show.fstat = TRUE)
  Pts
Predictors Estimates std. Error Statistic p
(Intercept) 29.66 3.12 9.49 <0.001
valor 0.71 0.17 4.15 <0.001
Observations 30
R2 / adjusted R2 0.381 / 0.359

Y analizamos los residuos

par(mfrow = c(2,2))
plot(regresion_valor)

Si bien existe correlación significativa, los residuos nos señalan que el modelo no es perfecto. Podríamos intentar algunas tranformaciones a ver si mejora. El valor, que está en dinero podemos transformarlo con el logaritmo. Y los puntos los podemos transformar con raíz cuadrada.

datos_transf <- data.frame(datos_2015$Pts, datos_2015$valor,
                           "sqrt.Pts" = sqrt(datos_2015$Pts), 
                           "log.valor" = log10(datos_2015$valor)) %>% 
  arrange(desc(datos_2015.Pts))

kable(datos_transf) %>% 
  kable_styling(bootstrap_options = "striped", full_width = F)
datos_2015.Pts datos_2015.valor sqrt.Pts log.valor
64 38.83 8.000000 1.5891674
61 37.85 7.810250 1.5780659
59 16.33 7.681146 1.2129862
57 19.40 7.549834 1.2878017
54 16.90 7.348469 1.2278867
51 20.15 7.141428 1.3042751
51 8.78 7.141428 0.9434945
50 16.93 7.071068 1.2286570
49 43.40 7.000000 1.6374897
46 13.95 6.782330 1.1445742
45 13.23 6.708204 1.1215598
44 10.33 6.633250 1.0141003
42 18.32 6.480741 1.2629255
41 6.58 6.403124 0.8182259
40 6.73 6.324555 0.8280151
40 32.56 6.324555 1.5126844
37 6.83 6.082763 0.8344207
36 10.00 6.000000 1.0000000
34 8.55 5.830952 0.9319661
33 14.48 5.744563 1.1607686
32 13.08 5.656854 1.1166077
32 6.90 5.656854 0.8388491
30 9.63 5.477226 0.9836263
30 6.20 5.477226 0.7923917
30 3.38 5.477226 0.5289167
29 3.75 5.385165 0.5740313
29 20.00 5.385165 1.3010300
27 14.28 5.196152 1.1547282
23 8.58 4.795832 0.9334873
14 4.43 3.741657 0.6464037

Y hago un gráfico, a ver si mejoró y ploteo los graficos de diagnostico.

regresion_transf <- lm(sqrt.Pts ~ log.valor, data = datos_transf)
plot(y = datos_transf$sqrt.Pts, x = datos_transf$log.valor, 
     ylab = "Raíz cuadrada de Puntos", xlab = "Logaritmo de valor")
abline(regresion_transf, lwd = 3, col = "red")

tab_model(regresion_transf, 
          show.se = TRUE,  show.df = TRUE,  show.stat = TRUE,  show.ci = FALSE,  show.fstat = TRUE)
  sqrt.Pts
Predictors Estimates std. Error Statistic p
(Intercept) 3.89 0.55 7.10 <0.001
log.valor 2.20 0.49 4.51 <0.001
Observations 30
R2 / adjusted R2 0.421 / 0.400 
#residudos
par(mfrow = c(2,2))
plot(regresion_transf)

Este modelo encaja mucho mejor!

Intervalo de confianza

Luego de las transformaciones, vimos que nuestra regresión es de un modelo lineal, con la siguiente fórmula.

\(y = \beta{_0} + \beta{_1}x\)

\(\sqrt{(Pts)} = intercept + pendiente*\log({valor})\)

\(\sqrt{(Pts)} = 3.88 + 2.22* \log(valor)\)

Y en base a esto, para calcular el intervalo de confianza podemos utilizar la técnica de Boostrap no paramétrico. Para ello, haremos una tabla con 10.000 replicaciones de \(\beta{_0}\) y de \(\beta{_1}\).

#funcion que reordena y recalcula la regresion
set.seed(999)
regresiones <- function(X,Y){
  boot <- sample.int(n = nrow(datos_transf), size = nrow(datos_transf), replace = TRUE)
  boot_fit <- lm(sqrt.Pts ~ log.valor, data=datos_transf[boot,])
  return(coef(boot_fit))
  }

# matriz de 10000 replicaciones de los coeficientes de la regresion
boot <- t(replicate(10000,regresiones(X,Y))) 
head(boot, 10) %>% 
  kable() %>% 
  kable_styling(bootstrap_options = "striped", full_width = F)
(Intercept) log.valor
3.724547 2.469323
4.185238 1.989383
4.045299 2.048250
4.070316 2.117026
4.079419 1.906992
3.921058 2.149061
4.475351 1.748147
3.705875 2.406653
4.439895 1.756008
3.405725 2.703779

El promedio de la intercepción nos da 3.892, y de la pendiente 2.201. Esto es bastante similar de nuestros valores originales ( que eran 3.888, 2.204).

Y podemos graficarlo

par(mfrow=c(1,2))
MASS::truehist(boot[,1], main = "Boostrap del intercept", xlab = "intercept")
abline(v =median(boot[,1]), col = "red", lty = 2)
MASS::truehist(boot[,2], main = "Boostrap de la pendiente", xlab = "pendiente")
abline(v = median(boot[,2]), col = "red", lty = 2)

Luego, un gráfico de la regresión

plot(x = boot[,2], y = boot[,1], pch=19, cex=0.1, 
     ylab = "Intercept", xlab = "pendiente",
     main = "Scaterplot del Boostrap de la regresión")

Y calcular el Intervalo de confianza para cada parámetro debemos analizar los quartiles de los parametros.

# Intervalos de Confianza para beta_0 

beta0 <- quantile(boot[,1], p= c(.025,.25,.5,.75,.975))
beta1 <- quantile(boot[,2], p= c(.025,.25,.5,.75,.975))

kable(list(beta0, beta1)) %>% 
  kable_styling(bootstrap_options = "striped", full_width = F)
x
2.5% 2.755488
25% 3.544417
50% 3.925619
75% 4.269352
97.5% 4.864679
x
2.5% 1.307157
25% 1.868165
50% 2.183405
75% 2.509989
97.5% 3.194552

De este modo, el intervalo del 95% de confianza para la intercepción va de 2.7554875 a 4.864679; y el de la pendiente va de 1.3071572 a 3.1945521.

Test de hipótesis

Para el modelo elegido, ¿se puede decir que existe algún coeficiente de las variables explicativas estadísticamente distinto de cero? Pruebe el test de hipótesis (Busque para contextualizar en Agresti, el test F, capítulo 11, en particular 11.4 -11.5)

???