CN101908102B - 基于茎区的核糖核酸二级结构的预测方法及装置 - Google Patents

Abstract

本发明公开了基于茎区的RNA二级结构的预测方法及装置。该预测方法包括：定义碱基对数据、茎区数据和祖先茎区数据；输入RNA序列；在所输入的RNA序列中进行查找，确定各种碱基对的第一分量在所输入的RNA序列所处的第一位置，并将第一位置存储在碱基对数据中；在所输入的RNA序列中进行查找，确定各种碱基对的第二分量在所输入的RNA序列所处的第二位置，并将第二位置存储在碱基对数据中；搜索在碱基对数据中的各元素的第一分量和第二分量之间的RNA序列中的所有茎区，将与所搜索到的与各茎区相对应的有关信息存储在茎区数据中，将与所搜索到的各祖先茎区有关的信息存储在祖先茎区数据中；根据茎区数据和祖先茎区数据确定最大茎区。本发明能降低时间与空间复杂度，提高查全率和查准率。

Description

基于茎区的核糖核酸二级结构的预测方法及装置

技术领域

本发明属于生物信息工程领域，涉及一种对核糖核酸(在下文中，简称为RNA)的二级结构进行预测的方法，尤其涉及基于茎区包含假结的RNA二级结构进行预测的方法及装置。

背景技术

RNA二级结构预测是RNA识别及其功能研究的根本途径与核心基础，是从海量的生物学数据中提取有用的知识，揭示其蕴含的生物学意义的手段和方法，为人类基因组计划的实施作出贡献。RNA二级结构预测问题是计算分子生物学中的一个经典问题，自2000年以来，有关RNA的研究已经多年连续被世界顶级刊物《科学》(《Science》)杂志列入世界十大科技进展。

RNA是合成蛋白质的模板，RNA二级结构预测用于蛋白质功能分析，是RNA三级结构预测的基础。RNA与DNA结构预测问题是计算机科学、数学和生物学交叉学科的基本课题之一。假结(pseudoknot)是RNA中最广泛的结构单元，是非常复杂和稳定的RNA结构，假结在RNA分子中具有构造、催化和调节功能，是目前RNA结构预测研究的关键点。

目前实验方法测定RNA二级结构成本很高且耗时太多，通常把待测结构RNA的一级序列作为输入，仅根据输入的序列预测其二级结构。这样做丢失了待测RNA的类别信息，进而无法利用同类别RNA二级结构的保守性。因此用计算机预测RNA二级结构成为一种重要的方法。由于茎区结构相对稳定，包含假结的RNA二级结构预测问题是NPC问题。现在急需一种可行的计算机高效方法来预测RNA二级结构，降低时间和空间复杂性，特别是寻找一种基于茎区包含假结的RNA二级结构预测方法来提高预测的准确性。

RNA一级结构如图7所示，RNA二级结构如图8所示，碱基之间可形成氢键配对，相邻碱基对构成茎区。

当前，国内外许多计算机科学家和生物学家提出了用于预测RNA二级结构的方法，例如，Nussinov方法、Mfold方法、Rivas方法、LyngsΦ方法等。然而，Nussinov方法和Mfold方法均不能预测假结结构。Nussinov方法忽略环的自由能量，只考虑临近碱基的相互作用，时间复杂性为O(n³)。Mfold方法只对不包括邻接和交叉内边的半环结构进行预测，其时间复杂性也为O(n³)。Rivas方法提出了对包含假结的RNA二级结构进行预测的动态规划算法，但其时间复杂性为O(n⁶)，空间复杂度为O(n⁴)。再如，LyngsΦ方法虽然对Rivas方法进行了一些改进，但当预测任意的平面假结或受限的非平面假结时，时间复杂性和空间复杂度分别为O(n⁵)和O(n⁴)，但只能预测一个平面假节点。在当前预测任意平面假结和部分非平面假结的方法中，较好的为Pknots方法，需O(n⁶)时间和O(n⁴)空间。由此可见，这些预测方法存在准确性不高、时间复杂性和空间复杂性偏高、仅能预测符合严格限制的RNA二级结构等问题。

由于对RNA二级结构预测方法的空间复杂性和时间复杂性直接影响预测成本，因此，如何使RNA二级结构、尤其是包含假结的RNA的二级结构的预测方法的时间复杂性和空间复杂性尽可能小，且确保预测结果的准确性，已经成为生物信息工程领域的一项重要的研究课题。

发明内容

本发明所要解决的主要技术问题是使得对RNA二级结构、尤其是使基于茎区包含假结的RNA二级结构进行预测方法的时间复杂性和空间复杂性尽可能小，提高预测准确性。

