Enriquecimento GO (topGO): principalmente para montagens de novo - lmigueel/Bioinformatica GitHub Wiki
1. Sobre
O pacote topGO fornece ferramentas para testar os termos GO enquanto leva em consideração a topologia do gráfico GO. Diferentes estatísticas de teste e diferentes métodos para eliminar semelhanças e dependências locais entre os termos GO podem ser implementados e aplicados.
Uma das principais vantagens do topGO é a estrutura de teste de conjunto de genes unificado que oferece. Além de fornecer um conjunto de funções fáceis de usar para realizar a análise de enriquecimento GO, também permite ao usuário implementar facilmente novos testes estatísticos ou novos algoritmos que lidam com a estrutura gráfica do GO. Este framework unificado também facilita a comparação entre diferentes metodologias de enriquecimento GO. Há uma série de estatísticas de teste e algoritmos lidando com o grafo GO estruturado pronto para usar no topGO. A Tabela 1 apresenta a tabela de compatibilidade entre as estatísticas do teste e os métodos do gráfico GO.
Tabela 1. Algoritmos suportados pelo GO
| algol | fisher | ks | t | globaltest | sum |
|---|---|---|---|---|---|
| classic | ✓ | ✓ | ✓ | ✓ | ✓ |
| elim | ✓ | ✓ | ✓ | ✓ | ✓ |
| weight | ✓ | - | - | - | - |
| weight01 | ✓ | ✓ | ✓ | ✓ | ✓ |
| lea | ✓ | ✓ | ✓ | ✓ | ✓ |
| parentchild | ✓ | - | - | - | - |
Os algoritmos de elim e weight foram introduzidos em Alexa et al. (2006). O algoritmo padrão usado pelo pacote topGO é uma mistura entre os algoritmos elim e weight e será referido como weight01. O algoritmo parentchild foi introduzido por Grossmann et al. (2007).
Usaremos principalmente o teste exato de Fisher, que é baseado na contagem de genes.
2. Instalação
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("topGO")
Acesse o MANUAL AQUI.
3. Exemplos
O script abaixo está completo. Ele fornece todas as métricas de corte (que possuem um p-valor associado), mas também a lista de genes presentes por processo GO, algo que originalmente não vem no script base. Eu tive que executar este procedimento pela necessidade e trago à vocês.
Eu utilizou muito o topGo após realizar uma montagem de novo e realizar uma anotação GO com o Pannzer2. Acesse sua anotação AQUI.
Você precisa de três arquivos:
- universeFile : arquivo contendo todos os genes do organismo e o GO associado a cada gene. Ele está no formado
gene GO GO GO. Você pode fazer um script que percorre o arquivo output GO.out do Pannzer2 e montar esse arquivo. Como, por exemplo:
TRINITY_DN2789_c0_g1_i4 GO:0055085 GO:0006810
TRINITY_DN28613_c1_g1_i2 GO:0042744 GO:0042743 GO:0072593
TRINITY_DN8344_c2_g1_i3 GO:0097502 GO:0070085
...
-
interestingGenesFile: lista de genes de interesse que você deseja saber quais processos GO estão enriquecidos
-
output_file: arquivo de saída
Como exemplo, os arquivos (1), (2) e (3) possuem nomes isoforma_GO.txt, transcritos_diferenciais.txt e resultado_topGO.txt, respectivamente. Segue o script abaixo com os passos descritos.
Você pode trocar as primeiras linhas por:
args <- commandArgs(trailingOnly = TRUE)
universeFile = args[1]
interestingGenesFile = args[2]
output_file = args[3]
sendo possível rodar em linha de comando: Rscript topGO.R isoforma_GO.txt transcritos_diferenciais.txt resultado_topGO.txt .
Saliento que o dataset final "AllRes_final" possui todos os p-values calculados, mas irei utilizar como corte o teste de estrutura pelo weight01 + fisher, pela preferência na Tabela 1.
# topGO.R
# eg. run:
# Rscript topGO.R annotations.txt interestinggenes.txt output_file
#args <- commandArgs(trailingOnly = TRUE)
universeFile = "isoforma_GO.txt"
interestingGenesFile = "transcritos_diferenciais.txt"
output_file = "resultado_topGO.txt"
# load topGO
library("topGO")
# read in the 'gene universe' file
geneID2GO <- readMappings(file = universeFile)
geneUniverse <- names(geneID2GO)
# read in the genes of interest
genesOfInterest <- read.table(interestingGenesFile,header=FALSE)
genesOfInterest <- as.character(genesOfInterest$V1)
geneList <- factor(as.integer(geneUniverse %in% genesOfInterest))
names(geneList) <- geneUniverse
# build the GOdata object in topGO
myGOdata <- new("topGOdata", description="My project", ontology="BP", allGenes=geneList, annot = annFUN.gene2GO, gene2GO = geneID2GO)
# run the Fisher's exact tests
resultClassic <- runTest(myGOdata, algorithm="classic", statistic="fisher")
resultElim <- runTest(myGOdata, algorithm="elim", statistic="fisher")
resultTopgo <- runTest(myGOdata, algorithm="weight01", statistic="fisher")
resultParentchild <- runTest(myGOdata, algorithm="parentchild", statistic="fisher")
# Eu gosto do weight01 + fisher. Vou cortar nele.
# see how many results we get where weight01 gives a P-value <= 0.001:
mysummary <- summary(attributes(resultTopgo)$score <= 0.001)
numsignif <- as.integer(mysummary[3](/lmigueel/Bioinformatica/wiki/3)) # how many terms is it true that P <= 0.001
# print out the top 'numsignif' results:
allRes <- GenTable(myGOdata, classicFisher = resultClassic, elimFisher = resultElim, topgoFisher = resultTopgo, parentchildFisher = resultParentchild, orderBy = "topgoFisher", ranksOf = "classicFisher", topNodes = numsignif)
# print a graph (to a pdf file) with the top 'numsignif' results:
output_file2 = paste(output_file,"Topgo", sep="_")
#printGraph(myGOdata, resultTopgo, firstSigNodes = numsignif, fn.prefix = output_file2, useInfo = "all", pdfSW = TRUE)
# print out the genes that are annotated with the significantly enriched GO terms:
myterms <- allRes$GO.ID
mygenes <- genesInTerm(myGOdata, myterms)
AllRes_final <- allRes
res <- NULL
for (i in 1:length(myterms)){
myterm <- myterms[i]
mygenesforterm <- mygenes[myterm][1](/lmigueel/Bioinformatica/wiki/1)
myfactor <- mygenesforterm %in% genesOfInterest # find the genes that are in the list of genes of interest
mygenesforterm2 <- mygenesforterm[myfactor == TRUE]
mygenesforterm2 <- paste(mygenesforterm2, collapse=',')
#print(paste("Term",myterm,"genes:",mygenesforterm2))
res[i] <- mygenesforterm2
}
AllRes_final$genes <- res
write.table(AllRes_final,file= output_file,sep="\t")
Citação
Alexa A, Rahnenfuhrer J (2020). topGO: Enrichment Analysis for Gene Ontology. R package version 2.42.0.
Alexa, A., Rahnenf ̈uhrer, J., and Lengauer, T. (2006). Improved scoring of functional groups from geneexpression data by decorrelating go graph structure.Bioinformatics (Oxford, England), 22:1600–1607.10.1093/bioinformatics/btl140.
Grossmann, S., Bauer, S., Robinson, P. N., and Vingron, M. (2007). Improved detection of overrepresentationof gene-ontology annotations with parent child analysis.Bioinformatics, 23:3024. 10.1093/bioinformat-ics/btm440.
https://www.bioconductor.org/packages/devel/bioc/vignettes/topGO/inst/doc/topGO.pdf