CN109872770B - 一种结合排挤度评价的多变异策略蛋白质结构预测方法 - Google Patents

Abstract

一种结合排挤度评价的多变异策略蛋白质结构预测方法，在进化算法框架下，首先，建立三种不同的变异策略，根据轮盘赌的方式选择一个变异策略对构象进行变异，并且对生成的变异构象进行一次3片段组装，生成变异构象；其次对变异构象进行交叉操作；最后用Rosetta能量函数score3、蒙特卡洛玻尔兹曼接收准则对索引为

的构象进行选择，用排挤度指标Exclusion、蒙特卡洛玻尔兹曼接收准则对索引为

的构象进行选择，结合排挤度评价的多变异策略蛋白质结构预测方法不仅可以增强种群的多样性，而且能够缓解能量函数不精确的问题，提高采样效率。本发明提供一种采样效率高、预测精度高的结合排挤度评价的多变异策略蛋白质结构预测方法。

Description

一种结合排挤度评价的多变异策略蛋白质结构预测方法

技术领域

本发明涉及生物信息学、计算机应用领域，尤其涉及的是一种结合排挤度评价的多变异策略蛋白质结构预测方法。

背景技术

1953年，Watson和Crick在《Nature》杂志上发表了DNA分子双螺旋结构模型，标志着分子生物学真正意义上的诞生；1958年，Crick提出分子生物学“中心法则”设想，深刻地揭示了生物遗传信息由DNA向蛋白质传递与表达的一般过程；1965年，Nirenberg、Khorana和Holley发现三联体遗传密码(即第一遗传密码)，DNA以三个核苷酸为一组的密码子转译为蛋白质的氨基酸序列；1990年，美国启动人类基因组计划，并于2003年正式宣布完成。十几年以来，人类基因组计划不断深化人类对自身和疾病的认识，对生物学、医学、数学和计算机科学带来了深远影响。然而至今为止，尚未呈现美国前总统克林顿当时所描绘的蓝图：“彻底改变我们对绝大多数疾病的诊断、预防和治疗手段”。究其原因，是基因图谱仅描绘了蛋白质的氨基酸序列(即蛋白质一级结构)，而蛋白质只有折叠成特定的三维结构(即蛋白质三级结构)才能产生其特定的生物学功能。相对于第一遗传密码，蛋白质序列一级结构与其三级结构之间的对应关系(即第二遗传密码或称折叠密码)仍为未解之谜。《Science》杂志在纪念创刊125周年之际，把“能否预测蛋白质折叠？”列为21世纪125个科学前沿问题之一，指出在今后的25年内，人类将有望发现这一重要基础性问题的答案。为了解决蛋白质折叠这一“世纪之问”，越来越多不同学科背景的研究者参与其中，特别是针对蛋白质折叠过程的终点—蛋白质结构预测，受到了科学界广泛的关注和研究。相对蛋白质折叠而言，蛋白质结构预测具有更强的实用性，只有获得蛋白质的三维结构，才能真正实现基因诊断，并最终达到基因治疗的目的。

测定蛋白质三维结构的实验方法主要包括X射线晶体衍射和多维核磁共振(NMR)。X射线晶体衍射是目前测定蛋白质结构最有效的方法，所达到的精度是其它方法所不能比拟的，主要缺点是蛋白质晶体难以培养且晶体结构测定的周期较长；NMR方法可以直接测定蛋白质在溶液中的构象，但是对样品的需要量大、纯度要求高，目前只能测定小分子蛋白质。实验测定结构方法存在的主要问题在于两个方面：一方面，对于现代药物设计的主要靶标—膜蛋白而言，极难获得其结构；另外，实验测定过程费时费钱费力，代价不菲，比如，使用NMR方法测定一个蛋白质结构通常需要15万美元以及半年的时间。

在理论探索和应用需求的双重推动下，根据Anfinsen法则，利用计算机设计适当的算法，以序列为起点，三维结构为目标的蛋白质结构预测自20世纪末蓬勃发展。计算机硬件和软件技术的快速发展，为从蛋白质三级结构预测的发展提供了结实的基础平台条件。蛋白质三级结构预测是生物信息学的一个重要任务。蛋白质构象优化问题现在面临最大的挑战是对及其复杂的蛋白质能量模型函数曲面进行搜索，而且模型越精细，越能提供细节知识，所需的计算资源也就越多。

