CN101957892A - 一种全基因组复制事件的检测方法和系统 - Google Patents

Abstract

本发明公开一种全基因组复制的检测方法和系统。该方法包括：获取与目标物种基因组上基因集的蛋白质序列以及基因在基因组上的位置信息；对基因集的蛋白质序列进行蛋白比蛋白两两比对获得基因集的同源基因对；利用动态规划获得基因组上具有共线性的同源区块；获得同源区块中同源基因对的核酸序列的同源位点对位排列；计算同源区块中同源基因对的遗传距离，获得年龄分布；根据年龄分布确定全基因组复制事件。本发明的方法和系统，全面考虑了全基因组上的所有基因的序列信息，排除了基因组上非保守和选择压力受到释放的区域的干扰，并利用动态规划算法实现最佳复制区域的寻找，确定的同源基因对和同源区块全面、准确，从而能够更准确地实现全基因组复制事件的检测。

Description

一种全基因组复制事件的检测方法和系统

技术领域

本发明涉及生物信息学技术领域，尤其涉及一种全基因组复制事件的检测方法和系统。

背景技术

全基因组复制事件(Whole Genome Duplication，WGD)，也称古多倍化事件，在植物界是一个频繁发生和普遍存在的现象，这些全基因组范围内的复制对植物的多样化和适应性产生了广泛而深远的影响。随着大量的植物物种全基因组被测序，以及对其基因组内功能元件的研究，推断基因组的演变历史及其对物种形成和多样化的影响已越来越受到人们的关注。基因组上的序列信息蕴含了一个物种所有进化历史所遗留下来的痕迹，因此，鉴定基因组上序列的进化过程和变化模型为人们了解该物种提供了很好的信息。在基因组测序技术成熟之前，一般的方法是通过估计相关物种基因组大小，染色体数目的变化，或通过染色体荧光定位杂交的方法来推断某植物物种是否经历过全基因组复制。但是，由于全基因组复制事件是一个特别复杂的动态过程，并且随后会伴随大量基因的丢失和染色体重排等复杂变化，该方法检测出来的结果只能提供一小部分信息或仅针对少数物种，或只能针对较近期发生的比较明显的全基因组复制事件。随着测序技术的逐渐成熟，出现基于EST(序列表达标签)数据鉴定全基因组复制事件的方法，但该方法也只是基于能够表达的在基因组上比较保守的基因区域，由于表达数据的获取不全面或者说表达数据只是基因组序列水平上的一个间接产物，该方法鉴定出来的植物基因组内全基因组复制事件的次数和时间往往模糊而不能确定。

发明内容

本发明要解决的一个技术问题是提供一种全基因组复制事件的检测方法，使得确定的全基因组复制事件更准确。

本发明提供一种全基因组复制事件的检测方法，包括：

获取与目标物种基因组上基因集的核酸序列和对应的基因集的蛋白质序列以及基因在基因组上的位置信息；

对基因集的蛋白质序列进行蛋白比蛋白两两比对获得基因集的同源基因对；

根据基因集中基因在基因组上的位置信息以及基因集的同源基因对，利用动态规划获得基因组上具有共线性的同源区块；

获得同源区块中同源基因对的核酸序列的同源位点对位排列；

计算同源区块中同源基因对的遗传距离，获得年龄分布；

根据年龄分布确定全基因组复制事件。

根据本发明方法的一个实施例，根据基因集中基因在基因组上的位置信息以及基因集的同源基因对利用动态规划获得基因组的具有共线性的同源区块的步骤包括：根据基因集的同源基因对，形成所有基因间的同源关系矩阵；利用基因集中每个基因在染色体上的位置信息，利用聚类算法，将相互连锁的共线性的同源基因块聚集在一起，形成具有共线性的同源区块。

根据本发明方法的一个实施例，获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质序列的步骤包括：获取目标物种基因组上基因集的核酸序列；根据密码子表将基因集的核酸序列转化为基因集的蛋白质序列。

根据本发明方法的一个实施例，获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质序列的步骤还包括：当基因集中的基因有多个转录本时，选择最长的转录本对应的蛋白质序列作为基因的核酸序列对应蛋白质序列。

