基于序列相似性和Z曲线方法重注释原核生物蛋白编码基因

刘硕，曾志，曾凡才，杜萌泽

研究报告

基于序列相似性和Z曲线方法重注释原核生物蛋白编码基因

刘硕1，曾志1，曾凡才2，杜萌泽2

1. 电子科技大学生命科学与技术学院，成都 611731 2. 西南医科大学基础医学院，分子生物与生物化学教研室，泸州 646000

随着测序技术的不断发展，产生了海量的基因组测序数据，极大地丰富了公共遗传数据资源。同时为了应对大量基因组数据的产生，基因组比较和注释算法、工具不断更新，使得联合多种注释工具得到更准确的蛋白编码基因的注释信息成为可能。目前公共数据库的原核生物基因组测序和装配有些是10多年前的，存在大量预测的功能未知的编码基因。为了提升美国国家生物信息中心(National Center for Biotechnology Information, NCBI)数据库中基因组的注释质量，本研究联合使用多种原核基因识别算法/软件和基因表达数据重注释1587个细菌和古细菌基因组。首先，利用Z曲线的33个变量从177个基因组原注释中识别获得3092个被过度注释为蛋白编码基因的序列；其次，通过同源比对为939个基因组中的4447个功能未知的蛋白编码基因注释上具体功能；最后，通过联合采用ZCURVE 3.0和Glimmer 3.02以及Prodigal这3种高精度的、广泛使用且基于算法不同而互补的基因识别软件来寻找漏注释基因。最终，从9个基因组中找到了2003个被漏注释的蛋白编码基因，这些基因属于多个蛋白质直系同源簇(clusters of orthologous groups of proteins, COG)。本研究使用新的工具并结合多组学数据重新注释早期测序的细菌和古细菌基因组，不仅为新测序菌株提供注释方法参考，而且这些重注释后得到的细菌基因序列也会对后续基础研究有所帮助。

细菌；重注释；Z曲线；假定ORFs；非蛋白编码ORFs

基因组分析是生命科学研究中非常关键的基础工作，只有获得准确的基因组注释才能为后续的分析和研究提供有力支撑并得到有意义的结果。目前，美国国家生物信息中心(National Center for Biotech­nology Information, NCBI)积累了大量10多年前完成测序的基因组。受限于早期有限的基因注释工具[1]，这些基因组注释的精确度有待于进一步提升。例如，在原核生物基因组中寻找小的基因时容易出现此类基因被漏注释而丢失的问题[2]。随着测序效率的提升，几何级数增长的测序量意味着这些注释错误将会通过同源搜索等方式被放大[3]。Breitwieser等[4]检查了所有的细菌和古细菌基因组后发现2250个原核生物基因组居然包含了人类基因组的信息。为了确保基因组信息的准确性[3,5]，相关数据库需要时常更新原始的基因组注释信息。当前随着各种组学数据的日益完善，通过联合多工具多组学的方式对基因组精确重注释也具备可行性。

自1995年第1个原核生物流感嗜血杆菌()基因组被注释以来，大量的原核生物基因组陆续完成测序和注释[5]。根据GenBank[6]的记录，自1999年4月至2019年10月，被测序序列的数量从3,525,418增长至216,763,706，增长倍数约为60倍。近年来基因识别软件，尤其是识别原核生物中基因的新方法和软件在被陆续开发和不断升级[7~14]，其中较广泛使用的有IPred[9]、Gene­MarkS[11]、Glimmer[12]、EasyGene[13]和ZCURVE[14]等。这些新的基因探测工具不仅可以用于注释基因组从而获得更精确的结果，而且还可赋予假定蛋白详尽的功能[15~19]。这些工具中基于Z-curve理论的软件ZCURVE 3.0[8]具有93.7%的准确率，与具有93.0%准确率的Glimmer 3.02[20]相当，而且Z-curve通过把单核苷酸、双核苷酸和三核苷酸的碱基组成转换成33个空间变量来表示DNA序列的特征，并通过对正样本和负样本的学习来对未知序列进行编码蛋白能力的预测[15,16]。该方法在原核生物的编码蛋白的预测上准确度高、应用广[21]。这两个软件联合使用可以获得互相补充的结果。更准确的细菌和古细菌基因组信息有助于更精确地推导生命的起源和演化，如Weiss等[22]对原核生物基因组的610万个蛋白编码基因的全部簇和进化树进行研究以推断最后的共同祖先(the last universal common ancestor, LUCA)，而不准确的注释会影响该研究结论的可靠性。

