Análise de Componentes Principais - ACP

---

library(factoextra)
library(FactoMineR)
library(ggplot2)
library(ggcorrplot)
library(psych)

Perfil de profissional de empregados de uma empresa

Variáveis

1. gênero 1 H 2 M
2. educação - em anos
3. ordenado em 100
4. sal inicial em 100
5. tempo serv em meses
6. exp anterior em meses

empresa <- rio::import(file = "Base Perfil dos Empregados de uma Empresa.sav")
names(empresa) <- c("genero", "educacao", "odenado", "salinicial","temposerv", "expanterior")
empresa$genero <- factor(c(1,2), labels = c("Homem", "Mulher"))
empresa

# retirar gênero da amostra
dados <- empresa[, -1]
# sumário estatístico 
describe(dados)

            vars  n     mean       sd  median  trimmed     mad   min    max range
educacao       1 30    13.67     2.73    15.0    13.83    3.71     8     19    11
odenado        2 30 38660.33 19687.22 30750.0 34242.08 8784.40 21150 100000 78850
salinicial     3 30 18437.00  7787.54 15750.0 16926.25 3002.26  9000  44100 35100
temposerv      4 30    78.37    10.02    76.5    77.79   11.12    64     97    33
expanterior    5 30   127.20   117.85    93.0   106.71   65.98     4    451   447
             skew kurtosis      se
educacao    -0.42    -0.27    0.50
odenado      1.84     2.62 3594.38
salinicial   1.92     3.20 1421.80
temposerv    0.29    -1.22    1.83
expanterior  1.40     1.19   21.52

# gráfico boxplot 
boxplot(dados) # dados originais

boxplot(scale(dados)) # dados padronizados

correlações

# Correlation matrix
matcor <- round(cor(dados), 2)
matcor

            educacao odenado salinicial temposerv expanterior
educacao        1.00    0.57       0.58     -0.18       -0.36
odenado         0.57    1.00       0.96     -0.01       -0.07
salinicial      0.58    0.96       1.00     -0.04       -0.05
temposerv      -0.18   -0.01      -0.04      1.00        0.35
expanterior    -0.36   -0.07      -0.05      0.35        1.00

# Plot
ggcorrplot(matcor, hc.order = TRUE, 
           type = "lower", 
           lab = TRUE, 
           lab_size = 3, 
           method="circle", 
           colors = c("tomato2", "white", "springgreen3"), 
           title="Correlograma", 
           ggtheme=theme_bw)

teste de Bartlett

Bartlett.sphericity.test

Bartlett.sphericity.test <- function(x)
{
  method <- "Bartlett's test of sphericity"
  data.name <- deparse(substitute(x))
  x <- subset(x, complete.cases(x)) # Omit missing values
  n <- nrow(x)
  p <- ncol(x)
  chisq <- (1-n+(2*p+5)/6)*log(det(cor(x)))
  df <- p*(p-1)/2
  p.value <- pchisq(chisq, df, lower.tail=FALSE)
  names(chisq) <- "X-squared"
  names(df) <- "df"
  return(structure(list(statistic=chisq, parameter=df, p.value=p.value,
                        method=method, data.name=data.name), class="htest"))
}
Bartlett.sphericity.test(dados)

    Bartlett's test of sphericity

data:  dados
X-squared = 88.724, df = 10, p-value = 9.589e-15

PCA

A análise de componentes principais será baseada na matriz de correlações amostral

# PCA com a matriz de cor
res.pca.cor <- PCA(dados, scale.unit = T, graph = FALSE)
# matriz de covariância
round(cor(dados),4)

            educacao odenado salinicial temposerv expanterior
educacao      1.0000  0.5709     0.5752   -0.1792     -0.3556
odenado       0.5709  1.0000     0.9625   -0.0099     -0.0658
salinicial    0.5752  0.9625     1.0000   -0.0405     -0.0539
temposerv    -0.1792 -0.0099    -0.0405    1.0000      0.3506
expanterior  -0.3556 -0.0658    -0.0539    0.3506      1.0000

