基因注释资源（1）

$\mathcal{T}ranslate$ from BioConductor:Genomic Annotation Resources

1.版本信息

R version: R version 3.5.0 (2018-04-23)
Bioconductor version: 3.7
Package version: 1.2.0

2. 引言

注释资源在Bioconductor项目中占很大比例[1]。此外，还有一系列可用的在线资源，可以使用特定的软件包访问。本演练将介绍这些资源中最受欢迎的资源，并提供一些有关如何使用它们的高级示例。

传统上，Bioconductor注释资源在分析结束时使用。在大量数据分析之后，注释将被解释性地用于了解最重要的结果。但是，它们越来越多地被用作起点，甚至作为中间步骤来帮助指导仍在进行中的研究。除此之外，作为注释的东西意味着什么也变得不像以前那么清晰。过去很明显，注释只是在进行了多项不同研究后建立的事物（例如基因产物的主要作用）。但是今天许多大型数据集都被社区处理，就像经典注释一样：作为额外比较的参考。

Bioconductor中注释正在进行的另一项变化是获得它们的方式。在过去，注释几乎完全作为单独的注释包存在[2],[3],[4]。今天的包仍然是一个巨大的注释来源。当前版本存储库包含八百多个注释包。此表总结了一些通常使用包访问的更重要的注释对象类：

Object Type	Example Package Name	Contents
TxDb	TxDb.Hsapiens.UCSC.hg19.knownGene	Transcriptome ranges for the known gene track of Homo sapiens, e.g., introns, exons, UTR regions.
OrgDb	org.Hs.eg.db	Gene-based information for Homo sapiens; useful for mapping between gene IDs, Names, Symbols, GO and KEGG identifiers, etc.
BSgenome	BSgenome.Hsapiens.UCSC.hg19	Full genome sequence for Homo sapiens.
Organism.dplyr	src_organism	Collection of multiple annotations for a common organism and genome build.
AnnotationHub	AnnotationHub	Provides a convenient interface to annotations from many different sources; objects are returned as fully parsed Bioconductor data objects or as the name of a file on disk.

注释资源在Bioconductor项目中占很大比例[1]。此外，还有一系列可用的在线资源，可以使用特定的软件包访问。本演练将介绍这些资源中最受欢迎的资源，并提供一些有关如何使用它们的高级示例。

传统上，Bioconductor注释资源在分析结束时使用。在大量数据分析之后，注释将被解释性地用于了解最重要的结果。但是，它们越来越多地被用作起点，甚至作为中间步骤来帮助指导仍在进行中的研究。除此之外，作为注释的东西意味着什么也变得不像以前那么清晰。过去很明显，注释只是在执行了多项不同研究后已经建立的那些事物（例如主要角色）但是，尽管注释包的流行，注释也越来越多地从像biomaRt这样的Web服务中被删除[5],[6],[7]或来自AnnotationHub [8]。这两者都代表了注释数据的巨大资源。

部分由于快速发展的景观，目前在单个文档中不可能涵盖Bioconductor中存在的每种可能的注释甚至每种注释。因此，我们将在此处讨论最流行的注释资源，并以旨在揭示用于访问它们的常用模式的方式对其进行描述。希望拥有此信息的用户能够对如何查找和使用以后不可避免地添加的其他资源进行有根据的猜测。将涵盖的主题包括以下内容：

3.安装

