Introduccion:
El análisis de componentes principales (ACP) permite resumir y visualizar la información de un conjunto de datos que contiene individuos descritos por múltiples variables cuantitativas interrelacionadas. Cada variable puede considerarse una dimensión diferente. Si tenemos más de 3 variables puede ser muy difÃcil visualizar un hiperespacio multidimensional. El análisis de componentes principales se utiliza para extraer la información importante de una tabla de datos multivariante y expresar esta información como un conjunto de pocas variables nuevas denominadas componentes principales. Estas nuevas variables corresponden a una combinación lineal de las originales. El objetivo del ACP es identificar las direcciones (o componentes principales) en las que la variación de los datos es máxima. En otras palabras, el ACP reduce la dimensionalidad de los datos multivariantes a uno o tres componentes principales, que pueden visualizarse gráficamente, con una pérdida mÃnima de información.
Instale los dos paquetes como se indica a continuación:
Cargarlos en R, escribiendo esto:
library(FactoMineR)
## Warning: package 'FactoMineR' was built under R version 4.3.2
library(factoextra)
## Warning: package 'factoextra' was built under R version 4.3.2
## Loading required package: ggplot2
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
Utilizaremos conjuntos de datos de demostración "decathlon2" del paquete "factoextra"
data("decathlon2")
los datos aquà utilizados describen el rendimiento de los atletas durante dos eventos deportivos (Descatar y OlymicG). Contiene 27 individuos (atletas) descritos por 13 variables.
Empezamos por subconjuntar los individuos activos y las variables activas para el análisis de componentes principales
decathlon2.active <- decathlon2[1:23, 1:10]
head(decathlon2.active[, 1:6], 4)
## X100m Long.jump Shot.put High.jump X400m X110m.hurdle
## SEBRLE 11.04 7.58 14.83 2.07 49.81 14.69
## CLAY 10.76 7.40 14.26 1.86 49.37 14.05
## BERNARD 11.02 7.23 14.25 1.92 48.93 14.99
## YURKOV 11.34 7.09 15.19 2.10 50.42 15.31
El código R siguiente, calcula el análisis de componentes principales en los individuos/varibles activos
library("FactoMineR")
res.pca <- PCA(decathlon2.active, graph = FALSE)
La salida de la funciónPCA() es una lista que incluye los siguientes componentes
print(res.pca)
## **Results for the Principal Component Analysis (PCA)**
## The analysis was performed on 23 individuals, described by 10 variables
## *The results are available in the following objects:
##
## name description
## 1 "$eig" "eigenvalues"
## 2 "$var" "results for the variables"
## 3 "$var$coord" "coord. for the variables"
## 4 "$var$cor" "correlations variables - dimensions"
## 5 "$var$cos2" "cos2 for the variables"
## 6 "$var$contrib" "contributions of the variables"
## 7 "$ind" "results for the individuals"
## 8 "$ind$coord" "coord. for the individuals"
## 9 "$ind$cos2" "cos2 for the individuals"
## 10 "$ind$contrib" "contributions of the individuals"
## 11 "$call" "summary statistics"
## 12 "$call$centre" "mean of the variables"
## 13 "$call$ecart.type" "standard error of the variables"
## 14 "$call$row.w" "weights for the individuals"
## 15 "$call$col.w" "weights for the variables"
Ahora, examinamos los valores propios para determinar el número de componentes principales a considerar. los "valores propios" y la proporción de varianzas retenidas por los componentes principales (PCs) pueden extraerse utilizando la función "get_eigenvalue()"[factoextra]"
library("factoextra")
eig.val <- get_eigenvalue(res.pca)
eig.val
## eigenvalue variance.percent cumulative.variance.percent
## Dim.1 4.1242133 41.242133 41.24213
## Dim.2 1.8385309 18.385309 59.62744
## Dim.3 1.2391403 12.391403 72.01885
## Dim.4 0.8194402 8.194402 80.21325
## Dim.5 0.7015528 7.015528 87.22878
## Dim.6 0.4228828 4.228828 91.45760
## Dim.7 0.3025817 3.025817 94.48342
## Dim.8 0.2744700 2.744700 97.22812
## Dim.9 0.1552169 1.552169 98.78029
## Dim.10 0.1219710 1.219710 100.00000
la suma de todos los valores propios da una varianza total de 10. La proporción de variación explicada por cada valor propio figura en la segunda columna. Por ejemplo, 4,124 dividido por 10 es igual a 0,4124, es decir, aproximadamente el 41,24% de la variación se explica por las proporciones sucesivas de variación explicadas para obtener el total acumulado. Por ejemplo, 41,242% más 18,385% es igual a 59,627%, y asà sucesivamente. Por lo tanto, aproximadamente el 59,62% de la variación se explica por los dos primeros valores propios juntos.
Desgraciadamente, no existe un método objetivo bien aceptado para decidir cuántos componentes principales son suficientes, lo que dependerá del campo de aplicación especÃfico y del conjunto de datos concreto. En la práctica, tendemos a fijarnos en los primeros componentes principales para encontrar patrones interesantes en los datos. El número de componentes se determina en el punto a partir del cual los valores propios restantes son todos relativamente pequeños y de tamaño comparable. el "scree plot" se puede producir utilizando la función "fviz_eig() o fviz_screeplot() [paquete factoextra].
fviz_eig(res.pca, addlabels = TRUE, ylim = c(0, 50))
Un método sencillo para extraer los resultados, para las variables, de una salida PCA es utilizar la función get_var() [paquete factoextra}. esta función proporciona una lista de matrices que contienen todos los resultados para la variabless activa 8coordenadas, correlación entre variables y ejes, conseno al cuadrado y contribuciones)
var <- get_pca_var(res.pca)
var
## Principal Component Analysis Results for variables
## ===================================================
## Name Description
## 1 "$coord" "Coordinates for the variables"
## 2 "$cor" "Correlations between variables and dimensions"
## 3 "$cos2" "Cos2 for the variables"
## 4 "$contrib" "contributions of the variables"
los componentes de get_pca_var() se pueden utilizar en el trazado de variables como se indica a continuación: varcoord:coordenadasdelasvariablesparacrearungráficodedispersiónvar�����:��������������������������������������á��������������ó����cos2: representa la calidad de la representación de las variables en el mapa de factores. se calcula como el cuadrado de las coordenadas: var.cos2*var.coord. var$contrib: contiene las contribuciones (en porcentaje) de las variables a los componentes principales.
head(var$coord, 4)
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## X100m -0.8506257 -0.17939806 0.3015564 0.0335732 -0.1944440
## Long.jump 0.7941806 0.28085695 -0.1905465 -0.1153896 0.2331567
## Shot.put 0.7339127 0.08540412 0.5175978 0.1284684 -0.2488129
## High.jump 0.6100840 -0.46521415 0.3300852 0.1445501 0.4027002
fviz_pca_var(res.pca, col.var = "black")
El gráfico anterior también se conoce como gráfico de correlación de variables. muestra las relaciones entre todas las variables.
head(var$cos2, 4)
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## X100m 0.7235641 0.032183664 0.09093628 0.00112716 0.03780845
## Long.jump 0.6307229 0.078880629 0.03630798 0.01331475 0.05436203
## Shot.put 0.5386279 0.007293864 0.26790749 0.01650412 0.06190783
## High.jump 0.3722025 0.216424207 0.10895622 0.02089474 0.16216747
library ("corrplot")
## Warning: package 'corrplot' was built under R version 4.3.2
## corrplot 0.92 loaded
corrplot (var$cos2, is.corr=FALSE)
Tambien es posible crear un diagramade barras de las variables cos2 usando lafuncion fviz_cos2()
fviz_cos2(res.pca, choice = "var", axes = 1:2)
Un Cos2 alto indica una buena representacion de la variable en el componente principal. en este caso la variable se posiciona cerca de la circunferencia del circulpo de correlacion. Un cos2 bajo ndica que la variable no esta perfectamente representada por las PC. en est caso la variable esta cerca del centro del circulo. Los valores de cos2 se utilizan para estimarla calidad de la represntacion.cuanto mas cerca este una variable del circulo de correlaciones, mejor sera su representacion.
fviz_pca_var(res.pca, col.var = "cos2",
gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
repel = TRUE
)
head(var$contrib, 4)
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## X100m 17.544293 1.7505098 7.338659 0.1375524 5.389252
## Long.jump 15.293168 4.2904162 2.930094 1.6248594 7.748815
## Shot.put 13.060137 0.3967224 21.620432 2.0140727 8.824401
## High.jump 9.024811 11.7715838 8.792888 2.5498795 23.115504
Se uso la funcion corrplot()para resaltar los elementos mas variables para cada elemento.
library ("corrplot")
corrplot (var$contrib, is.corr=FALSE)
El codigo siguiente muestra las 10 principales variables que contribuyen a los componentes principales
fviz_contrib(res.pca, choice = "var", axes = 1, top = 10)
La linea discontinua roja de los graficosanteriores indica la contribucion promedio esperada. Las variables mas importantesse puedenresaltaren el grafico siguiente
fviz_contrib(res.pca, choice = "var", axes = 2, top = 10)
fviz_pca_var(res.pca, col.var = "contrib",
gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07")
)
set.seed(123)
my.cont.var <- rnorm(10)
fviz_pca_var(res.pca, col.var = my.cont.var,
gradient.cols = c("blue", "yellow", "red"),
legend.title = "Cont .Var")
set.seed(123)
res.km <- kmeans(var$coord, centers = 3, nstart = 25)
grp <- as.factor(res.km$cluster)
fviz_pca_var(res.pca, col.var = grp,
palette = c("#0073C2FF", "#EFC000FF", "#868686FF"),
legend.title = "Cluster")
res.desc <- dimdesc(res.pca, axes = c(1,2), proba = 0.05)
res.desc$Dim.1
##
## Link between the variable and the continuous variables (R-square)
## =================================================================================
## correlation p.value
## Long.jump 0.7941806 6.059893e-06
## Discus 0.7432090 4.842563e-05
## Shot.put 0.7339127 6.723102e-05
## High.jump 0.6100840 1.993677e-03
## Javeline 0.4282266 4.149192e-02
## X400m -0.7016034 1.910387e-04
## X110m.hurdle -0.7641252 2.195812e-05
## X100m -0.8506257 2.727129e-07
res.desc$Dim.2
##
## Link between the variable and the continuous variables (R-square)
## =================================================================================
## correlation p.value
## Pole.vault 0.8074511 3.205016e-06
## X1500m 0.7844802 9.384747e-06
## High.jump -0.4652142 2.529390e-02
ind <- get_pca_ind(res.pca)
ind
## Principal Component Analysis Results for individuals
## ===================================================
## Name Description
## 1 "$coord" "Coordinates for the individuals"
## 2 "$cos2" "Cos2 for the individuals"
## 3 "$contrib" "contributions of the individuals"
fviz_pca_ind(res.pca)
fviz_pca_ind(res.pca, col.ind = "cos2",
gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
repel = TRUE
)
fviz_pca_ind(res.pca, pointsize = "cos2",
pointshape = 21, fill = "#E7B800",
repel = TRUE)
la generacion de un grafico de barrar de que represente (cos2)) de lo s individuos en el mapa de factores, puede usar la funcion fuiz_cos2().
fviz_contrib(res.pca, choice = "ind", axes = 1:2)
head(iris, 3)
## Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1 5.1 3.5 1.4 0.2 setosa
## 2 4.9 3.0 1.4 0.2 setosa
## 3 4.7 3.2 1.3 0.2 setosa
iris.pca <- PCA(iris[,-5], graph = FALSE)
Aquà describimos cómo colorear individuos por grupos. Además, mostraremos cómo añadir elipses de concentración y elipses de confianza por grupos. para ello, utilizaremos los datos del iris como conjuntos de datos de demostración.
fviz_pca_ind(iris.pca,
geom.ind = "point",
col.ind = iris$Species,
palette = c("#00AFBB", "#E7B800", "#FC4E07"),
addEllipses = TRUE,
legend.title = "Groups"
)
ind.p <- fviz_pca_ind(iris.pca, geom = "point", col.ind = iris$Species)
ggpubr::ggpar(ind.p,
title = "Principal Component Analysis",
subtitle = "Iris data set",
caption = "Source: factoextra",
xlab = "PC1", ylab = "PC?",
legend.title = "Species", legend.position = "top",
ggtheme = theme_gray(), palette = "jco"
)
fviz_pca_biplot(res.pca, repel = TRUE,
col.var = "#2E9FDF",
col.ind = "#696969"
)
fviz_pca_biplot(iris.pca,
col.ind = iris$Species, palette = "jco",
addEllipses = TRUE, label = "var",
col.var = "black", repel = TRUE,
legend.title = "Species")
fviz_pca_biplot(iris.pca,
geom.ind = "point",
pointshape = 21,
pointsize = 2.5,
fill.ind = iris$Species,
col.ind = "black",
col.var = factor(c("sepal", "sepal", "petal", "petal")),
legend.title = list(fill = "Species", color = "Clusters"),
repel = TRUE
)+
ggpubr::fill_palette("jco")+
ggpubr::color_palette ("npg")
Ahora un ejemplo diferente que es el de colorear los individuos por grupos (colores discretos) y variables por sus contribuciones a los componentes principales. Tambien, se puede cambiar la trasnparencia de las variables de las ariables por sus contribuciones utilizando el argumento alfa.var
fviz_pca_biplot(iris.pca,
geom.ind = "point",
fill.ind = iris$Species, col.ind = "black",
pointshape = 21, pointsize =,
palette = "jco",
addEllipses = TRUE,
alpha.var ="contrib", col.var = "contrib",
gradient.cols = "RdYlBu",
legend.title = list(fill = "Species", color = "Contrib",
alpha = "Contrib")
)
res.pca <- PCA(decathlon2, ind.sup = 24:27,
quanti.sup = 11:12, quali.sup = 13, graph=FALSE)
res.pca$quanti.sup
## $coord
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Rank -0.7014777 -0.24519443 -0.1834294 0.05575186 -0.07382647
## Points 0.9637075 0.07768262 0.1580225 -0.16623092 -0.03114711
##
## $cor
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Rank -0.7014777 -0.24519443 -0.1834294 0.05575186 -0.07382647
## Points 0.9637075 0.07768262 0.1580225 -0.16623092 -0.03114711
##
## $cos2
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Rank 0.4920710 0.060120310 0.03364635 0.00310827 0.0054503477
## Points 0.9287322 0.006034589 0.02497110 0.02763272 0.0009701427
fviz_pca_var (res.pca)
Estos son los resultados previstos para los individuos suplementarios (ind.sup) -Visualizar todo los individuos(activos y suplementarios). En grafico, puede agregar tambien las variables cualitativas(quali.sup), cuyas coordenadas son accesibles usando res.pcaquali.supp coord.
res.pca$ind.sup
## $coord
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## KARPOV 0.7947206 0.77951227 -1.6330203 1.7242283 -0.75070396
## WARNERS -0.3864645 -0.12159237 -1.7387332 -0.7063341 -0.03230011
## Nool -0.5591306 1.97748871 -0.4830358 -2.2784526 -0.25461493
## Drews -1.1092038 0.01741477 -3.0488182 -1.5343468 -0.32642192
##
## $cos2
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## KARPOV 0.05104677 4.911173e-02 0.21553730 0.24028620 0.0455487744
## WARNERS 0.02422707 2.398250e-03 0.49039677 0.08092862 0.0001692349
## Nool 0.02897149 3.623868e-01 0.02162236 0.48108780 0.0060077529
## Drews 0.09207094 2.269527e-05 0.69560547 0.17617609 0.0079736753
##
## $dist
## KARPOV WARNERS Nool Drews
## 3.517470 2.482899 3.284943 3.655527
p <- fviz_pca_ind(res.pca, col.ind.sup = "blue", repel = TRUE)
p <- fviz_add(p, res.pca$quali.sup$coord, color = "red")
p
res.pca$quali
## $coord
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Decastar -1.343451 0.1218097 -0.03789524 0.1808357 0.1343364
## OlympicG 1.231497 -0.1116589 0.03473730 -0.1657661 -0.1231417
##
## $cos2
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Decastar 0.9051233 0.007440939 0.0007201669 0.01639956 0.009050062
## OlympicG 0.9051233 0.007440939 0.0007201669 0.01639956 0.009050062
##
## $v.test
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Decastar -2.970766 0.4034256 -0.1528767 0.8971036 0.7202457
## OlympicG 2.970766 -0.4034256 0.1528767 -0.8971036 -0.7202457
##
## $dist
## Decastar OlympicG
## 1.412108 1.294433
##
## $eta2
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Competition 0.4011568 0.00739783 0.001062332 0.03658159 0.02357972
fviz_pca_ind(res.pca, habillage = 13,
addEllipses =TRUE, ellipse.type = "confidence",
palette = "jco", repel = TRUE)
fviz_pca_var(res.pca, select.var = list(cos2 = 0.6))
fviz_pca_var(res.pca, select.var = list(cos2= 5))
name<- list( name =c("Long jump", "high.jump", "X100m"))
fviz_pca_var(res.pca, select.var = name)
fviz_pca_biplot(res.pca, select.ind= list(contrib =5),
select.var = list(contrib = 5),
ggtheme = theme_minimal())
pdf("myplot.pdf")
dev.off()
## png
## 2
scree.plot <-fviz_eig(res.pca)
ind.plot <-fviz_pca_ind(res.pca)
var.plot <-fviz_pca_var(res.pca)
pdf("PCA.pdf")
print(scree.plot)
print(ind.plot)
print(var.plot)
png("pca-scree-plot.png")
print(scree.plot)
dev.off()
## png
## 2
png("pca_variables.png")
print(var.plot)
dev.off()
## png
## 2
png("pca-individuals.pnag")
print(ind.plot)
dev.off()
## png
## 2
library(ggpubr)
## Warning: package 'ggpubr' was built under R version 4.3.2
ggexport(plotlist= list(scree.plot, ind.plot, var.plot),
filename= "PCA.pdf")
## file saved to PCA.pdf
ggexport( plotlist = list(scree.plot, ind.plot, var.plot), nrow = 2, ncol = 2, filename = "PCA.pdf")
## file saved to PCA.pdf
ggexport(plotlist = list(scree.plot, ind.plot, var.plot), filename = "PCA.pdf")
## file saved to PCA.pdf
write.infile(res.pca, "pca.txt", sep = "t")
write.infile(res.pca, "pca.csv", sep = ";")
res.pca <-prcomp(iris [, -5],scale. =TRUE)
res.pca<- princomp(iris[, -5], cor =TRUE)