根据本发明方法的一个实施例，将基因集的蛋白质序列进行蛋白比蛋白两两比对获得基因集的同源基因对的步骤包括：通过blastp对基因集的蛋白质序列进行蛋白比蛋白两两比对，设置p-value＜＝1e-7；对比对输出的同源基因对通过e-value＜＝1e-10进行筛选；当同一基因对有多个比对结果仅以不同的e-value值给出时，选择e-value值最小的比对结果作为同源基因对。

根据本发明方法的一个实施例，获得同源区块中同源基因对的核酸序列的同源位点对位排列的步骤包括：对同源区块上的同源基因对的蛋白质序列进行全局序列比对获得同源基因对的蛋白质序列的同源位点对位排列；根据同源基因对的蛋白质序列的同源位点对位排列构建基因集的核酸序列的同源位点对位排列。

根据本发明方法的一个实施例，根据同源基因对的蛋白质序列的同源位点对位排列构建基因集的核酸序列的同源位点对位排列的步骤包括：根据同源基因对的蛋白质序列的同源位点对位排列构建基因集的核酸序列的同源位点对位排列；其中，利用三联密码子将蛋白质“反翻译”(CDS-back-translation)成相应核酸序列后，关注的同源位点中不包含起始密码子和终子密码子。

根据本发明方法的一个实施例，计算同源区块中同源基因对的遗传距离建立年龄分布图的步骤包括：每一个同源比对块上的每对同源基因对，计算遗传距离，对同一个同源区块上的所有基因对的遗传距离计算一个平均值；将遗传距离作为其分化年龄的一种衡量，建立年龄分布图。

通过本发明方法的实施例，全面考虑了全基因组上的所有基因的序列信息，排除了基因组上非保守和选择压力受到释放的区域的干扰，并利用动态规划算法实现最佳复制区域的寻找，确定的同源基因对和同源区块全面、准确，从而能够更准确地实现全基因组复制事件的检测。

本发明要解决的一个技术问题是提供一种全基因组复制事件的检测系统，使得确定的全基因组复制事件更准确。

本发明提供一种全基因组复制事件的检测系统，包括：

蛋白质序列获取模块，用于获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质序列以及基因在基因组上的位置信息；

蛋白质序列比对模块，用于将基因集的蛋白质序列进行蛋白比蛋白两两比对，获得基因集的同源基因对；

同源区块确定模块，用于根据基因集中基因在基因组上的位置信息以及基因集的同源基因对利用动态规划获得基因组的具有共线性的同源区块；

对位排列确定模块，用于获得同源区块中同源基因对的核酸序列的同源位点对位排列；

年龄分布获取模块，用于计算同源区块中同源基因对的遗传距离，获得年龄分布；

全基因组复制确定模块，用于根据年龄分布确定全基因组复制事件。

根据本发明系统的一个实施例，蛋白质序列获取模块包括：核酸序列获取单元，用于获得目标物种基因序列的基因集的核酸序列以及基因在基因组上的位置信息，发送基因集的核酸序列；蛋白质序列转化单元，用于接收基因集的核酸序列，根据密码子表将基因集中各个基因的核酸序列转化为蛋白质序列，获得基因集的蛋白质序列。

根据本发明系统的一个实施例，蛋白质序列比对模块包括：

同源基因对获取单元，用于通过BLASTP蛋白质组比对对基因集的蛋白质序列进行蛋白比蛋白两两比对，获得同源基因对，其中blastp蛋白质组比对的参数设置为p-value＜＝1e-7；

同源基因对过滤单元，用于同源基因对获取单元获得的同源基因对进行筛选，筛选的参数为设置为e-value＜＝1e-10。

根据本发明系统的一个实施例，对位排列确定模块包括：

蛋白质对位排列获取单元，用于对同源区块上的同源基因对的蛋白质序列进行全局序列比对获得同源基因对的蛋白质序列的同源位点对位排列；

核酸对位排列获取单元，用于根据同源基因对的蛋白质序列的同源位点对位排列构建基因集的核酸序列的同源位点对位排列。

通过本发明系统的实施例，蛋白质序列获取模块获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质，全面考虑了全基因组上的所有基因的序列信息，排除了基因组上非保守和选择压力受到释放的区域的干扰，同源区块确定模块利用动态规划算法和有效的聚类算法实现最佳复制区域的寻找，确定的同源基因对和同源区块全面、准确，能够更准确地实现全基因组复制事件的检测。

附图说明

图1示出本发明全基因组复制事件的检测方法的一个实施例的流程图；

