CN106096326A - 一种基于质心变异策略的差分进化蛋白质结构预测方法 - Google Patents

Abstract

一种基于质心变异策略的差分进化蛋白质结构预测方法，首先，根据各构象的能量值进行升序排列，并计算各构象与能量最低构象的平均能量误差值；然后，选取部分能量较低的构象计算质心构象；最后，根据平均能量误差值判断算法所达到的搜索状态，从而设计不同的质心变异策略生成测试构象，即如果平均能量误差值大于设定的阈值，则设计DE/rand‑to‑centroid/1策略进行变异，通过提取质心构象中的部分片段替换随机选取的构象中的对应片段生成测试构象，否则设计DE/centroid/2策略进行变异，通过提取随机选择的构象中的片段替换质心构象中的对应片段生成测试构象，从而提高算法搜索效率和预测精度。

Description

一种基于质心变异策略的差分进化蛋白质结构预测方法

技术领域

本发明涉及一种生物学信息学、智能优化、计算机应用领域，尤其涉及的是，一种基于质心变异策略的差分进化蛋白质结构预测方法。

背景技术

1953年，J.Watson和F.Crick在英国《Nature》杂志上发表了DNA分子双螺旋结构模型，标志着分子生物学真正意义上的诞生；五年后，F.Crick提出分子生物学“中心法则”的设想，揭示了生命遗传信息传递的一般规律。作为该法则的关键部分，从DNA到蛋白质氨基酸序列的三联遗传密码(简称“第一密码”)的破译工作早在1965年就已经全部完成；然而，从氨基酸序列到空间结构的折叠密码(简称“第二密码”)至今尚未破解。随着2003年人类基因组测序工作的完成，蛋白质氨基酸序列数量激增，蛋白质折叠密码的理论研究成为当前蛋白质工程领域迫切需要解决的一个关键问题。

结构基因组学利用实验手段来测定蛋白质的三维结构。X射线晶体学方法是至今为止研究蛋白质结构最有效的方法，所能达到的精度是任何其他方法所不能比拟的，它的缺点主要是蛋白质的晶体难以培养且晶体结构测定的周期较长。多维核磁共振方法可以直接测定蛋白质在溶液中的构象，但是由于对样品的需要量大、纯度要求高，目前只能测定小分子蛋白质。总体上，蛋白质结构实验测定方法极其费时费钱费力。

从头预测方法被誉为蛋白质结构预测领域的圣杯，鉴于其重要的生物学意义和问题的复杂性，2005年《Science》杂志将其列为当前科学界亟待解决的100个最具挑战性问题之一。蛋白质从头预测方法必须考虑以下两个因素：(1)蛋白质结构能量函数；(2)构象空间搜索方法。第一个因素本质上属于分子力学问题，主要是为了能够计算得到每个蛋白质结构对应的能量值。第二个因素本质上属于全局优化问题，通过选择一种合适的优化方法，对构象空间进行快速搜索，得到与某一全局最小能量对应的构象。其中，蛋白质构象空间优化属于一类非常难解的NP-Hard问题。群体进化类算法是研究蛋白质分子构象优化的重要方法，主要包括差分进化算法(DE)、遗传算法(GA)、粒子群算法(PSO)，这些算法不仅结构简单，易于实现，而且鲁棒性强，因此，经常被用于从头预测方法中的全局最小能量构象搜索。然而群体优化算法属于一类随机优化算法，现有蛋白质构象优化方面的文献主要研究如何从一个局部最小解跳到另一个局部最小解，没有提供一种机制有效利用群体进化过程的智能信息指导搜索，从而导致算法效率较低。此外，受选择压力和随机采样过程中遗传漂变的影响，群体中所有个体将不可避免收敛到某个吸收态。对于蛋白质构象这类优化问题，该吸收态并不一定就是全局最优解，从而影响预测精度。

因此，现有的基于群体的蛋白质结构预测方法在搜索效率和预测精度方面存在着缺陷，需要改进。

发明内容

