Me pide mi amigo Jesús Lagos que hagamos un vídeo hablando del análisis de componentes principales para un canal que tiene junto a Miguel Angel.
El caso es que llevo muchos años usándolo y lo estudié en la carrera, haciendo varios a mano, como no podía ser de otra manera, pero desde que empecé a usar software estadístico se me habían olvidado los detalles de la matemática subyacente.
Bien, dejando un poco todo eso, el análisis de componentes principales es una técnica de interdependencia que busca describir la relación entre las variables de un conjunto de datos e incluso obtener una reducción de su dimensión (si nos quedamos con k componentes principales siendo k menor que el número de variables del conjunto de datos).
Las componentes son una combinación lineal de las variables tales que, la primera recoge la dirección de la máxima varianza de los datos, la segunda es aquella que siendo ortogonal a la primera recoge la máxima varianza que no ha recogido la primera y así hasta llegar a tantas componentes como variables originales tienes. Visto de esta forma, no se trata más que de un cambio de base elegida de cierta forma y ahí es dónde interviene la diagonalización de matrices.
En el blog de Joaquín Amat, el cual recomiendo encarecidamente, se cuenta todo esto de manera mucho más clara, aquí os lo dejo
Una explicación más algebraica la podemos encontrar en el blog de Carlos
Las 2 primeras componentes serían algo así como
\[ \begin{aligned} PC1 = \phi_{11}X_1 + \phi_{21}X_2 + \ldots + \phi_{p1}X_p \\ PC2 = \phi_{12}X_1 + \phi_{22}X_2 + \ldots + \phi_{p2}X_p \\ \end{aligned} \]
Y para encontrar los \(\phi_{ij}\) es cuando se recurre a la diagonalización de la matriz de covarianzas.
Veamos un ejemplo simple con datos de fbref que me ha pasado Jesús
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──
## ✓ ggplot2 3.3.2 ✓ purrr 0.3.4
## ✓ tibble 3.0.4 ✓ dplyr 1.0.2
## ✓ tidyr 1.1.2 ✓ stringr 1.4.0
## ✓ readr 1.4.0 ✓ forcats 0.5.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
df <- readRDS("../../data/FBREF_team_kpi_90.rds")
# nos quedamos con temporada 2018/2019 de España
df_2018_2019 <- df %>%
filter(temporada == "2018/2019" & Country == "Spain")
glimpse(df_2018_2019)
## Rows: 20
## Columns: 73
## $ Squad <chr> "Alavés", "Athletic Club", "Atlético…
## $ CODE <chr> "ESP_ALA", "ESP_ATH", "ESP_ATM", "ES…
## $ Country <chr> "Spain", "Spain", "Spain", "Spain", …
## $ temporada <chr> "2018/2019", "2018/2019", "2018/2019…
## $ Creacion_Toques <dbl> 489.0789, 572.1316, 613.1053, 811.23…
## $ Creacion_Toques_AP <dbl> 61.34211, 55.07895, 53.78947, 62.578…
## $ Creacion_Toques_Z1 <dbl> 168.6842, 170.2105, 172.2368, 208.34…
## $ Creacion_Toques_Z2 <dbl> 232.7895, 286.2632, 301.1842, 415.97…
## $ Finalizacion_Toques_Z3 <dbl> 113.5526, 148.9474, 176.7632, 240.81…
## $ Finalizacion_Toques_AR <dbl> 17.68421, 20.42105, 23.71053, 31.500…
## $ Creacion_Conducciones <dbl> 270.6316, 357.8947, 416.1842, 617.02…
## $ Creacion_Conducciones_Dist <dbl> 1563.132, 1974.105, 2238.579, 3127.5…
## $ Creacion_Conduccioones_Dist_Prog <dbl> 866.8684, 1070.6579, 1173.0263, 1759…
## $ Defensa_Entradas <dbl> 16.73684, 16.94737, 20.63158, 15.052…
## $ Defensa_Ent_Z1 <dbl> 8.657895, 8.263158, 10.947368, 6.289…
## $ Defensa_Ent_Z2 <dbl> 6.157895, 6.289474, 7.815789, 6.5263…
## $ Defensa_Ent_Z3 <dbl> 1.921053, 2.394737, 1.868421, 2.2368…
## $ Defensa_Presiones <dbl> 154.3947, 154.2368, 184.0789, 163.81…
## $ Defensa_Presiones_Z1 <dbl> 51.65789, 45.26316, 58.31579, 43.342…
## $ Defensa_Presiones_Z2 <dbl> 72.50000, 73.18421, 86.71053, 77.921…
## $ Defensa_Presiones_Z3 <dbl> 30.23684, 35.78947, 39.05263, 42.552…
## $ Finalizacion_Pases_Clave <dbl> 7.157895, 7.184211, 8.921053, 10.789…
## $ Finalizacion_Pases_Z3 <dbl> 21.86842, 28.94737, 32.76316, 46.815…
## $ Finalizacion_Pases_Area <dbl> 5.947368, 7.947368, 9.000000, 12.605…
## $ Finalizacion_Centros_Area <dbl> 2.736842, 2.500000, 1.605263, 1.5789…
## $ Creacion_Toques_Juego <dbl> 438.6579, 523.5000, 568.7105, 768.36…
## $ Creacion_Regates <dbl> 15.10526, 15.81579, 18.23684, 24.000…
## $ Defensa_Regates_en_contra <dbl> 17.34211, 19.21053, 19.55263, 14.815…
## $ Defensa_Bloqueos <dbl> 15.94737, 15.15789, 15.31579, 15.421…
## $ Defensa_Bloqueos_Disp <dbl> 3.052632, 2.526316, 2.736842, 2.8947…
## $ Defensa_Bloqueos_Pases <dbl> 12.89474, 12.63158, 12.57895, 12.526…
## $ Defensa_Intercepciones <dbl> 10.263158, 11.868421, 12.184211, 10.…
## $ Defensa_Despejes <dbl> 30.21053, 25.94737, 21.84211, 15.368…
## $ Creacion_Pases_vivo <dbl> 314.9211, 407.7105, 455.0789, 658.31…
## $ Creacion_Pases_ABP <dbl> 51.36842, 50.00000, 45.81579, 45.210…
## $ Creacion_Pases_Falta <dbl> 15.39474, 14.92105, 11.92105, 14.947…
## $ Creacion_Pases_Largos_Z1 <dbl> 0.3421053, 0.7631579, 1.5263158, 3.0…
## $ Creacion_Pases_Presion <dbl> 55.68421, 72.84211, 88.23684, 109.15…
## $ Creacion_Pases_Amplitud <dbl> 11.39474, 14.57895, 11.65789, 15.342…
## $ Finalizacion_Pases_Centros <dbl> 19.18421, 19.13158, 15.52632, 10.368…
## $ Finalizacion_Pases_Corner <dbl> 3.815789, 4.052632, 5.657895, 5.3157…
## $ Creacion_Pases_Rasos <dbl> 190.3421, 267.6053, 330.6053, 559.34…
## $ Creacion_Pases_Media_Altura <dbl> 56.68421, 72.31579, 69.07895, 70.078…
## $ Creacion_Pases_Altos <dbl> 119.26316, 117.78947, 101.21053, 74.…
## $ Creacion_Pases_Completados <dbl> 258.2632, 349.1316, 398.5000, 617.47…
## $ Creacion_Pases_FdJ <dbl> 1.263158, 1.789474, 1.815789, 2.4473…
## $ Creacion_Pases_Interceptados <dbl> 6.657895, 7.552632, 9.736842, 12.236…
## $ Creacion_Pases_Bloqueados <dbl> 10.97368, 11.89474, 12.18421, 12.552…
## $ Creacion_Pases_Distancia <dbl> 5355.526, 7134.263, 7263.395, 11080.…
## $ Creacion_Pases_Dist_Prog <dbl> 2307.974, 2639.342, 2592.842, 3306.6…
## $ Creacion_Pases_Cortos <dbl> 13.97368, 14.15789, 16.55263, 19.710…
## $ Creacion_Pases_Cercanos <dbl> 225.4737, 301.6053, 356.4737, 534.10…
## $ Creacion_Pases_Largos <dbl> 126.8421, 141.9474, 127.8684, 149.71…
## $ Finalizacion_Asistencias <dbl> 0.7105263, 0.6842105, 0.9473684, 1.6…
## $ Finalizacion_xA <dbl> 0.6289474, 0.6631579, 0.8236842, 1.4…
## $ Creacion_Pases_Progresivos <dbl> 36.28947, 41.52632, 43.71053, 50.605…
## $ Creacion_Posesion <dbl> 1.123684, 1.318421, 1.300000, 1.6894…
## $ Finalizacion_Goles <dbl> 1.0263158, 1.0000000, 1.3684211, 2.3…
## $ Finalizacion_xG <dbl> 0.9815789, 1.0578947, 1.2078947, 1.9…
## $ Finalizacion_Disparos <dbl> 11.03, 9.68, 11.55, 14.53, 11.58, 11…
## $ Finalizacion_Disparos_Puerta <dbl> 3.11, 3.50, 4.21, 6.39, 4.05, 4.13, …
## $ Finalizacion_Acciones_Disparo <dbl> 15.13, 14.97, 18.53, 24.71, 18.76, 1…
## $ Finalizacion_Acciones_Gol <dbl> 1.53, 1.74, 2.29, 4.13, 2.08, 2.34, …
## $ Defensa_Goles_Contra <dbl> 1.32, 1.18, 0.76, 0.95, 1.37, 1.63, …
## $ Defensa_Disparos_Contra <dbl> 4.184211, 3.052632, 3.421053, 3.3947…
## $ Defensa_Goles_Contra_Corner <dbl> 0.05263158, 0.15789474, 0.07894737, …
## $ Porteria_Saque_Largo <dbl> 22.289474, 23.026316, 16.131579, 7.8…
## $ Porteria_Pases <dbl> 21.65789, 23.92105, 15.31579, 25.078…
## $ Porterioa_Saques_Puerta <dbl> 8.421053, 6.947368, 7.684211, 5.8421…
## $ Porteria_Centros <dbl> 10.710526, 9.447368, 9.605263, 8.394…
## $ Porteria_Salidas_Area <dbl> 0.82, 0.74, 0.32, 0.87, 1.11, 0.76, …
## $ Defensa_Duelos_Aereos <dbl> 37.44737, 34.28947, 26.05263, 16.552…
## $ Defensa_xGA <dbl> 1.336842, 1.152632, 1.018421, 1.0684…
Nos quedamos sólo con 4 variables para el ejemplo, escalamos los datos y obtenemos la matriz de correlaciones
df_numerico <- df_2018_2019 %>%
select(Squad,Creacion_Toques:Creacion_Toques_Z2)
df_centrado <- scale(df_numerico[, 2:5], center= TRUE, scale=TRUE)
matriz_cov <- cov(df_centrado)
matriz_cov
## Creacion_Toques Creacion_Toques_AP Creacion_Toques_Z1
## Creacion_Toques 1.0000000 0.1919393 0.6053904
## Creacion_Toques_AP 0.1919393 1.0000000 0.8552548
## Creacion_Toques_Z1 0.6053904 0.8552548 1.0000000
## Creacion_Toques_Z2 0.9941664 0.1501945 0.5673953
## Creacion_Toques_Z2
## Creacion_Toques 0.9941664
## Creacion_Toques_AP 0.1501945
## Creacion_Toques_Z1 0.5673953
## Creacion_Toques_Z2 1.0000000
Esta matriz simétrica es justo la materia prima del PCA, diagonalizándola es como encontramos las direcciones de mayor variabilidad y se hace el cambio de base
Con la función eigen obtenemos tanto los valores propios como los vectores propios
v_propios <- eigen(matriz_cov)
v_propios$values
## [1] 2.722832215 1.232848636 0.039628182 0.004690967
head(v_propios$vectors)
## [,1] [,2] [,3] [,4]
## [1,] -0.5353245 0.4191853 -0.1263221 0.72232542
## [2,] -0.3766687 -0.6970094 -0.6098767 0.01868295
## [3,] -0.5469743 -0.3635967 0.7514344 -0.06295168
## [4,] -0.5218883 0.4541575 -0.2178062 -0.68842867
Ahora vamos a calcular para cada club su valor en las 4 componentes principales para eso necesitamos multiplicar los vectores propios por los datos.
# transponemos los vectores propios y los datos para obtener las componentes
t_vec_propios <- t(v_propios$vectors)
t_df_centrado <- t(df_centrado)
# Producto matricial
pc_scores <- t_vec_propios %*% t_df_centrado
rownames(pc_scores) <- c("PC1", "PC2","PC3","PC4")
colnames(pc_scores) <- df_numerico$Squad
head(t(pc_scores), 10)
## PC1 PC2 PC3 PC4
## Alavés 1.4358255 -0.8720034 -0.21233168 -0.019668901
## Athletic Club 0.6999736 0.6767305 0.08296797 -0.081576648
## Atlético Madrid 0.3199010 1.1202263 0.14681784 0.046972256
## Barcelona -3.4646008 1.5976916 -0.32442720 0.056962318
## Betis -2.6909040 -0.1125068 -0.21107497 0.013592645
## Celta Vigo -0.7776694 -0.1446199 -0.02618855 -0.152061073
## Eibar 1.8179428 2.5941963 -0.20134263 -0.006140546
## Espanyol -0.9644433 -1.5216226 -0.17289972 0.086669116
## Getafe 3.2694678 0.5147319 -0.15659309 0.033695280
## Girona -0.6842965 -1.2661054 0.28673745 -0.032474417
Estas son las coordenadas de los clubes en esta nueva base, pero podemos recuperar la info original.
Reconstrucción de los datos originales
Como bien dice Joaquín en su post se pueden recuperar los datos originales cuando nos quedamos con todos los vectores propios
datos_recuperados <- t(v_propios$vectors %*% pc_scores)
head(datos_recuperados,10)
## [,1] [,2] [,3] [,4]
## Alavés -1.12154884 0.1960926 -0.6266173 -1.08557962
## Athletic Club -0.16044314 -0.7874700 -0.5614442 -0.01987703
## Atlético Madrid 0.31371452 -0.9899681 -0.4749212 0.27749175
## Barcelona 2.60654241 0.3903256 1.0667608 2.56518599
## Betis 1.42982745 1.2209814 1.3532979 1.38987128
## Celta Vigo 0.24915356 0.4068560 0.4678421 0.45056353
## Eibar 0.13525827 -2.3702619 -2.0885184 0.27749175
## Espanyol -0.03710733 1.5309276 0.9454039 -0.20973153
## Getafe -1.49033813 -1.4941473 -2.0952605 -1.46161740
## Girona -0.22409050 0.9647592 1.0521530 -0.25798184
head(df_centrado,10)
## Creacion_Toques Creacion_Toques_AP Creacion_Toques_Z1 Creacion_Toques_Z2
## [1,] -1.12154884 0.1960926 -0.6266173 -1.08557962
## [2,] -0.16044314 -0.7874700 -0.5614442 -0.01987703
## [3,] 0.31371452 -0.9899681 -0.4749212 0.27749175
## [4,] 2.60654241 0.3903256 1.0667608 2.56518599
## [5,] 1.42982745 1.2209814 1.3532979 1.38987128
## [6,] 0.24915356 0.4068560 0.4678421 0.45056353
## [7,] 0.13525827 -2.3702619 -2.0885184 0.27749175
## [8,] -0.03710733 1.5309276 0.9454039 -0.20973153
## [9,] -1.49033813 -1.4941473 -2.0952605 -1.46161740
## [10,] -0.22409050 0.9647592 1.0521530 -0.25798184
head(df_centrado - datos_recuperados, 10)
## Creacion_Toques Creacion_Toques_AP Creacion_Toques_Z1 Creacion_Toques_Z2
## [1,] -2.220446e-16 0.000000e+00 2.220446e-16 2.220446e-16
## [2,] -5.551115e-17 1.110223e-16 4.440892e-16 1.353084e-16
## [3,] -1.110223e-16 2.220446e-16 3.330669e-16 0.000000e+00
## [4,] -4.440892e-16 -2.220446e-16 -2.220446e-16 -8.881784e-16
## [5,] 2.220446e-16 -4.440892e-16 -8.881784e-16 -6.661338e-16
## [6,] -2.775558e-17 -1.665335e-16 -2.775558e-16 -1.665335e-16
## [7,] -1.110223e-16 4.440892e-16 1.332268e-15 2.220446e-16
## [8,] 4.163336e-17 -2.220446e-16 -6.661338e-16 -2.220446e-16
## [9,] 0.000000e+00 6.661338e-16 8.881784e-16 6.661338e-16
## [10,] 1.110223e-16 -4.440892e-16 -8.881784e-16 0.000000e+00
Afortunadamente hay funciones y librerías que hacen estas cosas y dan ayudas a la interpretación.
Reducción de dimensionalidad
Ahora si simplemente queremos quedarnos con el mejor resumen de las 4 variables originales nos quedamos solo con primer PC
sort(t(pc_scores)[,1])
## Barcelona Betis Real Madrid Real Sociedad Espanyol
## -3.46460082 -2.69090398 -2.19175378 -1.63031251 -0.96444334
## Celta Vigo Girona Sevilla Rayo Vallecano Villarreal
## -0.77766939 -0.68429649 -0.20760011 0.01532123 0.13996320
## Valencia Atlético Madrid Valladolid Athletic Club Levante
## 0.25217154 0.31990099 0.42063611 0.69997361 1.20774846
## Huesca Alavés Leganés Eibar Getafe
## 1.42757573 1.43582553 1.60505341 1.81794277 3.26946783
Dónde vemos que en un extremo está el Barcelona y en el otro el Getafe, ¿pero como pondera cada variable en esa componente principal? Pues simplemente viendo el primer vector propio (primera columna) se puede ver lo que pesa cada variable.
pesos_variables <- v_propios$vectors
rownames(pesos_variables) <- colnames(df_centrado)
colnames(pesos_variables) <- paste0("PC",1:4)
pesos_variables
## PC1 PC2 PC3 PC4
## Creacion_Toques -0.5353245 0.4191853 -0.1263221 0.72232542
## Creacion_Toques_AP -0.3766687 -0.6970094 -0.6098767 0.01868295
## Creacion_Toques_Z1 -0.5469743 -0.3635967 0.7514344 -0.06295168
## Creacion_Toques_Z2 -0.5218883 0.4541575 -0.2178062 -0.68842867
Asi vemos, que en la primera componente todas las variables tienen coordenada del mismo signo, es decir, se trata de una variable que pondrá en un extremo todos los equipos con alto valor en todas esas variables y en el otro los que tienen bajo valor en ese conjunto de variables.
La segunda PC ya si diferencia entre Creacion_toques y Creacion_toques_z2 con peso positivo vs Creacion_toques en AP y Z1.
Este comportamiento de un PCA es usual y se suele decir que la primera componente recoge el tamaño (equipos con alto valor en la mayoría de variables vs los que tienen bajo valor) y la segunda y siguientes componentes añaden matices.
Para facilitar la interpretación a veces se utiliza otro cambio de base haciendo rotaciones, la más usual es la rotación varimax, para facilitar la interpretación
PCA en R
En R ya hay multitud de funciones y librerías que hacen el pca directamente, por ejemplo princomp
y obtenemos lo mismo que antes salvo el signo en alguna componente, lo cual no varía su interpretación
res <- princomp(df_numerico[,2:5], cor=TRUE)
res$loadings
##
## Loadings:
## Comp.1 Comp.2 Comp.3 Comp.4
## Creacion_Toques 0.535 0.419 0.126 0.722
## Creacion_Toques_AP 0.377 -0.697 0.610
## Creacion_Toques_Z1 0.547 -0.364 -0.751
## Creacion_Toques_Z2 0.522 0.454 0.218 -0.688
##
## Comp.1 Comp.2 Comp.3 Comp.4
## SS loadings 1.00 1.00 1.00 1.00
## Proportion Var 0.25 0.25 0.25 0.25
## Cumulative Var 0.25 0.50 0.75 1.00
Y aquí los scores
res$scores
## Comp.1 Comp.2 Comp.3 Comp.4
## [1,] -1.47312591 -0.894656624 0.21784771 -0.020179866
## [2,] -0.71815778 0.694310797 -0.08512334 -0.083695875
## [3,] -0.32821149 1.149327966 -0.15063193 0.048192518
## [4,] 3.55460544 1.639197045 0.33285528 0.058442105
## [5,] 2.76080923 -0.115429590 0.21655835 0.013945759
## [6,] 0.79787196 -0.148376931 0.02686889 -0.156011369
## [7,] -1.86516992 2.661589241 0.20657318 -0.006300067
## [8,] 0.98949799 -1.561151804 0.17739137 0.088920636
## [9,] -3.35440322 0.528103818 0.16066112 0.034570627
## [10,] 0.70207339 -1.298996700 -0.29418642 -0.033318049
## [11,] -1.46466180 -0.249386278 0.02551176 -0.022539838
## [12,] -1.64675005 -0.024226568 -0.06234356 -0.061721168
## [13,] -1.23912377 -0.918827637 0.18189450 0.072835609
## [14,] -0.01571925 -0.996347013 -0.16727041 0.026577141
## [15,] 2.24869193 1.744009038 -0.32301661 -0.035840606
## [16,] 1.67266534 -1.476954610 0.14195906 -0.081997377
## [17,] 0.21299322 0.031819187 -0.29816043 0.146668682
## [18,] -0.25872255 0.082636708 -0.30022953 0.046260068
## [19,] -0.43156355 -0.844484815 0.01989843 0.023146309
## [20,] -0.14359922 -0.002155231 -0.02705741 -0.057955239
Esta entrada era solo para recordar que podemos obtener un PCA simplemente obteniendo los vectores propios de la matriz de correlaciones, y que aquello que vimos en los primeros cursos de la universidad sirve para algo.