图2示出本发明全基因组复制事件的检测方法的另一个实施例的流程图；

图3示出本发明全基因组复制事件的检测系统的一个实施例的框图；

图4示出本发明全基因组复制事件的检测系统的另一个实施例的框图。

图5示出本发明全基因组复制事件的检测方法的一个应用例的流程图；

图6示出根据图5的应用例获得的土豆全基因组复制时间的4DTV图。

具体实施方式

下面参照附图对本发明进行更全面的描述，其中说明本发明的示例性实施例。

图1示出本发明全基因组复制事件的检测方法的一个实施例的流程图。

如图1所示，在步骤102，获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质序列以及基因在基因组上的位置信息。获取目标物种基因组序列的基因注释信息，例如，基因集中各个基因的核酸序列和各个基因在基因组上的位置信息；根据密码子表将基因集中各个基因的核酸序列转化为该基因的蛋白质序列。基因集指物种基因组包括的所有基因的集合。

步骤104，对基因集的蛋白质序列进行蛋白比蛋白两两比对获得基因集的同源基因对。将基因集中每个基因的蛋白质序列和其他每个基因的蛋白质序列进行比对，获得基因集的同源基因对。

步骤106，根据基因集中基因的位置信息以及基因集的同源基因对，利用动态规划获得基因组上具有共线性的同源区块。通过步骤104的基因的同源比对，获得基因集的同源基因对，形成了一个所有基因间的同源关系矩阵。然后，利用每个基因在染色体上的位置信息，利用聚类算法，将相互连锁的共线性的同源基因块聚集在一起，即将单对的同源基因对聚集成多对的共线性的连锁基因块，即具有共线性的同源区块。由于基因组的进化过程中会涉及到大量的染色体重组，缺失等因素，所以要保证所鉴定出来的共线性区块确实是由于全基因组复制所遗留下来的信号，在这里，相邻基因对之间的距离不能太远(一般平均设为10k)，整个共线性基因区块的长度要足够大(例如，至少＞＝100K)来确定一个基因区块是否为同源区块。

步骤108，获得同源区块中同源基因对的核酸序列的同源位点对位排列。对于所有同源区块上的同源基因对，用其蛋白质序列进行全局序列比对(MUSCLE)，获得同源基因对的蛋白质序列的同源位点对位排列方式，然后将蛋白质序列比对反溯至核酸序列水平上去，构建对应的同源基因对的核酸序列的同源位点对位排列方式。

步骤110，计算同源区块中同源基因对的遗传距离，获得年龄分布。例如，对于同源区块上的每对同源基因对，计算其遗传距离(也即分歧度，ks或4dtv值)，对同一个同源区块上的所有同源基因对的遗传距离进行平均，以其平均值作为该同源区块的遗传距离。也可以通过例如加权平均的方法根据每个同源基因对的遗传距离获得同源区块的遗传距离。

步骤112，根据年龄分布确定全基因组复制事件。由于每一次全基因组复制(随后可能会发生基因丢失、染色体重排等事件)都是同源基因对的一次大量涌现，在分布统计上会对应地形成一个显著的峰，据此可对全基因组复制事件发生的次数和时间形成判断。

在上述实施例中，全面考虑了全基因组上的所有基因的序列信息，排除了基因组上非保守和选择压力受到释放的区域的干扰，并利用动态规划算法和有效的聚类算法实现最佳复制区域的寻找，确定的同源基因对和同源区块全面、准确，建立基因序列最有可能的进化模型，从而能够更准确地实现全基因组复制事件的检测。

图2示出本发明全基因组复制事件的检测方法的另一个实施例的流程图。

如图2所示，在步骤202，获取目标物种基因组序列的基因注释信息，包括基因集核酸序列；根据密码子表将基因集的核酸序列转化为基因集的蛋白质序列。当基因集中的基因有多个转录本时，选择最长的转录本对应的蛋白质序列作为该基因的核酸序列对应蛋白质序列。

在步骤204，将基因集的蛋白质序列进行蛋白比蛋白两两比对，对比对结果进行过滤后获得基因集的同源基因对。该过滤包括去除基因的蛋白质序列自身的比对结果，去除同一对基因中多重比对的结果，只保留一个基因对的多个比对结果中的最佳同源基因对。在blastp比对结果中，可以根据e-value值参数进行过滤，有时也需要结合相似度和覆盖度。通过对比对结果进行过滤，使得最终结果可靠性更高。