为了克服现有的蛋白质结构预测方法在搜索效率和预测精度方面的不足，本发明通过提取能量较低的构象信息，设计质心变异策略，同时基于片段组装技术，提出一种搜索效率高、预测精度高的基于质心变异策略的差分进化蛋白质结构预测方法。

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

一种基于质心变异策略的差分进化蛋白质结构预测方法，所述优化方法包括以下步骤：

1)选取蛋白质力场模型，即能量函数E(X)；

2)给定输入序列信息；

3)初始化：设置种群大小NP，交叉因子CR，最大迭代次数，由输入序列产生初始构象种群并初始化迭代次数G＝0，其中，N表示维数，表示第i个构象Cⁱ的第N维元素；

4)计算当前种群各构象的能量函数值E(Cⁱ),i＝1,2,…,N，并根据当前种群中各构象能量值对各构象进行升序排列；

5)找出当前种群中能量最低的构象C_best，并计算其他构象的能量与C_best的能量E(C_best)的平均能量误差如果迭代次数G＝0，则令δ_max＝δ；

6)针对种群中的每个构象个体Cⁱ，i∈{1,2,3,…,NP}，令C_target＝Cⁱ，C_target表示目标构象个体，提取当前种群中能量较低的构象信息，执行以下操作生成变异构象C_mutant：

6.1)选取排名前CT个构象其中CT＝rand(NP/3,NP/2)，rand(NP/3,NP/2)表示NP/3和NP/2之间的随机整数，表示第m个选取构象的第N维元素；

6.2)计算所选取的CT个构象的质心构象C_centroid＝(x_centroid,1,x_centroid,2,…,x_centroid,N)，其中，构象C_centroid的第j维元素j＝1,2,…,N；

6.3)设置序列长度L，在1和L之间随机生成4个整数randint1、randint2、randint3和randint4，其中randint1和randint2，randint3和randint4互不相同，令a＝min(randint1,randint2)，b＝max(randint1,randint2)，c＝min(randint3,randint4)，d＝max(randint3,randint4)，其中min表示取两个数的最小值，max表示取两个数的最大值；

6.4)如果δ＞0.5δ_max，则设计DE/rand-to-centroid/策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取质心构象C_centroid位置a到位置b的片段的氨基酸所对应的二面角替换构象C^rand1的相同位置所对应的二面角，同时提取构象C^rand2位置c到位置d的片段的氨基酸所对应的二面角替换构象C^rand1相同位置所对应的二面角，然后将所得C^rand1进行片段组装得到变异构象个体C_mutant；

6.5)如果δ≤0.5δ_max，则设计DE/centroid/2策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取构象C^rand1位置a到位置b的片段的氨基酸所对应的二面角替换质心构象C_centroid的相同位置所对应的二面角，同时使用C^rand2上位置c到位置d的片段的氨基酸所对应的二面角替换质心构象C_centroid相同位置所对应的二面角，然后将所得C_centroid进行片段组装得到变异构象个体C_mutant；；

7)对变异构象C_mutant执行交叉操作生成测试构象C_trial：

7.1)在0和1之间随机生成小数rand3；

7.2)若rand3≤CR，则在1和L之间随机生成整数rand4，利用变异构象C_mutant中的片段rand4替换目标构象C_target中对应的片段，从而生成测试构象C_trial，若rand3>CR，则C_trial直接等于变异构象C_mutant；

8)计算测试构象C_trial的能量值E(C_trial)，如果E(C_trial)-E(C_target)＜0，表明测试构象优于目标构象，则测试构象C_trial替换目标构象C_target；

9)判断是否满足终止条件，若满足则输出结果并退出，否则返回步骤4)。

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

本发明的技术构思为：首先，根据各构象的能量值进行升序排列，并计算各构象与能量最低构象的平均能量误差值；然后，选取部分能量较低的构象计算质心构象；最后，根据平均能量误差值判断算法所达到的搜索状态，从而设计不同的质心变异策略生成测试构象，即如果平均能量误差值大于设定的阈值，则设计DE/rand-to-centroid/1策略进行变异，通过提取质心构象中的部分片段替换随机选取的构象中的对应片段生成测试构象，否则设计DE/centroid/2策略进行变异，通过提取随机选择的构象中的片段替换质心构象中的对应片段生成测试构象，从而提高算法搜索效率和预测精度。