从头预测方法直接基于蛋白质物理或知识能量模型，利用优化算法在构象空间搜索全局最低能量构象解。构象空间优化(或称采样)方法是目前制约蛋白质结构从头预测精度最关键的因素之一。优化算法应用于从头预测采样过程必须首先解决以下三个方面的问题：(1)能量模型的复杂性。蛋白质能量模型考虑分子体系成键作用以及范德华力、静电、氢键、疏水等非成键作用，致使其形成的能量曲面极其粗糙，局部极小解数量随着序列长度的增加呈指数增长；能量模型的漏斗特性也必然会产生局部高能量障碍，导致算法极易陷入局部解。(2)能量模型高维特性。就目前而言，从头预测方法只能应对尺寸较小(<150残基)的目标蛋白，一般不超过100。对尺寸超过150残基以上的目标蛋白，除了CASP11中的仅有个例之外(如T0806)，现有优化方法均无能为力。这也就进一步说明了随着尺寸规模的增加，必然造成维数灾问题，完成如此浩瀚的构象搜索过程所涉及的计算量是目前最先进的计算机也难以承受的。(3)能量模型的不精确性。对于蛋白质这类复杂的生物大分子，除了考虑各种物理成键和知识推理的作用之外，还要考虑它与周围溶剂分子的相互作用，目前还无法给出精确的物理描述。考虑到计算代价问题，近十年来研究者陆续提出了一些列基于物理的力场简化模型(AMBER，CHARMM等)、基于知识的力场简化模型(Rosetta,QUARK等)。然而，我们还远远无法构建起能引导目标序列朝正确方向折叠的足够精确力场，导致数学上的最优解并不一定对应于目标蛋白的天然态结构；此外，模型的不精确性也必然会导致无法对算法性能进行客观地分析，从而阻碍了高性能算法在蛋白质结构从头预测领域中的应用。

随着氨基酸序列的增长，蛋白质分子体系自由度也增大，利用传统群体算法采样获得大规模蛋白质构象空间的全局最优解成为一项挑战性的工作；其次，粗粒度模型虽然减小了构象搜索空间，但是也导致了相互作用力之间的信息丢失，从而直接影响预测精度。

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

发明内容

为了克服现有的蛋白质结构预测方法采样效率低、种群多样性差、预测精度较低的不足，本发明在基本差分进化算法框架下，引入多变异策略来指导构象空间搜索，同时结合排挤度评价的策略选择构象，提出一种采样效率高、预测精度高的结合排挤度评价的多变异策略蛋白质结构预测方法。

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

一种结合排挤度评价的多变异策略蛋白质结构预测方法，所述预测方法包括以下步骤：

1)给定目标蛋白的序列信息；

2)根据目标蛋白序列从ROBETTA服务器(http://www.robetta.org/)上得到片段库文件，其中包括3片段库文件和9片段库文件；

3)设置参数：种群大小NP，算法的最大迭代代数G，交叉因子CR，温度因子β，置迭代代数g＝0；

4)种群初始化：随机片段组装生成NP个初始构象C_i，i＝{1,2，…,NP}；

5)对种群中的每个个体C_i进行如下操作：

5.1)将C_i设为目标个体

从以下变异策略池5.1.1)～5.1.3)中等概率的随机选择一个变异策略指导构象进行变异操作：

5.1.1)先从种群中选择一个能量最低的个体C_best，再从种群中随机选出两个互不相同的个体C_a、C_b，

分别从C_a和C_b中随机选择一个位置不同的3片段，分别替换C_best对应位置的片段生成变异构象

5.1.2)从种群中随机选出三个互不相同的个体C_a、C_b和C_c，

分别从C_b、C_c中随机选择一个位置不同的9片段，分别替换C_a对应位置的片段生成变异构象

5.1.3)从种群中随机选出两个互不相同的个体C_a和C_b，

分别从C_a中选择一个3片段，并且从C_b中随机选择一个位置不重叠的9片段，分别替换

对应位置的片段生成变异构象

5.2)对生成的变异构象

进行一次片三段组装生成新变异构象

5.3)生成随机数pCR，其中pCR∈(0,1)，若pCR＜CR，则从

中随机选择一个3片段，替换到

对应的位置的片段生成测试构象

否则直接把

记为

5.4)根据种群个体的索引对种群中索引为

的个体进行步骤5.4.1)，对索引为

的个体进行步骤5.4.2)～5.4.3)：

5.4.1)用Rosetta能量函数score3计算

的能量

如果

则

替换

如果

则按照概率

以蒙特卡洛准则接收构象，若接收构象，则

替换

5.4.2)从当前种群中选择一个能量最低的个体C_best，用公式(1)分别计算

与种群中能量最低个体之间的排挤度

和

排挤度Exclusion(C)计算如下，其中，L是蛋白质序列的残基数，

和

分别是个体C和个体C_best第k个残基CA原子的三维坐标：

5.4.3)如果

则

替换

否则按照概率

以玻尔兹曼准则接收构象，若接收构象，则

替换

6)g＝g+1，迭代运行步骤5)～6)，至g>G为止；

7)输出能量最低的构象为最终结果。