在本章中，我们使用了几个Bioconductor包。您可以使用biocLite()安装它们：

    ## try http:// if https:// URLs are not supported
    source("https://bioconductor.org/biocLite.R")
             biocLite(c(`"AnnotationHub", "Homo.sapiens",
               "Organism.dplyr",
               -  "TxDb.Hsapiens.UCSC.hg19.knownGene",
               "TxDb.Hsapiens.UCSC.hg38.knownGene",
               "BSgenome.Hsapiens.UCSC.hg19", "biomaRt",
               "TxDb.Athaliana.BioMart.plantsmart22"))

已使用的部分将详细介绍已安装软件包的使用情况。

4.AnnotationHub的使用

学习注释资源的列表顶部是相对较新的AnnotationHub包[8]。创建AnnotationHub是为了为最终用户提供方便的访问点，以便查找大量不同的注释对象以供Bioconductor使用。AnnotationHub中的资源很容易被发现，并作为熟悉的Bioconductor数据对象呈现给用户。因为它是最近添加的，AnnotationHub允许访问类似对象的广泛注释，其中一些注释甚至可能在几年前才被认为是注释。要开始使用AnnotationHub，用户只需要加载包，然后创建一个本地AnnotationHub对象，如下所示：

## snapshotDate(): 2018-04-23

ah <- AnnotationHub()

第一次调用AnnotationHub时，它将在您的系统上创建一个缓存目录，并下载集线器当前内容的最新元数据。从那时起，每当您下载其中一个集线器数据对象时，它也会将这些文件缓存在本地目录中，这样，如果再次请求该信息，您将能够快速访问它。

AnnotationHub对象的show方法将告诉您当前可以使用该对象访问的资源数量，并提供最常见的数据类型的高级概述。

ah

## AnnotationHub with 44923 records
## # snapshotDate(): 2018-04-23 
## # $dataprovider: BroadInstitute, Ensembl, UCSC, ftp://ftp.ncbi.nlm.nih.go...
## # $species: Homo sapiens, Mus musculus, Drosophila melanogaster, Bos taur...
## # $rdataclass: GRanges, BigWigFile, FaFile, TwoBitFile, Rle, OrgDb, Chain...
## # additional mcols(): taxonomyid, genome, description,
## #   coordinate_1_based, maintainer, rdatadateadded, preparerclass,
## #   tags, rdatapath, sourceurl, sourcetype 
## # retrieve records with, e.g., 'object[["AH2"]]' 
## 
##             title                                               
##   AH2     | Ailuropoda_melanoleuca.ailMel1.69.dna.toplevel.fa   
##   AH3     | Ailuropoda_melanoleuca.ailMel1.69.dna_rm.toplevel.fa
##   AH4     | Ailuropoda_melanoleuca.ailMel1.69.dna_sm.toplevel.fa
##   AH5     | Ailuropoda_melanoleuca.ailMel1.69.ncrna.fa          
##   AH6     | Ailuropoda_melanoleuca.ailMel1.69.pep.all.fa        
##   ...       ...                                                 
##   AH63652 | phastCons46wayPrimates.UCSC.hg19.chrUn_gl000247.rds 
##   AH63653 | phastCons46wayPrimates.UCSC.hg19.chrUn_gl000248.rds 
##   AH63654 | phastCons46wayPrimates.UCSC.hg19.chrUn_gl000249.rds 
##   AH63655 | phastCons46wayPrimates.UCSC.hg19.chrX.rds           
##   AH63656 | phastCons46wayPrimates.UCSC.hg19.chrY.rds

从上面的对象可以看出，有很多不同的资源可供使用。因此，通常当您获得AnnotationHub对象时，您要做的第一件事就是过滤它以删除不需要的资源。

幸运的是，AnnotationHub有几种不同类型的元数据可用于搜索和子集化。要查看不同的类别，您只需键入AnnotationHub对象的名称，然后从’$'运算符选项卡完成。要查看其中一个类别的所有可能内容，您可以将该值传递给唯一，如下所示：

unique(ah$dataprovider)

##  [1] "Ensembl"                              
##  [2] "UCSC"                                 
##  [3] "RefNet"                               
##  [4] "Inparanoid8"                          
##  [5] "NHLBI"                                
##  [6] "ChEA"                                 
##  [7] "Pazar"                                
##  [8] "NIH Pathway Interaction Database"     
##  [9] "Haemcode"                             
## [10] "BroadInstitute"                       
## [11] "PRIDE"                                
## [12] "Gencode"                              
## [13] "CRIBI"                                
## [14] "Genoscope"                            
## [15] "MISO, VAST-TOOLS, UCSC"               
## [16] "UWashington"                          
## [17] "Stanford"                             
## [18] "dbSNP"                                
## [19] "BioMart"                              
## [20] "GeneOntology"                         
## [21] "KEGG"                                 
## [22] "URGI"                                 
## [23] "ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/"

标记数据的最有价值的方法之一是根据将返回给您的R对象的类型。

unique(ah$rdataclass)

##  [1] "FaFile"           "GRanges"          "data.frame"      
##  [4] "Inparanoid8Db"    "TwoBitFile"       "ChainFile"       
##  [7] "SQLiteConnection" "biopax"           "BigWigFile"      
## [10] "AAStringSet"      "MSnSet"           "mzRpwiz"         
## [13] "mzRident"         "list"             "TxDb"            
## [16] "Rle"              "EnsDb"            "VcfFile"         
## [19] "igraph"           "OrgDb"

一旦确定了要用于查找感兴趣数据的元数据类型，就可以使用子集或查询方法将中心对象的大小减小到更易于管理的范围。例如，您只能选择字符串’GRanges’在元数据中的记录。正如您所看到的，GRanges是来自AnnotationHub的更流行的数据格式之一。

grs <- query(ah, "GRanges")
grs

## AnnotationHub with 19550 records
## # snapshotDate(): 2018-04-23 
## # $dataprovider: BroadInstitute, UCSC, Ensembl, Haemcode, Pazar, Gencode,...
## # $species: Homo sapiens, Mus musculus, Bos taurus, Pan troglodytes, Cani...
## # $rdataclass: GRanges
## # additional mcols(): taxonomyid, genome, description,
## #   coordinate_1_based, maintainer, rdatadateadded, preparerclass,
## #   tags, rdatapath, sourceurl, sourcetype 
## # retrieve records with, e.g., 'object[["AH5012"]]' 
## 
##             title                                            
##   AH5012  | Chromosome Band                                  
##   AH5013  | STS Markers                                      
##   AH5014  | FISH Clones                                      
##   AH5015  | Recomb Rate                                      
##   AH5016  | ENCODE Pilot                                     
##   ...       ...                                              
##   AH61310 | Vicugna_pacos.vicPac1.92.gtf                     
##   AH61311 | Xenopus_tropicalis.JGI_4.2.92.abinitio.gtf       
##   AH61312 | Xenopus_tropicalis.JGI_4.2.92.gtf                
##   AH61313 | Xiphophorus_maculatus.Xipmac4.4.2.92.abinitio.gtf
##   AH61314 | Xiphophorus_maculatus.Xipmac4.4.2.92.gtf

或者，您可以使用子集来仅选择特定字段的匹配项。

grs <- ah[ah$rdataclass == "GRanges",]

还提供了子集函数。

orgs <- subset(ah, ah$rdataclass == "OrgDb")
orgs

## AnnotationHub with 1691 records
## # snapshotDate(): 2018-04-23 
## # $dataprovider: ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/
## # $species: Escherichia coli, 'Caballeronia concitans', 'Chlorella vulgar...
## # $rdataclass: OrgDb
## # additional mcols(): taxonomyid, genome, description,
## #   coordinate_1_based, maintainer, rdatadateadded, preparerclass,
## #   tags, rdatapath, sourceurl, sourcetype 
## # retrieve records with, e.g., 'object[["AH61768"]]' 
## 
##             title                                      
##   AH61768 | org.Ag.eg.db.sqlite                        
##   AH61769 | org.At.tair.db.sqlite                      
##   AH61770 | org.Bt.eg.db.sqlite                        
##   AH61771 | org.Cf.eg.db.sqlite                        
##   AH61772 | org.Gg.eg.db.sqlite                        
##   ...       ...                                        
##   AH63468 | org.Salmonella_typhimurium_LT2.eg.sqlite   
##   AH63469 | org.Acinetobacter_baumannii.eg.sqlite      
##   AH63470 | org.Acinetobacter_genomosp._2.eg.sqlite    
##   AH63471 | org.Acinetobacter_genomospecies_2.eg.sqlite
##   AH63472 | org.Bacterium_anitratum.eg.sqlite

如果您真的需要访问所有元数据，可以使用mcols()将其作为DataFrame提取，如下所示：

meta <- mcols(ah)

此外，如果您是GUI的粉丝，您可以使用display方法在浏览器中查看数据，并将选定的行作为较小的AnnotationHub对象返回，如下所示：

sah <- display(ah)

调用此方法将生成一个基于Web的界面，如下图所示：

一旦你将AnnotationHub对象削减到合理的大小，并确定要检索哪些记录，那么你只需要使用’[[‘运算符来提取它们。使用’[['运算符，可以通过数字索引（1,2,3）或AnnotationHub ID提取。如果您选择使用前者，则只需提取您感兴趣的元素。因此，对于我们的链示例，您可能只想要这样的第一个：

res <- grs[[1]]

## downloading 0 resources

## loading from cache 
##     '/home/biocbuild//.AnnotationHub/5012'

head(res, n=3)

## UCSC track 'cytoBand'
## UCSCData object with 3 ranges and 1 metadata column:
##       seqnames          ranges strand |        name
##          <Rle>       <IRanges>  <Rle> | <character>
##   [1]     chr1       1-2300000      * |      p36.33
##   [2]     chr1 2300001-5400000      * |      p36.32
##   [3]     chr1 5400001-7200000      * |      p36.31
##   -------
##   seqinfo: 93 sequences (1 circular) from hg19 genome

4.1 AnnotationHub 练习

练习1：使用AnnotationHub提取来自智人（Homo sapiens）的UCSC数据，特别是来自hg19基因组的数据。在每个步骤过滤数据时，中心对象会发生什么？

练习2：现在你已经基本上缩小了UCSC基因组浏览器的hg19注释，让我们得到其中一个注释。找到oreganno轨道并将其保存到局部变量中。

5.OrgDb 对象

此时您可能想知道：这个OrgDb对象是什么？ OrgDb对象是一系列注释对象的一个​​成员，它们都通过一组共享方法表示隐藏数据。因此，如果仔细观察下面创建的狗对象，您会发现它包含Canis familiaris的数据（分类ID = 9615）。通过了解columns方法，您可以了解更多相关信息。

dog <- query(orgs, "Canis familiaris")[[1]]

## downloading 0 resources

## loading from cache 
##     '/home/biocbuild//.AnnotationHub/68517'

dog

## OrgDb object:
## | DBSCHEMAVERSION: 2.1
## | Db type: OrgDb
## | Supporting package: AnnotationDbi
## | DBSCHEMA: CANINE_DB
## | ORGANISM: Canis familiaris
## | SPECIES: Canine
## | EGSOURCEDATE: 2018-Apr4
## | EGSOURCENAME: Entrez Gene
## | EGSOURCEURL: ftp://ftp.ncbi.nlm.nih.gov/gene/DATA
## | CENTRALID: EG
## | TAXID: 9615
## | GOSOURCENAME: Gene Ontology
## | GOSOURCEURL: ftp://ftp.geneontology.org/pub/go/godatabase/archive/latest-lite/
## | GOSOURCEDATE: 2018-Mar28
## | GOEGSOURCEDATE: 2018-Apr4
## | GOEGSOURCENAME: Entrez Gene
## | GOEGSOURCEURL: ftp://ftp.ncbi.nlm.nih.gov/gene/DATA
## | KEGGSOURCENAME: KEGG GENOME
## | KEGGSOURCEURL: ftp://ftp.genome.jp/pub/kegg/genomes
## | KEGGSOURCEDATE: 2011-Mar15
## | GPSOURCENAME: UCSC Genome Bioinformatics (Canis familiaris)
## | GPSOURCEURL: 
## | GPSOURCEDATE: 2017-Apr6
## | UNIPROTGOSOURCENAME: UNIPROTGO
## | UNIPROTGOSOURCEDATE: Mon Apr  9 19:58:07 2018
## | UNIPROTGOSOURCEURL: ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/idmapping/idmapping_selected.tab.gz
## | ENSOURCEDATE: 2017-Dec04
## | ENSOURCENAME: Ensembl
## | ENSOURCEURL: ftp://ftp.ensembl.org/pub/current_fasta

## 
## Please see: help('select') for usage information

columns(dog)

##  [1] "ACCNUM"       "ALIAS"        "ENSEMBL"      "ENSEMBLPROT" 
##  [5] "ENSEMBLTRANS" "ENTREZID"     "ENZYME"       "EVIDENCE"    
##  [9] "EVIDENCEALL"  "GENENAME"     "GO"           "GOALL"       
## [13] "ONTOLOGY"     "ONTOLOGYALL"  "PATH"         "PMID"        
## [17] "REFSEQ"       "SYMBOL"       "UNIGENE"      "UNIPROT"

columns方法为您提供了一个数据类型向量，可以从您调用它的对象中检索这些数据类型。所以上面的调用表明可以从tetra对象中检索到几种不同的数据类型。

一个非常类似的方法是keytypes方法，它将列出所有也可以用作键的数据类型。

keytypes(dog)

##  [1] "ACCNUM"       "ALIAS"        "ENSEMBL"      "ENSEMBLPROT" 
##  [5] "ENSEMBLTRANS" "ENTREZID"     "ENZYME"       "EVIDENCE"    
##  [9] "EVIDENCEALL"  "GENENAME"     "GO"           "GOALL"       
## [13] "ONTOLOGY"     "ONTOLOGYALL"  "PATH"         "PMID"        
## [17] "REFSEQ"       "SYMBOL"       "UNIGENE"      "UNIPROT"

在许多情况下，列为列的大部分内容也将从keytypes调用返回，但由于这两件事不能保证相同，我们维护两个不同的方法。

现在您可以看到哪些类型的东西可以用作键，您可以调用keys方法来提取给定键类型的所有键。

head(keys(dog, keytype="ENTREZID"))

## [1] "399518" "399530" "399544" "399545" "399653" "403152"

如果您需要获取特定类型的所有ID，这很有用，但keys方法有一些额外的参数可以使它更灵活。例如，使用keys方法，您还可以提取包含“COX”的基因SYMBOLS，如下所示：

keys(dog, keytype="SYMBOL", pattern="COX")

##  [1] "COX5B"   "COX7A2L" "COX8A"   "COX15"   "COX5A"   "COX4I1"  "COX6A2" 
##  [8] "COX20"   "COX18"   "ACOX1"   "COX4I2"  "ACOX3"   "COX10"   "COX17"  
## [15] "COX11"   "ACOXL"   "COX7A1"  "COX1"    "COX2"    "COX3"    "COX19"  
## [22] "COX7B2"  "COX14"   "COX16"   "ACOX2"

或者，如果您确实需要其他键类型，则可以使用column参数为包含字符串“COX”的基因SYMBOLS提取ENTREZ GENE ID：

keys(dog, keytype="ENTREZID", pattern="COX", column="SYMBOL")

## 'select()' returned 1:1 mapping between keys and columns

##  [1] "474567"    "475739"    "476040"    "477792"    "478370"    "479623"   
##  [7] "479780"    "480099"    "482193"    "483322"    "485825"    "488790"   
## [13] "489515"    "503668"    "609555"    "611729"    "612614"    "804478"   
## [19] "804479"    "804480"    "100685945" "100687434" "100688544" "100846976"
## [25] "100855488"

但通常，您确实希望提取与特定键或一组键匹配的其他数据。为此，您可以使用两种方法。其中更强大的可能是选择。以下是查找SYMBOL基因和特定entrez基因ID的REFSEQ id的方法。

select(dog, keys="804478", columns=c("SYMBOL","REFSEQ"), keytype="ENTREZID")

## 'select()' returned 1:1 mapping between keys and columns

##   ENTREZID SYMBOL    REFSEQ
## 1   804478   COX1 NP_008473

当您调用它时，select将返回一个data.frame，它尝试为您请求的所有列填充匹配值。但是，如果要求select选择与键具有多对一关系的东西，则可能导致返回的数据对象的扩展。例如，观察当我们要求相同entrez基因ID的GO术语时会发生什么：

select(dog, keys="804478", columns="GO", keytype="ENTREZID")

## 'select()' returned 1:many mapping between keys and columns

##    ENTREZID         GO EVIDENCE ONTOLOGY
## 1    804478 GO:0004129      IEA       MF
## 2    804478 GO:0005506      IEA       MF
## 3    804478 GO:0005750      IEA       CC
## 4    804478 GO:0005751      IEA       CC
## 5    804478 GO:0006119      IEA       BP
## 6    804478 GO:0006123      IEA       BP
## 7    804478 GO:0009060      IEA       BP
## 8    804478 GO:0015990      IEA       BP
## 9    804478 GO:0016021      IEA       CC
## 10   804478 GO:0020037      IEA       MF
## 11   804478 GO:0045277      IEA       CC

因为有几个GO术语与基因“804478”相关联，所以在data.frame中最终会有很多行。如果您随后要求与原始密钥具有多对一关系的多个列，则可能会出现问题。如果你这样做，不仅结果会增加，而且也会变得非常难以使用。使用select时，更好的策略是选择性。

有时，您可能希望以比选择返回的data.frame对象更简单的方式查找匹配结果。当您只想查找每个键的一种值时尤其如此。对于这些情况，我们建议您查看mapIds方法。让我们看看如果请求与我们最近的select调用中相同的基本信息会发生什么，而是使用mapIds方法：

mapIds(dog, keys="804478", column="GO", keytype="ENTREZID")

## 'select()' returned 1:many mapping between keys and columns

##       804478 
## "GO:0004129"

如您所见，mapIds方法允许您简化返回的结果。默认情况下，mapIds仅返回每个键的第一个匹配元素。但是如果你真的需要在调用mapIds时返回的所有GO术语怎么办？那么你可以使用mapIds multiVals参数。这个参数有几个选项，我们已经看到默认情况下你只能返回’first’元素。但是你也可以返回’list’或’CharacterList’对象，或者你可以’过滤’或者返回’asNA’任何有多个匹配的键。您甚至可以定义自己的规则（作为函数）并将其作为参数传递给multiVals。让我们看一下返回列表时会发生什么：

mapIds(dog, keys="804478", column="GO", keytype="ENTREZID", multiVals="list")

## 'select()' returned 1:many mapping between keys and columns

## $`804478`
##  [1] "GO:0004129" "GO:0005506" "GO:0005750" "GO:0005751" "GO:0006119"
##  [6] "GO:0006123" "GO:0009060" "GO:0015990" "GO:0016021" "GO:0020037"
## [11] "GO:0045277"

现在您知道如何从OrgDb对象中提取信息，您可能会发现知道有一整套其他AnnotationDb派生对象也可以使用这些相同的五种方法（keytypes()，columns()，key() ，select()和mapIds()）。例如，存在ChipDb对象，InparanoidDb对象和TxDb对象，其分别包含关于微阵列探针，非精子同源伴侣或转录物范围信息的数据。还有更多专门的对象，如GODb或ReactomeDb对象，它们提供对GO和反应堆数据的访问。在下一节中，我们将查看这些对象中更受欢迎的类之一：TxDb对象。

5.1 OrgDb 练习

练习3：使用help（“SYMBOL”）查看不同列和键类型值的帮助页面。现在使用这些信息以及我们刚刚描述的查找entrez基因和染色体的基因符号“MSX2”。

练习4：在上一个练习中，我们不得不使用基因符号作为键。但在过去，这种行为有时是不可取的，因为一些基因符号被用作多个基因的官方符号。要了解这是否仍在发生，请利用entrez基因id唯一分配的事实，并从org.Hs.eg.db包中提取所有基因符号及其相关的entrez基因ID。然后检查符号是否冗余。

To be continue···