El análisis de componentes principales (ACP) es un método estadístico de extracción de datos que permite simplificar la complejidad de espacios muestrales con múltiples dimensiones pero manteniendo su información. El ACP consiste en reducir la dimensión de un gran conjunto de variables a un menor número a las que se le conocerán como componentes principales, esto lo convierte en un método muy útil de usar antes de hacer uso de otras técnicas estadísticas, pues condensa la información de los datos por medio de la búsqueda que realiza en otras dimensiones del dataset de los componentes principales para definir la estructura base, reduciéndolo a menos variables no relacionadas entre sí y manteniendo la mayor cantidad de información original, por lo que es primordial disponer del valor de las variables originales para calcular los componentes.

Primero, elegimos nuestra base de datos, en este caso se hará uso de los paquetes basicos de R que es “VADeaths”: Tasas de mortalidad en Virginia (Estados Unidos) por cada 1000 habitantes (1940)

data("VADeaths")
head(VADeaths)
##       Rural Male Rural Female Urban Male Urban Female
## 50-54       11.7          8.7       15.4          8.4
## 55-59       18.1         11.7       24.3         13.6
## 60-64       26.9         20.3       37.0         19.3
## 65-69       41.0         30.9       54.6         35.1
## 70-74       66.0         54.3       71.1         50.0
str(VADeaths)
##  num [1:5, 1:4] 11.7 18.1 26.9 41 66 8.7 11.7 20.3 30.9 54.3 ...
##  - attr(*, "dimnames")=List of 2
##   ..$ : chr [1:5] "50-54" "55-59" "60-64" "65-69" ...
##   ..$ : chr [1:4] "Rural Male" "Rural Female" "Urban Male" "Urban Female"

Se puede observar los tipo de datos de las variables que tiene la data.

Ahora debemos asegurarnos de que no hayan datos faltantes, en caso de haberlos deben ser eliminados:

missing(VADeaths)
## [1] FALSE

Al comprobarse que no hay datos perdidos o faltantes, vemos la variabilidad de las columnas:

library(PerformanceAnalytics)
## Loading required package: xts
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## Attaching package: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
## 
##     legend
chart.Correlation(VADeaths, histogram = F, pch = 19)

apply(VADeaths, 2, var)
##   Rural Male Rural Female   Urban Male Urban Female 
##      466.393      339.452      509.967      291.157

Procedemos con la comprobar si es necesario aplicar un ACP, pues es necesaria una base de datos con una alta correlación entre las variables a analizar para llevar a cabo este analisis:

Hipotesis para comprobar correlación

H0: Variables no correlacionadas entre sí

H1:Variables correlacionadas entre sí

library(psych)
cortest.bartlett(VADeaths, n=50)
## R was not square, finding R from data
## $chisq
## [1] 27.24037
## 
## $p.value
## [1] 0.0001305366
## 
## $df
## [1] 6

como el p-valor es menor a 0.05, no existe suficiente información para aceptar la hipótesis nula, así que la matriz es válida para hacer el análisis

Test KMO:

KMO(VADeaths)
## Kaiser-Meyer-Olkin factor adequacy
## Call: KMO(r = VADeaths)
## Overall MSA =  0.75
## MSA for each item = 
##   Rural Male Rural Female   Urban Male Urban Female 
##         0.68         0.69         0.81         0.86

Con este test podemos ver la relación de coeficientes de correlación entre variables si el valor es cercano a 1, mayor será la relación entre variables, dado los valores, se puede decir que es apto para la realización del análisis de componentes.

apply(X = VADeaths, MARGIN = 2, FUN = mean)
##   Rural Male Rural Female   Urban Male Urban Female 
##        32.74        25.18        40.48        25.28

la varianza es diferente entre las variables. Para Urban Male, la varianza es de ordenes de magnitud mayor al resto.

Se debe estandarizar la variables del data para que tengan media cero y desviación estándar de 1, para aplicar la técnica de PCA ya que urban Male tiene una varianza mayor al resto, si no se hace esta variable puede dominar o controlar los componentes principales.

apply(X = VADeaths, MARGIN = 2, FUN = var)
##   Rural Male Rural Female   Urban Male Urban Female 
##      466.393      339.452      509.967      291.157

Prcomp es una de las múltiples funciones de R que realiza PCA, esta se centra en las variables para que tengan media cero y el scale=TRUE hace que la desviación estándar sea uno.

pca <- prcomp(VADeaths, scale = TRUE)
names(pca)
## [1] "sdev"     "rotation" "center"   "scale"    "x"

Center y scale de pca contienen media y la desviación típica de las variables previa a la estandarización en escala original.