步骤206，获取每个基因在基因组上的位置信息，若基因在负链上，按5’-＞3’的方向给出位置。以基因为锚点利用动态规划算法找出所有基因组内部具有共线性的同源区块。对步骤204中所有的同源基因对利用基于图论的聚类算法，根据打分矩阵，将局部匹配的比对区域聚集在一起，形成最高分值的大比对区块。在确定同源区块的过程中考虑相邻基因之间的距离，基因缺失和插入所带来的影响，以及最小区块的设置。可以根据植物基因组的进化特点，进行了有效的模型设计。

步骤208，对于所有同源区块上的同源基因对，用其蛋白质序列例如通过MUSCLE进行全局序列比对，获得最佳同源位点对位排列方式。同源位点对位排列方式是指将可能是由同一个祖先传下来的同源区域或同源碱基位点对齐。在MUSLCE算法，可以根据打分获得一个最高分的同源位点对位排列结果，这个最高分值的同源位点对位排列结果作为最佳同源位点对位排列方式。

步骤210，将蛋白质序列比对反溯至核酸序列水平上去，构建对应的核酸序列对位排列方式。根据密码子将蛋白质“反翻译”回核酸序列(CDS-back-translation)，充分利用核酸序列的蛋白序列之间的对应密码子，去除起始密码子和终子密码子。对于残缺基因(有移码或中间有终子密码子的基因)在序列准备阶段给以抛弃。

步骤212，对于每个同源区块上的每对基因对，计算其遗传距离(也即分歧度，ks或4dtv值)，对同一同源区块上的所有基因对的遗传距离计算一个平均值，作为给同源区块的遗传距离。遗传距离的衡量，一个是通过ks(同义突变位点上的平均同义突变率)，它主要是基于分子钟假设，在完全中性的条件下基因的突变速率的衡量；而由于平等突变、回复突变等的影响，ks的衡量标准往往会受到影响而失真，所以需要对其进行矫正，第三位密码子四重简并位点的计算(4dtv，并考虑转换、颠换等参数)会更加接近真实的进化速率。对遗传距离的计算时，可以针对物种基因的特点进行参数训练，建立基因序列最有可能的进化模型，并对距离计算中可能会带来噪音的一些假设进行矫正。

步骤214，将同源区块的遗传距离作为其分化年龄的一种衡量，建立年龄分布图。

步骤216，每一次全基因组复制(伴随着随后的基因丢失，染色体重排等事件)都是同源基因对的一次大量涌现，在分布图上会对应地形成一个显著的峰，据此可对全基因组复制事件发生的次数和时间形成一个直观的认识和判断。

随后，可根据需要，逐步分析每次全基因组复制事件中产生的同源基因的进化趋势或命运，分析它们的功能分化或表达差异，分析它们分析所受到的选择压力和丢失、死亡等事件。

在上述实施例中，在鉴定同源区块时，首先以基因为锚点利用了高效的动态规划算法找出所有可能的最佳同源基因对，并利用基于图论的马尔可夫聚类算法对各对同源基因进行聚类，以得到所有可能的同源区块。处理过程中利用了一系列加权采样值来表示每对基因可能所在位置的后验概率分布，算法每一步都包括基因位置预测和位置更新两个阶段。位置预测阶段是利用m个加权采样值对后验概率分布进行描述的过程，位置更新阶段则是通过重要性采样操作对其进行及时不断更新，最终达到最佳同源区块的鉴定。此外，对遗传距离的计算时，针对物种基因的特点进行参数训练，建立基因序列最有可能的进化模型，并对距离计算中可能会带来噪音的一些假设进行矫正(如充分排除平行突变、回复突变，考虑转换巅换比及考虑基因所受到的选择压力等)，进一步保证了全基因组复制事件检测的准确性。

图3示出本发明全基因组复制事件的检测系统的一个实施例的框图。如图3所示，该实施例的检测系统包括蛋白质序列获取模块31、蛋白质序列比对模块32、同源区块确定模块33、对位排列确定模块34、年龄分布获取模块35和全基因组复制确定模块36。

