Convertir Ensembl IDs

Los Ensembl IDs son identificadores estables que se usan en las bases de datos para etiquetar inequivocamente, genes, transcripciones, exones o proteínas. A diferencia de los nombres de los genes que pueden cambiar como resultado de las mejoras en el conocimiento científico.

Cargar datos

Para el ejemplo vamos a utilizar un dataset de Drosophila melanogaster.

df.mcf <- read.csv('data/mcf_normalized_counts.csv', row.names = 1)


head(df.mcf)
SRX008026 SRX008174 SRX008201 SRX008239 SRX008008 SRX008168 SRX008211 SRX008255 SRX008261
FBgn0000008 577.30414 589.66391 536.2258 595.58471 560.6214 581.6257 516.9053 536.3239 567.9487
FBgn0000014 1147.76523 1197.37875 1121.8777 1217.76428 671.0268 695.3146 682.5334 668.7586 651.6000
FBgn0000015 395.03031 443.75218 452.2948 431.81641 229.7358 255.3881 236.6116 256.0892 209.1284
FBgn0000017 895.81677 1007.84392 952.1506 930.81982 780.1100 792.5268 789.9186 771.1942 838.7149
FBgn0000018 93.31425 96.27166 95.1218 97.28117 122.6359 154.8805 120.1259 118.5327 107.8662
FBgn0000024 1117.90467 1063.50098 1050.0701 1048.39706 562.2742 561.8537 607.9097 578.7615 517.3176

BioMart

El paquete biomaRt de Bioconductor, es una interfaz de conexion a las bases de datos de BioMart.

Seleccionamos la base de datos.

library(biomaRt)

listEnsembl()
biomart version
genes Ensembl Genes 113
mouse_strains Mouse strains 113
snps Ensembl Variation 113
regulation Ensembl Regulation 113
ensembl <- useEnsembl(biomart = 'genes')

ensembl
Object of class 'Mart':
  Using the ENSEMBL_MART_ENSEMBL BioMart database
  No dataset selected.

Seleccionamos el dataset.

searchDatasets(mart = ensembl, pattern = 'melanogaster')
dataset description version
55 dmelanogaster_gene_ensembl Drosophila melanogaster (Fruit fly) genes (BDGP6.46) BDGP6.46

Necesitamos la base de datos genes y el dataset dmelanogaster_gene_ensembl

ensembl <- useEnsembl(biomart = 'genes', dataset = 'dmelanogaster_gene_ensembl')

ensembl
Object of class 'Mart':
  Using the ENSEMBL_MART_ENSEMBL BioMart database
  Using the dmelanogaster_gene_ensembl dataset

Construimos una consulta y la lanzamos contra el servidor Ensembl BioMart.

genesBM <- getBM(filters = 'ensembl_gene_id', 
                  attributes = c("ensembl_gene_id","entrezgene_id",'external_gene_name'),
                  values = rownames(df.mcf),
                  mart = ensembl
                 )

head(genesBM)
ensembl_gene_id entrezgene_id external_gene_name
FBgn0000008 43852 a
FBgn0000014 42037 abd-A
FBgn0000015 47763 Abd-B
FBgn0000017 45821 Abl
FBgn0000018 44793 abo
FBgn0000024 41625 Ace
Note

En ocasiones el servidor y sus mirrors pueden estar caidos. Es cuestion de esperar un tiempo.

Observando los resultados de la consulta, se muestra que el nombre de los genes (external_gene_name) puede aparece duplicado, mientras que el esembl_gene_id es unico.

head(genesBM[duplicated(genesBM$external_gene_name), ], 10)
ensembl_gene_id entrezgene_id external_gene_name
522 FBgn0003356 43543 Jon99Cii
524 FBgn0003357 43544 Jon99Ciii
701 FBgn0004403 47219 RpS14a
703 FBgn0004404 47218 RpS14b
782 FBgn0004828 33736 His3.3B
1196 FBgn0013674 19893533 mt:CoI
1231 FBgn0013981 3771938 His4r
1232 FBgn0013981 318846 His4r
1233 FBgn0013981 3771893 His4r
1234 FBgn0013981 3771908 His4r

El dataset original (df.mcf) tiene 11342 genes, mientras que la consulta tiene 10474, luego hay genes que no han podido ser traducidos.