本研究对1587个测序和注释较早的细菌和古细菌菌株的基因组进行了重新注释，这些基因组包含功能已知的开放阅读框(open reading frame, ORF)和假定ORFs (功能未知的ORFs)，而假定ORFs中有些实际上是非蛋白编码基因的序列，此外基因组还包含被漏注释掉但实际上是可编码蛋白基因的那些序列。参考基因表达数据并联合3个注释准确率较高的原核生物识别软件(ZCURVE 3.0[8]、Glimmer 3.02[20]以及Prodigal[23])，移除掉一些之前被过度注释为基因的序列，同时给功能未知的ORFs注释了新的功能，并获得了之前漏注释的基因。

1 材料与方法

1.1 数据的选取

选用来自30个门(phylum)1587个20年内被陆续测序的细菌和古细菌的基因组来进行重注释(附表1)，基因组来自992个物种(附表2)。基因组的装配序列和注释信息均于2019年1月获取自NCBI数据库，它们的测序时间为1999年~2019年。

在使用33个Z-curve变量去除1587个基因组的过度注释基因时，选取它们的一些具有已知功能的ORFs作为正样本，并选取每个基因组的此类ORFs随机打乱1000次后产生的完全随机的序列作为负样本。

为判断计算识别的序列为真正的非蛋白编码序列，从GEO[24]和paxdb[25]下载了蛋白丰度或者RNA丰度数据。原基因组中第二类功能不确定的ORFs即假定ORFs如果被Z-curve 33变量方法识别为可能的非蛋白编码序列，同时其没有对应的蛋白或者mRNA被检测到，即认定其为过度注释为基因的序列，继而从基因组基因列表中移除。

本研究对1587个基因组中的9个模式物种基因组进行新基因的注释。它们分别是嗜酸氧化亚铁硫杆菌(ATCC 23270)、炭疽芽孢杆菌(str. Ames)、枯草芽孢杆菌(subsp. subtilis str. 168)、大肠杆菌(str. K-12 substr. MG1655)、流感嗜血杆菌(Rd KW20)、脑膜炎奈瑟球菌(MC58)、肠道沙门氏菌(subsp. enterica serovar Typhi str. CT18)、金黄色葡萄球菌(subsp. aureus NCTC 8325)和酿脓链球菌(SF370)。

1.2 重注释流程

整个重注释的流程分为3个部分：(1)识别过度注释为基因的序列；(2)未知功能ORFs的功能注释；(3)识别欠注释的基因(图1)，其具体的方法学细节部分见1.3、1.4、1.5和1.6。

1.3 寻找过度注释为基因的序列

Z-curve把DNA序列中出现在第一(1、4、7...)、第二(2、5、8...)和第三(3、6、9...)密码子位的A、G、C、T的碱基的频率分别用a1、g1、c1、t1；a2、g2、c2、t2；a3、g3、c3、t3来表示，根据Z曲线理论，一条DNA序列可以用三维平面中的点P(=1、2、3)来表示，而1、2、3所对应的x、y、z可以通过下边的公式来计算得出，据此得出Z-curve的单核苷酸的9个变量(9=3×3)。

图1 重注释流程图

当用12()或23()来表示DNA序列的第一第二密码子位或第二第三密码子位时，上边的公式可以衍生出以下公式，其中代表A、G、C或者T中的一种碱基。表示密码子相位组合，=12表示第一第二相位，同理=23表示第二第三相位，可以得到24个变量(24=3×4×2)。

得到以上Z-curve 33变量后，文章使用每个基因组对应的正样本和负样本进行训练。从正样本和负样本中随机抽取1/10的样本作为测试集，并进行十次交叉验证。最终平均预测准确率均大于99%，证明该预测方法稳定可靠。

将没有提供明确功能的ORFs作为测试集可寻找过度注释为基因的序列。预测为非蛋白编码基因的序列且同时没有在公共数据库中找到表达数据的假定ORFs将被归为过度注释为基因的序列。基因的表达数据来自不同的RNA测序集和GEO数据集(附表3)，这些数据可以帮助排除掉具有表达数据的基因。

1.4 给功能未知ORFs注释功能