根据本发明的一方面，提供一种基于茎区的核糖核酸二级结构的预测方法，包括：初始化步骤，定义碱基对数据、茎区数据和祖先茎区数据；输入步骤，输入核糖核酸序列；第一分量查找步骤，在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第一分量在所输入的核糖核酸序列所处的第一位置，并将所确定的所述第一位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第一分量中；第二分量查找步骤，在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第二分量在所输入的核糖核酸序列所处的第二位置，并将所确定的所述第二位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第二分量中；茎区搜索步骤，搜索在所述碱基对数据的各元素的第一分量和第二分量之间的核糖核酸序列中的所有茎区，将与所搜索到的与各茎区相对应的有关信息分别存储在所述茎区数据的各元素中，将与所搜索到的各祖先茎区有关的信息分别存储在所述祖先茎区数据的各元素中，其中所述所有茎区包含所搜索到的茎区的子茎区；最大茎区确定步骤，根据所述茎区数据和所述祖先茎区数据来确定在所输入的核糖核酸序列中的最大茎区。

根据本发明一方面的预测方法，所述祖先茎区数据包括与各祖先茎区一一对应的各元素，所述祖先茎区数据的各元素用于存储相应祖先茎区所在家族的家族能量值和相应祖先茎区的位置信息。

根据本发明一方面的预测方法，所述祖先茎区的位置信息是所述茎区数据中存储了与相应祖先茎区有关的信息的元素的位置或指向该位置的指针。

根据本发明的另一方面，提供一种基于茎区的核糖核酸二级结构的预测装置，该预测装置包括：初始化单元，其定义碱基对数据、茎区数据和祖先茎区数据；输入单元，其输入核糖核酸序列；第一分量查找单元，其在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第一分量在所输入的核糖核酸序列所处的第一位置，并将所确定的所述第一位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第一分量中；第二分量查找单元，其在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第二分量在所输入的核糖核酸序列所处的第二位置，并将所确定的所述第二位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第二分量中；茎区搜索单元，其搜索在所述碱基对数据的各元素的第一分量和第二分量之间的核糖核酸序列中的所有茎区，将与所搜索到的与各茎区相对应的有关信息分别存储在所述茎区数据的各元素中，将与所搜索到的各祖先茎区有关的信息分别存储在所述祖先茎区数据的各元素中，其中所述所有茎区包含所搜索到的茎区的子茎区；最大茎区确定单元，其根据所述茎区数据和所述祖先茎区数据来确定在所输入的核糖核酸序列中的最大茎区。

根据本发明的另一方面的预测装置中的祖先茎区数据包括与各祖先茎区一一对应的各元素，所述祖先茎区数据的各元素用于存储相应祖先茎区所在家族的家族能量值和所述相应祖先茎区的位置信息。

根据本发明一个技术方案的结构建模能够降低空间复杂度。

根据本发明的一个技术方案利用了基于茎区结构相对稳定的特征，采用了贪婪(Greedy)算法与动态规划算法相结合的方法，基于茎区结构来计算嵌套结构和包含交叉假接的非嵌套结构，能够预测RNA的嵌套结构、非嵌套结构和假结结构，以及能够降低时间复杂度。

根据本发明一个技术方案采用最小自由能量方法预测R NA二级结构，能够提高基于茎区包含假结的RNA二级结构进行预测的预测准确性。

本发明的一个技术方案采用了将已查找到的最大的茎区标记而非将其从RNA茎区中删除的方法，能够提高对包含假结RNA二级结构进行预测的预测准确性。

通过下面参考附图对示例性实施例的说明，本发明的其他特征将变得清楚。

附图说明

图1是示出根据本发明用于预测RNA的二级结构的预测方法的例子的流程图。

图2是示出根据本发明用于预测RNA的二级结构的预测方法的又一例子的流程图。

图3A是示出根据本发明的茎区搜索处理的例子的流程图。

图3B是示出根据本发明搜索处理中的茎区确定处理的例子的流程图。

图4示出与图1中用于预测RNA二级结构的预测方法对应的测量装置。

图5示出与图2中用于预测RNA二级结构的预测方法对应的测量装置

图6示出根据本发明第二实施例与Pknots技术方案的查全率与查准率对照表。

图7示出RNA一级结构的例子。

图8示出RNA二级结构的例子。

图9示出RNA二级结构中包含假结结构的例子。

此外，在图9的假结中，()和[]表示基对。

具体实施方式

本部分内容将根据附图详细说明本发明的各实施例。

在结合附图说明本发明的各实施例之前，说明几个概念。

RNA序列：将RNA序列表示为s₁s₂...s_n，其中s_i∈{A，C，G，U}，序列的长度为n，1≤i≤n。

碱基对：如果s_i·s_j∈{AU，CG，GU}，则s_i·s_j构成碱基对。典型地，存在AU、CG碱基对和GU错配碱基对，其中，s_i·s_j表示RNA序列S中两个相邻碱基。