在图3的实施例中，蛋白质序列获取模块31获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质序列以及基因在基因组上的位置信息。例如，蛋白质序列获取模块31获得目标物种基因序列的注释信息，其中包括基因集的核酸序列，然后根据密码子表将基因集中各个基因的核酸序列转化为该基因的蛋白质序列。当基因集中的基因有多个转录本时，选择最长的转录本对应的蛋白质序列作为基因的核酸序列对应蛋白质序列。蛋白质序列比对模块32将基因集的蛋白质序列进行蛋白比蛋白两两比对，获得基因集的同源基因对。蛋白质序列比对模块32将基因集中每个基因的蛋白质序列和分别和其他基因的蛋白质序列进行比对，从而获得基因集的同源基因对。同源区块确定模块33根据基因集基因在基因组上的位置信息以及基因集的同源基因对利用动态规划获得基因组的具有共线性的同源区块。例如，同源区块确定模块33根据蛋白质序列比对模块32获得的基因集的同源基因对，形成一个所有基因间的同源关系矩阵；然后利用每个基因在染色体上的位置信息，利用聚类算法，将相互连锁的共线性的同源基因块聚集在一起，即将单对的同源基因对聚集成多对的共线性的连锁基因块，即具有共线性的同源区块。对位排列确定模块34获得同源区块中同源基因对的核酸序列的同源位点对位排列。对位排列确定模块34对所有同源区块上的同源基因对，利用其蛋白质序列进行全局序列比对，获得同源基因对的蛋白质序列的同源位点对位排列方式，然后将蛋白质序列比对反溯至核酸序列水平上去，构建对应的彤云基因对的核酸序列的同源位点对位排列方式。年龄分布获取模块35计算同源区块中同源基因对的遗传距离，获得年龄分布。例如，年龄分布获取模块35对于同源区块上的每对同源基因对，计算其遗传距离，对同一个同源区块上的所有同源基因对的遗传距离进行平均，以其平均值作为该同源区块的遗传距离。全基因组复制确定模块36根据年龄分布确定全基因组复制事件。

在上述实施例中，蛋白质序列获取模块获取与目标物种基因组上基因集的核酸序列对应的基因集的蛋白质，全面考虑了全基因组上的所有基因的序列信息，排除了基因组上非保守和选择压力受到释放的区域的干扰，同源区块确定模块利用动态规划算法和有效的聚类算法实现最佳复制区域的寻找，确定的同源基因对和同源区块全面、准确，建立基因序列最有可能的进化模型，从而能够更准确地实现全基因组复制事件的检测。

图4示出本发明全基因组复制事件的检测系统的另一个实施例的框图。如图4所示，该实施例的检测系统包括蛋白质序列获取模块41、蛋白质序列比对模块42、同源区块确定模块33、对位排列确定模块44、年龄分布获取模块35和全基因组复制确定模块36。同源区块确定模块33、年龄分布获取模块35和全基因组复制确定模块36的功能和实现可以参见图3中对应模块的描述，为简洁起见在此不再详细描述。下面具体描述蛋白质序列获取模块41、蛋白质序列比对模块42、和对位排列确定模块44的实现。

在图4的实施例中，蛋白质序列获取模块41包括核酸序列获取单元411和蛋白质序列转化单元412。核酸序列获取单元411获得目标物种基因序列的基因集的核酸序列以及基因在基因组上的位置信息，将基因集的核酸序列发送给蛋白质序列转化单元412。蛋白质序列转化单元412接收基因集的核酸序列，根据密码子表将基因集中各个基因的核酸序列转化为蛋白质序列，从而获得基因集的蛋白质序列。蛋白质序列比对模块42包括同源基因对获取单元421和同源基因对过滤单元422。同源基因对获取单元421用于通过BLASTP蛋白质组比对对基因集的蛋白质序列进行蛋白比蛋白两两比对，获得同源基因对，其中blastp蛋白质组比对的参数设置为p-value＜＝1e-7。同源基因对过滤单元422对同源基因对获取单元421获得的同源基因对进行筛选，筛选的参数为设置为e-value＜＝1e-10。对位排列确定模块44包括蛋白质对位排列获取单元441和核酸对位排列获取单元442。蛋白质对位排列获取单元441对同源区块上的同源基因对的蛋白质序列进行全局序列比对获得同源基因对的蛋白质序列的同源位点对位排列。核酸对位排列获取单元442根据蛋白质对位排列获取单元411输出的同源基因对的蛋白质序列的同源位点对位排列来构建基因集的核酸序列的同源位点对位排列。