本发明的技术构思为：在进化算法框架下，首先，建立三种不同的变异策略，根据轮盘赌的方式选择一个变异策略对构象进行变异，并且对生成的变异构象进行一次3片段组装，生成变异构象；其次对变异构象进行交叉操作；最后用Rosetta能量函数score3、蒙特卡洛玻尔兹曼接收准则对索引为

的构象进行选择，用排挤度指标Exclusion、蒙特卡洛玻尔兹曼接收准则对索引为

的构象进行选择，结合排挤度评价的多变异策略蛋白质结构预测方法不仅可以增强种群的多样性，而且能够缓解能量函数不精确的问题，提高采样效率。

本发明的有益效果为：根据轮盘赌策略选用不同的变异策略指导构象变异，不仅能够提高种群的多样性，而且能够缓解传统进化算法采样效率低的问题；使用能量函数与排挤度指标分别对构象选择，改善选择效率，同时缓解了能量函数不精确导致的预测误差的问题，进而提高预测精度。

附图说明

图1是一种结合排挤度评价的多变异策略蛋白质结构预测方法对蛋白质1abv采样得到的构象分布图。

图2是一种结合排挤度评价的多变异策略蛋白质结构预测方法对蛋白质1abv采样时的构象更新示意图。

图3是一种结合排挤度评价的多变异策略蛋白质结构预测方法对蛋白质1abv结构预测得到的三维结构。

具体实施方式

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

参照图1～图3，一种结合排挤度评价的多变异策略蛋白质结构预测方法，所述预测方法包括以下步骤：

1)给定目标蛋白的序列信息；

2)根据目标蛋白序列从ROBETTA服务器(http://www.robetta.org/)上得到片段库文件，其中包括3片段库文件和9片段库文件；

3)设置参数：种群大小NP，算法的最大迭代代数G，交叉因子CR，温度因子β，置迭代代数g＝0；

4)种群初始化：随机片段组装生成NP个初始构象C_i，i＝{1,2，…,NP}；

5)对种群中的每个个体C_i进行如下操作：

5.1)将C_i设为目标个体

从以下变异策略池5.1.1)～5.1.3)中等概率的随机选择一个变异策略指导构象进行变异操作：

5.1.1)先从种群中选择一个能量最低的个体C_best，再从种群中随机选出两个互不相同的个体C_a、C_b，

分别从C_a和C_b中随机选择一个位置不同的3片段，分别替换C_best对应位置的片段生成变异构象

5.1.2)从种群中随机选出三个互不相同的个体C_a、C_b和C_c，

分别从C_b、C_c中随机选择一个位置不同的9片段，分别替换C_a对应位置的片段生成变异构象

5.1.3)从种群中随机选出两个互不相同的个体C_a和C_b，

分别从C_a中选择一个3片段，并且从C_b中随机选择一个位置不重叠的9片段，分别替换

对应位置的片段生成变异构象

5.2)对生成的变异构象

进行一次片三段组装生成新变异构象

5.3)生成随机数pCR，其中pCR∈(0,1)，若pCR＜CR，则从

中随机选择一个3片段，替换到

对应的位置的片段生成测试构象

否则直接把

记为

5.4)根据种群个体的索引对种群中索引为

的个体进行步骤5.4.1)，对索引为

的个体进行步骤5.4.2)～5.4.3)：

5.4.1)用Rosetta能量函数score3计算

的能量

如果

则

替换

如果

则按照概率

以蒙特卡洛准则接收构象，若接收构象，则

替换

5.4.2)从当前种群中选择一个能量最低的个体C_best，用公式(1)分别计算

与种群中能量最低个体之间的排挤度

和

排挤度Exclusion(C)计算如下，其中，L是蛋白质序列的残基数，

和

分别是个体C和个体C_best第k个残基CA原子的三维坐标：

5.4.3)如果

则

替换

否则按照概率

以玻尔兹曼准则接收构象，若接收构象，则

替换

6)g＝g+1，迭代运行步骤5)～6)，至g>G为止；

7)输出能量最低的构象为最终结果。

以序列长度为103的α蛋白1abv为实例，一种结合排挤度评价的多变异策略蛋白质结构预测方法，包括以下步骤：

1)给定目标蛋白的序列信息；

2)根据目标蛋白序列从ROBETTA服务器(http://www.robetta.org/)上得到片段库文件，其中包括3片段库文件和9片段库文件；

3)设置参数：种群大小NP＝100，算法的最大迭代代数G＝3000，交叉因子CR＝0.5，温度因子β＝2，置迭代代数g＝0；

4)种群初始化：随机片段组装生成NP个初始构象C_i，i＝{1,2，…,NP}；

5)对种群中的每个个体C_i进行如下操作：

5.1)将C_i设为目标个体