RNA二级结构S：RNA序列中的一组基对构成的集合构成RNA二级结构，以S表示。对于任意基对，如果s_i·s_j∈S、s_i′·s_j′∈S且i＝i′，则j＝j′，亦即，一个基不可同时与两个及两个以上的基构成基对。

茎区：将S中两个以上相邻基对s_i·s_j、s_(i+1)·s_(j-1)、s_(i+2)·s_(j-2)、s_(i+3)·s_(j-3)......∈S所构成封闭的结构，称为茎区。

假结：如果基对s_i·s_j与s_i′·s_j′∈S，如果i＜i′＜j＜j′，则序列s_i...s_i′...s_j...s_j′构成假结结构。

NPC问题：若判定问题P是非多项式的，若对每个其他判定问题P′是非多项式的，都有P′∝P，则称P是NP完全的(NPC)。

MFE原理：将最小折叠自由能量(Minimal Folding Free Energy)原理简称为MFE原理。

茎区家族：如果两个或两个以上茎区为包含或被包含关系，则称该两个或两个以上茎区属于同一茎区家族。

祖先茎区：将属于茎区家族中包含最多碱基对的茎区称为祖先茎区。

茎区的累加能量值，为茎区中的所有碱基对的能量值的和。家族能量值，为茎区家族中的所有茎区的能量和。

第一实施例

下面参考图1说明本发明的第一实施例。图1是示出根据本发明的用于预测RNA的二级结构的预测方法的例子的流程图。

参考图1，首先在步骤S100中，定义用于存储与碱基对的种类有关的数据结构，作为碱基对种类数据。例如，可定义包含36个元素，每个元素由三个值组成的矩阵M，作为碱基对种类数据，其中，M＝{m₀，m₁，m₂，......m₃₅}，以存储如下所示的各种(共36种)碱基对的信息：

{″UU″，″AA″，-0.9}，{″AC″，″GU″，-2.1}，{″UC″，″GA″，-2.3}，

{″AG″，″CU″，-1.7}，{″UA″，″UA″，-1.1}，{″CA″，″UG″，-1.8}，

{″CU″，″AG″，-1.7}，{″GA″，″UC″，-2.3}，{″GU″，″AC″，-2.1}，

{″CG″，″CG″，-2.0}，{″GC″，″GC″，-3.4}，{″GG″，″CC″，-2.9}，

{″AA″，″UU″，-0.9}，{″CC″，″GG″，-2.9}，{″AU″，″AU″，-0.9}，

{″GA″，″UU″，-0.5}，{″UG″，″CA″，-1.8}，{″GG″，″UU″，-0.5}，

{″GU″，″GU″，-0.5}，{″UG″，″UG″，-0.6}，{″UU″，″GG″，-0.5}，

{″AG″，″UU″，-0.9}，{″AU″，″GU″，-0.9}，{″CG″，″UG″，-1.6}，

{″CU″，″GG″，-1.9}，{″GG″，″UC″，-1.4}，{″GU″，″GC″，-2.3}，

{″UG″，″UA″，-0.9}，{″UU″，″GA″，-1.0}，{″UU″，″AG″，-0.5}，

{″GU″，″AU″，-0.7}，{″UG″，″CG″，-1.5}，{″GG″，″CU″，-1.5}，

{″UC″，″GG″，-1.3}，{″GC″，″GU″，-1.9}，{″UA″，″UG″，-0.7}。

在矩阵M中，每个元素m_i(0＜＝i＜＝35)包括三个值，即，m_i＝{m_(i，0)，m_(i，1)，m_(i，2)}，可将第一个和第二个值设置为二个字节的字符，分别用于存储碱基对的第一分量和第二分量，将第三个值设置为实数，用于存储对应碱基对的能量值。在利用程序来实现的情况下，例如，在利用C语言来实现的情况下，可以利用char[3]等来定义前两个分量，采用double等来定义第三个分量。

并且，还对矩阵M进行初始化。例如，将矩阵M第一个元素m₀(与碱基对{″UU″，″AA″，-0.9}相对应的元素)的第一个元素的第一分量设为“UU”、第二分量设为“AA”，第三分量的值设为-0.9，将M的第二个元素m₁(与碱基对{″AC″，″GU″，-2.1}相对应的元素)的第一分量m_(1，0)设为“AC”、第二分量m_(1，1)设为“GU”，第三分量m_(1，2)的值设为-2.1，依此类推，依次对M中的各个元素赋值，从而完成对M的初始化。

需要注意的是，上述通过定义矩阵M来存储与36种碱基对的属性值的方式仅是例子，本发明不限于上述的数据结构定义方式，例如，还可以通过定义一个长度为36×2×2的字符串CH来依次存储前述36个碱基对的第一和第二个分量，并且定义一个一维的实数数组RR来依次存储前述36个碱基对的能量值。此时，CH[0]＝‘U’，CH[1]＝‘U’，CH[2]＝‘A’，CH[3]＝‘A’，CH[4]＝‘A’，CH[5]＝‘C’，CH[6]＝G，CH[7]＝‘U’，依次类推，存储其它第三种至第36种碱基的第一及第二分量；RR[0]＝-0.9，RR[1]＝-2.1，依次类推，存储其它第3种至第36种碱基对的第三分量。

