PCA: Reduzir Dimensão sem Perder o que Importa (Nem Sempre)

Aula 10a e 10b. Depois de K-means e DBSCAN (agrupar sem rótulo), PCA ataca outro problema sem rótulo: como resumir um monte de variáveis em poucas, sem perder o que importa. Spoiler do post inteiro: "o que importa" é uma pergunta traiçoeira.

O que PCA calcula, de verdade

O professor volta pro dataset de câncer de mama que já apareceu nessa playlist, 30 variáveis, 569 pacientes:

from sklearn.decomposition import PCA
pca = PCA(n_components=30)
X_pca = pca.fit_transform(X)
print(pca.explained_variance_ratio_)

Saída: o primeiro componente sozinho explica 98.2% de toda a variância dos dados. Os 29 restantes, juntos, explicam os 1.8% que sobram.

Antes da fórmula, vale pegar a cena por trás dessas duas linhas. Imagina a nuvem de 569 pacientes espalhada num espaço de 30 dimensões, e você virando a cabeça, ângulo por ângulo, procurando o ponto de vista de onde essa nuvem parece mais espalhada possível (o oposto de olhar uma folha de papel de lado, que só mostra uma linha fina). É exatamente isso que pca.fit_transform(X) faz, só que testando matematicamente, não virando a cabeça na tentativa e erro: primeiro centraliza a nuvem (subtrai a média de cada variável, pra ela ficar "flutuando" em volta da origem), depois gira os eixos até achar essa direção de maior espalhamento, e chama ela de "componente 1", a segunda maior direção sobrando vira "componente 2", e assim por diante. PCA(n_components=30) só configura quantos ângulos novos eu quero guardar no final (aqui, todos os 30, de propósito, pra ver a variância se espalhar por todos antes de eu cortar qualquer coisa). explained_variance_ratio_ é a fatia de espalhamento que cada ângulo novo capturou, já em ordem decrescente, por isso o primeiro número da lista é sempre o maior.

Uma concentração absurda, e faz sentido: várias dessas 30 variáveis medem essencialmente a mesma coisa de jeitos diferentes (raio, perímetro e área do tumor são praticamente a mesma informação, só em unidades diferentes), então a variação real do dataset é bem mais "estreita" do que 30 números sugerem.

O Bishop define PCA como achar a direção que maximiza a variância dos dados projetados. Pra um único componente, isso vira um problema de autovalor: a direção ótima $\mathbf{u}_1$ satisfaz

$\mathbf{S}\mathbf{u}_1 = \lambda_1 \mathbf{u}_1$

onde $\mathbf{S}$ é a matriz de covariância dos dados. Ou seja, $\mathbf{u}_1$ precisa ser um autovetor de $\mathbf{S}$, e a variância capturada é exatamente o autovalor $\lambda_1$ correspondente. Quanto maior o autovalor, mais variância aquela direção captura, e é por isso que o "primeiro componente principal" é sempre o autovetor de maior autovalor. Os componentes seguintes repetem a receita, cada um maximizando a variância entre as direções que sobraram, ortogonais às anteriores.

Repara: em nenhum momento essa conta usa o rótulo y. PCA olha só pra X, procurando onde os dados variam mais. Guarda essa frase, ela é o fio condutor do post inteiro.

Reduzindo de verdade, com um porém

Com 9 componentes (30% das 30 variáveis originais), a variância capturada já é 99.99997%:

n_most_important = int(X.shape[1]*0.3)  # 9
scores = cross_val_score(model, X_pca[:, :n_most_important], y, cv=splitter)

O professor compara acurácia (espaço original de 30 variáveis vs os 9 primeiros componentes) em cinco classificadores diferentes:

Repara que o efeito não é o mesmo pra todo mundo. KNN nem percebe a diferença (PCA com quase toda a variância preservada é basicamente uma rotação + um corte quase sem perda, e distância euclidiana não liga pra rotação). Naive Bayes Gaussiano piora bastante, o que é meio contra-intuitivo à primeira vista: ele assume que as variáveis são independentes entre si, e os componentes do PCA são descorrelacionados por construção, então eu esperaria uma ajuda, não um tombo. Mas "descorrelacionado" não é a mesma coisa que "todo mundo importa igual": os 9 componentes que sobraram têm variância bem desigual entre si (o primeiro sozinho já carrega a fatia esmagadora da variância total), e Naive Bayes Gaussiano ajusta uma variância por variável por classe, então ele confia mais em componentes de variância maior. O problema é que "maior variância" e "melhor pra separar as classes" não são a mesma coisa, a mesma pegadinha que a próxima seção deixa ainda mais clara: descartar os componentes de menor variância pode estar jogando fora justamente o eixo onde a média de uma classe difere da outra, o que Naive Bayes Gaussiano sente na pele porque ele depende diretamente dessas variâncias por variável, ao contrário do KNN, que só olha a distância total. E árvore de decisão piora um pouquinho, pelo mesmo motivo que já vi no post de árvores: árvore corta um eixo de cada vez, e PCA gira os eixos, então a fronteira de decisão que ficava alinhada com uma variável original pode virar diagonal nos componentes, mais difícil de recortar com cortes retos.

O scikit-learn ainda deixa escolher a quantidade de componentes por fração de variância em vez de contagem fixa:

pca = PCA(n_components=0.9999)
X_pca = pca.fit_transform(X)
print(X_pca.shape)

Saída: (569, 5). Só 5 componentes já bastam pra 99.99% da variância, menos até que os 9 escolhidos na mão.

O aviso: PCA não sabe o que é "importante" pra você