# autovalores
round(res.pca.cor$eig,3)

# autovetores
round(res.pca.cor$svd$V,3)

       [,1]   [,2]   [,3]   [,4]   [,5]
[1,]  0.512 -0.158  0.277  0.798 -0.009
[2,]  0.583  0.257 -0.113 -0.292 -0.704
[3,]  0.585  0.247 -0.158 -0.264  0.709
[4,] -0.128  0.644  0.752 -0.051  0.029
[5,] -0.199  0.658 -0.566  0.454 -0.022

# A proporção de variação retida pelos componentes principais (CP) pode ser extraída da seguinte forma
res.pca.cor$eig

# A importância dos CP pode ser visualizada usando o scree plot :
fviz_screeplot(res.pca.cor, ncp=4)+ theme_minimal()

# A correlação entre uma variável e um CP é chamada de carga (loadings). 
round(res.pca.cor$var$cor,4)

              Dim.1   Dim.2   Dim.3   Dim.4   Dim.5
educacao     0.8078 -0.1859  0.2282  0.5107 -0.0017
odenado      0.9204  0.3019 -0.0932 -0.1868 -0.1348
salinicial   0.9228  0.2909 -0.1299 -0.1689  0.1357
temposerv   -0.2021  0.7574  0.6200 -0.0324  0.0056
expanterior -0.3142  0.7738 -0.4669  0.2908 -0.0043

O quadrado da correlação entre a variável e a CP representa a porcentagem de variância de uma das variáveis originais explicada por uma das CP

Conceito de Comunalidade

a comunalidade é a soma dos quadrados das correlações entre cada variável i e a componente principal j (ou o mesmo que o índice cos2). A soma limite-se ao número do componenentes retidos. Em nosso exemplo ilustrativo retemos todos os componentes que é igual ao número de variáveis.

# banseando-se na matriz de cor
round(res.pca.cor$var$cor^2,4)

             Dim.1  Dim.2  Dim.3  Dim.4  Dim.5
educacao    0.6525 0.0345 0.0521 0.2609 0.0000
odenado     0.8471 0.0911 0.0087 0.0349 0.0182
salinicial  0.8516 0.0846 0.0169 0.0285 0.0184
temposerv   0.0409 0.5736 0.3844 0.0011 0.0000
expanterior 0.0987 0.5987 0.2180 0.0846 0.0000

# ou  round(res.pca.cor$var$cos2,4)

Mapa Fatorial

Quando um subespaço projetivo bidimensional determinado por duas direções principais escolhidas (CP), sua imagem geométrica plana com os pontos projetados e o círculo de correlações é denominada MAPA FATORIAL

# Quanto mais próxima uma variável for do círculo de correlações, melhor sua representação no mapa fatorial (e mais importante é a variável para a interpretação desses componentes)
# As variáveis próximas ao centro do gráfico são menos importantes para os primeiros componentes.
# No gráfico abaixo os componentes são coloridas de acordo com os valores do coseno quadrado:
fviz_pca_var(res.pca.cor, col.var="cos2") +
scale_color_gradient2(low="white", mid="blue", 
                    high="red", midpoint=0.5) + theme_minimal()

# Coordenadas de variáveis
round(res.pca.cor$var$coord,2)

            Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
educacao     0.81 -0.19  0.23  0.51  0.00
odenado      0.92  0.30 -0.09 -0.19 -0.13
salinicial   0.92  0.29 -0.13 -0.17  0.14
temposerv   -0.20  0.76  0.62 -0.03  0.01
expanterior -0.31  0.77 -0.47  0.29  0.00

# Cos2: é uma medida que indica a qualidade da representação para variáveis no mapa fatorial
round(res.pca.cor$var$cos2,2)

            Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