此外，还定义用于存储与所查找到的各种碱基对所在位置及能量有关的信息的碱基对数据。该碱基对数据中的各元素具有两个分量，将各种碱基对的第一分量在所输入的核糖核酸序列所处的位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第一分量中，且将各种碱基对的第二分量在所输入的核糖核酸序列所处的位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第二分量中。例如，定义另一矩阵N，作为碱基对数据，N＝{n₀，n₁，n₂，......n₃₅}，矩阵N包含36个元素，每个元素包含三个值，即，n_i＝{n_(i，0)，n_(i，1)，n_(i，2)}，(0＜＝i＜＝35)。正如下文中所述，矩阵N中的36个元素将用于存储与上述36种碱基对有关的值。更具体地，36个元素中的每个元素的第一个值和第二个值将用于依次存储上述36种碱基对的第一分量和第二分量在所输入的碱基序列中所处的位置，第三个值用于依次存储上述36种碱基对对应的能量值。并且，将矩阵N中各元素的第三个值n_i2初始化为与36种碱基对相对应的能量值。注意，与矩阵M类似，本领域技术人员可以明白还可以定义其它数据结构形式来实现与上述N等同的功能。例如，可以定义三个一维数组，分别用于存储与上述矩阵N相对应的上述36种元素的36个第一个值、36个第二个值和36个第三个值，在此不再展开赘述。此外，在矩阵M中已经存储了上述36种碱基对对应的能量值的情况下，矩阵N中的各个元素也可以只包含上述第一个值和第二个值，而不包含用于依次存储上述36种碱基对对应的能量值的第三个值。

再定义茎区数据，茎区数据中的各元素分别用于存储与各茎区有关的信息，其中与各茎区有关的信息可以包括但不限于各茎区的起始位置、结束位置、所包含的碱基对数量、茎区的累加能量值。例如，可定义一个链表X，链表中的各个节点x包含五个变量，分别表示为x.next、x.start、x.end、x.energy和x.number，五个变量分别用于存储指向子节点的指针、所预测出的各个茎区的起始位置、所预测出的各个茎区的结束位置、所预测出的各个茎区中的碱基对的累加能量值以及所预测出的各个茎区中的碱基对的数量。然而，虽然优先为定义链接，但也可以定义一数组，该数组长度小于等于要进行二级结构确定的RNA序列的长度的四分之一，且利用数组中的各个元素来存储各个茎区的起始位置、结束位置及其累加能量值。也就是说，用于存储上述这些值的数据结构的形式多种多样，只要能用于存储这些值即可。

还定义祖先茎区数据，祖先茎区数据中的各元素分别存储与各祖先茎区有关的信息，其中与各祖先茎区有关的信息可以是在该祖先茎区所在茎区家族的家族能量值以及祖先茎区的位置信息。其中，祖先茎区的位置信息可以是茎区数据中存储了与相应祖先茎区有关的信息的元素的位置或指向该位置的指针等信息。例如，可定义用于存储与祖先茎区的信息链表E，作为祖先茎区数据，E中的节点e包含三个变量e.xpointer、e.energy、e.next。其中，e.xpointer是指针，用于指向茎区数据中存储了与相应祖先茎区有关的信息的节点；e.energy表示该祖先茎区所在的茎区家族的家族能量值；e.next指向链接的下一节点。

在本实施例中，还定义搜索初始位置F_b和搜索结束位置F_e。用于在搜索时对搜索范围进行界定。

然后，进入RNA序列输入和变量初始化步骤S101。

在步骤S101中，输入RNA序列。具体地，在从S100进入步骤S101时，所输入的RNA序列是要进行二级结构预测的整个RNA序列；与之对比，在从步骤S109返回到步骤S101时，所输入的RNA序列是经步骤S108处理(后面将说明)后得到的RNA序列。

在步骤S101中，还对所定义的数据进行初始化。例如，将矩阵N中的各个元素的第一分量n_i0及第二个分量n_i1均初始化为零，以及将F_b初始化为所输入的碱基序列的起始位置，将F_e初始化为所输入的碱基序列的长度减1。需要说明的是，在从S100进入步骤S101时，则相当于将F_b初始化为零(即要进行二级结构确定的RNA序列的初始位置)，以及将F_e初始化为要进行二级结构确定的RNA序列的长度减1；与之对比，在从步骤S109返回到步骤S101时，则相当于将F_b初始化为经步骤S108处理后所获得的碱基序列(即所输入的碱基序列)的起始位置，将F_e初始化为经步骤S108处理后所获得的碱基序列(即所输入的碱基序列)的长度减1。