本发明的有益效果表现在：根据能量较低的构象计算质心构象，并通过提取质心构象的进化信息设计质心变异策略生成测试构象，从而提高预测精度；其次，根据平均能量误差值判断算法所达到的搜索状态，从而设计适合对应状态的质心变异策略生成测试构象，达到提高算法搜索效率的效果。

附图说明

图1是本发明中蛋白质结构预测方法的流程图。

图2是本发明中的预测方法对蛋白质4ICB预测时的构象更新示意图。

图3是本发明中的预测方法对蛋白质4ICB预测时得到的构象分布图。

图4是本发明中的预测方法对蛋白质4ICB预测得到的三维结构。

具体实施方式

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

参照图1和图4，一种基于质心变异策略差分进化蛋白质结构预测方法，包括以下步骤：

1)选取蛋白质力场模型，即能量函数E(X)；

2)给定输入序列信息；

3)初始化：设置种群大小NP，交叉因子CR，最大迭代次数，由输入序列产生初始构象种群并初始化迭代次数G＝0，其中，N表示维数，表示第i个构象Cⁱ的第N维元素；

4)计算当前种群各构象的能量函数值E(Cⁱ),i＝1,2,…,N，并根据当前种群中各构象能量值对各构象进行升序排列；

5)记当前种群中能量最低的构象C_best，并计算其他构象的能量与C_best的能量E(C_best)的平均能量误差如果迭代次数G＝0，则令δ_max＝δ；

6)针对种群中的每个构象个体Cⁱ，i∈{1,2,3,…,NP}，令C_target＝Cⁱ，C_target表示目标构象个体，提取当前种群中能量较低的构象信息，执行以下操作生成变异构象C_mutant：

6.1)选取排名前CT个构象其中CT＝rand(NP/3,NP/2)，rand(NP/3,NP/2)表示NP/3和NP/2之间的随机整数，表示第m个选取构象的第N维元素；

6.2)计算所选取的CT个构象的质心构象C_centroid＝(x_centroid,1,x_centroid,2,…,x_centroid,N)，其中，构象C_centroid的第j维元素j＝1,2,…,N；

6.3)设置序列长度L，在1和L之间随机生成4个整数randint1、randint2、randint3和randint4，其中randint1和randint2，randint3和randint4互不相同，令a＝min(randint1,randint2)，b＝max(randint1,randint2)，c＝min(randint3,randint4)，d＝max(randint3,randint4)，其中min表示取两个数的最小值，max表示取两个数的最大值；

6.4)如果δ＞0.5δ_max，则设计DE/rand-to-centroid/策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取质心构象C_centroid位置a到位置b的片段的氨基酸所对应的二面角替换构象C^rand1的相同位置所对应的二面角，同时提取构象C^rand2位置c到位置d的片段的氨基酸所对应的二面角替换构象C^rand1相同位置所对应的二面角，然后将所得C^rand1进行片段组装得到变异构象个体C_mutant；

6.5)如果δ≤0.5δ_max，则设计DE/centroid/2策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取构象C^rand1位置a到位置b的片段的氨基酸所对应的二面角替换质心构象C_centroid的相同位置所对应的二面角，同时使用C^rand2上位置c到位置d的片段的氨基酸所对应的二面角替换质心构象C_centroid相同位置所对应的二面角，然后将所得C_centroid进行片段组装得到变异构象个体C_mutant；；

7)为了提高种群的多样性，对变异构象C_mutant执行交叉操作生成测试构象C_trial：

7.1)在0和1之间随机生成小数rand3；

7.2)若rand3≤CR，则在1和L之间随机生成整数rand4，利用变异构象C_mutant中的片段rand4替换目标构象C_target中对应的片段，从而生成测试构象C_trial，若rand3>CR，则C_trial直接等于变异构象C_mutant；

8)计算测试构象C_trial的能量值E(C_trial)，如果E(C_trial)-E(C_target)＜0，表明测试构象优于目标构象，则测试构象C_trial替换目标构象C_target；

