CN106055920B - 一种基于阶段性多策略副本交换的蛋白质结构预测方法 - Google Patents

Abstract

一种基于阶段性多策略副本交换的蛋白质结构预测方法，在差分进化算法(DE)的框架下，首先，在每个温度层中生成具有多样化折叠类型的初始构象种群；然后，对每个温度层，根据迭代次数将构象搜索分为两个阶段，第一个阶段从种群中随机选取一个构象作为目标个体，第二个阶段根据能量将种群分为两个部分，从能量低的前50％的种群中随机选取一个个体作为目标个体；之后，随机选取与目标个体不一样的三个构象个体经过变异、交叉和片段组装策略生成测试个体，并根据构象的能量判断是否接受测试个体，进而将相邻的温度层的对应个体进行副本交换；在阶段性策略的指导下，辅以副本交换策略，通过不断更新种群获得一系列亚稳态构象。

Description

一种基于阶段性多策略副本交换的蛋白质结构预测方法

技术领域

本发明涉及生物信息学、计算机应用领域，尤其涉及的是一种基于阶段性多策略副本交换的蛋白质结构预测方法。

背景技术

发展一种根据氨基酸序列就能预测蛋白质的结构的计算机算法是一个巨大的挑战。一方面，了解天然结构的相关知识是了解生物机制、发现抑制或激活蛋白质的药物的一个出发点。另一方面，我们对结构的了解多过序列，而且由于高通量测序技术的发展，这个差距在增大。因此，能根据序列准确预测其结构的方法是有相当大的价值的。

在1994开始举办的CASP(蛋白质结构预测评估)竞赛中，Moult和他的同事加快了基于计算机的蛋白质结构预测的进程。每两年的夏天举行一次CASP是一个全社会的无目的竞赛，来自全世界的150多个研究小组会得到100多个不同的目标序列(蛋白质，结构已知但并未向外界公开)。每个参与小组运用一些算法方案来预测这些目标蛋白质的三维结构。每次CASP举行后，组委会随后会公布用于实验中蛋白质真正的结构，评估所有小组的表现，出版整个评估结果。

蛋白质结构从头预测是不依赖模板仅从氨基酸序列信息得到天然结构。它的关键是正确定义能量函数、精确选用计算机搜索算法来寻找能量最低值。与同源模建、折叠识别法相比，从头预测方法具有一定的特点和优势，它适用于同源性小于25％的大多数蛋白质，仅从序列产生全新结构，所以对于蛋白质分子设计以及蛋白质折叠的研究等具有重要意义，同时该方法及其理论在其它的蛋白蛋结构分析中也有重要作用，如同源模建方法中侧链安装，模型优化，X射线衍射晶体结构的优化。但是到目前还没有一种十分完善的方法来预测蛋白质的三维结构，即使获得了很好的预测结果，但也只是针对某些蛋白质而言的，目前主要的技术瓶颈在于两个方面，第一方面在于采样方法，现有技术对构象空间采样能力不强，另一方面在于构象更新方法，现有技术对构象的更新精度仍然不足。

因此，现有的蛋白质结构预测方法存在采样效率、复杂度及预测精度方面存在不足，需要改进。

发明内容

为了克服现有蛋白质结构预测存在采样效率低、复杂度较高及预测精度较低的不足，本发明提出一种预测精度较高、复杂度较低的基于阶段性多策略副本交换的蛋白质结构预测方法。

本发明解决其技术问题所采用的技术方案是：

一种基于阶段性多策略副本交换的蛋白质结构预测方法，所述预测方法包括以下步骤：

1)给定查询序列信息；

2)初始化：设置群体规模NP，变异因子F，交叉概率CR，n个温度层T，迭代次数iteration，片段长度为L，能量函数选用Rosetta Score3，首先通过对查询序列随机的折叠和变换，在每个温度层生成一个规模为NP的初始构象种群，初始群体为P＝{xⁱ|i∈I}，计算初始种群构象个体的能量值f(xⁱ),i∈I，并设其中i为种群个体编号，I为种群个体编号集合，I＝{1,2,...,NP}，为对应温度层能量函数值中的最小值，T为温度层；

