基因功能注释

生信喵

背景
       预测得到一个物种的全部基因之后，接下来自然而然的问题就出现了。这些基因都有哪些功能呢?这就需要进行基因功能注释。基因功能注释就是将待查基因与已知数据库进行比对，如果比对上则认为二者为同源基因，执行相同的功能。宏基因组中通常包括很多新发现的基因，无法比对上已知数据库。所以，在宏基因组研究中，一部分或者大部分基因无法注释得到功能属于正常现象。

一、基因功能注释原理
       基因功能无法实现软件预测，只能与数据库进行比对。基因功能注释操作起来并不十分困难，主要就是与数据库比对。一般的数据库至少包含两部分内容：一是基因的序列，包含核酸和氨基酸，一般为 fasta 格式。另外就是基因的功能信息。这两部分可以在一个文件内，例如将功能信息写到每个基因的 ID 行，也可以放到不同文件中，每个基因 ID 对应相应功能。首先，将物种的基因集序列，一般我们采用氨基酸序列，与数据库进行 blast 同源比对。比对完了之后对 blast 结果进行过滤，因为我们需要在数据库中找到与物种基因集里面基因保持同源的序列。因为二者是同源关系的基因，所以就认为二者执行同样的功能。基因功能注释大体上就是这样的过程，原理并不复杂。这里我们其实注意到。对于基因功能注释，数据库的影响非常大。如果数据库中没有找到同源的基因，那么这个基因就无法注释出来。另外，如果数据库中的信息有错误，基因功能的注释也会出现错误，而且这个错误会逐渐累积，问题就会非常严重。所以，数据库的准确性非常重要。

二、基因功能注释数据库
2.1 NCBI nr

1. nucleotide sequence database, with entries from all traditional divisions of
   
