Gene Set Enrichment Analysis (GSEA) - lmigueel/Bioinformatica GitHub Wiki
Gene Set Enrichment Analysis (GSEA) (também análise de enriquecimento funcional) é um método para identificar classes de genes ou proteínas que estão super-representados em um grande conjunto de genes ou proteínas e podem ter uma associação com fenótipos de doenças.
Para isso, iremos utilizar um pacote em Python para estas análises: o GSEApy.
GSEApy é uma implementação em Python para GSEA e wrapper para Enrichr. GSEApy pode ser usado para dados de RNA-Seq, ChIP-seq, Microarray. Ele pode ser usado para enriquecimento GO conveniente e para produzir números de qualidade de publicação em Python. GSEApy tem seis subcomandos disponíveis: gsea, prerank, ssgsea, replot richr, biomart. Para mais informações, acesse: https://pypi.org/project/gseapy/.
1. Instalação
pip install gseapy
2. Script em python
Eu particularmente uso o GSEApy para obter o enriquecimento de vias KEGG e Gene Ontology (GO). Ele possui várias bibliotecas que são utilizados no módulo do Enrichr. Você pode acessá-las aqui: https://maayanlab.cloud/Enrichr/#stats
O script se baseia em passar um ou mais arquivos de entrada, e o enriquecimento ser executado. Os gráficos são bem auto explicativos, mas a lista em .txt do enriquecimento KEGG e GO também são fornecidas como reports (variável gene_sets). Ele apenas aceita cinco organismos para o Enrichr: human, mouse, fly, yeast, worm e fish. Abaixo eu tenho como entrada duas listas de genes (up_genes.txt e down_genes.txt).
import os
import pandas as pd
import numpy as np
from pandas import DataFrame
import gseapy
files = ['down_genes.txt','up_genes.txt']
for arq in files:
input_file = arq
taxid = 9606
aux = arq.replace(".txt","")
output_folder = aux + "_kegg_enrichment"
fdr_cutoff = 0.05
genes_list = pd.read_csv(input_file, sep='\t',header=None)
genes = genes_list.iloc[:,0]
my_genes = genes.tolist()
organism_validated = [9606,10090,7227,4932,6396,7955]
org = {9606:'Human',10090:'Mouse',7227:'Fly',4932:'Yeast',6396:'Worm',7955:'Fish'}
if(taxid in organism_validated):
l = my_genes
output = output_folder
gseapy.enrichr(gene_list=l,
description='pathway',
gene_sets=['KEGG_2019_Human','GO_Biological_Process_2018'],
organism = org[taxid],
outdir=output,
top_term=20,
cutoff = fdr_cutoff )
3. Output
Um exemplo de output do gráfico de enriquecimento KEGG e GO pode ser encontrado aqui: Exemplo - GSEApy output
Citation
Fang, Z. "GSEAPY: Gene Set Enrichment Analysis in Python." (2020).