然后进入第一分量查找步骤S102，在步骤S102中，在输入的RNA序列中进行查找，以确定矩阵M中各碱基对的第一分量在所输入的碱基序列中所处的位置，并将所确定的位置分别存储在N中与各碱基对相对应的元素的第一分量中。例如，如果在所输入的碱基序列中的第5个位置处查找出与上述36种碱基对中的第二种碱基对相对应的第一分量“AC”，亦即，所输入的碱基序列中的第5个和第6个元素为“A”和“C”，则使得N的第2个元素的第1个分量n_(1，0)设置为4(本实施例中，将第1个碱基的地址设为0)。

然后进入第二分量查找步骤S103，在步骤S103中，在输入的碱基序列中进行查找，以确定矩阵M中各碱基对的第二分量在所输入的碱基序列中所处的位置，并将所确定的位置分别存储在矩阵N中与各碱基对相对应的元素的第二分量中。例如，如果在所输入的碱基序列中的第100个位置处查找出上述36种碱基对中的第二种碱基对相对应的第二分量“GU”，亦即，所输入的碱基序列中的第100个和第101个元素为“G”和“U”，则将N的第2个元素的第2个分量n_(1，1)赋值为99(将第1个碱基的地址设为0)。此外，还将n_(1，2)赋值为与M中所存在的与碱基对{“AC”，“GU”}对应的能量值-2.1，即，n_(1，2)＝-2.1。

注意，在步骤S102及S103中，如果在所输入的碱基序列中未找到矩阵M中的部分碱基对的第一分量或第二分量，则使得矩阵N中与未找到的碱基对相对应的元素第一分量或第二分量的值以及相应的元素的能量值保持不变。例如，假设未找到与m₅相对应的碱基的第一分量和第二分量，则保持相应元素n₅中的n_(5，0)、n_(5，1)、n_(5，2)不变。

后面将结合图3A和图3B对茎区搜索处理进行说明。

茎区搜索处理用于针对矩阵N中各个元素n_i，搜索在n_i的第一个值所表示的位置和第二个值所表示的位置之间的RNA序列中存在的所有茎区。为了说明的简单，下面将“在n_i的第一个值所表示的位置和第二个值所表示的位置之间的RNA序列”简称为“与矩阵N的下一个元素n_i所对应的RNA子序列”。

在图3A中，初始化变量J_b、J_e，使得J_b＝n_(i，0)；J_e＝n_(i，1)。然后进入步骤S302，判断J_e-J_b是否大于等于4。如果步骤S302中的判断结果为否，则转至步骤S304，使得i＝i+1，如果加1后的i大于35，则处理结束，否则返回步骤S301，以判断与矩阵N的下一个元素n_i+1所对应的RNA子序列是否包含子茎区。如果步骤S302中的判断结果为是，则进入茎区确定处理步骤S303。茎区确定处理用于确定出在所输入的RNA序列中在位置F_b到位置F_e之间的RNA子序列中存在的所有子茎区。

下面结合图3B来说明茎区确定处理步骤S303。

如图3B所示，在步骤S311中初始化子茎区查找初始位置j_b和子茎区查找结束位置j_e，使得j_b＝n_(i，0)；j_e＝n_(i，1)。

然后，在步骤S312中，查找从位置j_b至位置j_e之间的茎区JQ。具体地，根据茎区定义，茎区是碱基对的连续区域，可知，可通过查找位置j_b和位置j_e之间的连续碱基对来找出茎区JQ。例如，可按照从两端位置(即位置j_b和位置j_e)向中间的方向、依次逐个碱基对地进行判断，判断在位置j_b和位置j_e之间是否存在连续的碱基对，如果发现某一段RNA子序列的碱基对是连续的，则将该段RNA子系统确定为茎区JQ。但是，茎区JQ的确定方法显然多种多样，例如，也可以在位置j_b和位置j_e之间的中间位置向两端依次判断碱基对是否连续。甚至可以找出从位置j_b和位置j_e之间多个茎区并从中选择碱基对数量最多的那个作为茎区JQ。