educacao     0.65  0.03  0.05  0.26  0.00
odenado      0.85  0.09  0.01  0.03  0.02
salinicial   0.85  0.08  0.02  0.03  0.02
temposerv    0.04  0.57  0.38  0.00  0.00
expanterior  0.10  0.60  0.22  0.08  0.00

Contribuições das variáveis para os componentes principais

As variáveis que são correlacionadas com PC1 e PC2 são as mais importantes para explicar a variabilidade no conjunto de dados. Variáveis que não se correlacionam com nenhum PC ou correlacionadas com as últimas dimensões são variáveis com baixa contribuição e podem ser removidas para simplificar a análise geral. As contribuições das variáveis na contabilização da variabilidade em uma determinada componente principal são (em porcentagem): (variável.cos2 * 100) / (cos2 total da componente)

# A contribuição das variáveis pode ser extraída da seguinte forma:
round(res.pca.cor$var$contrib,2)

            Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
educacao    26.20  2.50  7.65 63.64  0.01
odenado     34.01  6.59  1.28  8.51 49.61
salinicial  34.19  6.12  2.48  6.96 50.25
temposerv    1.64 41.49 56.53  0.26  0.08
expanterior  3.96 43.30 32.05 20.63  0.05

# veja que a soma é igual a 100%
sum (res.pca.cor$var$contrib[,1])  

[1] 100

# Quanto maior o valor da contribuição, mais a variável contribui para o componente.
# As variáveis mais importantes associadas a um determinado PC podem ser visualizadas, usando a função fviz_contrib () [factoextra package], da seguinte forma:
# Contribuições de variáveis no PC1
fviz_contrib(res.pca.cor, choice = "var", axes = 1)+ theme_minimal()

# Contribuições de variáveis no PC2
fviz_contrib(res.pca.cor, choice = "var", axes = 2)+ theme_minimal()

# Contribuição total nos PC1 e PC2
fviz_contrib(res.pca.cor, choice = "var", axes = 1:2)+ theme_minimal()

# Controle as cores das variáveis usando suas contribuições
# a cor representa a contribuição conjunta dim1-dim2
fviz_pca_var(res.pca.cor, col.var="contrib")+ theme_minimal()

# Alterar a cor 
fviz_pca_var(res.pca.cor, col.var="contrib") +
scale_color_gradient2(low="white", mid="blue", 
                  high="red", midpoint=50) + theme_minimal()

A função dimdesc () [em FactoMineR] pode ser usada para identificar as variáveis mais correlacionadas com uma determinada componente principal.

res.desc <- dimdesc(res.pca.cor, axes = c(1,2))
# Descrição da dimensão 1
res.desc$Dim.1

$quanti
           correlation      p.value
salinicial   0.9227974 4.062575e-13
odenado      0.9203827 6.158434e-13
educacao     0.8078005 6.798363e-08

# Descrição da dimensão 2
res.desc$Dim.2

$quanti
            correlation      p.value
expanterior   0.7737622 5.309827e-07
temposerv     0.7573869 1.264469e-06

# Descrição da dimensão 3
res.desc$Dim.3

NULL

Análise de pontos (escores, objetos, indivíduos)

# Gráfico de escores (indivíduos ou pontos objetos)
# As coordenadas dos escores nos componentes principais são:
round(res.pca.cor$ind$coord,2)

   Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
