PCA: reduzindo a dimensionalidade dos dados
Renata Brandt e Fábio Barros
July 19, 2015
PCA (Principal Component Analysis)
A Análise de Componentes Principais é indicada quando possuímos um conjunto de variáveis contínuas muito amplo. Após a análise, cada componente principal representa uma combinação linear das variáveis originais inseridas na matriz dados.
Os componentes principais agrupam elementos que variam num mesmo sentido, o que facilita a compreensão de como as variáveis interagem entre si. Eles têm a vantagem de serem independentes, podendo ser utilizados em testes de regressão linear, por exemplo.
Nós usaremos os dados de variáveis químicas de sedimento amostradas em 24 locais para exemplificar a construção de uma PCA. O conjunto de dados se encontra no arquivo varechem.csv.
# setwd("C:/Users/Fábio Barros/Desktop/CursoR")
# Para construir a PCA devemos definir que a primeira coluna
# do conjunto de dados referem-se aos nomes.
# Portanto, pedimos para abrir o arquivo com o comando row.names=1
solo <- read.csv("http://renatabrandt.github.io/EBC2015/data/varechem.csv", row.names=1)
is.data.frame(solo)## [1] TRUEUsando PRCOMP
Neste tópico, nós usaremos o pacote vegan - function prcomp - para construir a PCA. Existem dois tipos de PCA: a PCA de covariância e a PCA de correlação. Quando as variáveis foram medidas em unidades diferentes ou quando a variância de cada variável é muito diferente umas das outras devemos usar a PCA de correlação, que é construida usando os dados que foram padronizados para ter média 0 e desvio padrão 1.
Geralmente é exatamente esta situação que encontramos lidando com problemas biológicos!!!
Vamos observar as variâncias dos dados químicos do arquivo varechem?
round(apply(solo,2,var),2) ## N P K Ca Mg S Al Fe
## 30.56 223.39 4204.53 59332.17 1681.93 136.14 14962.28 3654.07
## Mn Zn Mo Baresoil Humdepth pH
## 1149.98 8.90 0.06 270.26 0.44 0.05Veja que a variância de cada variável é muito diferente. Portanto, temos que usar uma PCA de correlação para que a variável com maior variância não “domine” a análise.
Para construir a PCA de correlação usando prcomp devemos usar o argumento scale=T. Vamos agora rodar a análise no conjunto de dados solo?
install.packages(“vegan”)
require(vegan)## Loading required package: vegan## Warning: package 'vegan' was built under R version 3.1.3## Loading required package: permute## Warning: package 'permute' was built under R version 3.1.3## Loading required package: lattice## Warning: package 'lattice' was built under R version 3.1.3## This is vegan 2.3-0# para "rodar" a análise e visualizar a contribuição de
# cada variável no componente
prc.solo <- prcomp(solo, scale=TRUE); prc.solo## Standard deviations:
## [1] 2.2784976 1.7868529 1.2982707 1.0339221 0.9033240 0.8401250 0.6606103
## [8] 0.6073232 0.4132034 0.3866882 0.2919936 0.2643011 0.1872348 0.1536722
##
## Rotation:
## PC1 PC2 PC3 PC4 PC5
## N 0.09462709 0.12435859 -0.16204826 -0.79318367 0.351576656
## P -0.35589258 -0.20731069 0.09086948 -0.08242977 -0.219919629
## K -0.36318854 -0.16172669 -0.04259682 -0.21401371 -0.176886144
## Ca -0.35975485 -0.08457237 0.32644437 0.01826328 0.236121676
## Mg -0.35228219 -0.09776116 0.08467641 0.03784222 0.486392688
## S -0.33246862 -0.29668916 -0.21512714 0.01410909 -0.070039686
## Al 0.11552207 -0.48048232 -0.26540100 -0.07767669 -0.005266422
## Fe 0.20299687 -0.38301396 -0.15691339 -0.31769008 0.022858769
## Mn -0.28756460 0.19318776 0.08921960 -0.37885343 -0.489948050
## Zn -0.34603172 -0.16198093 0.04178581 0.05229606 0.150267479
## Mo -0.04154531 -0.22863241 -0.58319905 0.23286801 0.029395588
## Baresoil -0.16349721 0.31570370 -0.27274616 0.06312303 0.407367062
## Humdepth -0.23531857 0.33736134 -0.28352933 -0.02541797 0.046678588
## pH 0.16613845 -0.32216601 0.45164891 -0.06492250 0.259486172
## PC6 PC7 PC8 PC9 PC10
## N 0.31876043 -0.020445294 0.106079920 -0.029711941 0.25325314
## P 0.08346185 -0.202993610 -0.291937970 0.529497898 0.48691417
## K -0.30156737 -0.083967201 0.345810459 -0.347935071 0.19105757
## Ca -0.10026235 -0.056467690 -0.241640431 -0.132902952 0.15581868
## Mg -0.11335320 0.422483826 0.002896128 0.185900569 -0.18858761
## S -0.09229889 -0.001743922 0.271322793 -0.041707742 -0.01394397
## Al -0.20467312 -0.090888045 -0.037652261 -0.325530088 -0.02736515
## Fe -0.33950608 0.197580497 -0.192886023 0.424437201 -0.32652805
## Mn 0.07123733 -0.136543653 -0.274959999 -0.080205160 -0.54341919
## Zn 0.54306692 0.032492006 0.238361969 -0.005992104 -0.36679941
## Mo 0.41761742 -0.120157357 -0.352212047 -0.082184113 0.02092970
## Baresoil -0.32655495 -0.662017409 -0.006567573 0.150678212 -0.20591182
## Humdepth -0.16741045 0.425702015 -0.465061249 -0.315633320 0.12595591
## pH 0.05330910 -0.262135909 -0.374832351 -0.349371603 -0.07947962
## PC11 PC12 PC13 PC14
## N -0.005103292 0.03375847 0.12866179 -0.023577121
## P -0.212772662 -0.15369608 -0.16014570 -0.146803748
## K 0.118589614 0.28135388 -0.50806589 0.191151759
## Ca 0.669374213 -0.03834446 0.36929697 0.035843646
## Mg -0.239307771 0.42332112 -0.05980040 -0.345891242
## S -0.378600375 -0.07818118 0.63006839 0.345755661
## Al 0.018994101 -0.30010097 0.04257464 -0.654805487
## Fe 0.238879497 -0.09798831 -0.12533928 0.354078026
## Mn -0.035024496 0.18554028 0.13661147 -0.181857434
## Zn 0.104956310 -0.48566331 -0.29628189 0.070169749
## Mo 0.159957373 0.43567172 -0.03761246 0.117855258
## Baresoil -0.053374989 -0.10187147 -0.07709019 0.007445089
## Humdepth -0.168085863 -0.37001287 -0.13026298 0.156358945
## pH -0.406958380 0.04733921 -0.12013250 0.265236656# para acessar a estrutura do elemento "prc.solo"
str(prc.solo)## List of 5
## $ sdev : num [1:14] 2.278 1.787 1.298 1.034 0.903 ...
## $ rotation: num [1:14, 1:14] 0.0946 -0.3559 -0.3632 -0.3598 -0.3523 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:14] "N" "P" "K" "Ca" ...
## .. ..$ : chr [1:14] "PC1" "PC2" "PC3" "PC4" ...
## $ center : Named num [1:14] 22.4 45.1 162.9 569.7 87.5 ...
## ..- attr(*, "names")= chr [1:14] "N" "P" "K" "Ca" ...
## $ scale : Named num [1:14] 5.53 14.95 64.84 243.58 41.01 ...
## ..- attr(*, "names")= chr [1:14] "N" "P" "K" "Ca" ...
## $ x : num [1:24, 1:14] 1.539 4.255 -0.689 3.634 3.373 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:24] "local01" "local02" "local03" "local04" ...
## .. ..$ : chr [1:14] "PC1" "PC2" "PC3" "PC4" ...
## - attr(*, "class")= chr "prcomp"# para visualizar o número de componentes e as respectivas
# contribuições para a variância total.
screeplot(prc.solo)
# para mostrar a porcentagem de variância capturada por cada eixo
summary(prc.solo)## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6
## Standard deviation 2.2785 1.7869 1.2983 1.03392 0.90332 0.84013
## Proportion of Variance 0.3708 0.2281 0.1204 0.07636 0.05829 0.05042
## Cumulative Proportion 0.3708 0.5989 0.7193 0.79564 0.85392 0.90434
## PC7 PC8 PC9 PC10 PC11 PC12
## Standard deviation 0.66061 0.60732 0.4132 0.38669 0.29199 0.26430
## Proportion of Variance 0.03117 0.02635 0.0122 0.01068 0.00609 0.00499
## Cumulative Proportion 0.93551 0.96185 0.9740 0.98473 0.99082 0.99581
## PC13 PC14
## Standard deviation 0.1872 0.15367
## Proportion of Variance 0.0025 0.00169
## Cumulative Proportion 0.9983 1.00000# para fazer o plot da PCA
biplot(prc.solo)