Análise de componentes principais
Quando temos dados multivariados, a análise de componentes principais (PCA) é um recurso muito interessante e relativamente simples, em termos de conceito teórico e interpretação prática. Para exemplificar, vamos trabalhar com os dados climáticos de algumas cidades brasileiras. Os dados climáticos foram compilados a partir de estações automáticas do INMET.
No R, temos a facilidade de poder fazer o cálculo dos componentes principais e logo em seguida poder apresentá-los em gráficos elegantes e de fácil entendimento. O Objetivo deste post é apresentar uma rápida demonstração de como rodar um PCA e gerar os gráficos derivados.
if (!require("pacman")) install.packages("pacman")
pacman::p_load(readr, dplyr, ggplot2, ggrepel)
pacman::p_load_gh("vqv/ggbiplot")dados <- read_csv2(
"https://github.com/italocegatta/italocegatta.github.io_source/raw/master/content/dados/base_clima.csv"
)
print(dados, n=31)## # A tibble: 31 x 6
## Cidade Koppen Tmed PPT ETP DEF
## <chr> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 Bom Despacho Cwa 22.6 802. 1112. 238.
## 2 Niquelandia Aw 24.6 562. 1372. 656.
## 3 Arapoti Cfb 18.5 1367 670. 25.4
## 4 Rio Verde Aw 23.3 1245. 1196. 425.
## 5 Belo Oriente Aw 22.6 1054. 1105. 125.
## 6 Guanhaes Cwa 20.6 818. 866 188.
## 7 Eldorado do Sul Cfa 20.6 1787. 859. 10.7
## 8 Sao Gabriel Cfa 20.2 1783. 927. 2.8
## 9 Inhambupe As 24.2 715. 1318. 607.
## 10 Botucatu Cfb 21.9 1013. 1030 184.
## 11 Estrela do Sul Cwa 23.4 1134. 1208. 206.
## 12 Buri Cfa 20.0 1404. 805 3.6
## 13 Inocencia Am 24.5 1020. 1375. 389.
## 14 Chapadao do Sul Am 22.6 1027. 1098. 326.
## 15 Aracruz Aw 23.9 849. 1277. 336.
## 16 Tres Lagoas Aw 25.2 944. 1502. 538.
## 17 Tres Marias Aw 22.2 811. 1052. 405.
## 18 Peixe Aw 26.3 1208. 1630. 674.
## 19 Mogi Guacu Cwa 22.4 924. 1100. 210.
## 20 Brejinho de Nazare Aw 25.9 1507. 1563. 496.
## 21 Monte Dourado Am 27.4 2529. 1821. 530.
## 22 Otacilio Costa Cfb 16.9 2092. 548. 0
## 23 Telemaco Borba Cfa 18.5 1367 670. 28
## 24 Borebi Cfa 22.1 948. 1058. 201.
## 25 Coracao de Jesus As 23.9 413. 1275. 744.
## 26 Antonio Olinto Cfb 17.7 1740. 616. 0
## 27 Tres Barras Cfb 17.3 1123. 581 16.5
## 28 Urbano Santos Aw 27.0 1438. 1750. 935
## 29 Eunapolis Am 22.9 1419. 1128. 31.2
## 30 Itagimirim Aw 25.2 491. 1460. 870.
## 31 Bocaiuva Aw 23.9 413. 1276. 642.A análise de componentes principais nos mostra o quanto cada grupo de variáveis explicam a variabilidade total dados. No nosso caso, o primeiro componente responde por 72% da variabilidade e tem efeito quase que igual da temperatura (Tmed), evapotranspiração (ETP) e déficit hídrico (DEF). O segundo componente é majoritariamente o efeito da chuva (PPT). Juntos, os dois componentes explicam 95% dos dados.
pca <- select(dados, Tmed:DEF) %>%
princomp(cor = T)
summary(pca); loadings(pca)## Importance of components:
## Comp.1 Comp.2 Comp.3 Comp.4
## Standard deviation 1.7007490 0.9709348 0.40137904 0.0602784852
## Proportion of Variance 0.7231368 0.2356786 0.04027628 0.0009083739
## Cumulative Proportion 0.7231368 0.9588153 0.99909163 1.0000000000##
## Loadings:
## Comp.1 Comp.2 Comp.3 Comp.4
## Tmed 0.567 0.201 0.437 0.668
## PPT -0.241 0.933 -0.251
## ETP 0.561 0.284 0.261 -0.732
## DEF 0.553 -0.823
##
## 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.00A Figura 1 ajuda-nos a visualizar a disposição das cidades em função dos dois principais componentes. Se analisarmos por quadrantes, podemos agrupar as cidades de clima semelhante e ainda verificar a relação com as variáveis de clima. As setas indicam o efeito positivo ou negativo da variável. Por exemplo, o quadrante Q4 é caracterizado por valores altos de chuva e praticamente nenhum deficit hídrico. No oposto, temos o Q2 com baixa precipitação e alto déficit hídrico.
ggbiplot(pca) +
geom_point() +
geom_vline(xintercept = 0, size = 1.2, linetype = 6) +
geom_hline(yintercept = 0, size = 1.2, linetype = 6) +
geom_label_repel(aes(label = dados$Cidade), size = 3, nudge_x = .2) +
annotate(
"text",
x = c(-2, 2, 2, -2),
y = c(2, 2, -2, -2),
label = paste0("Q", 1:4), size = 6
) +
lims(x = c(-2,2), y = c(-2,2)) +
theme_bw()
Figura 1: Representação gráfica dos componentes principais.
Como também temos a informação do clima Koppen, podemos colorir o gráfico em função deste atributo (Figura 2).
ggbiplot(pca) +
geom_point(aes(color = dados$Koppen)) +
geom_vline(xintercept = 0, size = 1.2, linetype = 6) +
geom_hline(yintercept = 0, size = 1.2, linetype = 6) +
geom_label_repel(
aes(color = dados$Koppen, label = dados$Cidade),
size = 3, nudge_x = .2, show.legend = F
) +
lims(x = c(-2,2), y = c(-2,2)) +
scale_color_brewer("Clima Koppen", palette = "Dark2") +
theme_bw()+
theme(legend.position = "top")
Figura 2: Representação gráfica dos componentes principais com classificação Koppen.
Caso tenha alguma dúvida ou sugestão sobre o post, fique à vontade para fazer um comentário ou me contatar por E-mail.
sessioninfo::session_info(c("readr", "dplyr", "ggplot2", "ggrepel", "ggbiplot"))## - Session info ----------------------------------------------------------
## setting value
## version R version 3.5.3 (2019-03-11)
## os Windows 10 x64
## system x86_64, mingw32
## ui RTerm
## language (EN)
## collate Portuguese_Brazil.1252
## ctype Portuguese_Brazil.1252
## tz America/Sao_Paulo
## date 2019-08-25
##
## - Packages --------------------------------------------------------------
## package * version date lib source
## assertthat 0.2.1 2019-03-21 [1] CRAN (R 3.5.3)
## backports 1.1.4 2019-04-10 [1] CRAN (R 3.5.3)
## BH 1.69.0-1 2019-01-07 [1] CRAN (R 3.5.2)
## cli 1.1.0 2019-03-19 [1] CRAN (R 3.5.3)
## clipr 0.7.0 2019-07-23 [1] CRAN (R 3.5.3)
## colorspace 1.4-1 2019-03-18 [1] CRAN (R 3.5.3)
## crayon 1.3.4 2017-09-16 [1] CRAN (R 3.5.1)
## digest 0.6.20 2019-07-04 [1] CRAN (R 3.5.3)
## dplyr * 0.8.3 2019-07-04 [1] CRAN (R 3.5.3)
## ellipsis 0.2.0.1 2019-07-02 [1] CRAN (R 3.5.3)
## fansi 0.4.0 2018-10-05 [1] CRAN (R 3.5.1)
## ggbiplot * 0.55 2018-12-31 [1] Github (vqv/ggbiplot@7325e88)
## ggplot2 * 3.2.1 2019-08-10 [1] CRAN (R 3.5.3)
## ggrepel * 0.8.1 2019-05-07 [1] CRAN (R 3.5.3)
## glue 1.3.1 2019-03-12 [1] CRAN (R 3.5.3)
## gtable 0.3.0 2019-03-25 [1] CRAN (R 3.5.3)
## hms 0.5.0 2019-07-09 [1] CRAN (R 3.5.3)
## labeling 0.3 2014-08-23 [1] CRAN (R 3.5.0)
## lattice 0.20-38 2018-11-04 [2] CRAN (R 3.5.3)
## lazyeval 0.2.2 2019-03-15 [1] CRAN (R 3.5.3)
## magrittr 1.5 2014-11-22 [1] CRAN (R 3.5.1)
## MASS 7.3-51.1 2018-11-01 [2] CRAN (R 3.5.3)
## Matrix 1.2-17 2019-03-22 [1] CRAN (R 3.5.3)
## mgcv 1.8-28 2019-03-21 [1] CRAN (R 3.5.3)
## munsell 0.5.0 2018-06-12 [1] CRAN (R 3.5.1)
## nlme 3.1-137 2018-04-07 [2] CRAN (R 3.5.3)
## pillar 1.4.2 2019-06-29 [1] CRAN (R 3.5.3)
## pkgconfig 2.0.2 2018-08-16 [1] CRAN (R 3.5.1)
## plogr 0.2.0 2018-03-25 [1] CRAN (R 3.5.1)
## plyr * 1.8.4 2016-06-08 [1] CRAN (R 3.5.1)
## purrr 0.3.2 2019-03-15 [1] CRAN (R 3.5.3)
## R6 2.4.0 2019-02-14 [1] CRAN (R 3.5.2)
## RColorBrewer 1.1-2 2014-12-07 [1] CRAN (R 3.5.0)
## Rcpp 1.0.2 2019-07-25 [1] CRAN (R 3.5.3)
## readr * 1.3.1 2018-12-21 [1] CRAN (R 3.5.2)
## reshape2 1.4.3 2017-12-11 [1] CRAN (R 3.5.1)
## rlang 0.4.0 2019-06-25 [1] CRAN (R 3.5.3)
## scales * 1.0.0 2018-08-09 [1] CRAN (R 3.5.1)
## stringi 1.4.3 2019-03-12 [1] CRAN (R 3.5.3)
## stringr 1.4.0 2019-02-10 [1] CRAN (R 3.5.2)
## tibble 2.1.3 2019-06-06 [1] CRAN (R 3.5.3)
## tidyselect 0.2.5 2018-10-11 [1] CRAN (R 3.5.1)
## utf8 1.1.4 2018-05-24 [1] CRAN (R 3.5.1)
## vctrs 0.2.0 2019-07-05 [1] CRAN (R 3.5.3)
## viridisLite 0.3.0 2018-02-01 [1] CRAN (R 3.5.1)
## withr 2.1.2 2018-03-15 [1] CRAN (R 3.5.1)
## zeallot 0.1.0 2018-01-28 [1] CRAN (R 3.5.2)
##
## [1] C:/Users/Italo/Documents/R/win-library/3.5
## [2] C:/Program Files/R/R-3.5.3/library