PCA (Principal Component Analysis) em dados de RNA Seq - lmigueel/Bioinformatica GitHub Wiki

1. Sobre

A Análise de Componentes Principais ou Principal Component Analysis (PCA) é um procedimento matemático que utiliza uma transformação ortogonal (ortogonalização de vetores) para converter um conjunto de observações de variáveis possivelmente correlacionadas num conjunto de valores de variáveis linearmente não correlacionadas chamadas de componentes principais. O número de componentes principais é sempre menor ou igual ao número de variáveis originais. Os componentes principais são garantidamente independentes apenas se os dados forem normalmente distribuídos (conjuntamente).

PCA é um algoritmo determinístico.

A ideia é verificar quanto que as condições de interesse explicam a variância total. As condições que tiverem próximos no plot do PCA possuem viriâncias iguais. A tendência é que as amostras controle fiquem proximas, e as de tratado também.

Caso queira saber mais dos arquivos finais de análises de RNA-Seq que podem ser utilizados nestas análises de PCA, acesse AQUI.

2. DESeq2

library("DESeq")
countsTable <- read.delim("kallisto.isoform.counts.matrix", header=TRUE, stringsAsFactors=TRUE)
rownames(countsTable) <- countsTable[,1]

#supondo que temos 2 condições em triplicata (6 colunas)
countsTable <- countsTable[,2:7]
conds <- factor(c(names(countsTable)))

countsTable_novo <- apply(countsTable,2,as.integer)
countsTable_novo[is.na(countsTable_novo)] <- 0

cds<-newCountDataSet(countsTable_novo,conds)
cds<-estimateSizeFactors(cds)
sizeFactors(cds)
cds <- estimateDispersions(cds,method='blind')
vsd <- varianceStabilizingTransformation(cds)

pdf("PCA.pdf")
plotPCA(vsd)
dev.off()

3. Python

O importante aqui é você utilizar a matriz em normalizada pelo TMM, que é normalizada pelo edgeR. Esse output normalizado também pode ser gerado através do script abundance_estimates_to_matrix.pl que está no pacote TRINITY dentro da pasta /util/. Acesse aqui: script. Acesse o pacote Trinity aqui.

import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import matplotlib.pyplot as plt
import os
import matplotlib.patches as mpatches

foo = pd.read_csv('kallisto.isoform.TMM.EXPR.matrix', sep="\t",index_col=0)
foo = foo+1
foo = foo.apply(np.log10)

foo.columns

from sklearn.decomposition import PCA

pca = PCA(n_components=2)
bar = pca.fit_transform(foo.transpose())
bar = pd.DataFrame(bar, index=foo.columns)

targets = foo.columns
#o vetor de colors deve ter o tamanho das colunas. Neste caso, por exemplo, temos 6 colunas
colors = ['r', 'r','r','b','b','b']

legend_dict = { 'A' : 'red', 'B' : 'blue' }

pc1_var = round(pca.explained_variance_ratio_[0]*100,2)
pc2_var = round(pca.explained_variance_ratio_[1]*100,2)

fig = plt.figure(figsize=(18,18))
ax = fig.add_axes([0.1,0.1,0.8,0.8])
ax.set_xlabel('Principal Component 1 ('+str(pc1_var)+'%)', fontsize = 15)
ax.set_ylabel('Principal Component 2 ('+str(pc2_var)+'%)', fontsize = 15)

patchList = []
for key in legend_dict:
        data_key = mpatches.Patch(color=legend_dict[key], label=key)
        patchList.append(data_key)

plt.legend(handles=patchList,loc=1, prop={'size': 30})
for i,color in zip(bar.index,colors):
    ax.scatter(bar.loc[i][0], bar.loc[i][1],s=1000,c=color)
    ax.text(bar.loc[i][0], bar.loc[i][1], i,size=20,weight="bold",
           rotation=45,ha='center', va='center')


plt.tight_layout()
fig.savefig('PCA-geral.png', dpi=600)