对于图3至图4中各个装置或单元的功能，可以参考上文中关于本发明方法的实施例中对应部分的说明，为简洁起见，在此不再详述。

本领域的技术人员应当理解，对于图3至图4中的各个装置，可以通过单独的计算处理设备实现，或者将其集成为一个独立的设备实现。在图3至图4中用框示出以说明它们的功能。这些功能块可以用硬件、软件、固件、中间件、微代码、硬件描述语音或者它们的任意组合来实现。举例来说，一个或者两个功能块都可以利用运行在微处理器、数字信号处理器(DSP)或任何其他适当计算设备上的代码实现。代码可以表示过程、功能、子程序、程序、例行程序、子例行程序、模块或者指令、数据结构或程序语句的任意组合。代码可以位于计算机可读介质中。计算机可读介质可以包括一个或者多个存储设备，例如，包括RAM存储器、闪存存储器、ROM存储器、EPROM存储器、EEPROM存储器、寄存器、硬盘、移动硬盘、CD-ROM或本领域公知的其他任何形式的存储介质。计算机可读介质还可以包括编码数据信号的载波。

本领域技术人员将意识到硬件、固件和软件配置在这些情况下的可替换性，以及如何最好地实现每个特定应用地所述功能。

图5示出本发明全基因组复制事件的检测方法的一个应用例的流程图。该应用例通过本发明的检测方法鉴定土豆全基因组复制事件。土豆是茄科植物中非常重要的一个物种，其基因组大小约760M，功能基因有约4320个。每一次全基因组复制事件都会是同源基因数目的大爆发，其后伴随的大量基因的丢失和染色体重排虽然会抹去许多全基因组复制的痕迹，但功能基因在整体上会保留这些曾经发生过大爆发的证据。当把所有的今天依然存留下来的同源区块都鉴定出来后，计算每对同源区块之间分化的可能时间(遗传距离)，根据遗传距离累积所反映出来的“年龄”分布，就可判断时全基因复制事件是否发生和可能发生的时间。具体步骤如下：

步骤502，获得土豆基因集的核酸序列和对应蛋白质序列。土豆基因集包括4320个基因。如果一个基因有多个转录本，则选择其最长的转录本对应的蛋白质序列。

步骤504，blastp蛋白质组比对，比对时一般可设置p-value＜＝1e-7；比对输出结果中可以e-value＜＝1e-10来筛选，进一步保证比对结果的可信度。如果包括自身比对，则要过滤自身比对结果；若同一对基因对有多个比对结果只是以不同的e-value值给出，选择其最优比对结果(即e-value值最小)。土豆的蛋白质组比对结果过滤后有1,716,931条纪录。

步骤506，获取每个基因在染色体上的位置信息，若基因在负链上，按5’-＞3’的方向给出位置。一般可从基因注释结果中提取基因在染色体上的位置信息。

步骤508，以基因为锚点，利用以上基因对比上的结果和位置信息，实现动态规划算法获取最优局部匹配。对步骤504得到的比对结果利用基于图论的聚类算法，根据打分矩阵，将局部匹配的比对区域聚集在一起，形成两边延伸的大比对区块，即同源区块。将重点考虑相邻基因之间的距离，基因缺失和插入所带来的影响，以及最小区块的设置。土豆分析中相邻基因之间的最小距离是10k，每个区块中至少要有5个基因对。结果有284个大区块被鉴定出来，平均每个区块有9个基因，共2556个基因对涉及到3200个土豆基因。

步骤510，对2556个同源基因对的蛋白质序列进行全局比对，形成蛋白质的同源位点对位排列，随之转换为基于核酸序列的全局同源位点对位排列方式。蛋白质比对时选用MUSCLE全局比对软件包处理；CDS-back-translation中，要充分利用核酸序列的蛋白序列之间的对应密码子，并去除起始密码子和终子密码子。对于残缺基因(有移码或中间有终子密码子的基因)在序列准备阶段给以抛弃。

步骤512，基于核酸序列的比对结果，计算每个区块中每对同源基因的遗传距离。遗传距离的衡量，一个是通过ks(同义突变位点上的平均同义突变率)，它主要是基于分子钟假设，在完全中性的条件下基因的突变速率的衡量；而由于平等突变、回复突变等的影响，ks的衡量标准往往会受到影响而失真，所以需要对其进行矫正，第三位密码子四重简并位点的计算(4dtv，并考虑转换、颠换等参数)会更加接近真实的进化速率。