如果未找到茎区JQ，则结束处理。相反，如果找到了茎区JQ，则进行步骤S313，存储与茎区JQ有关的信息。例如，创建X的新节点x，使得x.start＝茎区JQ的开始位置、x.end＝茎区JQ的结束位置、x.energy＝茎区JQ中各碱基对的能量值的和以及X.number＝茎区JQ中包含的碱基对的数量，以存储与茎区JQ有关的信息。然后进入步骤S314，判断在步骤S313中找到的茎区JQ是否为祖先茎区，亦即，判断茎区JQ是否为在位置n_(i，0)和位置n_(i，1)找到的第一个茎区。例如，本实施例可通过判断j_b和j_e是否分别等于n_(i，0)和n_(i，1)来确定在步骤S313中找到的茎区JQ是否为祖先茎区。如果步骤S314中的判断结果为否，即，茎区JQ不是祖先茎区，则进入步骤S316，将茎区JQ的累加能量值加至祖先茎区数据中与茎区JQ的祖先节点相对应的元素的家族能量值中，例如在实施例中，可将x.energy加至e.energy。相反，如果步骤S314中的判断结果为是，即，茎区JQ是祖先茎区，则则进入步骤S317，创建祖先茎区数据的一个新元素，用于存储与茎区JQ相对应的与祖先茎区有关的信息。例如在实施例中，可以创建E的节点e，使得e.xpointer＝x、e.energy＝x.energy。

然后，从步骤S317进入步骤S315，使得j_b＝x.start+1且j_e＝x.end-1。在步骤S316中，判断j_e-j_b是否大于等于4，如果是，则返回步骤S312，以查找茎区JQ中的子茎区。如果步骤S316中为否，表示茎区JQ中不可能存在子茎区，然后茎区确定处理步骤结束。。

返回图1可知，根据上述结合图3A和图3B的说明，完成了图1中的步骤S104中的茎区搜索处理。

如图1所示，利用矩阵N中的所有元素进行茎区搜索处理之后，进入步骤S105，判断(F_e-1)-F_b＞＝4是否能满足(正如下文所述，如果存在已标记的碱基，则应判断是否满足(F_e-1)-F_b减去已标记的碱基数＞＝4，亦即在进行搜索时，将已标记的碱基视为不存在)。如果判断结果为是，则使得结束位置F_e减1，并返回步骤S103，否则进入步骤S106。

在步骤S106中，判断F_e-(F_b+1)＞＝4是否能满足(正如下文所述，如果存在已标记的碱基，则应判断是否满足F_e-(F_b+1)减去已标记的碱基数＞＝4，亦即已标记的碱基不参与新一轮的碱基对与茎区查找)，如果判断结果为是，则使得起始位置F_b加1，并返回步骤S102，否则进入步骤S107。

在步骤S107中，进行“最大茎区”确定处理，本发明根据MFE最小折叠自由能量原理，将在所搜索到的所有茎区中，与家族能量值最小的茎区家族对应的祖先茎区确定为“最大茎区”。也就是说，本发明可根据祖先茎区数据和茎区数据来确定“最大茎区”。更具体地，将祖先茎区中家族能量值最大的祖先茎区确定为“最大茎区”。在实施例中，可将E中e.energy最小(e.energy的绝对值最大)的元素所对应的祖先茎区确定为“最大茎区”。此外，本发明还可以将包含最多碱基对的祖先茎区确定为“最大茎区”。进而，本发明还可以根据祖先茎区数据中的与作为“最大茎区”的“祖先茎区”有关的信息来输出或标记整个“最大茎区”序列。例如，根据本发明的e.xpointer，可以在茎区数据中找到存储了与该“最大茎区”有关的信息的x节点，从而获得“最大茎区”的起始位置、结束位置、碱基对数量等信息。

此外，如果在确定了“最大茎区”后，还需要确定除“最大茎区”之外的“次最大茎区”(次大茎区)，则本实施例还可以在步骤S107之后，进入步骤S108。

在步骤S108中，删除被确定为最大茎区的碱基序列。然后进入步骤S109，以判断删除后所剩下的碱基序列中的碱基个数是否大于等于4，如是判断结果为是，则返回S102，将删除后剩下的碱基序列输入后，通过再次执行步骤S102-S107来确定当前所输入的碱基序列(被删除了上次执行S102-S107时被确定为“最大茎区”的序列后剩下的序列)中的“次最大茎区”，亦即经过步骤S108删除后所剩下的碱基序列中的最大茎区。反之，如果步骤S109的判断结果为否，则进入步骤S110，输出所有步骤S107中被确定为最大茎区的茎区。

如上所述，即使RNA序列中包含了非嵌套结构和假结结构，本发明的第一实施例能够较准确地确定RNA的二级结构。

第二实施例

接着将参考图2中所示的流程图说明本发明第二实施例。图2示出根据本发明的用于预测RNA的二级结构的预测方法的又一例子。根据本实施例的预测方法的步骤大体上可以与第一实施例中所述的预测装置的结构相同，因此，将省略对本实施例的预测方法相同的步骤的说明。图2中与第一实施例中相同的附图标记表示执行与第一实施例的大体相同的处理步骤，并将省略对这些大体相同步骤的说明。下面将仅对与第一实施例不同的步骤进行说明。