dim(df.mcf)
[1] 11342     9
dim(genesBM)
[1] 10474     3

Por ultimo, vamos a unir los resultados de la consulta al dataset original. Es evidente que no podemos añadirlos como columnas nuevas, ya que no tienen la misma longitud. Además nos tenemos que asegurar que añadimos los datos de manera que los ensembles coincidan.

En este caso vamos a usar la funcion merge() para unir los dataset. Lo vamos a realizar haciendo un left join, es decir vamos a conservar todas las observaciones del primer dataset (df.mcf) añadiendo las observaciones del segundo (genesBM) donde encuentre coincidencias. Si una observación del primer dataset no tiene coincidencia en el segundo las variables se añadiran con valores NA

samples <- names(df.mcf)
df.mcf$ensembl_gene_id <- rownames(df.mcf)

df <- merge(df.mcf, genesBM, 
            by =  'ensembl_gene_id', 
            all.x = TRUE)
df <- df[, c('ensembl_gene_id', 'entrezgene_id', 'external_gene_name', samples)]

head(df)
ensembl_gene_id entrezgene_id external_gene_name SRX008026 SRX008174 SRX008201 SRX008239 SRX008008 SRX008168 SRX008211 SRX008255 SRX008261
FBgn0000008 43852 a 577.30414 589.66391 536.2258 595.58471 560.6214 581.6257 516.9053 536.3239 567.9487
FBgn0000014 42037 abd-A 1147.76523 1197.37875 1121.8777 1217.76428 671.0268 695.3146 682.5334 668.7586 651.6000
FBgn0000015 47763 Abd-B 395.03031 443.75218 452.2948 431.81641 229.7358 255.3881 236.6116 256.0892 209.1284
FBgn0000017 45821 Abl 895.81677 1007.84392 952.1506 930.81982 780.1100 792.5268 789.9186 771.1942 838.7149
FBgn0000018 44793 abo 93.31425 96.27166 95.1218 97.28117 122.6359 154.8805 120.1259 118.5327 107.8662
FBgn0000024 41625 Ace 1117.90467 1063.50098 1050.0701 1048.39706 562.2742 561.8537 607.9097 578.7615 517.3176

Al realizar la consulta puede ocurrir que se generen duplicados.

head(df[duplicated(df$ensembl_gene_id), ], 10)
ensembl_gene_id entrezgene_id external_gene_name SRX008026 SRX008174 SRX008201 SRX008239 SRX008008 SRX008168 SRX008211 SRX008255 SRX008261
560 FBgn0003356 43543 Jon99Cii 922.5669 798.7539 839.3100 743.956 2329.09100 1962.368879 1923.834097 1667.50630 2474.319
562 FBgn0003357 43544 Jon99Ciii 3832.7272 3050.6082 3209.8945 2984.222 6375.41557 5092.932163 5578.572854 4558.38713 6456.564
758 FBgn0004403 47219 RpS14a 8884.1384 8935.2133 9724.8054 9143.730 8642.69277 8314.117015 9362.537932 8565.81671 9098.186
760 FBgn0004404 47218 RpS14b 735.3163 699.4738 814.1307 748.855 1044.88456 1008.370910 979.207894 1024.35666 1019.226
847 FBgn0004828 33736 His3.3B 3911.1111 3983.2399 4044.5417 4150.897 3700.56397 3830.820860 3873.149440 3760.85230 3823.748
1298 FBgn0013674 19893533 mt:CoI 0.0000 0.0000 0.0000 0.000 1.98333 1.647665 5.460267 2.19505 0.000
1336 FBgn0013981 3771938 His4r 4156.2165 4078.0073 4099.5631 4080.210 3043.75105 3193.174549 2952.184395 2960.39074 3125.919
1337 FBgn0013981 318846 His4r 4156.2165 4078.0073 4099.5631 4080.210 3043.75105 3193.174549 2952.184395 2960.39074 3125.919
1338 FBgn0013981 3771893 His4r 4156.2165 4078.0073 4099.5631 4080.210 3043.75105 3193.174549 2952.184395 2960.39074 3125.919
1339 FBgn0013981 3771908 His4r 4156.2165 4078.0073 4099.5631 4080.210 3043.75105 3193.174549 2952.184395 2960.39074 3125.919