步骤514，对每个区块中的多对遗传距离值计算一个平均值，作为该对区块分化时间的一个衡量。对所有区块分化时间作一个累积分布图。

步骤516，每一次全基因组复制事件都是重复区块的一个大涌现，根据分布图上的峰值可以判断全基因组复制事件发生的时间和发生的次数。

通过图5所示的应用例，土豆全基因组复制事件至少有两次被鉴定出来(如图6所示)。一次是在4dtv＝～0.33处，这是最近的一次全基因组复制事件；一次是在4dtv＝～0.75处，这是古老的一次全基因组复制事件。从图6中也可以看出，土豆古老的那次全基因组复制事件发生在与拟南芥分化之前，可能是与拟南芥甚至十字花科共有的一次全基因组复制事件；而土豆最近的那次全基因组复制事件发生在与拟南芥、葡萄分化之后，可能是土豆或说茄科所特有的。

上述应用例中，通过本发明的方法鉴定出来的两次全基因组复制事件非常清楚，比现有技术中只根据染色体数目进行的估计有足够的数据说服力，对土豆/茄科中确实有一次特有的全基因组复制事件给予了充分的回答。过去几年曾有不少研究根据EST数据对土豆、番茄等茄科植物的全基因组复制事件进行过推断，但都由于数据不全或方法不够有效而归于失败，或得到模凌两可的处于猜测阶段的结果。本发明实施例的方法得力于土全基因组数据及其注释基因的获得，更得力于高效的鉴定方法。本发明的方法可应用于其它所有的被子植物，已在木瓜、大豆、拟南芥、黄瓜等物种中做过测试，都取得显著的效果。

本发明的描述是为了示例和描述起见而给出的，而并不是无遗漏的或者将本发明限于所公开的形式。很多修改和变化对于本领域的普通技术人员而言是显然的。选择和描述实施例是为了更好说明本发明的原理和实际应用，并且使本领域的普通技术人员能够理解本发明从而设计适于特定用途的带有各种修改的各种实施例。

Claims (14)

1.一种全基因组复制事件的检测方法，其特征在于，包括：

获取与目标物种基因组上基因集的核酸序列对应的所述基因集的蛋白质序列以及基因在基因组上的位置信息；

对所述基因集的蛋白质序列进行蛋白比蛋白两两比对获得所述基因集的同源基因对；

根据所述基因集中基因在基因组上的位置信息以及所述基因集的同源基因对，利用动态规划算法获得所述基因组上具有共线性的同源区块；

获得所述同源区块中同源基因对的核酸序列的同源位点对位排列；

计算所述同源区块中同源基因对的遗传距离，获得年龄分布；

根据所述年龄分布确定全基因组复制事件。

2.根据权利要求1所述的检测方法，其特征在于，所述根据所述基因集中基因在基因组上的位置信息以及所述基因集的同源基因对利用动态规划获得所述基因组的具有共线性的同源区块的步骤包括：

根据所述基因集的同源基因对，形成所有基因间的同源关系矩阵；

利用基因集中每个基因在染色体上的位置信息，利用聚类算法，将相互连锁的共线性的同源基因块聚集在一起，形成具有共线性的同源区块。

3.根据权利要求1所述的检测方法，其特征在于，所述获取与目标物种基因组上基因集的核酸序列对应的所述基因集的蛋白质序列的步骤包括：

获取所述目标物种基因组上基因集的核酸序列；

根据密码子表将所述基因集的核酸序列转化为所述基因集的蛋白质序列。

4.根据权利要求3所述的检测方法，其特征在于，所述获取与目标物种基因组上基因集的核酸序列对应的所述基因集的蛋白质序列的步骤还包括：

当所述基因集中的基因有多个转录本时，选择最长的转录本对应的蛋白质序列作为所述基因的核酸序列对应蛋白质序列。

5.根据权利要求1所述的检测方法，其特征在于，所述将所述基因集的蛋白质序列进行蛋白比蛋白两两比对获得所述基因集的同源基因对的步骤包括：

通过blastp对所述基因集的蛋白质序列进行蛋白比蛋白两两比对，设置p-value＜＝1e-7；

