CN105205348B - 一种基于距离约束选择策略的群体构象空间优化方法 - Google Patents

Abstract

一种基于距离约束选择策略的群体构象空间优化方法，包括以下步骤：在差分进化算法框架下，采用Rosetta Score3粗粒度知识能量模型来有效降低构象空间搜索维数、提高算法的收敛速度；引入基于知识的片段组装技术可以有效提高预测精度；利用距离谱在蛋白质结构上的特性约束，使更多结构上相似的蛋白质构象进入种群，减少能量函数不精确带来的误差，以得到更为优良的局部构象；结合差分进化算法较强的全局搜索能力，可以对构象空间进行更为有效的采样，从而得到和天然态结构相似度较高的构象。本发明采用Rosetta Score3粗粒度知识能量模型来有效降低构象空间搜索维数、提高算法的收敛速度，提高预测精度。

Description

一种基于距离约束选择策略的群体构象空间优化方法

技术领域

本发明涉及生物信息学、计算机应用领域，尤其涉及的是一种基于距离约束选择策略的群体构象空间优化方法。

背景技术

蛋白质分子在生物细胞化学反应过程中起着至关重要的作用。它们的结构模型和生物活性状态对我们理解和治愈多种疾病有重要的意义。蛋白质只有折叠成特定的三维结构才能产生其特有的生物学功能。因此，要了解蛋白质的功能，就必须获得其三维空间结构。

蛋白质三级结构预测是生物信息学的一个重要任务。蛋白质构象优化问题现在面临最大的挑战是对极其复杂的蛋白质能量函数曲面进行搜索。蛋白质能量模型考虑了分子体系成键作用以及范德华力、静电、氢键、疏水等非成键作用，致使其形成的能量曲面极其粗糙，构象对应局部极小解数目随序列长度的增加呈指数增长。而蛋白质构象预测算法能够找到蛋白质稳定结构的机理是，大量的蛋白质亚稳定结构构成了低能量区域，所以能否找到蛋白质全局最稳定结构的关键是算法能够找到大量的蛋白质亚稳定结构，即增加算法的种群多样性。因此，针对更加精确的蛋白质力场模型，选取有效的构象空间优化算法，使新的蛋白质结构预测算法更具有普遍性和高效性成为生物信息学中蛋白质结构预测的焦点问题。

目前，蛋白质结构预测方法大致可以分为两类，基于模板的方法和不基于模板的方法。其中，不基于模板的从头预测(Ab-inito)方法应用最为广泛。它适用于同源性小于25％的大多数蛋白质，仅从序列产生全新结构，对蛋白质分子设计及蛋白质折叠的研究等具有重要意义。当前有以下几种比较成功的从头预测方法：张阳与Jeffrey Skolnick合作的TASSER(Threading/Assembly/Refinement)方法、David Baker及团队设计的Rosetta方法、Shehu等设计的FeLTr方法等。但是到目前还没有一种十分完善的方法来预测蛋白质的三维结构，即使获得了很好的预测结果，但也只是针对某些蛋白质而言的，目前主要的技术瓶颈在于两个方面，第一方面在于采样方法，现有技术对构象空间采样能力不强，另一方面在于构象更新方法，现有技术对构象的更新精度仍然不足。

因此，现有的构象空间搜索方法存在不足，需要改进。

发明内容

为了克服现有群体构象空间优化方法的构象空间搜索维数较高、收敛速度较慢、预测精度较低的不足，本发明在差分进化算法框架下，采用Rosetta Score3粗粒度知识能量模型来有效降低构象空间搜索维数、提高算法的收敛速度，提高预测精度，提供一种基于距离约束选择策略的群体构象空间优化方法。

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

一种基于距离约束选择策略的群体构象空间优化方法，所述搜索方法包括以下步骤：

1)给定输入序列信息；

2)根据序列信息从QUARK服务器(http://zhanglab.ccmb.med.umich.edu/QUARK/)上得到距离谱文件profile，rp_i为距离谱中记录下的残基对，D_i为该残基对之间的距离，其中i∈(1,N)，N为距离谱中残基对数量；

3)设置系统参数：种群大小NP，算法的迭代次数T，交叉因子CR；