假如一些假定/功能未知的ORFs没有被识别为非蛋白编码序列而是被预测为蛋白编码基因，研究就会通过同源搜索赋予其功能，使用的工具为BLAST[26]，赋予基因功能的同源搜索的参数设置为值(-value) < 1e-20，覆盖率(Coverage) > 80%以及一致性(Identity) > 70%。而假如被预测为非蛋白编码基因的序列其RNA-seq的数据显示为Count/TPM/ RPKM/CPM > 0，即存在，可以认为该序列具有表达数据，则需通过BLAST来注释其功能，根据经验参数，把赋予基因功能的同源搜索的参数同样设置为值 < 1e-20，覆盖率> 80%以及一致性 > 70%。

1.5 补充漏注释基因

基因组注释经常会丢掉一些真正的基因[2,27]，使用多种基因搜索工具可以尽可能减少该错误的发生。此处联合ZCURVE 3.0[8]、Glimmer 3.02[20]和Prodigal[23]来注释9种重要的模式生物的代表株。这3个软件识别出的不在原基因组基因列表中的ORFs将被全部加入候选基因序列集。对此集合中的序列用BLAST[26](https://blast.ncbi.nlm.nih.gov/Blast. cgi)进行同源比对后只保留那些与公共数据库中功能已知的基因有显著相似性的候选基因序列。由于此处目的为注释被遗漏的基因，研究根据经验参数把同源比对的参数设置为比1.4部分的阈值稍宽松的值，即< 1e-20，覆盖率 > 60%，一致性 > 60%，并用2个阈值对模式生物大肠杆菌的代表株str. K-12 substr. MG1655被注释后得到的246个新基因进行了2次同源比对，通过对比，发现宽松的阈值更能确保注释得到的新基因的完整性，详细结果参见2.4。

1.6 同源簇注释

结构和功能相似的蛋白编码序列属于同一个蛋白簇。在此，为注释得到的新基因进行同源簇COG注释可以帮助更好去理解基因组如何正常行使功能。研究同时联合使用了EggNOG[28]、WebMGA[29]和NCBI中的COG数据库[30]来注释这些基因序列的同源簇，其对应的COG编号可在附表7中找到。EggNOG包含2031个物种和19万个同源簇和对应功能注释数据。序列分析器WebMGA[29]和NCBI上的COG数据库[30]用于进一步完善同源簇功能注释。

2 结果与分析

2.1 1587个物种的注释现状

根据这1587个物种的具体测序时间，绘制了同一年份被测序的基因组数目占全部基因组比例的柱形图(图2，附表1)，其中2012年之前被测序的基因组占整体基因组的97%以上。根据CVTree[31]列出的门(phylum)的数量，与细菌和古细菌相关的门(phylum)为41个，选取的细菌和古细菌的基因组在门(phylum)水平的覆盖度为73%，研究选取的基因组尽可能包含了模式微生物的代表菌株如大肠杆菌和枯草芽孢杆菌等属于同一物种的注释时间不同的多种血清型，整体来说，选取的基因组在Refseq数据库中装配和注释的时间较早。

图2 1587个基因组的测序时间分布

2.2 原注释中3092个假定ORFs为非蛋白编码序列

在这1587个基因组的原始注释中，平均有1.26%的假定基因重注释后信息发生改变(图3A)。使用功能已知的ORFs作为正样本，选用对正样本序列随机打乱1000次后产生的序列作为负样本，进行训练后通过十重交叉验证后得到预测序列编码能力的平均准确率为0.9985 (图3B)。接着，研究依次对1587个细菌基因组中的功能未知ORFs进行预测，其中56,462个ORFs被预测为可能的非蛋白编码序列(附表4)。

由于计算方法本身的偏差，本研究使用实验数据进一步确定这些序列是否为非蛋白编码序列。当且仅当这些序列在公共数据库中查询不到表达数据，才认为它们确实不具有编码功能。研究查询了GEO等数据库并最终确定177个基因组中共3092个ORFs为过度注释为基因的序列(附表5)。需要重点提及的是其中43个基因组有20个以上的假定ORFs被识别为这类序列(表1)。

图3 假定ORFs的比率及不同准确率对应的基因组数量及大豆根瘤菌(B. japonicum USDA 110)3类序列两个密码子位GC含量

A：不同比率的重注释假定ORFs对应基因组的数目；B：不同识别非编码ORFs的准确率对应基因组的数目；C：大豆根瘤菌3类序列对应的第2和第3位密码子GC含量。

之前的研究已经报道过编码的ORFs在第二和第三密码子位的位置上有着不同的GC含量分布[16]，即大多数功能已知的ORFs的第三密码子位的GC含量比其第二密码子位的GC含量都高，而对于非编码ORFs，第二密码子位的GC含量则是接近第三密码子位的GC含量。在43个基因组中，大豆根瘤菌(USDA 110)基因组的假定ORFs被识别到了最多的非蛋白编码序列(147条)，其第二密码子位和第三密码子位的GC含量的分布和以上的论断是类似的(图3C)，这验证了本研究鉴定出的非编码ORFs具有可信性。

表1 识别出过度注释的ORFs多于20的菌株基因组的信息

2.3 4447个功能未知ORFs被注释上确定的功能

本研究用同源比对的方法为那些功能未知的ORFs注释上功能。最终，939个基因组中共计4447个ORFs被注释了确切功能(附表6)。其中在33个基因组中，有超过20个假定ORFs被注释上具体功能(表2)。

这些新注释的功能对生命活动非常重要。例如本研究发现1587个基因组中的601个功能未知的ORFs与一些质膜蛋白基因序列相似。质膜蛋白控制细胞内外物质的交换和信息的传递，参与信号通路的调控[32]。其中沙眼衣原体(D/UW-3/CX)中共有5个功能未知的蛋白被注释为质膜蛋白，沙眼衣原体可以引发炎症病变[33]，如它与子宫颈癌直接或间接相关。这些被注释上新功能的假定ORFs的详细信息请见附表6。

2.4 对新基因选取宽松阈值和严格阈值进行同源比对的结果

鉴于1.4和1.5部分进行同源比对的阈值不同，研究对模式生物大肠杆菌的代表株str. K-12 substr. MG1655注释后得到的246个新基因进行了第二次同源比对，选择的阈值为值 < 1e-20，覆盖率 > 80%以及一致性 > 70%，对两次比对做比较后发现有9个基因满足宽松的阈值(值 < 1e-20，覆盖率 > 60%，一致性 > 60%)，而当设置严格的阈值(值 < 1e-20，覆盖率 > 80%以及一致性 > 70%)时会被丢失掉。这9个基因的功能和对应值在表3中被列出，这些功能没有出现在str. K-12 substr. MG1655的现有注释中，因此认为它们可能对菌株发挥功能起着不可或缺的作用，而这恰恰是对该菌株注释新基因的意义所在。

表2 含有20个以上的假定ORFs被注释上准确功能的基因组的信息

2.5 9个基因组中新识别出2003个新基因

ZCURVE 3.0具有93.7%的准确率[8]，Glimmer 3.02具有93.0%的准确率[20]，这两个软件和假阳性率较低的Prodigal[23]被一起用来识别漏注释的基因，预测出的新的编码ORFs形成的并集被用来同源比对。最终，在9个基因组中识别到了2003个新基因(表4，附表7)。在最初的注释中基因数量相对较少的酿脓链球菌和流感嗜血杆菌分别得到了104和123个新基因，新基因的数量大约是原始基因数目的6%。脑膜炎奈瑟球菌和肠道沙门氏菌获得了很多新注释到的基因，它们占据了其原始基因数目的14%。这些被新识别的基因的名字和其起始位置的信息存放在附表7中。

注释同源蛋白簇对研究具有相似功能的基因帮助颇多。通过综合EggNOG[28]、WebMGA[29]和COG数据库[30]，共有1073个新识别的基因被注释到对应的COGs(图4)，其丰度就是某蛋白质直系同源簇所包含的新识别的基因占全部新识别的基因的比例。这1073条序列的详细COG注释可以在附表7中查询到，而所有1587个基因组重注释的具体信息则可以在附表8中被查到，该表列出的具体信息包含其基本信息即NC_ID和基因组大小(bp)，还包含其原始注释信息即蛋白质数目、功能已知的基因和假定基因，并列出了通过重注释得到的注释信息即假定ORFs中被注释上新功能的ORFs数目和非编码ORFs的数目，还提供了9个基因组新基因的数量，并提供了通过Z-curve 33变量方法识别非编码ORFs的准确率。

3 讨论

本研究结合3种原核基因组注释软件和表达数据重注释了1587个细菌和古细菌的基因组。首先，识别了过度注释为基因的序列，隶属于177个基因组的3092个假定ORFs被识别为非蛋白编码序列；接下来，还赋予了假定ORFs以功能，939个基因组的4447个基因被通过相似性搜索而注释获得准确功能，并且在9个基因组中识别了2003个属于多个蛋白质直系同源簇的新基因。

重注释可以保证进一步分析的数据的准确性。在真核生物中有研究表明多拷贝的基因会出现功能的分化[34]，原核生物也存在多拷贝的基因，它们是否也存在这种分化，可以通过开展重注释后获得准确的信息来进一步研究。而正如前文所提及的那样，伴随着测序数据的增加和技术的不断进步，对已测序菌株的基因组进行重新注释的必要性也在增大，同时根据Breitwieser的研究[4]，微生物基因组可能会被人类基因组所污染而导致基因的丢失，这就更加要求研究者们用多种方法来更新注释。本次进行重新注释的1587个菌株可归为992个物种，每个物种涵盖多种血清型，例如沙门氏菌根据抗原的不同可以分为不同的血清型[35]，其基因组层面可能会存在差异。Liu等[33]通过全基因组序列对6个沙门氏菌(Salmonella)菌株建立进化树显示亚种的相同血清型的C和A进化距离较远，却与不同血清型的进化距离较近。研究通过CVTree[31]对15个亚种构建进化树发现血清型之间会呈现不同的聚集，与Liu等[33]绘制的进化树一致(附图1)。因此，通过重注释不同血清型的菌株来使得基因组信息更完善，从而帮助深入探索血清型之间更准确的演化关系。而研究参考的数据集为1587个基因组中已经被注释为具有明确功能的ORFs，为了验证这些基因具有表达数据，针对大肠杆菌的功能已知的基因进行验证，研究下载了对应该物种的表达数据并寻找这些功能已知的基因中有多少具有表达数据，通过与GEO数据库中GSE56133[36]和GSE118058[37]的表达数据比较，发现2835个功能已知的基因中有2806个在这两个数据集中有表达数据，这进一步证实了现有明确功能注释的基因可以作为参考集来进行预测。

表3 大肠杆菌(E.coli str. K-12 substr. MG1655)满足宽松阈值的新基因

表4 9个菌株的名称、NC序列号、基因组大小、基因总数和新注释基因的数目

图4 特定同源簇对应的新基因占全部新基因的比例

A：RNA加工和修饰；B：染色质结构和动力学；Y：核结构；Z：细胞骨架；W：细胞外结构；D：细胞周期控制，细胞分裂，染色体分裂；O：翻译后修饰，蛋白反转，伴侣；Q：次生代谢产物的生物合成、运输和分解代谢；I：脂质转运与代谢；J：翻译，核糖体结构和生物发生；H：辅酶运输与代谢；F：核苷酸转运与代谢；V：防卫机制；N：细胞迁移；C：能量产生和转化；U：细胞内传输，分泌和囊泡转运；P：无机离子转运与代谢；T：信号转导机制；K：转录；E：氨基酸转运和代谢；M：细胞壁和细胞膜的生物发生；G：碳水化合物运输和代谢；S：功能未知；R：只能预测大致功能；L：复制重组和修复。

得益于众多基因识别软件的升级如ZCURVE 3.0[8]和不断更新的公共数据库如NCBI的GEO[24]，使得寻找新的编码蛋白ORFs和用表达数据验证成为可能。在识别新的基因方面，本研究联合ZCURVE 3.0[8]、Glimmer 3.02[20]和Prodigal[23]对原核基因组重新注释，其中Glimmer 3.02是基于隐马尔可夫模型，识别准确率高达93.0%，但是它对序列的局部特征有着强烈的依赖性，而Prodigal是通过对现有的细菌和古细菌基因组注释信息进行训练，因此对已知基因/保守基因的识别效果优于其余的软件，但预测未报道的基因时会有些许偏差。而ZCURVE 3.0是基于DNA序列的全局统计特征，它将Fisher线性判别替换为支持向量机(SVM)[38]以提高灵敏度，并且鉴于在预测过程中由于负样本的随机性容易产生假阳性的基因，ZCURVE 3.0依据了ORFs核酸分布的花瓣模型[39]在训练集中产生5类负样本并逐次分类，然后保留在多次分类中均被预测为基因的那些ORFs，并且算法还把33个包含零阶和一阶的Z曲线变量增加为额外包含二阶和三阶Z曲线变量的765个变量，可以最大化地优化程序的预测效果。选取三个程序混合预测更能保证预测的准确率。

总之，本研究基于序列的全局特征和局部特征，结合同源比对对多个物种的多种菌株的基因组进行了全面的重新注释。未来，伴随着基因组多组学数据的增加以及相应的注释工具和数据库的完善[40~43]，公共数据库里基因组的注释将会更加准确。

致谢

感谢电子科技大学生命科学与技术学院的郭锋彪教授在研究开展中给予的帮助。

附录：

附图和附表详见文章电子版www.chinagene.cn。

[1] Mørk S, Holmes I. Evaluating bacterial gene-finding HMM structures as probabilistic logic programs., 2012, 28(5): 636–642.

[2] Warren AS, Archuleta J, Feng WC, Setubal JC. Missing genes in the annotation of prokaryotic genomes., 2010, 11(1): 131.

[3] Salzberg SL. Next-generation genome annotation: we still struggle to get it right., 2019, 20(1): 92.

[4] Breitwieser FP, Pertea M, Zimin AV, Salzberg SL. Human contamination in bacterial genomes has created thousands of spurious proteins., 2019, 29(6): 954–960.

[5] Fleischmann RD, Adams MD, White O, Clayton RA, Kirkness EF, Kerlavage AR, Bult CJ, Tomb JF, Dougherty BA, Merrick JM, Mckenney K, Sutton G, FitzHugh W, Fields C, Gocayne JD, Scott J, Shirley R, Liu LL, Glodek A, Kelley JM, Weidman JF, Phillips CA, Spriggs T, Hedblom E, Cotton MD, Utterback TR, Hanna MC, Nguyen DT, Saudek DM, Brandon RC, Fine LD, Fritchman JL, Fuhrmann JL, Geoghagen NSM, Gnehm CL, McDonald LA, Small KV, Fraser CM, Smith HO, Venter JC. Whole-genome random sequencing and assembly ofRd., 1995, 269(5223): 496–512.

[6] Benson DA, Karsch-Mizrachi I, Lipman DJ, Ostell J, Sayers EW. Genbank., 2009, 37(Data­base issue): D26–D31.

[7] Yu JF, Xiao K, Jiang DK, Guo J, Wang JH, Sun X. An Integrative method for identifying the over-annotated protein-coding genes in microbial genomes., 2011, 18(6): 435–449.

[8] Hua ZG, Lin Y, Yuan YZ, Yang DC, Wei W, Guo FB. ZCURVE 3.0: identify prokaryotic genes with higher accuracy as well as automatically and accurately select essential genes., 2015, 43(W1): W85– W90.

[9] Zickmann F, Renard BY. IPred-integrating ab initio and evidence based gene predictions to improve prediction accuracy., 2015, 16(1): 134.

[10] Keilwagen J, Wenk M, Erickson JL, Schattat MH, Grau J, Hartung F. Using intron position conservation for homology-based gene prediction., 2016, 44(9): e89.

[11] Besemer J, Lomsadze A, Borodovsky M. GeneMarkS:a self-training method for prediction of gene starts in microbial genomes. Implications for finding sequence motifs in regulatory regions., 2001, 29(12): 2607–2618.

[12] Kelley DR, Liu B, Delcher AL, Pop M, Salzberg SL. Gene prediction with Glimmer for metagenomic sequences augmented by classification and clustering., 2012, 40(1): e9.

[13] Larsen TS, Krogh A. EasyGene-a prokaryotic gene finder that ranks ORFs by statistical significance., 2003, 4(1): 21.

[14] Guo FB, Ou HY, Zhang CT. ZCURVE: a new system for recognizing protein-coding genes in bacterial and archaeal genomes., 2003, 31(6): 1780–1789.

[15] Du MZ, Guo FB, Chen YY. Gene re-annotation in genome of the extremophileby using bioinformatics methods., 2011, 29(2): 391–401.

[16] Guo FB, Xiong LF, Teng JL, Yuen KY, Lau SK, Woo PC. Re-annotation of protein-coding genes in 10 complete genomes of Neisseriaceae family by combining similarity- based and composition-based methods., 2013, 20(3): 273–286.

[17] Lei Y, Kang SK, Gao JX, Jia XS, Chen LL. Improved annotation of a plant pathogen genomepv. oryzae PXO99A., 2013, 31(3): 342–350.

[18] Mao Y, Yang X, Liu Y, Yan Y, Du Z, Han Y, Song Y, Zhou L, Cui Y, Yang R. Reannotation ofstrain 91001 Based on Omics Data., 2016, 95(3): 562–570.

[19] Pfeiffer F, Bagyan I, Alfaro‐Espinoza G, Zamora‐Lagos MA, Habermann B, Marin‐Sanguino A, Oesterhelt D, Kunte HJ. Revision and reannotation of theDSM 2581T genome., 2017, 6(4): e00465.

[20] Delcher AL, Bratke KA, Powers EC, Salzberg SL. Identifying bacterial genes and endosymbiont DNA with Glimmer., 2007, 23(6): 673–679.

[21] Zhang R, Zhang CT. A Brief Review:The Z-curve theory and its application in genome analysis., 2014, 15(2): 78–94.

[22] Weiss MC, Sousa FL, Mrnjavac N, Neukirchen S, Roettger M, Nelson-Sathi S, Martin WF. The physiology and habitat of the last universal common ancestor., 2016, 1(9): 16116.

[23] Hyatt D, Chen GL, Locascio PF, Land ML, Larimer FW, Hauser LJ. Prodigal: prokaryotic gene recognition and translation initiation site identification., 2010, 11(1): 119.

[24] Barrett TT, Troup DB, Wilhite SE, Ledoux P, Rudnev D, Evangelista C, Kim IF, Soboleva A, Tomashevsky M, Marshall KA, Phillippy KH, Sherman PM, Muertter RN, Edgar R. NCBI GEO: archive for high-throughput functional genomic data., 2009, 37(Database issue): D885–D890.

[25] Wang M, Weiss M, Simonovic M, Haertinger G, Schrimpf SP, Hengartner MO, von Mering C. PaxDb, a database of protein abundance averages across all three domains of life., 2012, 11(8): 492–500.

[26] McGinnis S, Madden TL. BLAST: at the core of a powerful and diverse set of sequence analysis tools., 2004, 32(Suppl.2): W20–W25.

[27] Wood DE, Lin H, Levy-Moonshine A, Swaminathan R, Chang YC, Anton BP, Osmani L, Steffen M, Kasif S, Salzberg SL. Thousands of missed genes found in bacterial genomes and their analysis with COMBREX., 2012, 7(1): 37.

[28] Huerta-Cepas J, Szklarczyk D, Forslund K, Cook H, Heller D, Walker MC, Rattei T, Mende DR, Sunagawa S, Kuhn M, Jensen LJ, Mering CV, Bork P. eggNOG 4.5: a hierarchical orthology framework with improved functional annotations for eukaryotic, prokaryotic and viral sequences., 2016, 44(Database issue): D286–D293.

[29] Wu ST, Zhu ZW, Fu LM, Niu BF, Li WZ. WebMGA: a customizable web server for fast metagenomic sequence analysis., 2011, 12(1): 444.

[30] Tatusov RL, Galperin MY, Natale DA, Koonin EV. The COG database: a tool for genome-scale analysis of protein functions and evolution., 2000, 28(1): 33–36.

[31] Qi J, Luo H, Hao BL. CVTree: a phylogenetic tree reconstruction tool based on whole genomes., 2004, 32: W45–W47.

[32] Hockenbery D, Nuñez G, Milliman C, Schreiber RD, Korsmeyer SJ. Bcl-2 is an inner mitochondrial membrane protein that blocks programmed cell death., 1990, 348(6299): 334–336.

[33] Liu WQ, Feng Y, Wang Y, Zou QH, Chen F, Guo JT, Peng YH, Jin Y, Li YG, Hu SN, Johnson RN, Liu GR, Liu SL.C: Genetic divergence fromand pathogenic convergence with., 2009, 4(2): e4510.

[34] Vankuren NW, Long M. Gene duplicates resolving sexual conflict rapidly evolved essential gametogenesis functions., 2018, 2(4): 705–712.

[35] Minor LL, Bockemühl J. 1987 supplement (no.31) to the schema of Kauffmann-White., 1988, 139(3): 331–335.

[36] Dwyer DJ, Belenky PA, Yang JH, Macdonald IC, Martell JD, Takahashi N, Chan CT, lobritz MA, Braff D, Schwarz EG, Ye JD, Pati M, Vercruysse M, Ralifo PS, Allison KR, Khalil AS, Ting AY, Walker GC, Collins JJ. Antibiotics induce redox-related physiological alterations as part of their lethality., 2014, 111(20): E2100–E2109.

[37] Hadjeras L, Poljak L, Bouvier M, Morin-Ogier Q, Canal l, Cocaign-Bousquet M, Girbal L, Carpousis AJ. Detachment of the RNA degradosome from the inner membrane ofresults in a global slowdown of mRNA degradation, proteolysis of RNase E and increased turnover of ribosome-free transcripts., 2019, 111(6): 1715–1731.

[38] Kim S, Yu Z, Kil RM, Lee M. Deep learning of support vector machines with class probability output networks., 2015, 64: 19–28.

[39] Guo FB. The distribution patterns of bases of protein- coding genes, non-coding ORFs, and intergenic sequences inPA01 genome and its implication., 2007, 25(2): 127–133.

[40] Uyar B, Yusuf D, Wurmus R, Rajewsky N, Ohler U, Akalin A. RCAS: an RNA centric annotation system for transcriptome-wide regions of interest., 2017, 45(10): e91.

[41] Huang Y, Liu Q, Chi LJ, Shi CM, Wu Z, Hu M, Shi H, Chen H. Application of BIG-Annotator in the genome sequencing data functional annotation and genetic diagnosis., 2018, 40(11): 1015–1023.黄莹, 刘琪, 池连江, 石承民, 吴祯, 胡敏, 石宏, 陈华. BIG-Annotator: 基因组测序数据高效功能注释及其在遗传诊断中的应用. 遗传, 2018, 40(11): 1015–1023.

[42] Bick JT, Zeng SQ, Robinson MD, Ulbrich SE, Bauersachs S. Mammalian Annotation Database for improved annotation and functional classification of Omics datasets from less well-annotated organisms., 2019, 2019: 1–16.

[43] Ravindran SP, Lüneburg J, Gottschlich L, Tams V, Cordellier M. Daphnia stressor database: Taking advantage of a decade of Daphnia ‘-omics’ data for gene annotation., 2019, 9(1): 11135.

附图1亚种不同血清型的系统发育树

Supplementary Fig. 1 The phylogenetic tree of different serological types ofsub-species

Comprehensive re-annotation of protein-coding genes for prokaryotic genomes by Z-curve and similarity-based methods

Shuo Liu1, Zhi Zeng1, Fancai Zeng2, Mengze Du2

The development of sequencing technology has generated huge genomicsequencing information and largely enriched public genetic resources. To analyze such big data, the algorithms and tools for comparison and annotation of genomes are updated continually, enabling genome annotation with higher accuracyvarious annotation tools. Many prokaryotic genomes in public database were sequenced and assembled more than a decade ago, and they contained multiple genes with unknown functions. To improve the current annotation for those genomes in NCBI, we re-annotate 1587 bacterial and archaeal genomes using multiple prokaryotic gene recognition algorithms/softwares and gene expression data.The 33 Z-curve variableswere appliedto recognize sequences that were over-annotated to genes of 1587 bacterial and archaeal genomes deposited in public databases, anda total of 3092sequences belonging to 177 genomes were recognized as sequences over-annotated as protein-coding genes.Next, 4447 protein-coding genes with unknown functions from 939 genomes were annotated with definite functions by similarity search. Finally, we recognized 2003 missed protein-coding genesthat belong to known COG (clusters of orthologous groups of proteins) of nine genomes using three methods (ZCURVE 3.0, Glimmer3.02 and Prodigal), which are accurate and frequently used for gene finding. Their algorithms are different and complementary. This is a comprehensive study for re-annotation of bacterial and archaeal genomes with new tools combining multi-omics data, which should provide a reference for annotation of newly sequenced strains, and also benefit further fundamental researches with the bacterial gene sequences obtained after re-annotation.

bacteria; re-annotation; Z-curve; hypothetical ORFs; non-coding ORFs

2020-02-20;

2020-05-11

电子科技大学理科实力提升计划项目(编号：Y0301902610100202)资助[Supported by Science Strength Improvement Plan of University of Electronic Science and Technology of China (No. Y0301902610100202)]

刘硕，在读博士研究生，专业方向：微生物基因组学。E-mail: liushuo20022020@gmail.com

杜萌泽，博士，讲师，研究方向：生物信息学。E-mail: du_mengze@foxmail.com

10.16288/j.yczz.20-022

2020/5/28 11:15:04

URI: http://kns.cnki.net/kcms/detail/11.1913.R.20200528.0949.001.html

(责任编委: 包其郁)