如图2所示，在第二实施例中，执行完步骤S107之后，执行步骤S208，标记出(而不是删除)最大茎区中的碱基，然后进入步骤S209。在步骤S209中，判断所输入的碱基序列中未被标记的碱基数是否大于等于4，如果判断结果为否，则进入步骤S110，输出所有茎区。反之，如果判断结果为是，则返回步骤S101。需要注意的是，在本实施例中，步骤S101将RNA序列连同标记一起输入，步骤S102至步骤S107只对未被标记的茎区进行处理，亦即，在步骤S102至步骤S107中，将未被标记的茎区视为不存在。

需重点说明的是，对已预测的茎区进行标记而不是像传统的预测方法那样删除已经预测的茎区，可以提高预测方法的查准率、查全率和准确度。因为删除已经预测的茎区可能会因为改变了原有RNA序列的结构而导致预测的查准度下降。

综上可知，根据本发明的技术方案利用了基于茎区结构相对稳定的特征，采用了贪婪(Greedy)算法、动态规划算法，来计算嵌套结构和包含交叉假接的非嵌套结构，能够预测RNA的嵌套结构、非嵌套结构和假结结构，能够降低时间复杂性、空间复杂度。

图6示出根据本发明第二实施例与Pknots技术方案的查全率与查准率对照表。在图6中，查全率(Sensitivity)Sn＝TP/RP，查准率(Specificity)Sp＝TP/(TP+FP)，其中，TP表示RNA结构中正确预测出的碱基对数目，RP表示RNA结构中真实碱基对的数目，FP表示RNA结构中错误预测的碱基对数目。从图6可以看出，本发明在查准率、查全率、时间复杂度和空间复杂度等方面，均具有优势。

此外，发明人对PseudoBase数据库中的包含假接的RNA分子片段进行了实验，实验结果表明，本发明的预测方法对约22％的分子片段的预测准确率达到87.5％以上；对约51％分子片段的准确率达到76.3％以上。此外，本发明的预测方法可预测包含多达6356个核苷酸的R NA分子。

其它实施例

本发明还可以通过硬件装置实现。

下面详细说明与根据本发明各实施例对RNA的二级结构进行预测的预测方法相对应的各个预测装置的实施例。注意，根据本实施例对RNA的二级结构进行预测的预测装置可以但不限通过程序软件与计算机硬件设备相结合的方式实现。例如，还可以通过将代码嵌入至带有中央处理器和存储器的处理设备的方式来实现。

图4和图5分别示出与图1和图2中用于预测RNA的二级结构的预测方法对应的测量装置。

图4中的输入及初始化单元401与第一分量查找单元402连接，用于实现与图1的步骤S101相对应的功能。第一分量查找单元402与第二分量查找单元403和茎区搜索处理单元404连接，且用于实现与图1的步骤S102相对应的功能。第二分量查找单元403与第一分量查找单元402和茎区搜索处理单元404连接，且用于实现与图1的步骤S103相对应的功能。茎区搜索处理单元404与第一分量查找单元402、第二分量查找单元403和最大茎区确定单元405连接，且用于实现与图1的步骤S104-S106相对应的功能。最大茎区确定单元405与茎区搜索处理单元404和最大茎区删除单元406连接，用于实现与图1的步骤S107相对应的功能。最大茎区删除单元406与最大茎区确定单元405和输出单元407相连，用于实现与图1的步骤S108-S109相对应的功能。此外，输出单元407用于实现与图1的步骤S110相对应的功能。

在图5中，与图4功能大体相同的单元用相同的附图标记表示。将省略对这些大体相同的单元的说明。下面将仅对与第一实施例不同的单元进行说明。

相比图4，图5中的测量装置不包括最大茎区删除单元，取而代之的是最大茎区标记单元506。最大茎区标记单元506与最大茎区确定单元405和输出单元407连接，且用于实现与图2中的步骤S208-S209对应的功能。此外，需要重点说明的是，图5中的测量装置的第一分量查找单元402、第二分量查找单元403、茎区搜索处理单元404和最大茎区确定单元405只处理未被标记的碱基，而不是处理所有输入的碱基。

下面说明本发明各实施例的时间复杂度和空间复杂度。

首先说明时间复杂度。在查找第一个茎区的过程中使用两重循环，第一个循环寻找茎区的前两个碱基(i＝1...n)(n为碱基个数)，第二个循环寻找茎区的后两个碱基(j＝1...n)，则找到第一个茎区的时间复杂度为O(kn²)(k为常数)，假设最坏情况下找到第一个茎区的层数为l₁，...找到第m个茎区的层数为l_m，则l₁+l₂+...+l_m＜＝n/2，因此，查找最大茎区的时间复杂度应为O(l₁n²+l₂n²+...+l_mn²)＜＝O(n.n²/2)＝O(n³/2)。在标记最大茎区的过程中，仅有指针的移动，没有字符的交换，与n无关，所以该算法时间复杂度为O(n³)。