1  -1.31  0.85  1.64  1.06 -0.19
2   3.96  2.51  0.91  0.13  0.08
3  -1.46  1.94  0.37  0.41 -0.15
4  -2.28  2.74 -0.92 -0.32 -0.03
5  -1.11  1.11  0.43  0.00 -0.01
6   1.13  0.91  0.78  1.21  0.26
7  -2.25  2.52 -1.39 -0.22 -0.05
8  -1.01 -0.33  1.24 -0.63  0.00
9   0.09 -0.10  0.89  0.23 -0.25
10  1.36  0.22  0.63 -0.42  0.09
11  1.55  0.29  0.62 -0.06  0.40
12 -0.86 -0.49  0.72 -0.57  0.20
13 -1.74 -0.52  0.39 -1.77 -0.14
14 -0.11 -0.42  0.63 -1.01 -0.08
15  0.22 -0.60  0.38  0.23  0.11
16  0.42 -0.91  0.55 -0.02 -0.41
17 -1.03 -0.62  0.01 -0.24  0.21
18 -0.66 -0.92  0.31 -0.64  0.21
19 -0.99 -0.57 -0.29 -0.15  0.12
20  0.39 -0.55 -0.27  0.83  0.04
21  2.83 -0.22 -0.45 -0.16 -0.52
22 -1.04  0.54 -1.86  0.64  0.00
23 -0.19 -0.91 -0.51  0.27  0.10
24 -0.88 -0.85 -0.70 -0.14  0.18
25  0.08 -1.26 -0.28  0.42 -0.03
26  0.11 -1.15 -0.49  0.50 -0.11
27 -0.40 -1.27 -0.59  0.83  0.04
28  4.41  0.70 -1.59 -1.05  0.10
29  0.96 -0.92 -0.71  0.07 -0.21
30 -0.19 -1.74 -0.45  0.55  0.04

fviz_pca_ind(res.pca.cor)+ theme_minimal()

#Cos2: qualidade da representação para escores nos componentes principais
# O coseno quadrado mostra a importância de um componente para uma determinada observação.
round(res.pca.cor$ind$cos2,3)

   Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
1  0.273 0.115 0.428 0.179 0.006
2  0.686 0.277 0.036 0.001 0.000
3  0.341 0.606 0.022 0.028 0.004
4  0.381 0.550 0.061 0.007 0.000
5  0.466 0.463 0.071 0.000 0.000
6  0.302 0.196 0.144 0.342 0.016
7  0.378 0.474 0.144 0.004 0.000
8  0.333 0.036 0.502 0.129 0.000
9  0.009 0.011 0.854 0.059 0.067
10 0.743 0.020 0.162 0.072 0.003
11 0.793 0.028 0.126 0.001 0.052
12 0.400 0.130 0.277 0.172 0.022
13 0.459 0.041 0.023 0.474 0.003
14 0.007 0.108 0.246 0.635 0.004
15 0.080 0.582 0.232 0.084 0.021
16 0.120 0.563 0.203 0.000 0.114
17 0.687 0.247 0.000 0.038 0.028
18 0.241 0.463 0.051 0.222 0.023
19 0.688 0.227 0.057 0.017 0.011
20 0.124 0.247 0.061 0.567 0.001
21 0.937 0.006 0.023 0.003 0.031
22 0.207 0.055 0.659 0.079 0.000
23 0.030 0.684 0.214 0.063 0.009
24 0.379 0.354 0.241 0.009 0.016
25 0.003 0.858 0.043 0.095 0.000
26 0.006 0.720 0.133 0.134 0.006
27 0.057 0.573 0.124 0.245 0.001
28 0.824 0.021 0.108 0.047 0.000
29 0.398 0.363 0.218 0.002 0.020
30 0.010 0.849 0.056 0.085 0.000

fviz_pca_ind(res.pca.cor, col.ind="cos2") +
scale_color_gradient2(low="white", mid="blue", 
   high="red", midpoint=0.50) + theme_minimal()

# Contribuição dos escores para os componentes principais 
round(res.pca.cor$ind$contrib,2)

   Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