Essa seção (o professor batizou de "PCA fail" no próprio notebook) é o motivo pelo qual eu escrevi "guarda essa frase" lá em cima. Ele constrói um dataset artificial de propósito: 12 grupos de pontos, arranjados em duas fileiras (uma em $y\approx0.1$, outra em $y\approx0.4$), espalhados ao longo do eixo $x$ de 0 a 1.1. A classe alterna: fileira de baixo é uma classe, fileira de cima é outra, não importa o valor de $x$.

model = Perceptron()
model.fit(X, y)
print(accuracy_score(y, model.predict(X)))

Saída: 1.0. Cem por cento, porque separar por fileira (valor de y) é fácil, é só uma reta horizontal.

Agora aplica PCA com 1 componente antes de treinar:

pca = PCA(n_components=1)
X_pca = pca.fit_transform(X)
model.fit(X_pca, y)
print(accuracy_score(y, model.predict(X_pca)))

Saída: 0.47. Praticamente a cara de uma moeda sendo jogada.

O que aconteceu: o eixo $x$ (0 a 1.1) tem muito mais espalhamento que o eixo $y$ (0.1 a 0.4), então a direção de maior variância é quase totalmente ao longo de $x$. PCA, fazendo exatamente o que promete, escolhe essa direção como o único componente. Só que a informação que separa as classes mora no eixo $y$, o de menor variância, e PCA jogou ela fora inteira. O algoritmo não tem como saber que $y$ importava mais pra classificar, porque ele nunca olhou pro rótulo.

Interativo: pra onde o primeiro componente aponta

Reconstrução minha do mesmo dataset de 100 pontos (dois grupos, mesma disposição em fileiras), com os dois componentes principais desenhados por cima (vermelho sólido = 1º componente, verde tracejado = 2º, cada seta proporcional à raiz quadrada da variância que aquele componente explica):

Sem normalizar, o componente 1 (vermelho) fica quase deitado, seguindo o eixo $x$, exatamente a direção que não separa as classes. Clica em "Normalizado (z-score)": agora as duas variáveis competem em pé de igualdade, e o componente 1 gira pra uma diagonal, ainda não perfeitamente alinhado com a separação real (por isso o resultado com normalização, 0.81 de acurácia no notebook, é melhor mas não perfeito), mas bem menos cego ao que importa do que a versão sem normalizar.

PCA como compressão de verdade

aula10b troca de dataset pro MNIST (70 mil imagens, 784 pixels cada) e usa PCA pra achar quantos componentes bastam pra reter 95% da variância:

pca = PCA()
pca.fit(X_train)
cumsum = np.cumsum(pca.explained_variance_ratio_)
d = np.argmax(cumsum >= 0.95) + 1
print(d)

Saída: d = 154. De 784 pixels pra 154 componentes, uma redução de quase 80% no tamanho, ainda guardando 95% da variância.

A linha que acha esse 154 é densa, mas a cena por trás é simples: imagina enchendo um pote com "fatias de variância explicada", começando pela maior fatia (componente 1), depois a segunda maior, a terceira, e assim por diante, cada vez despejando mais uma em cima da anterior. cumsum é literalmente esse pote enchendo, guardando "quanto já foi despejado até aqui" a cada passo. Assim que o pote passa de 95% cheio, eu paro de despejar, e esse é o número de componentes que eu preciso guardar. np.argmax(cumsum >= 0.95) é só o jeito do numpy de achar esse ponto de parada exato, e o + 1 no final é porque contar posição de lista começa do zero, mas "quantos componentes eu preciso" começa do um.

Comparando acurácia, pixel cru contra os 154 componentes:

O mesmo padrão do dataset de câncer: modelo linear não se importa (ou até melhora um pouquinho), floresta de árvores piora, porque de novo perde o alinhamento de eixo que ela precisa pra cortar bem.

E porque PCA guarda a direção de cada componente, dá pra fazer o caminho de volta (com perda) e reconstruir uma aproximação da imagem original:

X_train_recovered = pca.inverse_transform(X_train_reduced)

inverse_transform desfaz a rotação que fit_transform fez: pega os 154 números de cada imagem e multiplica de volta pelos mesmos 154 eixos que o PCA aprendeu, depois soma de volta a média de cada pixel que tinha sido subtraída na centralização. É como tentar redesenhar o dígito original só com as 154 pinceladas mais importantes que eu guardei, sem lembrar dos detalhes finos que descartei junto com os outros 630 eixos.

Os dígitos reconstruídos ficam visivelmente mais "borrados" que o original, mas ainda dá pra reconhecer qual número é qual, com só 154 números guardados por imagem em vez de 784. É literalmente compressão com perda, no mesmo espírito de um JPEG, só que a "base" usada pra descrever a imagem foi aprendida a partir do próprio dataset em vez de ser genérica.

Fechando

Aplicação Prática

Sem normalizar, o componente 1 do dataset "PCA fail" captura 83.5% da variância total, quase todo mundo. Normalizado, cai pra 54.8%, quase empatado com o componente 2 (45.2%). Isso por si só já é um sinal de alerta útil: quando o primeiro componente domina demais a variância (como no câncer de mama, 98.2%, ou aqui sem normalizar, 83.5%), vale desconfiar que uma variável de escala grande está dominando a conta sozinha, o mesmo problema de escala que eu já vi lá atrás com KNN, só que agora afetando o próprio PCA em vez de uma distância.

# variância explicada pelo componente 1, com e sem normalizar
razao_sem_normalizar = 0.835  # x domina, y quase não conta
razao_normalizado = 0.548     # x e y competem em pé de igualdade

A régua prática que fica: sempre que eu for usar PCA em variáveis que não estão na mesma unidade (preço em dólar ao lado de idade em anos, por exemplo), normalizar antes não é opcional, é o mesmo cuidado de sempre, agora aplicado a um lugar novo.