4)种群初始化：由输入序列产生NP个种群个体P_init，k其中k∈(1,NP)，对每个个体全部位置片段组装；

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

5.1)设i＝1，其中i∈{1,2,3,…,NP}；令P_target＝P_i，其中i为序号，P_target表示目标个体；

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

5.3)针对个体P_j做变异操作，其中：j＝rand1；

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

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

5.4)针对测试个体P_trail做交叉操作；

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

5.5)针对目标个体P_target和测试个体P_trail进行选择操作；

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

5.5.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入5.6)；否则继续执行5.5.3)；

5.5.3)计算目标个体P_target和测试个体P_trail在残基对rp_i的距离相对于距离谱中D_i之间的距离差▽targetD_i和▽trailD_i；

5.5.4)判断▽targetD_i和▽trailD_i的值，若大于6，则令其等于6，去除一些距离差值较大的异常数值；

5.5.5)计算距离谱中所有位置的目标个体和测试个体的距离差值总和

5.5.6)若▽D_trail＞▽D_target，则进入5.6)；

若▽D_trail＜▽D_target，则产生一个(0,1)之间的随机数rand6，若rand6小于0.2，

则用P_trail替换P_target，否则进行5.6)；

5.6)i＝i+1；判断i是否大于等于NP，是则进入5.7)，否则返回到5.1)；

5.7)迭代的运行5.1)～5.6)，至满足终止条件。

本发明的技术构思为：在差分进化算法框架下，采用Rosetta Score3粗粒度知识能量模型来有效降低构象空间搜索维数、提高算法的收敛速度；引入基于知识的片段组装技术可以有效提高预测精度；利用距离谱在蛋白质结构上的特性约束，使更多结构上相似的蛋白质构象进入种群，减少能量函数不精确带来的误差，以得到更为优良的局部构象；结合差分进化算法较强的全局搜索能力，可以对构象空间进行更为有效的采样，从而得到和天然态结构相似度较高的构象。

本发明的有益效果为：有效降低构象空间搜索维数、提高算法的收敛速度；提高预测精度。

附图说明

图1是蛋白质1GYZ构象系综中构象更新示意图。

具体实施方式

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

参照图1，一种基于距离约束选择策略的群体构象空间优化方法，包括以下步骤：

1)给定输入序列信息；

2)根据序列信息从QUARK服务器(http://zhanglab.ccmb.med.umich.edu/QUARK/)上得到距离谱文件profile，rp_i为距离谱中记录下的残基对，D_i为该残基对之间的距离，其中i∈(1,N)，N为距离谱中残基对数量；

3)设置系统参数：种群大小NP，算法的迭代次数T，交叉因子CR；

4)种群初始化：由输入序列产生NP个种群个体P_init，k其中k∈(1,NP)，对每个个体全部位置片段组装；

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

5.1)设i＝1，其中i∈{1,2,3,…,NP}；令P_target＝P_i，其中i为序号，P_target表示目标个体；

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

5.3)针对个体P_j做变异操作，其中：j＝rand1；

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

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

5.4)针对测试个体P_trail做交叉操作；

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

5.5)针对目标个体P_target和测试个体P_trail进行选择操作；

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

5.5.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入5.6)；否则继续执行5.5.3)；

5.5.3)计算目标个体P_target和测试个体P_trail在残基对rp_i的距离相对于距离谱中D_i之间的距离差▽targetD_i和▽trailD_i；

5.5.4)判断▽targetD_i和▽trailD_i的值，若大于6，则令其等于6，去除一些距离差值较大的异常数值；

5.5.5)计算距离谱中所有位置的目标个体和测试个体的距离差值总和

5.5.6)若▽D_trail＞▽D_target，则进入5.6)；

若▽D_trail＜▽D_target，则产生一个(0,1)之间的随机数rand6，若rand6小于0.2，

则用P_trail替换P_target，否则进行5.6)；

5.6)i＝i+1；判断i是否大于等于NP，是则进入5.7)，否则返回到5.1)；

5.7)迭代的运行5.1)～5.6)，至满足终止条件。