从以下变异策略池5.1.1)～5.1.3)中等概率的随机选择一个变异策略指导构象进行变异操作：

5.1.1)先从种群中选择一个能量最低的个体C_best，再从种群中随机选出两个互不相同的个体C_a、C_b，

分别从C_a和C_b中随机选择一个位置不同的3片段，分别替换C_best对应位置的片段生成变异构象

5.1.2)从种群中随机选出三个互不相同的个体C_a、C_b和C_c，

分别从C_b、C_c中随机选择一个位置不同的9片段，分别替换C_a对应位置的片段生成变异构象

5.1.3)从种群中随机选出两个互不相同的个体C_a和C_b，

分别从C_a中选择一个3片段，并且从C_b中随机选择一个位置不重叠的9片段，分别替换

对应位置的片段生成变异构象

5.2)对生成的变异构象

进行一次片三段组装生成新变异构象

5.3)生成随机数pCR，其中pCR∈(0,1)，若pCR＜CR，则从

中随机选择一个3片段，替换到

对应的位置的片段生成测试构象

否则直接把

记为

5.4)根据种群个体的索引对种群中索引为

的个体进行步骤5.4.1)，对索引为

的个体进行步骤5.4.2)～5.4.3)：

5.4.1)用Rosetta能量函数score3计算

的能量

如果

则

替换

如果

则按照概率

以蒙特卡洛准则接收构象，若接收构象，则

替换

5.4.2)从当前种群中选择一个能量最低的个体C_best，用公式(1)分别计算

与种群中能量最低个体之间的排挤度

和

排挤度Exclusion(C)计算如下，其中，L是蛋白质序列的残基数，

和

分别是个体C和个体C_best第k个残基CA原子的三维坐标：

5.4.3)如果

则

替换

否则按照概率

以玻尔兹曼准则接收构象，若接收构象，则

替换

6)g＝g+1，迭代运行步骤5)～6)，至g>G为止；

7)输出能量最低的构象为最终结果。

以序列长度为103的α蛋白1abv为实施例，运用以上方法得到了该蛋白质的近天然态构象，运行1000代所得到的结构与天然态结构之间的平均均方根偏差为

最小均方根偏差为

预测得到的三维结构如图3所示。

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

Claims (1)

1.一种结合排挤度评价的多变异策略蛋白质结构预测方法，其特征在于：所述方法包括以下步骤：

1)给定目标蛋白的序列信息；

2)根据目标蛋白序列从ROBETTA服务器上得到片段库文件，其中包括3片段库文件和9片段库文件；

3)设置参数：种群大小NP，算法的最大迭代代数G，交叉因子CR，温度因子β，置迭代代数g＝0；

4)种群初始化：随机片段组装生成NP个初始构象C_i，i＝{1,2，…,NP}；

5)对种群中的每个个体C_i进行如下操作：

5.1)将C_i设为目标个体

从以下变异策略池5.1.1)～5.1.3)中等概率的随机选择一个变异策略指导构象进行变异操作：

5.1.1)先从种群中选择一个能量最低的个体C_best，再从种群中随机选出两个互不相同的个体C_a、C_b，

分别从C_a和C_b中随机选择一个位置不同的3片段，分别替换C_best对应位置的片段生成变异构象

5.1.2)从种群中随机选出三个互不相同的个体C_a、C_b和C_c，

分别从C_b、C_c中随机选择一个位置不同的9片段，分别替换C_a对应位置的片段生成变异构象

5.1.3)从种群中随机选出两个互不相同的个体C_a和C_b，

分别从C_a中选择一个3片段，并且从C_b中随机选择一个位置不重叠的9片段，分别替换

对应位置的片段生成变异构象

5.2)对生成的变异构象

进行一次3片段组装生成新变异构象

5.3)生成随机数pCR，其中pCR∈(0,1)，若pCR＜CR，则从

中随机选择一个3片段，替换到

对应的位置的片段生成测试构象

否则直接把

记为

5.4)根据种群个体的索引对种群中索引为

的个体进行步骤5.4.1)，对索引为

的个体进行步骤5.4.2)～5.4.3)：

5.4.1)用Rosetta能量函数score3计算

的能量

如果

则

替换

如果

则按照概率

以蒙特卡洛准则接收构象，若接收构象，则

替换

5.4.2)从当前种群中选择一个能量最低的个体C_best，用公式(1)分别计算

与种群中能量最低个体的排挤度

和

排挤度Exclusion(C)计算如下，其中，L是蛋白质序列的残基数，

和

分别是个体C和个体C_best第k个残基CA原子的三维坐标：

5.4.3)如果

则

替换

否则按照概率

以玻尔兹曼准则接收构象，若接收构象，则

替换

6)g＝g+1，迭代运行步骤5)～6)，至g>G为止；

7)输出能量最低的构象为最终结果。

Images