对比对输出的同源基因对通过e-value＜＝1e-10进行筛选；

当同一基因对有多个比对结果仅以不同的e-value值给出时，选择e-value值最小的比对结果作为同源基因对。

6.根据权利要求1所述的检测方法，其特征在于，所述获得所述同源区块中同源基因对的核酸序列的同源位点对位排列的步骤包括：

对所述同源区块上的同源基因对的蛋白质序列进行全局序列比对获得所述同源基因对的蛋白质序列的同源位点对位排列；

根据所述同源基因对的蛋白质序列的同源位点对位排列构建所述基因集的核酸序列的同源位点对位排列。

7.根据权利要求6所述的检测方法，其特征在于，所述根据所述同源基因对的蛋白质序列的同源位点对位排列构建所述基因集的核酸序列的同源位点对位排列的步骤包括：

利用核酸序列的蛋白序列之间的对应密码子根据所述同源基因对的蛋白质序列的同源位点对位排列构建所述基因集的核酸序列的同源位点对位排列。

8.根据权利要求1所述的检测方法，其特征在于，所述计算所述同源区块中同源基因对的遗传距离建立年龄分布图的步骤包括：

每一个同源比对块上的每对同源基因对，计算遗传距离，对同一个同源区块上的所有基因对的遗传距离计算一个平均值；

将遗传距离作为其分化年龄的一种衡量，建立年龄分布图。

9.一种全基因组复制的检测系统，其特征在于，包括：

蛋白质序列获取模块，用于获取与目标物种基因组上基因集的核酸序列对应的所述基因集的蛋白质序列以及基因在基因组上的位置信息；

蛋白质序列比对模块，用于将所述基因集的蛋白质序列进行蛋白比蛋白两两比对，获得所述基因集的同源基因对；

同源区块确定模块，用于根据所述基因在基因组上的位置信息以及所述基因集的同源基因对利用动态规划获得所述基因组的具有共线性的同源区块；

对位排列确定模块，用于获得所述同源区块中同源基因对的核酸序列的同源位点对位排列；

年龄分布获取模块，用于计算所述同源区块中同源基因对的遗传距离，获得年龄分布；

全基因组复制确定模块，用于根据所述年龄分布确定全基因组复制事件。

10.根据权利要求9所述的全基因组复制的检测系统，其特征在于，所述蛋白质序列获取模块包括：

核酸序列获取单元，用于获得目标物种基因序列的基因集的核酸序列以及基因在基因组上的位置信息，发送所述基因集的核酸序列；

蛋白质序列转化单元，用于接收所述基因集的核酸序列，根据密码子表将所述基因集中各个基因的核酸序列转化为蛋白质序列，获得所述基因集的蛋白质序列。

11.根据权利要求9所述的全基因组复制的检测系统，其特征在于，所述蛋白质序列比对模块包括：

同源基因对获取单元，用于通过BLASTP蛋白质组比对对所述基因集的蛋白质序列进行蛋白比蛋白两两比对，获得同源基因对，其中所述blastp蛋白质组比对的参数设置为p-value＜＝1e-7；

同源基因对过滤单元，用于同源基因对获取单元获得的同源基因对进行筛选，筛选的参数为设置为e-value＜＝1e-10。

12.根据权利要求9所述的全基因组复制的检测系统，其特征在于，所述对位排列确定模块包括：

蛋白质对位排列获取单元，用于对所述同源区块上的同源基因对的蛋白质序列进行全局序列比对获得所述同源基因对的蛋白质序列的同源位点对位排列；

核酸对位排列获取单元，用于根据所述同源基因对的蛋白质序列的同源位点对位排列构建所述基因集的核酸序列的同源位点对位排列。

13.根据权利要求9所述的全基因组复制的检测系统，其特征在于，所述同源区块确定模块根据所述基因集的同源基因对，形成所有基因间的同源关系矩阵；利用基因集中每个基因在染色体上的位置信息，利用聚类算法，将相互连锁的共线性的同源基因块聚集在一起，形成具有共线性的同源区块。

14.根据权利要求9所述的全基因组复制的检测系统，其特征在于，所述年龄分布获取模块对于同源区块上的每对同源基因对，计算其遗传距离，对同一个同源区块上的所有同源基因对的遗传距离进行平均，以其平均值作为所述同源区块的遗传距离。

Images