3)开始迭代，依次对每个温度层种群中的每个个体进行以下操作：

3.1)令g＝1，g为迭代计数器；

3.2)随机选取种群中的一个构象个体xⁱ，i∈{1,2,3,…,NP}，令P_target＝xⁱ，P_target表示目标构象个体，执行以下操作得到变异构象P_mutant；

3.3)随机生成正整数rand1,rand2,rand3∈{1,2,3,......NP}，且rand1≠rand2≠rand3≠i，令P_rand1＝x^rand1，P_rand2＝x^rand2，P_rand3＝x^rand3；再生成4个随机整数randrange1,randrange2,randrange3,randrange4；其中randrange1≠randrange2，randrange3≠randrange4∈{1,2,…,Length},Length为序列长度；

3.4)令a＝min(randrange1,randrange2)，b＝max(randrange1,randrange2)；令c＝min(randrange3,randrange4)，d＝max(randrange3,randrange4)；其中min表示取两个数的最小值，max表示取两个数的最大值；

3.5)如果g≤iteration/2，用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant，跳转至3.7)；

3.6)如果iteration/2<g<iteration，则执行如下操作：

3.6.1)根据能量对整个种群中的构象进行降序排列，然后从能量低的前NP/2个构象个体中随机选出一个构象记为P_pbest；

3.6.2)用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_pbest上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant；

3.7)针对变异个体P_mutant做交叉操作，过程如下：

3.7.1)生成随机数rand4，rand5，其中rand4∈(0,1)，rand5∈(1,Length)；

3.7.2)根据执行交叉过程：若随机数rand4≤CR，个体P_mutant的片段rand5替换为个体P_target中对应的片段，否则直接继承变异个体得到测试个体P_trail；

3.8)针对每个温度层中的目标个体P_target和测试个体P_trail进行选择操作：

3.8.1)计算P_target和P_trail的能量：E(P_target)和E(P_trail)；

3.8.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入3.8)；

4)在每个温度层都进行变异、交叉和选择更新完种群后，将相邻温度层中对应的个体进行副本交换：

4.1)首先选择种群中第i个个体，i∈I(i＝1,2，…，NP)；

4.2)选择两个相邻的温度层T_j和T_j+1，j＝1,2，…，n；

4.3)产生一个随机数randx∈(0,1)和判别数：

其中，表示在T_j个温度层中第i个个体，其中i＝1,2,···,NP，为种群个体的能量值，K为玻尔兹曼常数，T_j(j＝1,2,…,n)为第j个温度层的开尔文温度，j为索引编号，e为自然常数；

4.4)若randx＜judge，将相邻温度层中的和个体进行交换；

5)判断是否满足终止条件，如不满足则转至3)；若满足则输出结果。

进一步，所述步骤5)中，对种群中的每个构象个体都执行完步骤3)-4)以后，迭代次数g＝g+1，终止条件为迭代次数g达到最大迭代次数iteration。

本发明的技术构思为：在差分进化算法框架下，首先，对输入的查询序列进行随机的折叠、扰动，在每个温度层中生成具有多样化折叠类型的初始构象种群；然后，对每个温度层，根据迭代次数将构象搜索分为两个阶段，第一个阶段从种群中随机选取一个构象作为目标个体，第二个阶段根据能量将种群分为两个部分，从能量低的前50％的种群中随机选取一个个体作为目标个体；之后，随机选取与目标个体不一样的三个构象个体经过变异、交叉和片段组装策略生成测试个体；在种群更新时，以Rosetta Score3为优化目标函数，基于Anfinsen提出的蛋白质天然态结构是自由能最低的结构，根据构象的能量判断是否接受测试个体，进而将相邻的温度层的对应个体进行副本交换；在阶段性策略的指导下，辅以副本交换策略，通过不断更新种群获得一系列亚稳态构象。

本发明的有益效果为：预测精度较高、复杂度较低。