pca$center
##   Rural Male Rural Female   Urban Male Urban Female 
##        32.74        25.18        40.48        25.28
pca$scale
##   Rural Male Rural Female   Urban Male Urban Female 
##     21.59613     18.42422     22.58245     17.06332
pca$rotation
##                    PC1        PC2         PC3        PC4
## Rural Male   0.5014626 -0.3026103  0.09885977  0.8044806
## Rural Female 0.4993178 -0.6102618  0.23254315 -0.5693733
## Urban Male   0.4982092  0.6997763  0.50026270 -0.1088027
## Urban Female 0.5010037  0.2152226 -0.82818335 -0.1295643

Las siguientes son las coordenadas que corresponden a los scores de los componentes:

dim(pca$rotation)
## [1] 4 4
# número de distintos componentes
head(pca$x)
##              PC1          PC2          PC3          PC4
## 50-54 -1.9841069 -0.149398478 -0.040622594 -0.025465751
## 55-59 -1.4051659  0.002933054 -0.028688104  0.037866254
## 60-64 -0.5202149  0.060207074  0.124826129 -0.004563758
## 65-69  0.9466577  0.256202962 -0.053818854 -0.011668846
## 70-74  2.9628300 -0.169944612 -0.001696576  0.003832101

La función prcomp hace cálculos de los valores de los componentes principales para cada observación y esos resultados se almacenan en la matriz x

library(ggplot2)
## 
## Attaching package: 'ggplot2'
## The following objects are masked from 'package:psych':
## 
##     %+%, alpha
pca$sdev^2
## [1] 3.9640793316 0.0301186321 0.0052385281 0.0005635081
scree(VADeaths, main = "Grafico de Sedimentacion")
## Warning in fa.stats(r = r, f = f, phi = phi, n.obs = n.obs, np.obs = np.obs, :
## The estimated weights for the factor scores are probably incorrect. Try a
## different factor score estimation method.

Para la varianza explicada por cada componente

prop_varianza <- pca$sdev^2 / sum(pca$sdev^2)
prop_varianza
## [1] 0.991019833 0.007529658 0.001309632 0.000140877

A continuación, con la función biplot se obtiene la representación bidimensional de las dos primeras componentes y se indica 0 para que la flechas estén en la misma escala de los componentes

biplot(x = pca, scale = 0, cex = 0.6, col = c("blue4", "brown3"))

library(factoextra)
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
fviz_pca_var(pca, col.var = "cos2",
             geom.var = "arrow",
             labelsize = 2, repel = FALSE)

#valores de explicación de los eigenvalues.
#Marca el valor de la contribución de las variables a los componentes.
fviz_screeplot(pca, addlabels = TRUE, ylim = c(0, 100))

Se guardan en un vector los componentes

PC1<-pca[[2]][,1]
PC2<-pca[[2]][,2]
Individuos<-pca$x[,1:2]
head(Individuos)
##              PC1          PC2
## 50-54 -1.9841069 -0.149398478
## 55-59 -1.4051659  0.002933054
## 60-64 -0.5202149  0.060207074
## 65-69  0.9466577  0.256202962
## 70-74  2.9628300 -0.169944612
Componentes<-cbind(PC1,PC2)
head(Componentes)
##                    PC1        PC2
## Rural Male   0.5014626 -0.3026103
## Rural Female 0.4993178 -0.6102618
## Urban Male   0.4982092  0.6997763
## Urban Female 0.5010037  0.2152226
prop_varianza_acum <- cumsum(prop_varianza)
prop_varianza_acum
## [1] 0.9910198 0.9985495 0.9998591 1.0000000
ggplot(data = data.frame(prop_varianza_acum, pc = 1:4),
       aes(x = pc, y = prop_varianza_acum, group = 1)) +
  geom_point() +
  geom_line() +
  theme_bw() +
  labs(x = "Componente principal",
       y = "Prop. varianza explicada acumulada")

La primera componente explica el 99.1% de la varianza observada de los datos, la segunda componente explica el 0.8% de la varianza observada de los datos mientras que las dos últimas componentes (3 y 4) están muy por debajo del 1% o del 0.7%

Al emplearse solo las dos primeras componentes se explica el 99.85% de la varianza observada Si se emplean las tres primeras componentes, se explicaría el 99.99% de la varianza observada.

pc1<- apply(pca$rotation[,1]*VADeaths, 1, sum)
pc2<- apply(pca$rotation[,2]*VADeaths, 2, sum)
VADeaths$PC1 <- pc1
## Warning in VADeaths$PC1 <- pc1: Realizando coercion de LHD a una lista
VADeaths$PC2 <- pc2
VADeaths[1:4] <-  NULL
head(VADeaths)
## [[1]]
## [1] 66
## 
## [[2]]
## [1] 8.7
## 
## [[3]]
## [1] 11.7
## 
## [[4]]
## [1] 20.3
## 
## [[5]]
## [1] 30.9
## 
## [[6]]
## [1] 54.3