Generando un informe

Guillermo Ayala Gallego

2025-04-16

Justificación

  • Hemos visto cómo obtener subconjuntos de genes con expresión diferencial.
  • Ya tenemos la lista.
  • Queremos explorar este grupo de genes.
  • En la práctica estadística clásica se generaba un informe y el investigador interpreta los resultados y redacta un informe.
  • En este contexto el problema es mayor.
  • Básicamente se explora en la gran cantidad de genes y por ello los resultados han de ser validados y, posiblemente, nuevas hipótesis generadas.
  • El modo más útil posiblemente sea generar un archivo en html.
  • Nos proponemos el siguiente objetivo:
    • Generar un fichero en que cada fila contenga distintos identificadores de aquellos genes que declaramos con expresión diferencial significativa así como los p-valores originales y los ajustados.

Datos y análisis previo

  • Leemos los datos.
pacman::p_load(Biobase)
data(gse20986,package="tamidata")
  • Hacemos análisis de expresión diferencial.
eset = gse20986[, c(2, 3, 5, 10:12)]
tissue = factor(rep(1:2, each = 3), levels = 1:2,
                labels = c("retina", "huvec"))
tt = genefilter::rowttests(eset, tissue)
padj = p.adjust(tt[,"p.value"],method="BH")
sig = which(padj < 0.05)
  • Consideramos los p-valores y los p-valores ajustados.
pvalor = tt[sig,"p.value"]
pajustado = padj[sig]

Obteniendo identificadores con annotate

  • Utilizamos annotate::lookUp()
pacman::p_load(annotate)
  • ¿Qué paquete de anotación hemos de usar?
annotation(eset)
[1] "hgu133plus2"
  • Cargamos el paquete correspondiente.
pacman::p_load(hgu133plus2.db)
  • Obtenemos los identificadores de nuestros genes.
ID = featureNames(eset)[sig]
head(ID)
[1] "1552487_a_at" "1552626_a_at" "1552701_a_at" "1552703_s_at" "1552730_at"  
[6] "1552760_at"
  • Podemos determinar los nombres abreviados de nuestros genes.
lookUp(ID, "hgu133plus2.db", "SYMBOL")
getSYMBOL(ID,"hgu133plus2.db")
  • Sus nombres con
lookUp(ID, "hgu133plus2.db", "GENENAME")
  • Los identificadores ENSEMBL
lookUp(ID, "hgu133plus2.db", "ENSEMBL")
  • O en Gene Ontology,
lookUp(ID, "hgu133plus2.db", "GO")
  • O los identificadores ENTREZ.
lookUp(ID, "hgu133plus2.db", "ENTREZID")

Generando informe con ReportingTools

pacman::p_load(ReportingTools)
  • Obtenemos distintos identificadores.
ID = featureNames(eset)[sig]
Name = as.character(lookUp(ID,"hgu133plus2.db", "GENENAME"))
entrezid =  as.character(lookUp(ID, "hgu133plus2.db", "ENTREZID"))
ID[ID == "NA"] = NA
Name[Name == "NA"] = NA
entrezid = ifelse(entrezid == "NA",NA,
  paste0("<a href='http://www.ncbi.nlm.nih.gov/gene/?term=",
         entrezid,"'>",entrezid,"</a>"))
  • Posiblemente queramos añadir a estos descriptores los p-valores originales así como los ajustados.
df = data.frame(ID = ID,Name = Name,entrezid = entrezid,
    pvalor = pvalor,pajustado = pajustado,stringsAsFactors=FALSE)
  • Utilizando ReportingTools::HTMLReport() fijamos el nombre del fichero en que guardamos la información así como el directorio en donde estará dicho fichero.
foutput = "gse20986.DE"
htmlRep1 = HTMLReport(shortName = foutput,title = foutput,
                       reportDirectory = "./reports")
  • Guardamos la información y cerramos el fichero con las funciones ReportingTools::publish() y ReportingTools::finish().
publish(df,htmlRep1)
finish(htmlRep1)
[1] "./reports/gse20986.DE.html"
browseURL("./reports/gse20986.DE.html")

Utilizando AnnotationDbi::select

a = AnnotationDbi::select(hgu133plus2.db,keys=ID,
                          columns=c("ENTREZID","ENSEMBL"),
                          keytype="PROBEID")
b = BiocGenerics::match(ID,a[,"PROBEID"])
df0 =  a[b,]
names(df0)
[1] "PROBEID"  "ENTREZID" "ENSEMBL" 
df = data.frame(ID = ID,Name = Name,entrezid = df0$ENTREZID,
    pvalor = pvalor,pajustado = pajustado,stringsAsFactors=F)

Y lo demás igual que en la opción previa.

Utilizando DT

  • Tomamos data.frame previo.
ff = DT::datatable(df,escape=FALSE)
  • Generamos el informe en html.
ff = DT::datatable(df,escape=FALSE)
  • Guardamos informe generado.
DT::saveWidget(ff, "reports/gse20986_DE_DT.html")

Utilizando kableExtra

pacman::p_load(kableExtra)
df %>% kable(escape=FALSE) %>%
    kable_styling() %>%
    save_kable("reports/gse20986_DE_kE.html")

Generación de enlaces

tami::entrezid2url
function (id) 
ifelse(id == "NA", NA, paste("<a href='http://www.ncbi.nlm.nih.gov/gene/?term=", 
    id, "'>", id, "</a>", sep = ""))
<bytecode: 0x55a09a5274f8>
<environment: namespace:tami>
tami::ensembl2url
function (id, site = "http://www.ensembl.org") 
paste("<a href='", site, "/id/", id, "'>", id, "</a>", sep = "")
<bytecode: 0x55a09a7704f8>
<environment: namespace:tami>
tami::go2url
function (id) 
ifelse(id == "NA", NA, paste("<a href='http://amigo.geneontology.org/amigo/term/", 
    id, "'>", id, "</a>", sep = ""))
<bytecode: 0x55a09aad6f68>
<environment: namespace:tami>
tami::kegg2url
function (id) 
ifelse(id == "NA", NA, paste("<a href='http://www.genome.jp/dbget-bin/www_bget?", 
    id, "'>", id, "</a>", sep = ""))
<bytecode: 0x55a09ad3f2d0>
<environment: namespace:tami>
tami::WormBase2url
function (id) 
ifelse(id == "NA", NA, paste("<a href='http://www.wormbase.org/species/c_elegans/gene/", 
    id, "'>", id, "</a>", sep = ""))
<bytecode: 0x55a09b154d70>
<environment: namespace:tami>