发明内容
本公开的目的是提供一种基于动态规划和量子退火的RNA结构预测方法及系统,以通过结合动态规划和量子退火的优势,提高对于RNA二级结构预测的准确性。
为了实现上述目的,本公开实施例第一方面提供一种基于动态规划和量子退火的RNA结构预测方法,应用于基于动态规划和量子退火的RNA结构预测系统,所述系统包括隐马尔可夫模块、动态规划模块、量子退火模块和优化器,所述方法包括:
获取目标氨基酸序列的单字母表示序列;
将所述单字母表示序列输入隐马尔可夫模块,得到对应于所述单字母表示序列的目标密码子组合;
设置初始为空的密码子组合集合,将所述目标密码子组合输入动态规划模块,以使得所述动态规划模块从所述目标密码子组合中选取密码子组合片段,并判断所述密码子组合片段是否包含终止密码子;
在从所述目标密码子组合中选取的所述密码子组合片段不包含终止密码子时,基于当前的密码子组合集合计算出具有最优RNA二级结构的密码子组合初选片段,调用所述优化器对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,保存具有最高打分的密码子子结构的密码子组合片段至所述密码子组合集合,从所述目标密码子组合中重新选取密码子组合片段,重新选取的该密码子组合片段包含上一次从所述目标密码子组合中选取的密码子组合片段,并返回执行所述判断所述密码子组合片段是否包含终止密码子的步骤;
在从所述目标密码子组合中选取的所述密码子组合片段包含终止密码子时,从最后得到的密码子组合集合中输出具有最优RNA二级结构的完整密码子组合作为第一候选密码子组合;
确定所述目标密码子组合的哈密顿量,并将所述哈密顿量输入量子退火模块,基于所述量子退火模块的输出确定具有最优RNA二级结构的第二候选密码子组合;
从所述第一候选密码子组合以及所述第二候选密码子组合中选出具有最优RNA二级结构的密码子组合。
可选地,所述优化器包括JTVAE优化约束条件模块、热力学优化约束条件模块、CAI值和MFE值优化约束条件模块,所述调用所述优化器对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,保存具有最高打分的密码子子结构的密码子组合片段至所述密码子组合集合,包括:
调用所述优化器的JTVAE优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
调用所述优化器的热力学优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
调用所述优化器的CAI值和MFE值优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
针对所述密码子组合初选片段中具有同一预定义RNA子结构的密码子子结构,计算所述JTVAE优化约束条件模块、所述热力学优化约束条件模块、所述CAI值和MFE值优化约束条件模块打分分数的均值,保存具有最高所述均值的密码子子结构的密码子组合片段至所述密码子组合集合。
可选地,所述预定义RNA子结构包括发卡环、内环、多分支环以及碱基对。
可选地,所述调用所述优化器的JTVAE优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,包括:
调用所述优化器的JTVAE优化约束条件模块以对所述密码子组合初选片段进行分割,得到多个具有预定义RNA子结构的密码子子结构;
将多个所述密码子子结构输入图神经网络,得到所述密码子子结构的树形结构表示数据;
将所述树形结构表示数据输入编码器,得到所述树形结构表示数据对应的潜在表示向量;
基于重参数化技巧将潜在表示向量与预设噪声结合,得到最终的潜在表示向量;
将所述最终的潜在表示向量输入解码器,得到所述密码子子结构的打分分数。
可选地,所述调用所述优化器的CAI值和MFE值优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,包括:
计算所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构的CAI值和MFE值;
基于所述CAI值和MFE值计算所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构的打分分数。
可选地,所述确定所述目标密码子组合的哈密顿量,包括:
基于如下公式确定所述目标密码子组合的哈密顿量:
;
;
;
其中,H为所述目标密码子组合的哈密顿量,a、b、和/>均为可调参数,/>是第i个茎的碱基对长度,u是最长茎的长度,/>为对应第i个茎的发卡环惩罚项,/>和/>分别是对应于第i个茎和第j个茎的二进制变量,/>是对假结的惩罚函数,e是核苷酸的线性维度,/>是假结中单股核苷酸的数量,/>是长度为n的第k个内嵌茎的惩罚常数,/>是第i个茎和第j个茎的交叉惩罚函数。
本公开实施例第二方面,提供一种基于动态规划和量子退火的RNA结构预测系统,所述系统包括隐马尔可夫模块、动态规划模块、量子退火模块和优化器,其中:
所述隐马尔可夫模块用于对输入的目标氨基酸的单字母表示序列进行计算,得到对应于所述单字母表示序列的目标密码子组合;
所述动态规划模块用于在设置初始为空的密码子组合集合后,从所述目标密码子组合中选取密码子组合片段,并判断所述密码子组合片段是否包含终止密码子;在从所述目标密码子组合中选取的所述密码子组合片段不包含终止密码子时,基于当前的密码子组合集合计算出具有最优RNA二级结构的密码子组合初选片段,调用所述优化器对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,保存具有最高打分的密码子子结构的密码子组合片段至密码子组合集合,从所述目标密码子组合中重新选取密码子组合片段,重新选取的该密码子组合片段包含上一次从所述目标密码子组合中选取的密码子组合片段,并返回执行所述判断所述密码子组合片段是否包含终止密码子的步骤;在从所述目标密码子组合中选取的所述密码子组合片段包含终止密码子时,从最后得到的密码子组合集合中输出具有最优RNA二级结构的完整密码子组合作为第一候选密码子组合;
所述量子退火模块用于在确定所述目标密码子组合的哈密顿量后,对输入的所述哈密顿量进行计算,基于所述量子退火模块的输出确定具有最优RNA二级结构的第二候选密码子组合;以使得所述系统从所述第一候选密码子组合以及所述第二候选密码子组合中选出具有最优RNA二级结构的密码子组合。
可选地,所述优化器包括JTVAE优化约束条件模块、热力学优化约束条件模块、CAI值和MFE值优化约束条件模块,其中:
所述优化器的JTVAE优化约束条件模块用于对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
所述优化器的热力学优化约束条件模块用于对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
所述优化器的CAI值和MFE值优化约束条件模块用于对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
所述优化器还用于针对所述密码子组合初选片段中具有同一预定义RNA子结构的密码子子结构,计算所述JTVAE优化约束条件模块、所述热力学优化约束条件模块、所述CAI值和MFE值优化约束条件模块打分分数的均值,保存具有最高所述均值的密码子子结构的密码子组合片段至所述密码子组合集合。
可选地,所述预定义RNA子结构包括发卡环、内环、多分支环以及碱基对。
可选地,所述JTVAE优化约束条件模块还用于:
对所述密码子组合初选片段进行分割,得到多个具有预定义RNA子结构的密码子子结构;
将多个所述密码子子结构输入图神经网络,得到所述密码子子结构的树形结构表示数据;
将所述树形结构表示数据输入编码器,得到所述树形结构表示数据对应的潜在表示向量;
基于重参数化技巧将潜在表示向量与预设噪声结合,得到最终的潜在表示向量;
将所述最终的潜在表示向量输入解码器,得到所述密码子子结构的打分分数。
通过上述技术方案,利用包含隐马尔可夫模型的隐马尔可夫模块将目标氨基酸序列的单字母表示序列转化为对应的目标密码子组合,以得到目标氨基酸序列最可能的RNA一级结构。然后分别利用动态规划和量子退火计算目标密码子组合的二级结构,一方面利用动态规划避免重复计算,基于递归提高计算效率和准确性,另一方面利用量子退火克服解非全局最优的问题,将两者的优点结合,以找到具有最优RNA二级结构的密码子组合,提高计算的准确性。
本公开的其他特征和优点将在随后的具体实施方式部分予以详细说明。
具体实施方式
以下结合附图对本公开的具体实施方式进行详细说明。应当理解的是,此处所描述的具体实施方式仅用于说明和解释本公开,并不用于限制本公开。以下首先对一些基本概念进行解释:
RNA(Ribonucleic Acid,核糖核酸),即缩合核糖核苷酸经磷脂键,由此形成的长链状分子。在动植物、微生物、一些病毒以及噬箘体内都有RNA分子。一个核糖核苷酸分子由磷酸,核糖和碱基构成,它是一种以磷酸二酯键连接几十种核糖核苷酸的形式存在的核酸,而以核苷酸作为其基本构成单元。RNA之中,主要包含的碱基为四种,即A(Adenine,腺嘌呤)、G(Guanine,鸟嘌呤)、C(Cytosine,胞嘧啶)以及U(Uracil,尿嘧啶)。
RNA一级结构是指核苷酸序列单链,它由核酸的四个基本组成单位通过磷酸二酯键连接而成,形成一条线性的多聚体。
RNA二级结构,即处于自然条件之下,RNA分子开展的回折,处于各种碱基互补对茎区之间,运用最顶部单链、环区形式,由此交替产生的茎环结构。
茎区结构,可明确两方面构成了部分双螺旋区域,其中之一为依托于碱基互补对形式实现;其二,即在双螺旋区域外隔离,不存在碱基对的单链,由此构成了环结构。
假结,往往两大形式被包括于RNA分子结构内。其一,即为平面假结,即处于中环结构上、环结构外部的碱基配对,可运用平面图形予以体现的一种假结形式;与之相反的,即为非平面假结。
RNA分子的三维折叠在很大程度上是由分子内碱基间氢键的模式决定的。从序列中预测碱基配对网络,也称为RNA二级结构预测或RNA折叠,是一个不确定性多项式时间(NP)完全计算问题。
RNA分子折叠成复杂的二级结构,这决定了它们的分子性质,如热稳定性和致密性。此外,RNA折叠还影响RNA在蛋白质翻译、转录调控和其他重要细胞过程中的功能。因此,RNA二级结构的测定方法对于基础研究、应用生物技术和合理的药物发现具有重要的意义。为此目的开发的实验方法非常耗时和昂贵,因此在实践中的使用受到限制。
图1是根据一示例性实施例示出的一种基于动态规划和量子退火的RNA结构预测系统的框图,如图1所示,所述系统包括隐马尔可夫模块、动态规划模块、量子退火模块和优化器,其中:
所述隐马尔可夫模块用于对输入的目标氨基酸的单字母表示序列进行计算,得到对应于所述单字母表示序列的目标密码子组合;
所述动态规划模块用于在设置初始为空的密码子组合集合后,从所述目标密码子组合中选取密码子组合片段,并判断所述密码子组合片段是否包含终止密码子;在从所述目标密码子组合中选取的所述密码子组合片段不包含终止密码子时,基于当前的密码子组合集合计算出具有最优RNA二级结构的密码子组合初选片段,调用所述优化器对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,保存具有最高打分的密码子子结构的密码子组合片段至密码子组合集合,从所述目标密码子组合中重新选取密码子组合片段,重新选取的该密码子组合片段包含上一次从所述目标密码子组合中选取的密码子组合片段,并返回执行所述判断所述密码子组合片段是否包含终止密码子的步骤;在从所述目标密码子组合中选取的所述密码子组合片段包含终止密码子时,从最后得到的密码子组合集合中输出具有最优RNA二级结构的完整密码子组合作为第一候选密码子组合;
所述量子退火模块用于在确定所述目标密码子组合的哈密顿量后,对输入的所述哈密顿量进行计算,基于所述量子退火模块的输出确定具有最优RNA二级结构的第二候选密码子组合;以使得所述系统从所述第一候选密码子组合以及所述第二候选密码子组合中选出具有最优RNA二级结构的密码子组合。
可选地,所述优化器包括JTVAE优化约束条件模块、热力学优化约束条件模块、CAI值和MFE值优化约束条件模块,其中:
所述优化器的JTVAE优化约束条件模块用于对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
所述优化器的热力学优化约束条件模块用于对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
所述优化器的CAI值和MFE值优化约束条件模块用于对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
所述优化器还用于针对所述密码子组合初选片段中具有同一预定义RNA子结构的密码子子结构,计算所述JTVAE优化约束条件模块、所述热力学优化约束条件模块、所述CAI值和MFE值优化约束条件模块打分分数的均值,保存具有最高所述均值的密码子子结构的密码子组合片段至所述密码子组合集合。
可选地,所述预定义RNA子结构包括发卡环、内环、多分支环以及碱基对。
可选地,所述JTVAE优化约束条件模块还用于:
对所述密码子组合初选片段进行分割,得到多个具有预定义RNA子结构的密码子子结构;
将多个所述密码子子结构输入图神经网络,得到所述密码子子结构的树形结构表示数据;
将所述树形结构表示数据输入编码器,得到所述树形结构表示数据对应的潜在表示向量;
基于重参数化技巧将潜在表示向量与预设噪声结合,得到最终的潜在表示向量;
将所述最终的潜在表示向量输入解码器,得到所述密码子子结构的打分分数。
可选地,所述CAI值和MFE值优化约束条件模块还用于:
计算所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构的CAI值和MFE值;
基于所述CAI值和MFE值计算所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构的打分分数。
可选地,所述目标密码子组合的哈密顿量为:
;
;
;
其中,H为所述目标密码子组合的哈密顿量,a、b、和/>均为可调参数,/>是第i个茎的碱基对长度,u是最长茎的长度,/>为对应第i个茎的发卡环惩罚项,/>和/>分别是对应于第i个茎和第j个茎的二进制变量,/>是对假结的惩罚函数,e是核苷酸的线性维度,是假结中单股核苷酸的数量,/>是长度为n的第k个内嵌茎的惩罚常数,/>是第i个茎和第j个茎的交叉惩罚函数。
关于上述实施例中的系统,其中各个模块执行操作的具体方式通过下述基于动态规划和量子退火的RNA结构预测方法进行详细描述。
图2是根据一示例性实施例示出的一种基于动态规划和量子退火的RNA结构预测方法的流程图,该方法可以应用于如图1所示的基于动态规划和量子退火的RNA结构预测系统,参见图2,该基于动态规划和量子退火的RNA结构预测方法包括:
S201,获取目标氨基酸序列的单字母表示序列。
S202,将所述单字母表示序列输入隐马尔可夫模块,得到对应于所述单字母表示序列的目标密码子组合。
S203,设置初始为空的密码子组合集合,将所述目标密码子组合输入动态规划模块,以使得所述动态规划模块从所述目标密码子组合中选取密码子组合片段,并判断所述密码子组合片段是否包含终止密码子;
在从所述目标密码子组合中选取的所述密码子组合片段不包含终止密码子时,基于当前的密码子组合集合计算出具有最优RNA二级结构的密码子组合初选片段,调用所述优化器对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,保存具有最高打分的密码子子结构的密码子组合片段至所述密码子组合集合,从所述目标密码子组合中重新选取密码子组合片段,重新选取的该密码子组合片段包含上一次从所述目标密码子组合中选取的密码子组合片段,并返回执行所述判断所述密码子组合片段是否包含终止密码子的步骤;
在从所述目标密码子组合中选取的所述密码子组合片段包含终止密码子时,从最后得到的密码子组合集合中输出具有最优RNA二级结构的完整密码子组合作为第一候选密码子组合。
S204,确定所述目标密码子组合的哈密顿量,并将所述哈密顿量输入量子退火模块,基于所述量子退火模块的输出确定具有最优RNA二级结构的第二候选密码子组合。
S205,从所述第一候选密码子组合以及所述第二候选密码子组合中选出具有最优RNA二级结构的密码子组合。
具体来讲,在步骤S201中,目标氨基酸序列可以是目标蛋白质的氨基酸序列,氨基酸可以通过单字母进行表示,例如甘氨酸用字母G表示,丙氨酸用字母A表示,进而目标氨基酸序列可以用对应的单字母表示序列进行表示,如以下为一示例性氨基酸序列的单字母表示序列:
“MGVAARPPALRHWFSHSIPLAIFALLLLYLSVRSLGARSGCGPRAQPCVPGETAPFQVRQESGTLEAPERKQPPCLGPRGMLGRMMRRFHASLKPEGDVGLSPYLAGWRALVEFLTPLGSVFAFATREAFTKVTDLEARVHGPDAEHYWSLVAMAAWERRAGLLEQPGAAPRDPTRSSGSRTLLLLHRALRWSQLCLHRVATGALGGPDAGVQCSDAYRAALGPHHPWLVRQTARLAFLAFPGRRRLLELACPGATEAEARAALVRAAGTLEDVYNRTQSLLAERGLLQLA”
具体来讲,可以通过有线或无线通信输入的方式来获取目标氨基酸序列的单字母表示序列。获取单字母表示序列后,可以进入执行步骤S202。
在步骤S202中,隐马尔可夫模块采用隐马尔可夫模型,是一个能够通过可观测数据很好地捕捉真实空间统计性质的随机模型,通过事先建立的隐马尔可夫模型,将单字母表示序列输入隐马尔可夫模块后,可以利用其中的隐马尔可夫模型将单字母表示序列转化为对应的目标密码子组合。
密码子是生物细胞用来将遗传物质中的信息翻译成蛋白质的一套规则,每个密码子由3个相邻的碱基组成,每个氨基酸可能有多个对应的密码子,例如甘氨酸对应的密码子有GGU、GGC、GGA以及GGG。密码子组合则是多个密码子排成的序列组合。由于单个氨基酸可能有多个对应的密码子,因此采用隐马尔可夫模型,根据可观测的单字母表示序列来预测目标氨基酸对应RNA的密码子组合,以提高密码子组合一级结构预测的准确性,进而为准确预测RNA二级结构奠定基础。
隐马尔可夫模型通常包括:
状态集合S={s1,s2,…,sN}S’={s1’,s2’,…,sN’}:这是模型中所有可能的隐藏状态。
观测集合O={o1,o2,…,oM}O’={o1’,o2’,…,oM’}:这是所有可能的观测结果。
状态转移概率矩阵A=[aij]A’=[aij’],其中aij=P(qt+1=sj∣qt=si)aij’=P(qt+1’=sj’∣qt’=si’)表示在时间t处于状态sisi’的情况下,在时间(t+1)转移到状态 sjsj’的概率。
观测概率矩阵 B=[bj(k)]B’=[bj’(k)],其中bj(k)=P(ok at time t∣qt=sj)bj’(k)=P(ok’ at time t∣qt’=sj’) 表示在时间 t 处于状态 sjsj’ 的情况下,观测到okok’的概率。
初始状态概率分布π=[πi]π’=[πi’],其中πi=P(q1=si)πi’=P(q1’=si’) 表示模型在时间t=1处于状态sisi’的概率。
在隐马尔可夫模型中,通常有三个基本问题:
评估问题(Forward-Backward Algorithm):给定模型参数和观测序列,计算观测序列的概率。这通常通过前向算法(Forward Algorithm)来计算。
解码问题(Viterbi Algorithm):给定模型参数和观测序列,找到最有可能产生这些观测的隐藏状态序列。这通常通过维特比算法(Viterbi Algorithm)来解决。
学习问题(Baum-Welch Algorithm):调整模型参数以最大化给定观测序列的概率。这通常通过Baum-Welch算法(一种特殊的EM算法)来实现。
举例来讲,可以先确定模型参数如通过密码子偏好型确定转移概率,确定模型参数后,将单字母表示序列作为观测序列,通过上述维特比算法求解隐藏状态序列以得到对应的目标密码子组合。并在计算出目标密码子组合后,进入执行步骤S203和步骤S204。
在步骤S203中,预先设置一个空的集合作为密码子组合集合,用于存储动态规划计算过程产生的具有最高打分的密码子组合片段。将目标密码子组合输入动态规划模块,动态规划模块利用动态规划的思想来对目标密码子组合进行计算,通过先定出一小段密码子组合片段的最优二级结构,再用相同的方法将密码子组合片段扩展,找到相应的最好二级结构,这种方法不断迭代进行,直至全长序列。密码子组合片段为密码子组合中的部分连续碱基构成的片段。
具体来讲,每次迭代计算前,首先从目标密码子组合中选取密码子组合片段,例如初始可以选择包含4个碱基的密码子组合片段,并根据密码子组合片段判断是否满足迭代停止条件。
若选取的密码子组合片段中不包含终止密码子,说明对目标密码子组合的计算没有完成,进而对选取的密码子组合片段进行计算,计算过程利用动态规划算法的状态转移方程,例如可以利用基于动态规划的碱基最大配对方法或极小自由能方法,在一种可能的实施方式中,状态转移方程如下:
OPT(i,j)=max(OPT(i,j-1),max(1+OPT(i,t-1)+OPT(t+1,j-1)))
其中,i+1<t<j-2,OPT(i,j)表示表示碱基i和碱基j对应的密码子组合片段上的二级结构中碱基对的最大数目。通过遍历t来对OPT(i,j)进行计算。
计算过程可以使用之前已计算出的密码子组合片段的结果,进而在后续计算时直接通过查表获取相应结果,以计算出具有最优RNA二级结构的密码子组合初选片段,其最优可以通过如最大碱基配对数量或最小自由能进行评价。
在计算出密码子组合初选片段后,调用优化器基于相应的优化约束条件对密码子组合初选片段进行优化。具体可以先选出密码子组合初选片段中具有预定义RNA子结构的密码子子结构,然后对其分别进行打分,将包含最高打分的密码子子结构的密码子组合片段保存到密码子组合集合用于下一次迭代计算。保存后,重新选取新的密码子组合片段,新选取的密码子组合片段比上一次的密码子组合片段更长,且包含上一次选取的密码子组合片段,以便于利用动态规划进行计算。
若选取的密码子组合片段中包含终止密码子,说明对目标密码子组合的计算完成,此时从最后得到的密码子组合集合中找到包含所有密码子的完整密码子组合,该完整密码子组合包含了通过动态规划计算得到的最优RNA二级结构,进而将该密码子组合作为第一候选密码子组合。
传统的算法可能会在复杂度或准确性上遇到困难,尤其是在面对大型mRNA序列时。动态规划,作为一种强大的优化工具,可以有效地捕捉到目标氨基酸序列的对应密码子组合的最优性。引入束剪枝技术进一步加强了动态规划算法的效率,允许算法在广泛的搜索空间中找到前K最优的解。
在步骤S204中,首先确定目标密码子组合的哈密顿量,然后将该哈密顿量输入量子退火模块,以利用量子退火算法找到对应哈密顿量的基态,基态包含表征最优RNA二级结构的数据,进而根据该数据得到具有最优RNA二级结构的第二候选密码子组合。其中量子退火算法可以采用在量子计算机上实现的量子退火算法,也可以采用在经典计算机上实现的量子模拟退火算法。
可选地,步骤S204中,确定所述目标密码子组合的哈密顿量,包括:
基于如下公式确定所述目标密码子组合的哈密顿量:
;
;
;
其中,H为所述目标密码子组合的哈密顿量,a、b、和/>均为可调参数,/>是第i个茎的碱基对长度,u是最长茎的长度,/>为对应第i个茎的发卡环惩罚项,/>和/>分别是对应于第i个茎和第j个茎的二进制变量,/>是对假结的惩罚函数,e是核苷酸的线性维度,是假结中单股核苷酸的数量,/>是长度为n的第k个内嵌茎的惩罚常数,/>是第i个茎和第j个茎的交叉惩罚函数。
具体来讲,可以根据密码子组合进行预先计算确定该哈密顿量中的参数,或者通过类似于机器学习的训练方式来调整该哈密顿量的参数。通过训练来调整哈密顿量的参数的步骤可以包括:
获取用于训练的密码子组合样本及对应的RNA二级结构数据标签;
基于预设的初始参数以及对应哈密顿量H的公式确定所述密码子组合样本的初始哈密顿量;
利用量子退火算法对所述初始哈密顿量进行求解,基于求解结果得到RNA二阶结构数据结果;
计算输出的RNA二阶结构数据结果与RNA二级结构数据标签的马修斯相关系数评分;
基于同步扰动随机逼近算法和所述马修斯相关系数评分对所述哈密顿量的参数进行更新,得到更新后的哈密顿量;
若马修斯相关系数评分不满足预设要求且对哈密顿量参数的迭代更新次数小于预设次数,则利用量子退火算法对更新后的哈密顿量进行求解,得到输出的更新后的哈密顿量的RNA二阶结构数据结果,并返回执行所述计算输出的RNA二阶结构数据结果与RNA二级结构数据标签的马修斯相关系数评分的步骤;
若马修斯相关系数评分满足预设要求或对哈密顿量参数的迭代更新次数大于等于预设次数,将最后更新后的哈密顿量确定为所述目标密码子组合的哈密顿量。
具体来讲,用于训练的RNA的密码子组合样本的RNA二级结构数据是已知的,进而在获取RNA的密码子组合样本及对应的RNA二级结构数据标签后,可以给哈密顿量公式中的可调参数随机赋值,并基于密码子组合样本的结构确定其它参数,以得到初始哈密顿量。然后利用量子退火算法求解初始哈密顿量,根据求解结果得到输出的RNA二阶结构数据结果,其中变量为1表示存在相应的茎,变量为0表示不存在相应的茎。然后与机器学习的训练方法类似,通过计算马修斯相关系数评分来对计算结果与标签进行比较,并通过同步扰动随机逼近算法计算梯度,以通过梯度下降对哈密顿量的参数进行更新,直至满足相应条件时,将最后更新的哈密顿量作为目标密码子组合的哈密顿量。若不满足相应条件,则迭代执行相应更新步骤。在确定哈密顿量公式中的其它参数时,例如可以通过构建碱基序列的矩阵,在可配对位置设置元素1,其它位置设置元素0,并确定垂直于对角线的连续3个及以上元素1的碱基对为潜在茎,进而确定相关参数。
通过该哈密顿量,旨在确定一组潜在的茎,使连续碱基对的数量和茎的平均长度都最大化。其第一个求和的部分比较第i个茎的长度与潜在最长茎的长度,从而强制执行平均茎长度应该最大化的概念。此外还考虑奖励较大茎的形成,而不是一些在总能量上相同的较小茎。在考虑单个茎时,引入发卡环惩罚项,以禁止大小为1和2的环,而大小为3或更大的环则具有递减的惩罚,直到大小为7的环,此后所有环都被分配相同的惩罚。
哈密顿量的第二部分为较长的茎提供二次奖励,并引入对假结和茎交叉的惩罚,以减少涉及这些类型结构的碱基对的能量贡献。考虑到假结通常在其组成的串内嵌茎(in-line stems)之间展示短的单股区域,并且避免包含内嵌茎,同时考虑到计算的现实性,建立该哈密顿量的假结惩罚如上公式所示。如此以提高对RNA二级结构预测的求解准确度。
当然,在其它的实施例中,也可以基于如下公式确定所述目标密码子组合的哈密顿量:
;
;
;
其中,、/>、/>、/>为可调常数,/>是包含对数反密码子使用频率的向量,/>为第i个密码子中鸟嘌呤和胞嘧啶的个数,/>为第j个密码子中鸟嘌呤和胞嘧啶的个数,/>代表目标GC含量,N表示系统所包含的密码子数量,r(/>,/>)返回密码子/>和密码子/>之间的重复序列核苷酸的最大数目,/>表示在位置i找到密码子C的频率,/>、/>为二进制变量。
该哈密顿量同时考虑了多个优化目标,一是密码子使用偏好,二是GC含量与目标值之间的差异,三是连续核苷酸的数量。为了对不转译为查询序列的密码子组合施加能量惩罚,还添加了两个额外的约束条件。第一个约束条件对每个量子比特的一体项(one-bodyterm)添加了一个小的线性偏移。通过偏移势能,增加了将更多的密码子包含在序列中的能量有利性。类似地,另一个约束条件对映射到氨基酸序列中相同位置的密码子施加了显著的能量惩罚。相较于无效的密码子组合,这两个势能的结合对有效密码子组合的能量评分进行了优化。
在传统算法中,搜索最优mRNA序列的过程容易受到局部最优解的困扰,导致最终设计的序列并不是全局最优。而量子退火的独特之处在于,它利用了量子物理的特性,允许同时考虑多个可能的解决方案,并通过逐渐降低系统的能量来找到问题的全局最优解。在mRNA序列设计中,这意味着量子退火可以更有效地探索密码子的各种可能组合,从而更有可能找到那些能够提高mRNA的稳定性和翻译效率的序列。此外,量子退火的并行搜索特性也意味着在处理大型mRNA序列时,它可以提供更快的计算速度和更高的准确性。
步骤S203和步骤S204执行完之后,进入执行步骤S205。在步骤S205中,从基于动态规划模块计算出的第一候选密码子组合和基于量子退火模块计算出的第二候选密码子组合中选出具有最优RNA二级结构的密码子组合作为如图1所示的优化结果。举例来讲,可以计算第一候选密码子组合以及第二候选密码子组合的密码子适应指数(CAI,CodonAdaptation Index)和最小自由能(MFE,Minimum Free Energy)等参数,CAI值越大,MFE值越小,说明对应的密码子组合的RNA二级结构越优。
通过上述技术方案,利用包含隐马尔可夫模型的隐马尔可夫模块将目标氨基酸序列的单字母表示序列转化为对应的目标密码子组合,以得到目标氨基酸序列最可能的RNA一级结构。然后分别利用动态规划和量子退火对计算目标密码子组合的二级结构,一方面利用动态规划避免重复计算,基于递归提高计算效率和准确性,另一方面利用量子退火克服解非全局最优的问题,将两者的优点结合,以找到具有最优RNA二级结构的密码子组合,提高计算的准确性。
可选地,步骤S203中,调用所述优化器对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,保存具有最高打分的密码子子结构的密码子组合片段至所述密码子组合集合,包括:
S2031,调用所述优化器的JTVAE优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
S2032,调用所述优化器的热力学优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
S2033,调用所述优化器的CAI值和MFE值优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分;
S2034,针对所述密码子组合初选片段中具有同一预定义RNA子结构的密码子子结构,计算所述JTVAE优化约束条件模块、所述热力学优化约束条件模块、所述CAI值和MFE值优化约束条件模块打分分数的均值,保存具有最高所述均值的密码子子结构的密码子组合片段至所述密码子组合集合。
在一种可能的实施方式中,预定义RNA子结构包括发卡环、内环、多分支环以及碱基对。
可选地,步骤S2031中,调用所述优化器的JTVAE优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,包括:
S20311,调用所述优化器的JTVAE优化约束条件模块以对所述密码子组合初选片段进行分割,得到多个具有预定义RNA子结构的密码子子结构;
S20312,将多个所述密码子子结构输入图神经网络,得到所述密码子子结构的树形结构表示数据;
S20313,将所述树形结构表示数据输入编码器,得到所述树形结构表示数据对应的潜在表示向量;
S20314,基于重参数化技巧将潜在表示向量与预设噪声结合,得到最终的潜在表示向量;
S20315,将所述最终的潜在表示向量输入解码器,得到所述密码子子结构的打分分数。
JTVAE优化约束条件模块将结合树变分自编码器JTVAE(Junction TreeVariational Autoencoder)与图神经网络GNN(Graph Neural Network)结合。以进一步丰富基于动态规划和量子退火的RNA结构预测方法的能力,使其可以从多个层次上深入理解和优化RNA序列。这种深度学习结构不仅捕捉到了RNA序列的拓扑结构和复杂特性,而且还为其提供了更为细致和深入的分析,从而为研究人员带来了更为准确和高效的RNA设计结果。
具体来讲,JTVAE包括解码器和编码器。步骤S20311中,将密码子组合初选片段输入JTVAE优化约束条件模块后,首先按照预定义RNA子结构对密码子组合初选片段进行分割,得到多个密码子子结构如可能得到发卡环、内环、多分支环以及碱基对中的一种或多种结构。然后在步骤S20312中,将密码子子结构输入图神经网络,利用图神经网络捕捉和理解其拓扑结构,以得到树形结构表示数据用于表示密码子子结构。
图神经网络GNN是一种专门为处理图形数据而设计的深度学习框架。通过使用GNN,我们可以有效地捕捉和理解RNA的拓扑结构,从而为后续的计算和优化提供有力的支持。这个结构即树形结构表示数据由多个节点和边组成。每个节点可以代表RNA的某一部分,如一个碱基或一个密码子,而边则代表这些部分之间的相互作用或连接。
在GNN的工作过程中,信息会沿着边在节点之间传递。这意味着,每当信息从一个节点传递到另一个节点时,它都会经过一定的变换和整合,从而更好地捕捉到整个RNA拓扑结构中的特征和模式。
为了使这一过程更加高效,GNN中的每个节点都配备了一个或多个“邻居”节点。这使得节点可以不仅仅从直接连接的邻居获取信息,还可以从更远的节点获取信息。这种方法确保了在整个网络中的信息流动是高效和全面的。
在描述RNA的拓扑结构时,GNN的这一特性尤为重要。由于RNA的结构通常是复杂且多变的,因此捕捉其细微的拓扑特征对于理解其功能和特性至关重要。GNN的节点表示结构单元,边表示他们之间的关系和链接。每个节点还包含其化学性质、空间位置等信息。
JTVAE是一种为了解决分子结构和其他树形数据的表示和生成问题而设计的深度学习模型。这种模型结合了图结构的表示能力和变分自编码器(VAE)的生成能力,以确保在特定的应用中得到有效和准确的结果。
在步骤S20313中,将树形结构表示数据输入JTVAE中的编码器,以将其转化为一个低维的连续向量空间中的点。这个向量通常被称为“潜在变量”或“潜在表示”。在此过程中,编码器通过学习数据的内在模式和结构来产生两个输出:一个是潜在变量的均值μ,另一个是其标准差σ。这两个输出描述了一个正态分布,从中可以采样得到实际的潜在变量。
在步骤S20314中,潜在变量通过一个重参数化技巧与N(0,1)的噪声结合,产生最终的潜在表示向量z。这种技巧是为了使模型能够在训练期间进行反向传播。
在步骤S20315中,解码器接收这个最终的潜在表示向量z,并试图从中重建原始的树形数据。这一过程可以看作是将潜在空间中的点映射回到原始的数据空间。解码器的目标是确保重建的数据尽可能接近原始数据。在解码过程中,可以得到密码子子结构的打分分数。
JTVAE模型可以事先进行训练,为了训练JTVAE模型,通常采用两个损失函数。第一个是重建损失,它衡量原始数据和重建数据之间的差异。第二个是KL散度,它衡量编码器输出的正态分布与标准正态分布之间的差异。通过优化这两个损失,JTVAE能够学习有效的数据表示并生成相似的数据样本。
步骤S2032中,热力学优化约束条件模块利用相关的热力学优化约束条件对具有预定义RNA子结构的密码子子结构进行打分,如可以利用结构稳定性指标、结构动力学指标等对其进行打分。
可选地,步骤S2033中,调用所述优化器的CAI值和MFE值优化约束条件模块对所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构进行打分,包括:
S20331,计算所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构的CAI值和MFE值;
S20332,基于所述CAI值和MFE值计算所述密码子组合初选片段中具有预定义RNA子结构的密码子子结构的打分分数。
在步骤S20331中,首先计算密码子子结构的CAI值和MFE值,然后在步骤S20332中,根据计算出的这两个数值进行打分,例如可以对CAI值和MFE值进行加权求和作为打分分数。
在步骤S2034中,针对JTVAE优化约束条件模块、所述热力学优化约束条件模块、所述CAI值和MFE值优化约束条件模块计算的打分分数,计算其分数均值,将分数均值最高的密码子子结构的密码子组合片段保存到密码子组合集合中。
当然,在其它的实施例中,还可以考虑利用GC含量对密码子组合进行优化,以确保mRNA序列的整体稳定性和表达效率。
以上结合附图详细描述了本公开的优选实施方式,但是,本公开并不限于上述实施方式中的具体细节,在本公开的技术构思范围内,可以对本公开的技术方案进行多种简单变型,这些简单变型均属于本公开的保护范围。
另外需要说明的是,在上述具体实施方式中所描述的各个具体技术特征,在不矛盾的情况下,可以通过任何合适的方式进行组合,为了避免不必要的重复,本公开对各种可能的组合方式不再另行说明。
此外,本公开的各种不同的实施方式之间也可以进行任意组合,只要其不违背本公开的思想,其同样应当视为本公开所公开的内容。