附图说明

图1是测试序列在种群更新过程中RMSD和能量值的关系示意图。

图2是1ENH算法预测所得蛋白质三维结构示意图。

具体实施方式

下面结合附图对本发明作进一步描述。

参照图1～图2，一种基于阶段性多策略副本交换的蛋白质结构预测方法，包括以下步骤：

1)给定查询序列信息；

2)初始化：设置群体规模NP，变异因子F，交叉概率CR，n个温度层T，迭代次数iteration，片段长度为L，能量函数选用Rosetta Score3，首先通过对查询序列随机的折叠和变换，在每个温度层生成一个规模为NP的初始构象种群，初始群体为P＝{xⁱ|i∈I}，计算初始种群构象个体的能量值f(xⁱ),i∈I，并设其中i为种群个体编号，I为种群个体编号集合，I＝{1,2,...,NP}，为对应温度层能量函数值中的最小值，T为温度层；

3)开始迭代，依次对每个温度层种群中的每个个体进行以下操作：

3.1)令g＝1，g为迭代计数器；

3.2)随机选取种群中的一个构象个体xⁱ，i∈{1,2,3,…,NP}，令P_target＝xⁱ，P_target表示目标构象个体，执行以下操作得到变异构象P_mutant；

3.3)随机生成正整数rand1,rand2,rand3∈{1,2,3,......NP}，且rand1≠rand2≠rand3≠i，令P_rand1＝x^rand1，P_rand2＝x^rand2，P_rand3＝x^rand3；再生成4个随机整数randrange1,randrange2,randrange3,randrange4；其中randrange1≠randrange2，randrange3≠randrange4∈{1,2,…,Length},Length为序列长度；

3.4)令a＝min(randrange1,randrange2)，b＝max(randrange1,randrange2)；令c＝min(randrange3,randrange4)，d＝max(randrange3,randrange4)；其中min表示取两个数的最小值，max表示取两个数的最大值；

3.5)如果g≤iteration/2，用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant，跳转至3.7)；

3.6)如果iteration/2<g<iteration，则执行如下操作：

3.6.1)根据能量对整个种群中的构象进行降序排列，然后从能量低的前NP/2个构象个体中随机选出一个构象记为P_pbest；

3.6.2)用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_pbest上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant；

3.7)针对变异个体P_mutant做交叉操作，过程如下：

3.7.1)生成随机数rand4，rand5，其中rand4∈(0,1)，rand5∈(1,Length)；

3.7.2)根据执行交叉过程：若随机数rand4≤CR，个体P_mutant的片段rand5替换为个体P_target中对应的片段，否则直接继承变异个体得到测试个体P_trail；

3.8)针对每个温度层中的目标个体P_target和测试个体P_trail进行选择操作：

3.8.1)计算P_target和P_trail的能量：E(P_target)和E(P_trail)；

3.8.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入3.8)；

4)在每个温度层都进行变异、交叉和选择更新完种群后，将相邻温度层中对应的个体进行副本交换：

4.1)首先选择种群中第i个个体，i∈I(i＝1,2，…，NP)；

4.2)选择两个相邻的温度层T_j和T_j+1，j＝1,2，…，n；

4.3)产生一个随机数randx∈(0,1)和判别数：

其中，表示在T_j个温度层中第i个个体，其中i＝1,2,···,NP，为种群个体的能量值，K为玻尔兹曼常数，T_j(j＝1,2,…,n)为第j个温度层的开尔文温度，j为索引编号，e为自然常数；

4.4)若randx＜judge，将相邻温度层中的和个体进行交换；

5)判断是否满足终止条件，如不满足则转至3)；若满足则输出结果。

进一步，所述步骤5)中，对种群中的每个构象个体都执行完步骤3)-4)以后，迭代次数g＝g+1，终止条件为迭代次数g达到最大迭代次数iteration。

本实施例以序列长度为54的蛋白质1ENH为实施例，一种基于阶段性多策略副本交换的蛋白质结构预测方法，其中包含以下步骤：