9)判断是否满足终止条件，若满足则输出结果并退出，否则返回步骤4)。

所述步骤9)中，对种群中的每个构象个体都执行完步骤6)-8)以后，迭代次数G＝G+1，终止条件为迭代次数G达到步骤3)中预设的最大迭代次数

本实施例序列长度为76的α折叠蛋白质4ICB为实施例，一种基于质心变异策略的差分进化蛋白质结构预测方法，其中包含以下步骤：

1)选取Rosetta score3力场模型，即能量函数E(X)；

2)输入蛋白质4ICB的序列信息；

3)初始化：设置种群大小NP＝50，交叉因子CR＝0.5，最大迭代次数为10000，由输入序列产生初始构象种群 并初始化迭代次数G＝0，其中，N表示维数，表示第i个构象Cⁱ的第N维元素；

4)计算当前种群各构象的能量函数值E(Cⁱ),i＝1,2,…,N，并根据当前种群中各构象能量值对各构象进行升序排列；

5)记当前种群中能量最低的构象C_best，并计算其他构象的能量与C_best的能量E(C_best)的平均能量误差如果迭代次数G＝0，则令δ_max＝δ；

6)针对种群中的每个构象个体Cⁱ，i∈{1,2,3,…,NP}，令C_target＝Cⁱ，C_target表示目标构象个体，提取当前种群中能量较低的构象信息，执行以下操作生成变异构象C_mutant：

6.1)选取排名前CT个构象其中CT＝rand(NP/3,NP/2)，rand(NP/3,NP/2)表示NP/3和NP/2之间的随机整数，表示第m个选取构象的第N维元素；

6.2)计算所选取的CT个构象的质心构象C_centroid＝(x_centroid,1,x_centroid,2,…,x_centroid,N)，其中，构象C_centroid的第j维元素j＝1,2,…,N；

6.3)设置序列长度L＝76，在1和L之间随机生成4个整数randint1、randint2、randint3和randint4，其中randint1和randint2，randint3和randint4互不相同，令a＝min(randint1,randint2)，b＝max(randint1,randint2)，c＝min(randint3,randint4)，d＝max(randint3,randint4)，其中min表示取两个数的最小值，max表示取两个数的最大值；

6.4)如果δ＞0.5δ_max，则设计DE/rand-to-centroid/策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取质心构象C_centroid位置a到位置b的片段的氨基酸所对应的二面角替换构象C^rand1的相同位置所对应的二面角，同时提取构象C^rand2位置c到位置d的片段的氨基酸所对应的二面角替换构象C^rand1相同位置所对应的二面角，然后将所得C^rand1进行片段组装得到变异构象个体C_mutant；

6.5)如果δ≤0.5δ_max，则设计DE/centroid/2策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取构象C^rand1位置a到位置b的片段的氨基酸所对应的二面角替换质心构象C_centroid的相同位置所对应的二面角，同时使用C^rand2上位置c到位置d的片段的氨基酸所对应的二面角替换质心构象C_centroid相同位置所对应的二面角，然后将所得C_centroid进行片段组装得到变异构象个体C_mutant；；

7)为了提高种群的多样性，对变异构象C_mutant执行交叉操作生成测试构象C_trial：

7.1)在0和1之间随机生成小数rand3；

7.2)若rand3≤CR，则在1和L之间随机生成整数rand4，利用变异构象C_mutant中的片段rand4替换目标构象C_target中对应的片段，从而生成测试构象C_trial，若rand3>CR，则C_trial直接等于变异构象C_mutant；

8)计算测试构象C_trial的能量值E(C_trial)，如果E(C_trial)-E(C_target)＜0，表明测试构象优于目标构象，则测试构象C_trial替换目标构象C_target；

9)对种群中的每个构象个体都执行完步骤6)-8)以后，迭代次数G＝G+1，若迭代次数G达到最大迭代次数10000，则输出结果并退出，否则返回步骤4)。

以序列长度为76的α折叠蛋白质4ICB为实施例，运用以上方法得到了该蛋白质的近天然态构象，最小均方根偏差为平均均方根偏差为预测得到的三维结构如图4所示。