Tambien se puede dar el caso de datos faltantes de los genes que no se han ‘traducido’.

head(df[!complete.cases(df),])
ensembl_gene_id entrezgene_id external_gene_name SRX008026 SRX008174 SRX008201 SRX008239 SRX008008 SRX008168 SRX008211 SRX008255 SRX008261
21 FBgn0000054 NA NA 457.86190 497.9050 494.26034 468.90923 582.10747 604.69301 526.00573 566.32290 633.98924
35 FBgn0000099 NA NA 83.36073 84.2377 85.79613 82.58401 74.37489 72.49725 91.00445 82.68022 92.45676
49 FBgn0000142 NA NA 903.90400 812.2921 934.43182 930.11996 537.81309 486.06114 587.88875 561.20111 576.75410
58 FBgn0000171 NA NA 8949.45833 8890.0860 8856.58579 9277.40418 6353.92949 6338.56682 5791.52327 6530.27370 7053.13026
66 FBgn0000210 NA NA 314.78006 371.5484 361.83587 328.93633 99.16652 102.15522 96.46472 109.02082 112.26893
75 FBgn0000242 NA NA 255.05894 221.1240 224.74857 235.85435 79.00266 80.73558 89.18436 93.65547 118.87298
dim(df[!complete.cases(df),])
[1] 909  12

AnnotationDbi

Bioconductor proporciona varias fuentes de anotaciones, entre ellos hay varios paquetes. Los paquetes a nivel de organismo (p. ej.org.Mm.eg.db) usan un identificador de gen central (p. ej. Entrez Gene id) y contiene asignaciones entre este y otros tipos de identificadores. El nombre de los paquetes es del tipo org.<Ab>.<id>.db (p. ej. org.Mm.eg.db), donde <Ab> es una abreviatura de 2 letras del organismo (p. ej. Mm, mus musculus), y <id> es una abreviatura en minuscula que describe el tipo de indentificador central (p. ej. eg, entrez gene).

En nuestro caso nos interesa el paquete org.Dm.eg.db.

library(org.Dm.eg.db)

keytypes(org.Dm.eg.db)
 [1] "ACCNUM"       "ALIAS"        "ENSEMBL"      "ENSEMBLPROT"  "ENSEMBLTRANS"
 [6] "ENTREZID"     "ENZYME"       "EVIDENCE"     "EVIDENCEALL"  "FLYBASE"     
[11] "FLYBASECG"    "FLYBASEPROT"  "GENENAME"     "GENETYPE"     "GO"          
[16] "GOALL"        "MAP"          "ONTOLOGY"     "ONTOLOGYALL"  "PATH"        
[21] "PMID"         "REFSEQ"       "SYMBOL"       "UNIPROT"     
genes <- mapIds(org.Dm.eg.db, 
                keys = rownames(df.mcf), 
                keytype = 'ENSEMBL', 
                column = 'SYMBOL')
head(genes)
FBgn0000008 FBgn0000014 FBgn0000015 FBgn0000017 FBgn0000018 FBgn0000024 
        "a"     "abd-A"     "Abd-B"       "Abl"       "abo"       "Ace" 

En este caso el resultado es un vector con los simbolos de los genes. Si lo prefereimos podemos convertirlo en un dataframe facilmente.

genes <- as.data.frame(genes)

head(genes)
genes
FBgn0000008 a
FBgn0000014 abd-A
FBgn0000015 Abd-B
FBgn0000017 Abl
FBgn0000018 abo
FBgn0000024 Ace

Biological id Translator

En el paquete para análisis de enriquecimiento, clusterProfiler, podemos encontrar la función bitr() que permite cambiar de identificador.

library(clusterProfiler)

genes <- bitr(rownames(df.mcf), 
              fromType = 'ENSEMBL', 
              toType = 'SYMBOL', 
              OrgDb = org.Dm.eg.db)
head(genes)
ENSEMBL SYMBOL
FBgn0000008 a
FBgn0000014 abd-A
FBgn0000015 Abd-B
FBgn0000017 Abl
FBgn0000018 abo
FBgn0000024 Ace

Referencias

Open source software for Bioinformatics