1)给定查询序列信息；

2)初始化：设置群体规模NP＝50，变异因子F＝0.5，交叉概率CR＝0.5，8个温度层T，迭代次数iteration＝10000，片段长度为L＝3，能量函数选用RosettaScore3，首先通过对查询序列随机的折叠和变换，生成一个规模为NP的初始构象种群，初始群体为P＝{xⁱ|i∈I}，计算初始种群构象个体的能量值f(xⁱ),i∈I，并设其中i为种群个体编号，I为种群个体编号集合，I＝{1,2,...,NP}，为对应温度层能量函数值中的最小值，T为温度层；

3)开始迭代，依次对每个温度层种群中的每个个体进行以下操作：

3.1)令g＝1，g为迭代计数器；

3.2)随机选取种群中的一个构象个体xⁱ，i∈{1,2,3,…,NP}，令P_target＝xⁱ，P_target表示目标构象个体，执行以下操作得到变异构象P_mutant；

3.3)随机生成正整数rand1,rand2,rand3∈{1,2,3,......NP}，且rand1≠rand2≠rand3≠i，令P_rand1＝x^rand1，P_rand2＝x^rand2，P_rand3＝x^rand3；再生成4个随机整数randrange1,randrange2,randrange3,randrange4；其中randrange1≠randrange2，randrange3≠randrange4∈{1,2,…,Length},Length为序列长度；

3.4)令a＝min(randrange1,randrange2)，b＝max(randrange1,randrange2)；令c＝min(randrange3,randrange4)，d＝max(randrange3,randrange4)；其中min表示取两个数的最小值，max表示取两个数的最大值；

3.5)如果g≤iteration/2，用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant，跳转至3.7)；

3.6)如果iteration/2<g<iteration，则执行如下操作：

3.6.3)根据能量对整个种群中的构象进行降序排列，然后从能量低的前NP/2个构象个体中随机选出一个构象记为P_pbest；

3.6.4)用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_pbest上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant；

3.7)针对变异个体P_mutant做交叉操作，过程如下：

3.7.1)生成随机数rand4，rand5，其中rand4∈(0,1)，rand5∈(1,Length)；

3.7.2)根据执行交叉过程：若随机数rand4≤CR，个体P_mutant的片段rand5替换为个体P_target中对应的片段，否则直接继承变异个体得到测试个体P_trail；

3.8)针对每个温度层中的目标个体P_target和测试个体P_trail进行选择操作：

3.8.1)计算P_target和P_trail的能量：E(P_target)和E(P_trail)；

3.8.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入3.8)；

4)在每个温度层都进行变异、交叉和选择更新完种群后，将相邻温度层中对应的个体进行副本交换：

4.1)首先选择种群中第i个个体，i∈I(i＝1,2，…，NP)；

4.2)选择两个相邻的温度层T_j和T_j+1，j＝1,2，…，n；

4.3)产生一个随机数randx∈(0,1)和判别数：

其中，表示在T_j个温度层中第i个个体，其中i＝1,2,…,NP，为种群个体的能量值，K为玻尔兹曼常数，T_j(j＝1,2,…,n)为第j个温度层的开尔文温度，j为索引编号，e为自然常数；

4.4)若randx＜judge，将相邻温度层中的和个体进行交换；

5)如果迭代次数g等于最大迭代次数iteration，则输出结果，否则转至3)。

以序列长度为54的蛋白质1ENH为实施例，运用以上方法得到了该蛋白质的近天然态构象，构象系综中构象更新图如图1所示，算法预测所得蛋白质三维结构展示如图2所示。

以上阐述的是本发明给出的一个实施例表现出来的优良效果，显然本发明不仅适合上述实施例，在不偏离本发明基本精神及不超出本发明实质内容所涉及内容的前提下可对其做种种变化加以实施。

Claims (2)

1.一种基于阶段性多策略副本交换的蛋白质结构预测方法，其特征在于：所述预测方法包括以下步骤：

1)给定查询序列信息；

