4 Análise de Componentes Principais
A Análise de Componentes Principais (ACP) é uma técnica para estudar e reorganizar a variabilidade de um conjunto de variáveis. A ACP constrói um novo sistema de coordenadas para os dados, possivelmente simplificando as análises. Os novos eixos são combinações lineares das variáveis originais, são mutuamente ortogonais e aparecem em ordem decrescente de variância.
Definição 4.1 Dado um vetor aleatório p-dimensional \(\boldsymbol{x}= (X_1, \dots, X_p)^T\), com matriz de covariância \(\boldsymbol{\Sigma}\), os componentes principais (CPs) são um novo sistema de coordenadas \(Y_1, \dots, Y_p\), tal que \[ \begin{aligned} Y_1 &= \boldsymbol{u}_1^{T} \boldsymbol{x}= u_{11}X_1 + u_{12}X_2 + \dots + u_{1p}X_p \\ Y_2 &= \boldsymbol{u}_2^{T} \boldsymbol{x}= u_{21}X_1 + u_{22}X_2 + \dots + u_{2p}X_p \\ &\vdots \\ Y_p &= \boldsymbol{u}_p^{T} \boldsymbol{x}= u_{p1}X_1 + u_{p2}X_2 + \dots + u_{pp}X_p \end{aligned} \]
Onde \(\boldsymbol{u}_i^T\) são vetores que definem as direções das novas coordenadas. Os CPs tem duas propriedades importantes:
- São ortogonais entre sí, ou seja, para \(i \neq j\), \(\boldsymbol{u}_i^T \boldsymbol{u}_j = 0\) e consequentemente \(\operatorname{Var}\left(Y_i, Y_j\right) = 0\)
- Suas variâncias são decrescentes, isto é, \(\operatorname{Var}\left(Y_1\right) \ge \operatorname{Var}\left(Y_2\right) \ge \dots \ge \operatorname{Var}\left(Y_p\right)\).
Essa mudança de coordenadas pode ter dois propósitos relacionados: descrever a estrutura de dependência entre as variáveis e/ou reduzir a dimensionalidade. Frequentemente, a ACP é uma ferramenta intermediária em uma análise. Isso porque muitas vezes a análise dos dados observados no novo sistema de coordenadas é mais conveniente. Uma característica importantíssima para tal conveniência é a independência.
Importante
A variabilidade que a ACP reorganiza é medida em torno da média. Por isso, antes de construir os componentes, deslocamos a nuvem de dados para que seu centro coincida com a origem. No contexto populacional, o vetor centrado é \(\boldsymbol{x}^c=\boldsymbol{x}-\boldsymbol{\mu}\); em uma amostra, cada observação é centralizada por \(\boldsymbol{x}_i^c=\boldsymbol{x}_i-\bar{\boldsymbol{x}}\). Esse deslocamento não altera a forma da nuvem nem as distâncias entre as observações, apenas seu ponto de referência. Assim, as projeções nos novos eixos descrevem desvios em relação ao centro dos dados.
Para simplificar a notação, daqui em diante \(\boldsymbol{x}\) representa o vetor de variáveis já centralizado.
Exemplo 4.1 Considere dados fictícios de 160 estudantes matriculados em uma mesma disciplina. Para cada estudante, foram registradas a nota na prova teórica (\(X_1\)) e a nota em um trabalho prático (\(X_2\)), ambas na escala de 0 a 10. Aplicando a centralização descrita acima, subtraímos de cada nota a média da respectiva avaliação. Seguindo a convenção adotada no capítulo, continuaremos chamando as variáveis centralizadas de \(X_1\) e \(X_2\). Assim, valores positivos indicam desempenho acima da média da turma e valores negativos, desempenho abaixo da média. A Figura 4.1 (esquerda) mostra essas duas variáveis, que apresentam forte associação positiva; cada ponto representa um estudante.
Os eixos \(X_1\) e \(X_2\) são perfeitamente válidos, mas não acompanham a orientação principal da nuvem. A dispersão é muito maior ao longo de uma direção diagonal: estudantes com nota acima da média na prova tendem também a ficar acima da média no trabalho, e o mesmo ocorre com desempenhos abaixo da média. Isso sugere que podemos definir um novo eixo, isto é, uma nova dimensão, para descrever cada estudante combinando as notas na prova e no trabalho. A Figura 4.1 (direita) ilustra essa nova dimensão, que vamos chamar de \(Y_1\).
4.1 Variância como Medida de Informação
No Exemplo 4.1 Note que o eixo \(Y_1\) acompanha a direção de maior alongamento da nuvem de pontos. Em estatística, a variância é frequentemente usada como uma medida de informação. Uma variável com alta variância indica que seus valores são bem espalhados, o que nos ajuda a diferenciar as observações. Se a variância fosse zero, todos os pontos seriam idênticos, não nos fornecendo nenhuma informação sobre suas diferenças.
Definição 4.2 A variância total de um conjunto de dados com \(p\) variáveis é a soma das variâncias de cada variável individual. Matematicamente, se \(\boldsymbol{x}= (X_1, \dots, X_p)^T\) é o vetor de variáveis aleatórias com matriz de covariâncias \(\boldsymbol{\Sigma}\), a variância total é definida como:
\[ \text{Variância Total} = \sum_{j=1}^{p} \operatorname{Var}\left(X_j\right) = \sum_{j=1}^{p} \sigma_{jj} = \operatorname{tr}\left(\boldsymbol{\Sigma}\right) \]
onde \(\sigma_{jj}\) é a variância da \(j\)-ésima variável e \(\operatorname{tr}\left(\boldsymbol{\Sigma}\right)\) é o traço da matriz de covariâncias (a soma dos elementos da diagonal principal). Essa medida representa a dispersão total na nuvem de pontos, somando a variabilidade em cada uma das direções dos eixos originais.
4.2 Interpretação e solução geométrica
Seguindo o Exemplo 4.1. Propositalmente, a Figura 4.1 mostra que a busca pelo eixo mais informativo \(Y_1\) pode ser traduzido como a busca um um ângulo ideal \(\theta^*\) que define uma rotação do eixo original \(X_1\). Para um ângulo \(\theta\) qualquer, a rotação pode ser expressa pela matriz de rotação
\[ \begin{bmatrix} \boldsymbol{u}_1(\theta) & \boldsymbol{u}_2(\theta) \end{bmatrix} = \begin{bmatrix} \cos\theta & -\sin\theta\\ \sin\theta & \phantom{-}\cos\theta \end{bmatrix} \]
cujas colunas \(\boldsymbol{u}_1(\theta),\, \boldsymbol{u}_2(\theta)\) formam uma base ortonormal, isso garante que as direções \(Y_1\) e \(Y_2\) sejam perpendiculares entre si. As projeções dos dados originais nesse sistema são
\[ \begin{aligned} Y_1 &= \boldsymbol{u}_1(\theta)^T\boldsymbol{x} =\cos(\theta)X_1+\sin(\theta)X_2 \\ Y_2 &= \boldsymbol{u}_2(\theta)^T\boldsymbol{x} =-\sin(\theta)X_1+\cos(\theta)X_2. \end{aligned} \]
Ainda não sabemos como obter o ângulo ideal \(\theta^*\). Podemos, entretanto, formular precisamente o que significa alinhar o primeiro eixo com a maior dispersão. A variância da coordenada \(Y_1\) é
\[ \begin{aligned} \operatorname{Var}\left(Y_1\right) &= \operatorname{Var}\left(\boldsymbol{u}_1(\theta)^T\boldsymbol{x}\right) \\ &= \boldsymbol{u}_1(\theta)^T\boldsymbol{\Sigma}\boldsymbol{u}_1(\theta) \\ &= \begin{bmatrix} \cos\theta \sin\theta \end{bmatrix} \begin{bmatrix} \sigma_{11} & \sigma_{12} \\ \sigma_{12} & \sigma_{22} \end{bmatrix} \begin{bmatrix} \cos\theta \\ \sin\theta \end{bmatrix} \\ &= \sigma_{11}\cos^2\theta + 2\sigma_{12}\sin\theta\cos\theta + \sigma_{22}\sin^2\theta. \end{aligned} \]
em que \(\boldsymbol{\Sigma}\) é a matriz de covariâncias de \((X_1, X_2)^T\). Portanto, cada rotação produz uma variância diferente para o primeiro eixo. A ACP escolhe a rotação que maximiza essa função. Nesse exemplo simples, podemos observar isso gráficamente. A Figura 4.2 ilustra a solução.
Nota
Em duas dimensões, existe também uma solução trigonométrica para esse problema de maximização. Essa solução contudo, não se estende de maneira conveniente a dezenas ou centenas de variáveis, por isso não prosseguiremos com ela.
As variâncias das notas nos eixos originais são \(\operatorname{Var}\left(X_1\right) = 1.08\) e \(\operatorname{Var}\left(X_2\right) = 1.0\). Consequentemente, a variância total é 2.08.
Uma vez fixado o ângulo ideal \(\theta^*\), O eixo \(Y_1\) é exatamente a direção de maior variabilidade, e consequentemente, o mais informativo. Para calcular a variância dos dados nesse novo eixo consideramos as projeções ortogonais das observações sobre ele. Essas projeções são definidas exatamente pela matriz de rotação. Para um ponto \(\boldsymbol{x}_i = (x_1, x_2)\) arbitrário, a projeção sobre \(Y_1\) é dada por \(\boldsymbol{y}_i = \cos(\theta)x_1+\sin(\theta)x_2\). A Figura 4.1 (direita) ilustra a projeção. A variância dos dados projetados sobre \(Y_1\) é
\[ \begin{aligned} \operatorname{Var}\left(Y_1\right) &= \frac{1}{n} \sum_{i=1}^{n} (\boldsymbol{y}_i - \bar{Y}_1)^2 \\ &= 1.92 \end{aligned} \]
Note que \(Y_1\) possui uma variância superior a \(X_1\) e \(X_2\). Com relação a variância total, temos o resíduo \(\text{Variância Total} - \operatorname{Var}\left(Y_1\right) = 0.16\). É aqui que introduzimos o segundo novo eixo \(Y_2\). Ele é perpendicular a \(Y_1\) e captura exatamente essa variância residual \(\operatorname{Var}\left(Y_2\right) = 0.16\).
Fica claro como os novos eixos \(Y_1\) e \(Y_2\) definem um par de coordenadas que redistribui a informação dos dados. O primeiro eixo \(Y_1\) captura a maior parte da variabilidade, enquanto o segundo \(Y_2\) captura a variância residual.
Nota
Não é possível, e nem necessário, que o leitor verifique o calculo dessas variâncias, pois no exemplo apresentamos apenas a nuvem de pontos, mas não os dados originais.
4.3 A Solução dos componentes principais
Para encontrar o primeiro componente principal, partimos da combinação linear \(Y_1 = \boldsymbol{u}_1^T\boldsymbol{x}\), em que \(\boldsymbol{u}_1\) é, por enquanto, apenas uma direção candidata. Sua variância é \(\operatorname{Var}\left(Y_1\right) = \operatorname{Var}\left(\boldsymbol{u}_1^T\boldsymbol{x}\right) = \boldsymbol{u}_1^T \boldsymbol{\Sigma}\boldsymbol{u}_1\), onde \(\boldsymbol{\Sigma}\) é a matriz de covariâncias de \(\boldsymbol{x}\).
Para evitar que a variância seja aumentada simplesmente inflando os coeficientes em \(\boldsymbol{u}_1\), impomos a restrição de que seu comprimento seja unitário, \(\boldsymbol{u}_1^T\boldsymbol{u}_1 = 1\). Formalmente, o problema de maximização para o primeiro componente principal se torna:
\[ \begin{aligned} \max_{\boldsymbol{u}_1} \quad & \boldsymbol{u}_1^T \boldsymbol{\Sigma}\boldsymbol{u}_1 \\ \text{sujeito a} \quad & \boldsymbol{u}_1^T \boldsymbol{u}_1 = 1 \end{aligned} \]
Utilizando o método dos multiplicadores de Lagrange, a função a ser maximizada é:
\[ L(\boldsymbol{u}_1, \lambda_1) = \boldsymbol{u}_1^T \boldsymbol{\Sigma}\boldsymbol{u}_1 - \lambda_1 (\boldsymbol{u}_1^T \boldsymbol{u}_1 - 1) \]
Derivando em relação a \(\boldsymbol{u}_1\) e igualando a zero, obtemos:
\[ \frac{\partial L}{\partial \boldsymbol{u}_1} = 2 \boldsymbol{\Sigma}\boldsymbol{u}_1 - 2 \lambda_1 \boldsymbol{u}_1 = 0 \implies \boldsymbol{\Sigma}\boldsymbol{u}_1 = \lambda_1 \boldsymbol{u}_1 \]
Note que chegamos à equação fundamental de autovalores e autovetores. Esta equação mostra que a direção candidata \(\boldsymbol{u}_1\) que resolve o problema deve ser um autovetor da matriz de covariâncias \(\boldsymbol{\Sigma}\). Para encontrar a variância, pré-multiplicamos a equação por \(\boldsymbol{u}_1^T\):
\[ \boldsymbol{u}_1^T \boldsymbol{\Sigma}\boldsymbol{u}_1 = \lambda_1 \boldsymbol{u}_1^T \boldsymbol{u}_1 \]
Como \(\operatorname{Var}\left(Y_1\right) = \boldsymbol{u}_1^T \boldsymbol{\Sigma}\boldsymbol{u}_1\) e a restrição é \(\boldsymbol{u}_1^T \boldsymbol{u}_1 = 1\), temos:
\[ \operatorname{Var}\left(Y_1\right) = \lambda_1 \]
Para maximizar a variância de \(Y_1\), devemos escolher o maior autovalor possível. Portanto, \(\lambda_1\) é o maior autovalor de \(\boldsymbol{\Sigma}\). Denotamos por \(\boldsymbol{e}_1\) um autovetor unitário associado a \(\lambda_1\); na solução do problema, a direção candidata se torna \(\boldsymbol{u}_1 = \boldsymbol{e}_1\).
Nota
A matriz de covariâncias \(\boldsymbol{\Sigma}\) é, por construção, simétrica e positiva semidefinida. Conforme discutido em Teorema 2.1, o Teorema Espectral garante que seus autovalores são reais e não negativos e que existe uma base ortonormal de autovetores.
Para obter o segundo componente, consideramos uma nova direção candidata \(\boldsymbol{u}_2\) e escrevemos \(Y_2 = \boldsymbol{u}_2^T\boldsymbol{x}\). Como \(\boldsymbol{e}_1\) já determina a direção de maior variância, exigimos que \(\boldsymbol{u}_2\) tenha comprimento unitário e seja ortogonal a \(\boldsymbol{e}_1\). Essa ortogonalidade também garante que os componentes sejam não correlacionados: \(\boldsymbol{\Sigma}\), \(\operatorname{Cov}\left(Y_1,Y_2\right) = \boldsymbol{e}_1^T\boldsymbol{\Sigma}\boldsymbol{u}_2 = \lambda_1\boldsymbol{e}_1^T\boldsymbol{u}_2 = 0\). Assim, o segundo problema de otimização é
\[ \begin{aligned} \max_{\boldsymbol{u}_2} \quad & \boldsymbol{u}_2^T \boldsymbol{\Sigma}\boldsymbol{u}_2 \\ \text{sujeito a} \quad & \begin{cases} \boldsymbol{u}_2^T \boldsymbol{u}_2 = 1 \\ \boldsymbol{e}_1^T \boldsymbol{u}_2 = 0 \end{cases} \end{aligned} \]
A restrição de ortogonalidade introduz um segundo multiplicador de Lagrange, \(\phi\). A função Lagrangiana e sua derivada são
\[ \begin{aligned} L(\boldsymbol{u}_2, \lambda_2, \phi) &= \boldsymbol{u}_2^T\boldsymbol{\Sigma}\boldsymbol{u}_2 - \lambda_2(\boldsymbol{u}_2^T\boldsymbol{u}_2 - 1) - \phi\boldsymbol{e}_1^T\boldsymbol{u}_2, \\ \frac{\partial L}{\partial \boldsymbol{u}_2} &= 2\boldsymbol{\Sigma}\boldsymbol{u}_2 - 2\lambda_2\boldsymbol{u}_2 - \phi\boldsymbol{e}_1 = \boldsymbol{0}. \end{aligned} \]
Pré-multiplicamos a derivada por \(\boldsymbol{e}_1^T\). Como \(\boldsymbol{e}_1^T\boldsymbol{\Sigma}= \lambda_1\boldsymbol{e}_1^T\), \(\boldsymbol{e}_1^T\boldsymbol{u}_2 = 0\) e \(\boldsymbol{e}_1^T\boldsymbol{e}_1 = 1\), obtemos
\[ \begin{aligned} 0 &= 2\boldsymbol{e}_1^T\boldsymbol{\Sigma}\boldsymbol{u}_2 - 2\lambda_2\boldsymbol{e}_1^T\boldsymbol{u}_2 - \phi\boldsymbol{e}_1^T\boldsymbol{e}_1 \\ &= 2\lambda_1\boldsymbol{e}_1^T\boldsymbol{u}_2 - 2\lambda_2\boldsymbol{e}_1^T\boldsymbol{u}_2 - \phi \\ &= -\phi. \end{aligned} \]
Logo, \(\phi=0\), e a derivada se reduz a \(\boldsymbol{\Sigma}\boldsymbol{u}_2 = \lambda_2\boldsymbol{u}_2\). A segunda direção também deve, portanto, ser um autovetor de \(\boldsymbol{\Sigma}\). Entre os autovetores ortogonais a \(\boldsymbol{e}_1\), a maior variância disponível é o segundo maior autovalor, \(\lambda_2\). Denotando o autovetor unitário correspondente por \(\boldsymbol{e}_2\), a solução é \(\boldsymbol{u}_2=\boldsymbol{e}_2\).
Este processo continua: o \(k\)-ésimo componente principal (\(Y_k\)) é definido pelo autovetor \(\boldsymbol{e}_k\) associado ao \(k\)-ésimo maior autovalor \(\lambda_k\), garantindo que \(\operatorname{Var}\left(Y_k\right) = \lambda_k\) e que todos os componentes sejam mutuamente não correlacionados.
A decomposição espectral reúne todos esses resultados de uma só vez. Na expressão \(\boldsymbol{\Sigma}= \boldsymbol{E}\Lambda\boldsymbol{E}^T\), as colunas de \(\boldsymbol{E}\) são as direções dos componentes principais e os elementos da diagonal de \(\Lambda\) são suas variâncias. Ordenando os autovalores do maior para o menor, obtemos simultaneamente todos os componentes principais na ordem desejada.
Proposição 4.1 (Conservação da variância total) Se \(Y_1,\dots,Y_p\) são todos os componentes principais de \(\boldsymbol{x}\), então a soma de suas variâncias é igual à variância total das variáveis originais:
\[ \sum_{k=1}^p \operatorname{Var}\left(Y_k\right)=\sum_{j=1}^p \operatorname{Var}\left(X_j\right). \]
Comprovação. Pela Definição 4.2, a variância total das variáveis originais é \(\operatorname{tr}\left(\boldsymbol{\Sigma}\right)\). Usando a decomposição espectral \(\boldsymbol{\Sigma}=\boldsymbol{E}\Lambda\boldsymbol{E}^T\), as propriedades do traço, a ortogonalidade de \(\boldsymbol{E}\) e o fato de que \(\operatorname{Var}\left(Y_k\right)=\lambda_k\), obtemos
\[ \begin{aligned} \sum_{j=1}^p \operatorname{Var}\left(X_j\right) &= \operatorname{tr}\left(\boldsymbol{\Sigma}\right) \\ &= \operatorname{tr}\left(\boldsymbol{E}\Lambda\boldsymbol{E}^T\right) \\ &= \operatorname{tr}\left(\Lambda\boldsymbol{E}^T\boldsymbol{E}\right) \\ &= \operatorname{tr}\left(\Lambda\right) \\ &= \sum_{k=1}^p \lambda_k \\ &= \sum_{k=1}^p \operatorname{Var}\left(Y_k\right). \end{aligned} \]
Portanto, a ACP não cria nem elimina variabilidade; ela apenas a redistribui entre os componentes principais.
Nota
A demonstração acima, utilizando multiplicadores de Lagrange, é uma maneira moderna e elegante de conduzir a derivação do problema de maximização. Uma abordagem clássica restringe a norma de \(\boldsymbol{u}_1\) através do quociente,
\[ \operatorname{Var}\left(Y_1\right) = \max_{\boldsymbol{u}_1} \frac{\boldsymbol{u}_1^T \boldsymbol{\Sigma}\boldsymbol{u}_1}{\boldsymbol{u}_1^T \boldsymbol{u}_1} \]
Este é um problema clássico na álgebra linear. Um teorema fundamental afirma que para qualquer matriz simétrica \(A\), o máximo da forma quadrática \(\boldsymbol{x}^T A \boldsymbol{x}\), sujeito à restrição \(\boldsymbol{x}^T \boldsymbol{x}= 1\), é o maior autovalor de \(A\). O vetor \(\boldsymbol{x}\) que atinge esse máximo é o autovetor correspondente. Como a matriz de covariâncias \(\boldsymbol{\Sigma}\) é simétrica, este teorema se aplica diretamente ao nosso problema.
4.4 Componentes Principais Amostrais
A derivação anterior foi formulada em termos populacionais, supondo conhecidos a matriz de covariâncias \(\boldsymbol{\Sigma}\) e seus autovalores e autovetores. O Exemplo 4.1, por outro lado, já utilizou uma amostra para construir a interpretação geométrica. Resta explicitar como a solução teórica é aplicada nesse contexto.
Em uma amostra, centralizamos as observações e substituímos \(\boldsymbol{\Sigma}\) pela matriz de covariâncias amostral \(\boldsymbol{S}\). Os componentes principais são obtidos da mesma maneira, e as quantidades resultantes são estimativas de seus análogos populacionais:
- O \(k\)-ésimo autovalor amostral, \(\hat{\lambda}_k\), é uma estimativa de \(\lambda_k\).
- O \(k\)-ésimo autovetor amostral, \(\hat{\boldsymbol{e}}_k\), é uma estimativa de \(\boldsymbol{e}_k\).
- O \(k\)-ésimo componente principal amostral, \(\hat{Y}_k = \hat{\boldsymbol{e}}_k^T \boldsymbol{x}\), é uma estimativa de \(Y_k\).
A teoria e a interpretação permanecem as mesmas. Para simplificar a notação, ao longo deste capítulo, omitimos o acento circunflexo (\(\hat{\phantom{a}}\)), mas é importante lembrar que, na aplicação prática, estamos sempre lidando com estimativas amostrais.
Outra questão importante é a escala das variáveis. Retome os dados do Exemplo 4.1 e imagine que, além das duas notas, registramos o tempo de estudo no semestre, medido em horas. As notas são expressas em pontos e o tempo de estudo em horas, de modo que suas variâncias têm unidades diferentes e não são diretamente comparáveis. Se por acaso o tempo fosse registrado em minutos ou até segundos, sua variância seria muito maior.
Importante
Para evitar que as escalas tenham esse efeito, padronizamos as variáveis dividindo cada uma por seu desvio padrão. Todas passam a ter variância igual a 1. Devemos notar que essa padronização é exatamente aquela que transforma a matriz de covariâncias \(\boldsymbol{S}\) em uma matriz de correlações \(\boldsymbol{R}\). Ou seja, para evitar efeitos de escala, aplicamos a decomposição sobre a matriz de covariâncias \(\boldsymbol{S}\).
Exemplo 4.2 (Um Cálculo Manual Completo) Até aqui, a decomposição espectral resolveu o problema de ACP em termos abstratos, com \(\boldsymbol{\Sigma}\) e seus autovalores e autovetores. Mas como isso se traduz quando temos dados reais em mãos? Vamos resolver um problema pequeno inteiramente à mão, com uma amostra cujos números foram escolhidos para deixar as contas simples.
Suponha que observamos \(n = 5\) pares de valores já centralizados de duas variáveis, \(X_1\) e \(X_2\), organizados na matriz de dados \(\boldsymbol{X}\) (uma observação por linha, como em Equação 1.1):
\[ \boldsymbol{X}= \begin{bmatrix} -2 & -2 \\ -1 & 1 \\ 0 & 0 \\ 1 & -1 \\ 2 & 2 \end{bmatrix} \]
Matriz de covariâncias amostral. Como os dados já estão centralizados, a soma de cada coluna de \(\boldsymbol{X}\) é zero. As variâncias e a covariância amostrais são:
\[ s_{11} = \frac{\sum x_{i1}^2}{n-1} = \frac{10}{4} = 2.5, \qquad s_{22} = \frac{\sum x_{i2}^2}{n-1} = \frac{10}{4} = 2.5, \qquad s_{12} = \frac{\sum x_{i1}x_{i2}}{n-1} = \frac{6}{4} = 1.5 \]
\[ \boldsymbol{S}= \begin{bmatrix} 2.5 & 1.5 \\ & 2.5 \end{bmatrix} \]
Nota
Neste exemplo, vamos assumir que \(X_1\) e \(X_2\) tem escalas comparáveis, e por isso, nao padronizaremos os dados, ou seja, vamos aplicar a decomposição diretamente em \(\boldsymbol{\Sigma}\).
Autovalores e Autovetores. Resolvemos \(\det(\boldsymbol{S}- \lambda \boldsymbol{I}) = 0\):
\[ (2.5-\lambda)^2 - 1.5^2 = 0 \quad \Longrightarrow \quad 2.5 - \lambda = \pm 1.5 \quad \Longrightarrow \quad \lambda_1 = 4, \quad \lambda_2 = 1 \]
Note que \(\lambda_1 + \lambda_2 = 5 = \operatorname{tr}\left(\boldsymbol{S}\right)\), como esperado.
Para cada autovalor, substituímos \(\lambda_k\) no sistema \((\boldsymbol{S}-\lambda_k\boldsymbol{I})\boldsymbol{e}_k=\boldsymbol{0}\). Os sistemas resultantes e uma solução unitária para cada um são
\[ \begin{aligned} \lambda_1=4:\quad &\begin{bmatrix} -1.5 & 1.5\\ 1.5 & -1.5 \end{bmatrix} \begin{bmatrix}e_{11}\\e_{21}\end{bmatrix} =\begin{bmatrix}0\\0\end{bmatrix}, \\ &e_{21}=e_{11} \quad\Longrightarrow\quad \boldsymbol{e}_1=\frac{1}{\sqrt{2}} \begin{bmatrix}1\\1\end{bmatrix}; \\[0.8em] \lambda_2=1:\quad &\begin{bmatrix} 1.5 & 1.5\\ 1.5 & 1.5 \end{bmatrix} \begin{bmatrix}e_{12}\\e_{22}\end{bmatrix} =\begin{bmatrix}0\\0\end{bmatrix}, \\ &e_{22}=-e_{12} \quad\Longrightarrow\quad \boldsymbol{e}_2=\frac{1}{\sqrt{2}} \begin{bmatrix}1\\-1\end{bmatrix}. \end{aligned} \]
Os vetores com sinais opostos também seriam soluções válidas, pois representam os mesmos eixos com orientação invertida.
Expressões dos componentes principais. Substituindo os autovetores na equação \(Y_k = \boldsymbol{e}_k^T \boldsymbol{x}\), obtemos as combinações lineares explícitas dos dois componentes:
\[ \begin{aligned} Y_1 &= \frac{1}{\sqrt{2}} X_1 + \frac{1}{\sqrt{2}} X_2 \approx 0.707 X_1 + 0.707 X_2 \\ Y_2 &= \frac{1}{\sqrt{2}} X_1 - \frac{1}{\sqrt{2}} X_2 \approx 0.707 X_1 - 0.707 X_2 \end{aligned} \]
em que \(X_1\) e \(X_2\) representam as variáveis centralizadas.
Banco de dados transformado. Reunindo os autovetores em \(\boldsymbol{E}=[\boldsymbol{e}_1\ \boldsymbol{e}_2]\), obtemos a matriz de dados expressa no sistema de coordenadas dos componentes principais:
\[ \begin{aligned} \boldsymbol{T} &=\boldsymbol{X}\boldsymbol{E} \\ &=\boldsymbol{X}\frac{1}{\sqrt{2}} \begin{bmatrix} 1 & 1\\ 1 & -1 \end{bmatrix} \\ &=\frac{1}{\sqrt{2}} \begin{bmatrix} -4 & 0\\ 0 & -2\\ 0 & 0\\ 0 & 2\\ 4 & 0 \end{bmatrix}. \end{aligned} \]
As linhas de \(\boldsymbol{T}\) continuam representando as cinco observações, enquanto suas colunas agora correspondem a \(Y_1\) e \(Y_2\). O elemento \(t_{ik}\) é chamado escore da observação \(i\) no componente \(k\).
4.5 Interpretando os Componentes Principais
Os coeficientes \(e_{jk}\) do autovetor \(\boldsymbol{e}_k\), frequentemente chamados de cargas (loadings), indicam o peso de cada variável original \(X_j\) na formação do componente \(Y_k\). A magnitude de \(e_{jk}\) reflete a importância relativa da variável \(j\) no componente \(k\). Infelizmente, a interpretação direta desses coeficientes não é simples.
Uma medida mais interpretável é a correlação entre os componentes principais e as variáveis originais, \(\operatorname{Corr}\left(Y_k, X_j\right)\). Ela nos diz o quão “alinhado” um componente está com cada variável original, numa escala padronizada de -1 a 1.
Proposição 4.2 (Correlação entre componentes principais e variáveis originais) Se \(Y_k=\boldsymbol{e}_k^T\boldsymbol{x}\) é o \(k\)-ésimo componente principal, associado ao autovalor \(\lambda_k>0\), e \(\sigma_{jj}=\operatorname{Var}\left(X_j\right)>0\), então
\[ \operatorname{Corr}\left(Y_k,X_j\right) = \frac{e_{jk}\sqrt{\lambda_k}}{\sqrt{\sigma_{jj}}}, \]
em que \(e_{jk}\) é a carga da variável \(X_j\) no autovetor \(\boldsymbol{e}_k\).
Comprovação. Escrevendo o componente como \(Y_k=\sum_{\ell=1}^p e_{\ell k}X_\ell\), podemos calcular diretamente sua covariância com \(X_j\). Na derivação abaixo, a terceira igualdade usa a simetria da matriz de covariâncias, enquanto a quarta corresponde à \(j\)-ésima equação do sistema \(\boldsymbol{\Sigma}\boldsymbol{e}_k=\lambda_k\boldsymbol{e}_k\). Para passar da covariância à correlação, usamos \(\operatorname{Var}\left(Y_k\right)=\lambda_k\) e \(\operatorname{Var}\left(X_j\right)=\sigma_{jj}\):
\[ \begin{aligned} \operatorname{Cov}\left(Y_k,X_j\right) &=\sum_{\ell=1}^p e_{\ell k}\operatorname{Cov}\left(X_\ell,X_j\right) \\ &=\sum_{\ell=1}^p e_{\ell k}\sigma_{\ell j} \\ &=\sum_{\ell=1}^p \sigma_{j\ell}e_{\ell k} \\ &=\lambda_k e_{jk}, \\[0.5em] \operatorname{Corr}\left(Y_k,X_j\right) &=\frac{\operatorname{Cov}\left(Y_k,X_j\right)} {\sqrt{\operatorname{Var}\left(Y_k\right)}\sqrt{\operatorname{Var}\left(X_j\right)}} \\ &=\frac{e_{jk}\lambda_k} {\sqrt{\lambda_k}\sqrt{\sigma_{jj}}} \\ &=\frac{e_{jk}\sqrt{\lambda_k}}{\sqrt{\sigma_{jj}}}. \end{aligned} \]
Note que a correlação reescala a carga \(e_{jk}\) pelo desvio padrão do componente, \(\sqrt{\lambda_k}\), e pelo desvio padrão da variável original, \(\sqrt{\sigma_{jj}}\). Isso facilita a comparação entre componentes. Quando a ACP é realizada sobre a matriz de correlação (ou seja, com dados padronizados), as variâncias \(\sigma_{jj}\) são todas iguais a 1. Nesse caso, a fórmula simplifica para \(\operatorname{Corr}\left(Y_k,X_j\right)=e_{jk}\sqrt{\lambda_k}\).
Para facilitar ainda mais a interpretação, existem representações visuais úteis para a ACP. A mais comum é o biplot. Ele combina um gráfico de dispersão das observações projetadas nos componentes principais com os autovetores, permitindo visualizar tanto as direções dos componentes quanto as cargas das observações. Vamos ilustrar e interpretar um biplot no Exemplo 4.3, mais tarde nesse capítulo.
4.6 Escolhendo o Número de Componentes
A ACP expressa os dados em um sistema de eixos rotacionado. A Proposição 4.1 mostra que esse processo não perde nenhuma informação. Agora, quando o interesse é na redução de dimensionalidade, precisamos aceitar descartar parte da informação.
Na prática, escolhemos reter \(q < p\) componentes principais, descartando \(Y_{q+1}, \dots, Y_p\). As propriedades dos componentes principais garante que esse descarte é tal que a informação preservada é maximizada.
Para que essa redução seja adequada, é preciso escolher \(q\) cuidadosamente. Essa escolha é uma balança: poucos componentes tornam a análise mais simples e fácil de interpretar, enquanto mais componentes preservam mais da variância original. Existem na literatura diversos critérios para escolha de \(q\), cada um com suas particularidades. Os mais comuns são:
O Critério da Variância Explicada Acumulada é o mais comum deles. Calculamos a proporção da variância total explicada por cada componente e acumulamos essa proporção.
\[ \text{Proporção da Variância por } CP_k = \frac{\lambda_k}{\sum_{j=1}^{p} \lambda_j} \]
Em seguida, escolhemos o menor número de componentes \(q\) cuja variância explicada acumulada atinja um limiar satisfatório, geralmente entre 70% e 90%. A escolha do limiar depende do contexto da análise.
Já o Critério do Autovalor (ou Critério de Kaiser) sugere reter apenas os componentes cujos autovalores (\(\lambda_k\)) são maiores que 1. A intuição por trás dessa regra é mais clara quando a ACP é aplicada sobre a matriz de correlação. Nesse caso, as variáveis originais são padronizadas para ter variância 1. Um componente com autovalor (variância) menor que 1 está, portanto, explicando menos variabilidade do que uma única variável original. Reter tal componente não traria uma “economia” de informação, tornando-o um candidato à exclusão.
Uma terceira alternativa, mais visual, é o Gráfico do cotovelo (ou Scree Plot): A Figura 4.3 ilustra esse método. O gráfico mostra a variância explicada acumulada em função do número de componentes principais, com um “cotovelo” indicando o ponto de corte. A ideia é reter os componentes que aparecem antes do cotovelo, pois eles são os que contribuem mais significativamente para a variância total.
4.7 Exemplo Prático
Exemplo 4.3 (Medidas corporais dos pinguins) Vamos exemplificar o uso da ACP no conjunto de dados Palmer Penguins, ele contem medidas corporais de pinguins: comprimento e profundidade do bico, comprimento da nadadeira e massa corporal. As três primeiras estão em milímetros; a última, em gramas. Por isso, uma ACP baseada diretamente na matriz de covariâncias daria peso excessivo à massa corporal. Começamos removendo os valores ausentes e padronizando cada variável.
Código
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
nomes = {
"bill_length_mm": "Comp. do bico",
"bill_depth_mm": "Prof. do bico",
"flipper_length_mm": "Comp. da nadadeira",
"body_mass_g": "Massa corporal",
}
pinguins = pd.read_csv("dados/penguins.csv")
pinguins = pinguins.dropna(subset=list(nomes) + ["species"])
X = pinguins[list(nomes)].rename(columns=nomes)
Z = (X - X.mean()) / X.std()
R = Z.cov()
R.round(2)| Comp. do bico | Prof. do bico | Comp. da nadadeira | Massa corporal | |
|---|---|---|---|---|
| Comp. do bico | 1.00 | -0.24 | 0.66 | 0.60 |
| Prof. do bico | -0.24 | 1.00 | -0.58 | -0.47 |
| Comp. da nadadeira | 0.66 | -0.58 | 1.00 | 0.87 |
| Massa corporal | 0.60 | -0.47 | 0.87 | 1.00 |
O conjunto de dados tem 342 observações completas. Note que a variância total é \(\operatorname{tr}\left(\boldsymbol{R}\right)=1+1+1+1=4\), pois todas as variáveis padronizadas têm variância 1. A matriz também antecipa parte da estrutura que a ACP encontrará. Comprimento da nadadeira e massa corporal, por exemplo, têm correlação 0.87 isso indica que essas variáveis carregam informação parecida.
Agora fazemos a decomposição espectral de \(\boldsymbol{R}\). As colunas da matriz E são os autovetores e definem as combinações lineares; a matriz T contém os escores das observações. A Tabela 4.2 apresenta a variância explicada, enquanto a Figura 4.4 mostra os mesmos autovalores graficamente.
Nota
Existem bibliotecas que facilitam a aplicação da ACP, como sklearn.decomposition.PCA. Agora, como os calculos são extremamente simples, vamos apenas usar a decomposição espectral fornecida pelo NumPy.
Código
autovalores, E = np.linalg.eigh(R)
# A função eigh() devolve os autovalores em ordem crescente, então precisamos inverter a ordem
ordem = autovalores.argsort()[::-1]
autovalores = autovalores[ordem]
E = E[:, ordem]
# O sinal de um autovetor é arbitrário. Vamos garantir que todos os autovetores tenham sinal positivo.
if E[2, 0] < 0:
E[:, 0] *= -1
if E[1, 1] < 0:
E[:, 1] *= -1
# Os escores T são
T = Z.to_numpy() @ E
# Variância explicada pelos componentes principais
proporcoes = autovalores / autovalores.sum()
acumuladas = np.cumsum(proporcoes)
pd.DataFrame({
"Componente": ["CP1", "CP2", "CP3", "CP4"],
"Autovalor": autovalores,
"Proporção (%)": 100 * proporcoes,
"Acumulada (%)": 100 * acumuladas,
}).round(2).set_index("Componente").rename_axis(None)
# Gráfico do cotovelo
componentes = np.arange(1, 5)
fig, ax = plt.subplots(figsize=(5, 3.5))
ax.plot(componentes, autovalores, "o-", linewidth=2)
ax.axhline(1, color="gray", linestyle="--", label="Autovalor = 1")
ax.set(
xlabel="Componente principal",
ylabel="Autovalor",
xticks=componentes,
ylim=(0, 3),
)
ax.legend()
plt.tight_layout()
plt.show()| Autovalor | Proporção (%) | Acumulada (%) | |
|---|---|---|---|
| CP1 | 2.75 | 68.84 | 68.84 |
| CP2 | 0.77 | 19.31 | 88.16 |
| CP3 | 0.37 | 9.13 | 97.29 |
| CP4 | 0.11 | 2.71 | 100.00 |
O primeiro componente tem variância 2.75 e concentra 68.8% da variância total. Os dois primeiros, juntos, concentram 88.2%. Se o interesse do problema fosse a redução da dimensionalidade, seria razoável manter dois componentes, isso porque eles preservam 88,2% da variância e permitem analisar os dados em um espaço bidimensional. Aqui, vamos focar na interpretação, vamos interpretar CP1 e CP2.
Código
correlacoes = E[:, :2] * np.sqrt(autovalores[:2])
pd.DataFrame({
"Variável": Z.columns,
"CP1": correlacoes[:, 0],
"CP2": correlacoes[:, 1],
}).round(3).set_index("Variável").rename_axis(None)
fig, ax = plt.subplots(figsize=(5, 3.5))
ajuste_vertical = {
"Comp. da nadadeira": -0.10,
"Massa corporal": 0.10,
}
for nome, (x, y) in zip(Z.columns, correlacoes):
ax.annotate(
"",
xy=(x, y),
xytext=(0, 0),
arrowprops={"arrowstyle": "->", "color": "firebrick"},
)
ax.text(
1.08 * x,
1.08 * y + ajuste_vertical.get(nome, 0),
nome,
color="firebrick",
fontsize=9,
ha="center",
)
ax.add_patch(plt.Circle((0, 0), 1, fill=False, color="gray", linestyle=":"))
ax.axhline(0, color="gray", linewidth=0.8)
ax.axvline(0, color="gray", linewidth=0.8)
ax.set(
xlabel="CP1",
ylabel="CP2",
xlim=(-1.2, 1.2),
ylim=(-1.2, 1.2),
)
ax.set_aspect("equal")
plt.tight_layout()
plt.show()| CP1 | CP2 | |
|---|---|---|
| Comp. do bico | 0.755 | 0.525 |
| Prof. do bico | -0.664 | 0.701 |
| Comp. da nadadeira | 0.956 | 0.002 |
| Massa corporal | 0.910 | 0.074 |
As correlações tornam a interpretação direta. O primeiro componente cresce com o comprimento da nadadeira, a massa corporal e o comprimento do bico, mas diminui com a profundidade do bico. Ele separa, portanto, pinguins maiores e mais pesados, com nadadeiras longas e bicos mais compridos e menos profundos, de pinguins com o perfil oposto.
O segundo componente está ligado principalmente às duas medidas do bico. Comprimento da nadadeira e massa corporal têm correlações próximas de zero com esse eixo. Assim, \(Y_2\) descreve uma variação conjunta no comprimento e na profundidade do bico.
O círculo de correlações, mostrado na Figura 4.5, é outra ferramente interpretativa visual. Ele permite visualizar as correlações entre as variáveis em um plano bidimensional. Note que, variáveis correlacionadas, como comprimento da nadadeira e massa corporal, tem vetores parecidos na representação. Enquanto isso, os vetores de profundidade do bico e comprimento da nadadeira são opostos, indicando uma correlação negativa. Variáveis pouco correlacionadas, como profundidade do bico e comprimento do bico, formam um angulo quase perpendicular entre eles.
O banco de dados original contém também as espécies dos pinguins. Elas são três: Adelie, Chinstrap e Gentoo. Vamos construir um biplot, colorindo os pontos de acordo com a espécie. No biplot da Figura 4.6, os pontos mostram os escores brutos nos dois primeiros componentes e as setas indicam a direção em que cada medida aumenta. A espécie não entrou no cálculo da ACP; ela aparece apenas como cor para ajudar na leitura.
Os Gentoo aparecem principalmente à direita, na direção das setas de massa corporal e comprimento da nadadeira, indicando que esses pinguins são mais pesados e maiores. Os Adelie se concentram à esquerda, ilustrando que esses pinguins são mais leves e menores.
Os Chinstrap ocupam a parte superior, associada a valores altos nas duas medidas do bico, mostrando que esses pinguin tem os bicos maiores e mais profundos.
Como os dois eixos preservam 88,2% da variância, podemos concluir que essas interpretações são robustas. Em casos onde os componentes descartam mais informação, as interpretações devem ser mais cautelosas.
4.8 Exercícios
Exercício 4.1 Seja \(\boldsymbol{x}\) um vetor aleatório de dimensão \(p\) com matriz de covariâncias \(\boldsymbol{\Sigma}= \operatorname{Cov}\left(\boldsymbol{x}\right)\), real e simétrica. Os componentes principais de \(\boldsymbol{x}\) são definidos pela transformação linear \(\boldsymbol{y} = \boldsymbol{E}^T\boldsymbol{x}\), em que \(\boldsymbol{E}\) é a matriz ortogonal \(p \times p\) cujas colunas são os autovetores normalizados de \(\boldsymbol{\Sigma}\).
Mostre que os componentes principais \(Y_1, \dots, Y_p\) são não correlacionados entre si e têm variâncias \(\lambda_1, \dots, \lambda_p\), respectivamente. Ou seja, mostre que a matriz de covariâncias de \(\boldsymbol{y}\) é a matriz diagonal contendo os autovalores de \(\boldsymbol{\Sigma}\):
\[ \operatorname{Cov}\left(\boldsymbol{y}\right) = \begin{bmatrix} \lambda_1 & 0 & \dots & 0 \\ 0 & \lambda_2 & \dots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \dots & \lambda_p \end{bmatrix} \]
Exercício 4.2 Um vendedor de ferramentas automotivas coletou dados a respeito de um determinado produto de seu portfólio: \(X_1\): Qualidade das peças; \(X_2\): Facilidade de utilização; \(X_3\): Qualidade da embalagem. Após analisar os dados, ele calculou a seguinte matriz de covariância amostral:
\[ \boldsymbol{S}= \begin{bmatrix} 4 & 2 & 0 \\ 2 & 1 & 0 \\ 0 & 0 & 2 \end{bmatrix} \]
a) Apenas olhando para essa matriz, o que podemos esperar de uma ACP nesse caso?
b) Realize a decomposição da matriz para obter os componentes principais. Depois, compare os resultados com as presunções feitas no item (a).
Exercício 4.3 Dois grupos realizaram uma ACP sobre a mesma matriz de dados. Para um dos componentes, o primeiro grupo obteve o autovetor \(\boldsymbol{e}_k\) e os escores \(t_{ik}\). O segundo encontrou \(-\boldsymbol{e}_k\) e escores com os sinais invertidos. Os gráficos produzidos pelos grupos aparecem refletidos, e cada grupo acredita que o outro cometeu um erro.
a) Mostre que, se \(\boldsymbol{e}_k\) é um autovetor de \(\boldsymbol{\Sigma}\) associado a \(\lambda_k\), então \(-\boldsymbol{e}_k\) também é.
b) Compare a variância explicada, as cargas correlacionais e os escores obtidos pelos dois grupos. O que muda e o que permanece igual?
c) Afinal, os grupos encontraram componentes diferentes?
Exercício 4.4 (Correlação e Variância Explicada) Considere um vetor aleatório bidimensional \(\boldsymbol{x}=(X_1,X_2)^T\), em que \(X_1\) e \(X_2\) são padronizadas e \(\operatorname{Corr}\left(X_1,X_2\right)=\rho\), com \(0<\rho\leq 1\).
Mostre que a proporção da variância total explicada pelo primeiro componente principal varia linearmente com \(\rho\). Determine essa função e interprete seu comportamento quando \(\rho\to 0^+\) e quando \(\rho=1\).
Exercício 4.5 (Autovalores Repetidos) Dois pesquisadores aplicaram ACP à matriz de covariâncias
\[ \boldsymbol{\Sigma} = \begin{bmatrix} 4 & 0 & 0\\ 0 & 4 & 0\\ 0 & 0 & 1 \end{bmatrix}. \]
O primeiro escolheu \((1,0,0)^T\) e \((0,1,0)^T\) como os autovetores associados ao autovalor 4. O segundo escolheu, para algum ângulo \(\phi\),
\[ \boldsymbol{v}_1= \begin{bmatrix}\cos\phi\\ \sin\phi\\ 0\end{bmatrix} \qquad\text{e}\qquad \boldsymbol{v}_2= \begin{bmatrix}-\sin\phi\\ \cos\phi\\ 0\end{bmatrix}. \]
a) Verifique que \(\boldsymbol{v}_1\) e \(\boldsymbol{v}_2\) são autovetores unitários e ortogonais associados ao autovalor 4.
b) Algum dos pesquisadores está errado? Explique como a repetição dos autovalores afeta a ACP.
Exercício 4.6 (Projeção e Reconstrução) Uma equipe pretende comprimir uma matriz de dados centralizados \(\boldsymbol{X}\), de dimensão \(n \times p\), mantendo apenas os primeiros \(q\) componentes principais. Para isso, os primeiros \(q\) autovetores são reunidos na matriz \(\boldsymbol{E}_q\), e os escores correspondentes são dados por \(\boldsymbol{T}_q=\boldsymbol{X}\boldsymbol{E}_q\).
a) Mostre como reconstruir uma aproximação \(\widehat{\boldsymbol{X}}\) dos dados originais a partir de \(\boldsymbol{T}_q\) e \(\boldsymbol{E}_q\). Em seguida, escreva essa reconstrução diretamente em função de \(\boldsymbol{X}\).
b) Defina \(\boldsymbol{P}_q=\boldsymbol{E}_q\boldsymbol{E}_q^T\) e mostre que essa matriz é simétrica e idempotente, no sentido da Definição 2.4. Use esse resultado para interpretar geometricamente a reconstrução.
c) Mostre que o resíduo \(\boldsymbol{X}-\widehat{\boldsymbol{X}}\) é ortogonal às direções mantidas. O que essa propriedade diz sobre a informação que foi descartada?
d) Mostre que o erro quadrático de reconstrução satisfaz
\[ \frac{1}{n-1} \sum_{i=1}^n\sum_{j=1}^p \left(x_{ij}-\widehat{x}_{ij}\right)^2 =\sum_{k=q+1}^p\lambda_k. \]
Use essa relação para explicar por que os primeiros \(q\) componentes fornecem a melhor aproximação de dimensão \(q\).
Exercício 4.7 (Composição Química de Vinhos) O conjunto Wine, da UCI Machine Learning Repository, reúne 178 vinhos produzidos na mesma região da Itália a partir de três cultivares. Para cada vinho, foram registradas 13 medidas de sua composição química.
O que uma ACP permite concluir a respeito desses vinhos? Faça uma análise completa usando as medidas químicas como variáveis ativas. Justifique as decisões tomadas, interprete os componentes e use gráficos, incluindo um biplot, para sustentar suas conclusões. O cultivar não deve participar do cálculo dos componentes, mas pode ser usado depois para ajudar a interpretar os escores.