以上阐述的是本发明给出的一个实施例表现出来的优良优化效果，显然本发明不仅适合上述实施例，而且可以应用到实际工程中的各个领域，同时在不偏离本发明基本精神及不超出本发明实质内容所涉及内容的前提下可对其做种种变化加以实施。

Claims (2)

1.一种基于质心变异策略的差分进化蛋白质结构预测方法，其特征在于：所述蛋白质结构预测方法包括以下步骤：

1)选取蛋白质力场模型，即能量函数E(X)；

2)给定输入序列信息；

3)初始化：设置种群大小NP，交叉因子CR，最大迭代次数，由输入序列产生初始构象种群并初始化迭代次数G＝0，其中，N表示维数，表示第i个构象Cⁱ的第N维元素；

4)计算当前种群各构象的能量函数值E(Cⁱ),i＝1,2,…,N，并根据当前种群中各构象能量值对各构象进行升序排列；

5)找出当前种群中能量最低的构象C_best，并计算其他构象的能量与C_best的能量E(C_best)的平均能量误差如果迭代次数G＝0，则令δ_max＝δ；

6)针对种群中的每个构象个体Cⁱ，i∈{1,2,3,…,NP}，令C_target＝Cⁱ，C_target表示目标构象个体，提取当前种群中能量较低的构象信息，执行以下操作生成变异构象C_mutant：

6.1)选取排名前CT个构象其中CT＝rand(NP/3,NP/2)，rand(NP/3,NP/2)表示NP/3和NP/2之间的随机整数，表示第m个选取构象的第N维元素；

6.2)计算所选取的CT个构象的质心构象C_centroid＝(x_centroid,1,x_centroid,2,…,x_centroid,N)，其中，构象C_centroid的第j维元素

6.3)设置序列长度L，在1和L之间随机生成4个整数randint1、randint2、randint3和randint4，其中randint1和randint2，randint3和randint4互不相同，令a＝min(randint1,randint2)，b＝max(randint1,randint2)，c＝min(randint3,randint4)，d＝max(randint3,randint4)，其中min表示取两个数的最小值，max表示取两个数的最大值；

6.4)如果δ＞0.5δ_max，则设计DE/rand-to-centroid/策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取质心构象C_centroid位置a到位置b的片段的氨基酸所对应的二面角替换构象C^rand1的相同位置所对应的二面角，同时提取构象C^rand2位置c到位置d的片段的氨基酸所对应的二面角替换构象C^rand1相同位置所对应的二面角，然后将所得C^rand1进行片段组装得到变异构象个体C_mutant；

6.5)如果δ≤0.5δ_max，则设计DE/centroid/2策略进行变异：从当前种群中随机选取两个不同的构象C^rand1和C^rand2，其中rand1≠rand2∈[1,NP]，提取构象C^rand1位置a到位置b的片段的氨基酸所对应的二面角替换质心构象C_centroid的相同位置所对应的二面角，同时使用C^rand2上位置c到位置d的片段的氨基酸所对应的二面角替换质心构象C_centroid相同位置所对应的二面角，然后将所得C_centroid进行片段组装得到变异构象个体C_mutant；；

7)对变异构象C_mutant执行交叉操作生成测试构象C_trial：

7.1)在0和1之间随机生成小数rand3；

7.2)若rand3≤CR，则在1和L之间随机生成整数rand4，利用变异构象C_mutant中的片段rand4替换目标构象C_target中对应的片段，从而生成测试构象C_trial，若rand3>CR，则C_trial直接等于变异构象C_mutant；

8)计算测试构象C_trial的能量值E(C_trial)，如果E(C_trial)-E(C_target)＜0，表明测试构象优于目标构象，则测试构象C_trial替换目标构象C_target；

9)判断是否满足终止条件，若满足则输出结果并退出，否则返回步骤4)。

2.如权利要求1所述的一种基于质心变异策略的双层差分进化蛋白质结构预测方法，其特征在于：所述步骤9)中，对种群中的每个构象个体都执行完步骤6)-8)以后，迭代次数G＝G+1，终止条件为迭代次数G达到步骤3)中预设的最大迭代次数。