对于空间复杂度，在查找最大茎区过程中只使用一维数据结构(数组)存储碱基序列，与碱基个数n有关，所以其空间复杂度为O(n)。

综上且结合图6可知，本发明相比传统技术在时间复杂度、空间复杂度、查全率和查准率上均有优势。

此外，本发明的各实施例的单元和处理流程可以相互结合，以组成更多的技术方案。

虽然已经参考示例性实施例对本发明作了说明，但是应当理解，本发明并不限于所公开的示例性实施例。权利要求书的范围符合最宽的解释，以包括所有的修改以及等同结构和功能。

Claims (10)

1.一种基于茎区的核糖核酸二级结构的预测方法，其特征在于，包括：

初始化步骤，定义碱基对数据、茎区数据和祖先茎区数据；

输入步骤，输入核糖核酸序列；

第一分量查找步骤，在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第一分量在所输入的核糖核酸序列所处的第一位置，并将所确定的所述第一位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第一分量中；

第二分量查找步骤，在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第二分量在所输入的核糖核酸序列所处的第二位置，并将所确定的所述第二位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第二分量中；

茎区搜索步骤，搜索在所述碱基对数据的各元素的第一分量和第二分量之间的核糖核酸序列中的所有茎区，将与所搜索到的与各茎区相对应的有关信息分别存储在所述茎区数据的各元素中，将与所搜索到的各祖先茎区有关的信息分别存储在所述祖先茎区数据的各元素中，其中所述所有茎区包含所搜索到的茎区的子茎区；

最大茎区确定步骤，根据所述茎区数据和所述祖先茎区数据来确定在所输入的核糖核酸序列中的最大茎区。

2.根据权利要求1所述的预测方法，其特征在于，

所述祖先茎区数据包括与各祖先茎区一一对应的各元素，所述祖先茎区数据的各元素用于存储相应祖先茎区所在家族的家族能量值和相应祖先茎区的位置信息。

3.根据权利要求2所述的预测方法，其特征在于，所述祖先茎区的位置信息是所述茎区数据中存储了与相应祖先茎区有关的信息的元素的位置或指向该位置的指针。

4.根据权利要求2所述的预测方法，其特征在于，所述最大茎区确定步骤将所述祖先茎区数据中家族能量值最小的祖先茎区确定为在所输入的核糖核酸序列中的最大茎区。

5.根据权利要求2所述的预测方法，其特征在于，在完成所述最大茎区确定步骤后，删除所述最大茎区，并返回所述输入步骤，以输入经删除所述最大茎区后剩余的核糖核酸序列，以确定在经删除所述最大茎区后剩余的核糖核酸序列中的最大茎区。

6.根据权利要求2所述的预测方法，其特征在于，在完成所述最大茎区确定步骤后，标识所述最大茎区，并返回所述输入步骤，以输入经标记所述最大茎区后的核糖核酸序列，并在后续步骤中不对所输入的核糖核酸序列中的已标记碱基对进行处理。

7.根据权利要求2所述的预测方法，其特征在于，所述祖先茎区数据的各元素存储了与相应祖先茎区有关的信息的元素的位置。

8.根据权利要求5或6所述的预测方法，其特征在于，还包括：

茎区输出步骤，根据述茎区数据和所述祖先茎区数据，输出所有在所述茎区搜索步骤中搜索到的各茎区的累加能量值、碱基对数量、起始位置、结束位置。

9.一种基于茎区的核糖核酸二级结构的预测装置，其特征在于，包括：

初始化单元，其定义碱基对数据、茎区数据和祖先茎区数据；

输入单元，其输入核糖核酸序列；

第一分量查找单元，其在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第一分量在所输入的核糖核酸序列所处的第一位置，并将所确定的所述第一位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第一分量中；

第二分量查找单元，其在所输入的核糖核酸序列中进行查找，以确定各种碱基对的第二分量在所输入的核糖核酸序列所处的第二位置，并将所确定的所述第二位置分别存储在所述碱基对数据中与各种碱基对相对应的元素的第二分量中；

茎区搜索单元，其搜索在所述碱基对数据的各元素的第一分量和第二分量之间的核糖核酸序列中的所有茎区，将与所搜索到的与各茎区相对应的有关信息分别存储在所述茎区数据的各元素中，将与所搜索到的各祖先茎区有关的信息分别存储在所述祖先茎区数据的各元素中，其中所述所有茎区包含所搜索到的茎区的子茎区；

最大茎区确定单元，其根据所述茎区数据和所述祖先茎区数据来确定在所输入的核糖核酸序列中的最大茎区。

10.根据权利要求9所述的预测装置，其特征在于，

所述祖先茎区数据包括与各祖先茎区一一对应的各元素，所述祖先茎区数据的各元素用于存储相应祖先茎区所在家族的家族能量值和所述相应祖先茎区的位置信息。

Images