本实施例以序列长度为60的蛋白质1GYZ为实施例，一种基于距离约束选择策略的群体构象空间优化方法，其中包含以下步骤：

1)给定输入序列信息；

2)根据序列信息从QUARK服务器(http://zhanglab.ccmb.med.umich.edu/QUARK/)上得到距离谱文件profile，rp_i为距离谱中记录下的残基对，D_i为该残基对之间的距离，其中i∈(1,N)，N＝62为距离谱中残基对数量；

3)设置系统参数：种群大小NP＝30，算法的迭代次数T＝10000，交叉因子CR＝0.5；4)种群初始化：由输入序列产生NP个种群个体P_init，k其中k∈(1,NP)，对每个个体全部位置片段组装；

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

5.1)设i＝1，其中i∈{1,2,3,…,NP}；令P_target＝P_i，其中i为序号，P_target表示目标个体；

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

5.3)针对个体P_j做变异操作，其中：j＝rand1；

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

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

5.4)针对测试个体P_trail做交叉操作；

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

5.5)针对目标个体P_target和测试个体P_trail进行选择操作；

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

5.5.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入5.6)；否则继续执行5.5.3)；

5.5.3)计算目标个体P_target和测试个体P_trail在残基对rp_i的距离相对于距离谱中D_i之间的距离差▽targetD_i和▽trailD_i；

5.5.4)判断▽targetD_i和▽trailD_i的值，若大于6，则令其等于6，去除一些距离差值较大的异常数值；

5.5.5)计算距离谱中所有位置的目标个体和测试个体的距离差值总和

5.5.6)若▽D_trail＞▽D_target，则进入5.6)；

若▽D_trail＜▽D_target，则产生一个(0,1)之间的随机数rand6，若rand6小于0.2，

则用P_trail替换P_target，否则进行5.6)；

5.6)i＝i+1；判断i是否大于等于NP，是则进入5.7)，否则返回到5.1)；

5.7)迭代的运行5.1)～5.6)，至终止条件。

以序列长度为60的蛋白质1GYZ为实施例，运用以上方法得到了该蛋白质的近天然态构象，构象系综中构象更新图如图1所示。

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

Claims (1)

1.一种基于距离约束的群体构象空间优化方法，其特征在于：所述构象空间搜索方法包括以下步骤：

1)给定输入序列信息；

2)根据序列信息从QUARK服务器上得到距离谱文件profile，rp_n为距离谱中记录下的残基对，D_n为该残基对之间的距离，其中n∈(1,N)，N为距离谱中残基对数量；

3)设置系统参数：种群大小NP，算法的迭代次数T，交叉因子CR；

4)种群初始化：由输入序列产生NP个种群个体P_init，k其中k∈(1,NP)，对每个个体全部位置片段组装；

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

5.1)设i＝1，其中i∈{1,2,3,…,NP}；令P_target＝P_i，其中i为序号，P_target表示目标个体；

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

5.3)针对个体P_j做变异操作，其中：j＝rand1；

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

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

5.4)针对测试个体P_trail做交叉操作；

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

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

5.5)针对目标个体P_target和测试个体P_trail进行选择操作；

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

5.5.2)若E(P_target)>E(P_trail)则用P_trail替换P_target，进入5.6)；否则继续执行5.5.3)；

5.5.3)计算目标个体P_target和测试个体P_trail在残基对rp_n的距离相对于距离谱中D_n之间的距离差▽targetD_n和▽trailD_n；

5.5.4)判断▽targetD_n和▽trailD_n的值，若大于6，则令其等于6，去除一些距离差值较大的异常数值；

5.5.5)计算距离谱中所有位置的目标个体和测试个体的距离差值总和

5.5.6)若▽D_trail＞▽D_target，则进入5.6)；

若▽D_trail＜▽D_target，则产生一个(0,1)之间的随机数rand6，若rand6小于0.2，则用P_trail替换P_target，否则进行5.6)；

5.6)i＝i+1；判断i是否大于等于NP，是则进入5.7)，否则返回到5.1)；

5.7)迭代的运行5.1)～5.6)，至满足终止条件。