2)初始化：设置群体规模NP，变异因子F，交叉概率CR，n个温度层T，迭代次数iteration，片段长度为L，能量函数选用Rosetta Score3，首先通过对查询序列随机的折叠和变换，在每个温度层生成一个规模为NP的初始构象种群，初始群体为P＝{xⁱ|i∈I}，计算初始种群构象个体的能量值f(xⁱ),i∈I，并设其中i为种群个体编号，I为种群个体编号集合，I＝{1,2,...,NP}，为对应温度层能量函数值中的最小值，T为温度层；

3)开始迭代，依次对每个温度层种群中的每个个体进行以下操作：

3.1)令g＝1，g为迭代计数器；

3.2)随机选取种群中的一个构象个体xⁱ，i∈{1,2,3,… ,NP}，令P_target＝xⁱ，P_target表示目标构象个体，执行以下操作得到变异构象P_mutant；

3.3)随机生成正整数rand1,rand2,rand3，rand1,rand2和rand3均∈{1,2,3,......NP}，且rand1≠rand2≠rand3≠i，令P_rand1＝x^rand1，P_rand2＝x^rand2，P_rand3＝x^rand3；再生成4个随机整数randrange1,randrange2,randrange3,randrange4；其中randrange1≠randrange2，randrange3≠randrange4,randrange1,randrange2,randrange3和randrange4均∈{1,2,… ,Length}, Length为序列长度；

3.4)令a＝min(randrange1,randrange2)，b＝max(randrange1,randrange2)；令c＝min(randrange3,randrange4)，d＝max(randrange3,randrange4)；其中min表示取两个数的最小值，max表示取两个数的最大值；

3.5)如果g≤iteration/2，用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant，跳转至3.7)；

3.6)如果iteration/2<g<iteration，则执行如下操作：

3.6.1)根据能量对整个种群中的构象进行降序排列，然后从能量低的前NP/2个构象个体中随机选出一个构象记为P_pbest；

3.6.2)用P_rand2上位置a到位置b的片段的氨基酸所对应的二面角phi、psi、omega替换P_rand1的相同位置所对应的二面角phi、psi、omega；再使用P_rand3上位置c到位置d的片段的氨基酸所对应的二面角phi、psi、omega替换P_pbest上相同位置所对应的二面角phi、psi、omega，再将所得P_rand1进行片段组装得到变异个体P_mutant；

3.7)针对变异个体P_mutant做交叉操作，过程如下：

3.7.1)生成随机数rand4，rand5，其中rand4∈(0,1)，rand5∈(1,Length)；

3.7.2)根据执行交叉过程：若随机数rand4≤CR，个体P_mutant的片段rand5替换为个体P_target中对应的片段，否则直接继承变异个体得到测试个体P_trail；

3.8)针对每个温度层中的目标个体P_target和测试个体P_trail进行选择操作：

3.8.1)计算P_target和P_trail的能量：E(P_target)和E(P_trail)；

3.8.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入3.8)；

4)在每个温度层都进行变异、交叉和选择更新完种群后，将相邻温度层中对应的个体进行副本交换：

4.1)首先选择种群中第i个个体，i∈I(i＝1,2，… ，NP)；

4.2)选择两个相邻的温度层T_j和T_j+1，j＝1,2，… ，n；

4.3)产生一个随机数randx∈(0,1)和判别数：

其中，表示在T_j个温度层中第i个个体，其中i＝1,2,···,NP，为种群个体的能量值，K为玻尔兹曼常数，T_j(j＝1,2,…,n)为第j个温度层的开尔文温度，j为索引编号，e为自然常数；

4.4)若randx＜judge，将相邻温度层中的和个体进行交换；

5)判断是否满足终止条件，如不满足则转至3)；若满足则输出结果。

2.如权利要求1所述的一种基于阶段性多策略副本交换的蛋白质结构预测方法，其特征在于：所述步骤5)中，对种群中的每个构象个体都执行完步骤3)-4)以后，迭代次数g＝g+1，终止条件为迭代次数g达到最大迭代次数iteration。