1   2.30  1.75 13.21  9.15  3.44
2  20.96 15.24  4.02  0.14  0.58
3   2.84  9.09  0.66  1.39  2.03
4   6.96 18.09  4.11  0.83  0.08
5   1.65  2.95  0.92  0.00  0.00
6   1.72  2.00  3.01 11.83  6.01
7   6.77 15.32  9.47  0.39  0.22
8   1.35  0.27  7.49  3.19  0.00
9   0.01  0.03  3.87  0.44  5.65
10  2.46  0.12  1.96  1.45  0.70
11  3.23  0.21  1.88  0.03 14.31
12  1.00  0.59  2.54  2.62  3.66
13  4.05  0.66  0.76 25.43  1.76
14  0.02  0.42  1.93  8.27  0.60
15  0.07  0.87  0.71  0.43  1.20
16  0.24  2.00  1.47  0.00 15.36
17  1.42  0.92  0.00  0.47  3.86
18  0.59  2.05  0.46  3.31  3.86
19  1.31  0.78  0.40  0.19  1.40
20  0.20  0.72  0.36  5.57  0.14
21 10.75  0.12  0.97  0.21 24.58
22  1.45  0.70 16.95  3.35  0.00
23  0.05  1.97  1.26  0.61  1.00
24  1.04  1.75  2.41  0.15  3.06
25  0.01  3.81  0.38  1.42  0.07
26  0.02  3.18  1.20  2.00  1.05
27  0.21  3.90  1.72  5.63  0.16
28 26.05  1.19 12.45  8.95  0.91
29  1.23  2.02  2.46  0.04  4.15
30  0.05  7.30  0.97  2.48  0.15

# Contribuições de escores para PC1
fviz_contrib(res.pca.cor, choice = "ind", axes = 1)+ theme_minimal()

# Contribuições de escores para PC2
fviz_contrib(res.pca.cor, choice = "ind", axes = 2)+ theme_minimal()

# Contribuição total em PC1 e PC2
fviz_contrib(res.pca.cor, choice = "ind", axes = 1:2)+ theme_minimal()

# Contribuições dos escores para PC1  (apenas os "top")
fviz_contrib(res.pca.cor, choice = "ind", axes = 1:2, top = 5)+ theme_minimal()

# Mundando a cor
fviz_pca_ind(res.pca.cor, col.ind="contrib") +
scale_color_gradient2(low="white", mid="blue", 
                  high="red", midpoint=50) + theme_minimal()

O gráfico de dispersão dos escores dos dois primeiros componentes baseados na matriz de correlações juntamente com os respectivos autovetores

Este gráfico é chamado de biplot. É uma representação bidimensional de dados multivariados.

fviz_pca_biplot(res.pca.cor) + theme_minimal()

fviz_pca_biplot(res.pca.cor, habillage=empresa$genero) + theme_minimal()

sumário

# sumário
facto_summarize(res.pca.cor, "var") # para variáveis

facto_summarize(res.pca.cor, "ind") # para escores

Extra

simulação para reter o número de CP

Parallel analysis (see Hayton, Allen, and Scarpello, 2004 for more details)

Rotação varimax

A rotação pode facilitar a interpretação dos componentes

fa.parallel(dados, fa="pc", show.legend=FALSE,
            main="Scree plot with parallel analysis")

Parallel analysis suggests that the number of factors =  NA  and the number of components =  2 

# sugere duas CP
pc <- principal(r=dados, nfactors=2, rotate="none", scores=T)
pc$loadings

Loadings:
            PC1    PC2   
educacao     0.808 -0.186
odenado      0.920  0.302
salinicial   0.923  0.291
temposerv   -0.202  0.757
expanterior -0.314  0.774

                 PC1   PC2
SS loadings    2.491 1.383
Proportion Var 0.498 0.277
Cumulative Var 0.498 0.775

pc <- principal(r=dados, nfactors=2, rotate="varimax", scores=T)
pc$loadings

Loadings:
            RC1    RC2   
educacao     0.728 -0.397
odenado      0.968       
salinicial   0.967       
temposerv           0.784
expanterior         0.830

                 RC1   RC2
SS loadings    2.410 1.463
Proportion Var 0.482 0.293
Cumulative Var 0.482 0.775