2. GenBank, EMBL, and DDBJ excluding bulk divisions (gss, sts, pat, est, and htg
   
3. divisions. wgs entries are also excluded. Not non-redundant.

复制代码

      NR 数据库是 NCBI 中的一个非冗余的蛋白数据库。它包含从 GeneBank 核酸序列翻译而来的非冗余序列，并且还收录了其他蛋白数据库的非冗余序列，包括 RefSeq、PDB、SwissProt、PIR 和 PRF。
       nr 数据库地址：https://ftp.ncbi.nlm.nih.gov/blast/db/FASTA/nr.gz

2.2 uniport
       https://www.uniprot.org/
       全球蛋白质资源（Universal Protein Resource,UniProt）是全球有关蛋白质方面信息最全面的资源库。UniProt 由 UniprotKB、UniRef 和 UniParc 组成是蛋白质序列以及功能信息的集中资源，且其具有最小的冗余。UniProt 是对 PIR、TrEMBLe 以及 Swissprot 的信息进行组合而构成的。
      
       uniprot 数据库结构
       Uniprot 数据库包含蛋白质序列，功能信息，研究论文索引等信息，整合了包括 EBI（ European Bioinformatics Institute）、SIB（the Swiss Institute of Bioinformatics）、PIR（Protein Information Resource）三大数据库的资源。
       EBI（ European Bioinformatics Institute）：欧洲生物信息学研究所（EMBL-EBI）是欧洲生命科学旗舰实验室 EMBL 的一部分。位于英国剑桥欣克斯顿的惠康基因组校内，是世界上基因组学领域最强的存在之一。
       SIB（the Swiss Institute of Bioinformatics）：瑞士日内瓦的 SIB 维护着 ExPASy（专家蛋白质分析系统）服务器，这里包含有蛋白质组学工具和数据库的主要资源。
       PIR（Protein Information Resource）：PIR 是美国国家生物医学研究基金会（NBRF）于 1984 年创立，旨在协助研究人员识别和解释蛋白质序列信息。
       UniProt 提供了完全分类的、有丰富且准确注释信息的基于知识的蛋白质序列信息，且有广泛的交叉引用以及多种查询界面。

2.3 cog
       COG（Clusters of Orthologous Groups of proteins，直系同源蛋白簇）构成每个 COG 的蛋白都是被假定为来自于一个祖先蛋白，因此是 orthologs 或者是 paralogs。
       通过把所有完整基因组的编码蛋白一个一个的互相比较确定的。在考虑来自一个给定基因组的蛋白时，这种比较将给出每个其他基因组的一个最相似的蛋白（因此需要用完整的基因组来定义 COG），这些基因的每一个都轮番的被考虑。如果在这些蛋白（或子集）之间一个相互的最佳匹配关系被发现，那么那些相互的最佳匹配将形成一个 COG。这样，一个 COG 中的成员将与这个 COG 中的其他成员比起被比较的基因组中的其他蛋白更相像。
       NCBI COG 构建了不同系统分类分支的 COG 簇，比如 arCOG，KOG, POG 等。
       主页：https://www.ncbi.nlm.nih.gov/COG/
       COG 单字母描述，详见 http://www.sbg.bio.ic.ac.uk/~phunkee/html/old/COG_classes.html

2.4 kegg 库
       http://www.genome.jp/kegg/
       KEGG（Kyoto Encyclopedia of Genes and Genomes，京都基因与基因组百科全书，是一个关于基因组、酶促途径以及生物化学物质的在线数据库。虽然决定生物体基因分类的基因组测序工程有了飞速的发展，但对单个基因功能的研究仍然相差甚远。同时活细胞的生物学功能是许多分子相互作用的结果，不能仅仅归功于单个基因或单个分子。KEGG 将基因组信息和高一级的功能信息有机地结合起来，通过对细胞内已知生物学过程的计算机化处理和将现有的基因功能解释标准化，对基因的功能进行系统化的分析。KEGG 的另一个任务是一个将基因组中的一系列基因用一个细胞内的分子相互作用的网络连接起来的过程，如一个通路或是一个复合物，通过它们来展现更高一级的生物学功能。

2.5 go 功能数据库
       http://geneontology.org/
       GO(gene ontology)是基因本体联合会(Gene Onotology Consortium)所建立的数据库，旨在建立一个适用于各种物种的，堆积因和蛋白质功能进行限定和描述的，并能随着研究不断深入而更新的语言词汇标准．GO 是多种生物本体语言中的一种，提供了三层结构的系统定义方式，用于描述基因产物的功能。涵盖了基因的细胞组分（cellular component）、分子功能（molecular function）、生物学过程（biological process）。

三、使用案例
       这里我们使用 eggnog-mapper 工具来进行基因功能注释。eggnog-mapper 是一个非常方便 的 基 因 功 能 注 释 流 程 。 可 以 自 动 化 完 成 基 因 功 能 注 释 工 作 ， 其 内 置 了COG/KOG/KEGG/GO/BiGG 等数据库，也可以自行创建注释数据库。
       eggnog-mapper 在不同的系统分类水平都进行了构建直系同源簇，当前使用的 eggNOG v5.0 版本数据库，包含 5,090 代表性的基因组，其中包括 4445 个细菌，168 古细菌，477个真核生物以及 2502 个病毒。当前 v2 版本的 eggnog-mapper，不再支持 hmmer 比对。

3.1 软件安装
       网址：https://github.com/eggnogdb/eggnog-mapper
       在线版本：http://eggnog-mapper.embl.de/
       文档：https://github.com/eggnogdb/eggnog-mapper/wiki/eggNOG-mapper-v2

1. #安装 eggnog-emapper
   
2. mamba install -c bioconda -y eggnog-mapper
   
3. mkdir eggnog_database;
   
4. download_eggnog_data.py -y --data_dir eggnog_database
   
5. #镜像数据下载

复制代码

1. ftp://download.nmdc.cn/tools/eggnog/eggnog.db.gz
   
2. ftp://download.nmdc.cn/tools/eggnog/eggnog_proteins.dmnd.gz

复制代码

3.2 软件使用

1. #基因功能注释
   
2. emapper.py -i MGH78578.pep --output annotation -m diamond --data_dir
   
3. eggnog_database

复制代码

      选项参数：

1. -i: 输入文件，最好基因的氨基酸文件
   
2. -o: 输出结果前缀
   
3. -m: 使用 HMMER 策略还是 DIAMOND 策略，默认使用 HMMER，新版本只支持 diamond
   
4. --cpu：使用的线程数
   
5. --translate：如使用的核酸序列，选择 HMMER 策略时需要先翻译成氨基酸序列
   
6. --usemem： 将 emapper.db 读入内存
   
7. --output_dir：输出结果文件夹
   
8. --report_orthologs： 列出所有进行功能转移的直系同源基因
   
9. --no_annot: 只汇总鉴定的最佳 seed 序列以及相应的 E 值和 Bitscore 值
   
10. --no_search:可直接基于--no_annot 的结果进行后续功能注释
    
11. --target_orthologs: one2one,many2one,one2many,many2many,all 可选。
    
12. --data_dir: 数据库目录
    
13. --tax_scope: 指定选择的直系同源基因的物种分类范围，默认为自动判断。

复制代码

生信喵

四、结果解读
       最终会生成两个文件，分别是 mg.emapper.annotations 和 mg.emapper.seed_orthologs
       第一列：查询序列名称；
       第二列：eggNOG 种子序列；
       第三列：eggNOG 种子序列 evalue；
       第四列：eggNOG 种子序列 bit score；
       第五列：预测基因名称；
       第六列：GO_terms, 预测的 GO，分号分隔；
       第七列：KEGG_KO: 预测的 KO,分号分隔；
       第八列：BiGG_Reactions: BiGG 代谢反应预测，分号分隔；
       第九列：eggNOG Taxonomic Scope 信息；
       第十列：匹配的 OGs;
       第十一列：best_OG|evalue|score: Best matching Orthologous Groups (only in HMM mode)
       第十二列：COG 功能分类；
       第十三列：eggNOG 功能描述；

五、在线分析
       官方网址：http://eggnog-mapper.embl.de/
       输入文件：基因的氨基酸序列文件，最多可以包含 10 万条序列，基本上满足目前绝大部分物种了。
       简单，但是需要排队，且部分数据因单位保密不可以外传。