一种基于残基信息动态选择策略的蛋白质结构预测方法与流程

本发明涉及一种生物学信息学、智能优化、计算机应用领域，尤其涉及的是，一种基于残基信息动态选择策略的蛋白质结构预测方法。

背景技术：

生物细胞中包含许多由20多种氨基酸所形成的长链折叠而成的蛋白质，蛋白质结构预测问题是当今计算生物学领域中的研究热点，不仅具有非常重要的理论机制，而且对新蛋白的设计、蛋白质之间相互作用建模、药物标靶蛋白的设计具有十分重要的指导意义。目前，最常用的蛋白质结构实验测定方法包括x-晶体衍射和核磁共振，但是上述两种蛋白质结构测定方法存在一定的缺陷，无法满足所有蛋白质结构预测的需求。例如，对于一些不易结晶的蛋白，无法使用x-晶体衍射方法进行测定；使用核磁共振测定一个蛋白的结构费用较大(大约15万美元)，并且极其费时(大约需要半年)。因此，如何以计算机为工具，运用适当的算法，从氨基酸序列出发直接预测蛋白质的三维结构，成为当前生物信息学中一种重要的研究课题。

1965年，anfinsen等提出了蛋白质的一级结构完全决定了其三维空间结构的著名论断，这一论断使得根据蛋白质得氨基酸序列从理论上预测其相应的空间结构成为了现实。从热力学角度来看，也就是说蛋白质的天然态结构通常对应着能量最低的构象。因此，从氨基酸序列出发，以能量函数来衡量构象的质量，采用适当的算法搜索能量模型中最小能量构象，从而预测蛋白质的天然态结构，已成为计算生物信息学中重要的研究课题之一。目前，最常用的预测方法可以划分为以下三类：针对高相似序列的同源建模方法；针对较低相似性序列的折叠识别方法；以及不依赖模板的从头预测方法。

从头预测方法直接基于anfinsen假说建立蛋白质物理或知识能量模型，然后设计适当优化算法求解最小能量构象。可以看出，从头预测方法必须考虑以下两个因素：(1)蛋白质结构能量函数；(2)构象空间搜索方法。第一个因素本质上属于分子力学问题，主要是为了能够计算得到每个蛋白质结构对应的能量值。第二个因素本质上属于全局优化问题，通过选择一种合适的优化方法，对构象空间进行快速搜索，得到能量最低构象。其中，能量函数值的不精确是导致蛋白质结构预测精度的关键因素，能量函数的不精确会导致搜索过程中结构合理的构象丢失，从而影响预测精度。

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

技术实现要素：

为了克服现有的蛋白质结构预测方法在预测精度方面的不足，本发明提供一种预测精度较高的基于残基信息动态选择策略的蛋白质结构预测方法。

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

一种基于残基信息动态选择策略的蛋白质结构预测方法，所述预测方法包括以下步骤：

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

2)根据序列信息从quark服务器(https://zhanglab.ccmb.med.umich.edu/quark/)上得到距离谱文件；

3)根据序列信息从raptorx(http://raptorx.uchicago.edu/contactmap/)服务器上得到残基接触信息文件；

4)设置参数：种群大小np，算法的最大迭代代数g，交叉因子cr，选择因子sr，温度因子β，接触信息因子βcon，置迭代代数g＝0；

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

6)对种群中的每个构象个体ci，i＝{1,2，…,np}看作目标构象个体进行如下操作生成变异构象

6.1)在1到np范围内随机生成三个正整数n1,n2,n3,且n1≠n2≠n3≠i；

6.2)在构象cn1上随机选择一个loop区的片段替换cn3相同位置的所对应的片段，然后在cn2上随机选择一个不同于cn1上选则的loop区的片段替换cn3相同位置的所对应的片段，然后对cn3进行一次片段组装生成变异构象

7)对每个变异构象i＝{1,2，…,np}执行交叉操作生成交叉操作生成测试构象过程如下：

7.1)生成随机数rand1，其中rand1∈(0,1)；

7.2)若随机数rand1≤cr，则从目标构象中随机选择一个loop区的片段替换变异构象相同位置所对应的片段，否则变异构象保持不变；

8)对每个目标构象和测试构象进行如下操作；

8.1)用rosettascore3能量函数分别计算与的能量：和

8.2)若则构象替换否则产生一个随机数rand2，其中rand2∈(0,1)，若rand2<sr，则进行步骤8.3)；否则，进行步骤8.6)；

8.3)根据距离谱中的残基对分别计算和的残基间距离和然后根据公式(1)和(2)分别计算和的距离误差系数dtrial和dtarget，其中，t表示距离谱中残基对的数量，和分别表示和构象中第t对残基m,n对应cα原子之间的距离，表示距离谱中第t对残基m,n对应的距离谱的峰值距离；

8.4)若dtrial<dtarget，则构象替换构象接收次数count1加1，否则进行步骤8.5)；

8.5)计算目标构象和测试构象的距离误差系数差▽d＝dtrial-dtarget，按照概率以蒙特卡洛准则接受构象若接受，则count1加1，其中β为温度因子；

8.6)根据残基接触信息，由公式(3)分别计算和的接触能量contrial和contarget，其中，n表示残基接触信息文件中残基对的总数，k是残基接触信息文件中的第k对残基对，是第p个残基和第q个残基的接触概率，是第p个残基和第q个残基的cβ原子之间的距离；

8.7)contrial>contarget，则构象替换构象接收次数count2加1，否则进行步骤8.8)；

8.8)计算目标构象和测试构象的距离误差系数差▽con＝contarget-contrial，按照概率以蒙特卡洛准则接受构象若接受，则count2加1，其中βcon为接触信息因子；

9)g＝g+1，sr＝count1/(count1+count2)迭代运行步骤5)～8)，至g>g为止；

10)输出能量最低的构象为最终预测结果。

本发明的技术构思为：在rosetta的基本框架下，引入距离谱和残基接触信息来弥补能量函数精确性的不足。首先，初始化构象，并针对loop区进行交叉、变异来生成新构象；然后，分别利用距离谱和残基接触信息设计两种不同的选择策略，同时根据各选择策略的历史成功信息计算其选择概率；最后，根据选择概率动态自适应地选择策略指导种群更新；迭代上述过程直至满足条件后输出结果。

本发明的有益效果为：在新生成的构象中利用了能量函数、距离谱和残基接触信息动态自适应地对构象进行筛选，从而缓解了能量函数不精确引起的预测误差，提高预测精度。

附图说明

图1是基于残基信息动态选择策略的蛋白质结构预测方法对蛋白质3gwl采样得到的构象分布图。

图2是基于残基信息动态选择策略的蛋白质结构预测方法对蛋白质3gwl采样得到的均方根偏差和能量的散点图。

图3是基于残基信息动态选择策略的蛋白质结构预测方法对3gwl蛋白结构预测得到的三维结构。

具体实施方式

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

参照图1～图3，一种基于残基信息动态选择策略的蛋白质结构预测方法，所述方法包括以下步骤：

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

2)根据序列信息从quark服务器(https://zhanglab.ccmb.med.umich.edu/quark/)上得到距离谱文件；

3)根据序列信息从raptorx(http://raptorx.uchicago.edu/contactmap/)服务器上得到残基接触信息文件；

4)设置参数：种群大小np，算法的最大迭代代数g，交叉因子cr，选择因子sr，温度因子β，接触信息因子βcon，置迭代代数g＝0；

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

6)对种群中的每个构象个体ci，i＝{1,2，…,np}看作目标构象个体进行如下操作生成变异构象

6.1)在1到np范围内随机生成三个正整数n1,n2,n3,且n1≠n2≠n3≠i；

6.2)在构象cn1上随机选择一个loop区的片段替换cn3相同位置的所对应的片段，然后在cn2上随机选择一个不同于cn1上选则的loop区的片段替换cn3相同位置的所对应的片段，然后对cn3进行一次片段组装生成变异构象

7)对每个变异构象i＝{1,2，…,np}执行交叉操作生成交叉操作生成测试构象过程如下：

7.1)生成随机数rand1，其中rand1∈(0,1)；

7.2)若随机数rand1≤cr，则从目标构象中随机选择一个loop区的片段替换变异构象相同位置所对应的片段，否则变异构象保持不变；

8)对每个目标构象和测试构象进行如下操作；

8.1)用rosettascore3能量函数分别计算与的能量：和

8.2)若则构象替换否则产生一个随机数rand2，其中rand2∈(0,1)，若rand2<sr，则进行步骤8.3)；否则，进行步骤8.6)；

8.3)根据距离谱中的残基对分别计算和的残基间距离和然后根据公式(1)和(2)分别计算和的距离误差系数dtrial和dtarget，其中，t表示距离谱中残基对的数量，和分别表示和构象中第t对残基m,n对应cα原子之间的距离，表示距离谱中第t对残基m,n对应的距离谱的峰值距离；

8.4)若dtrial<dtarget，则构象替换构象接收次数count1加1，否则进行步骤8.5)；

8.5)计算目标构象和测试构象的距离误差系数差▽d＝dtrial-dtarget，按照概率以蒙特卡洛准则接受构象若接受，则count1加1，其中β为温度因子；

8.6)根据残基接触信息，由公式(3)分别计算和的接触能量contrial和contarget，其中，n表示残基接触信息文件中残基对的总数，k是残基接触信息文件中的第k对残基对，是第p个残基和第q个残基的接触概率，是第p个残基和第q个残基的cβ原子之间的距离；

8.7)contrial>contarget，则构象替换构象接收次数count2加1，否则进行步骤8.8)；

8.8)计算目标构象和测试构象的距离误差系数差▽con＝contarget-contrial，按照概率以蒙特卡洛准则接受构象若接受，则count2加1，其中βcon为接触信息因子；

9)g＝g+1，sr＝count1/(count1+count2)迭代运行步骤5)～8)，至g>g为止；

10)输出能量最低的构象为最终预测结果。

以序列长度为106的蛋白质3gwl为实例，一种基于残基信息动态选择策略的蛋白质结构预测方法，包括以下步骤：

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

2)根据序列信息从quark服务器(https://zhanglab.ccmb.med.umich.edu/quark/)上得到距离谱文件；

3)根据序列信息从raptorx(http://raptorx.uchicago.edu/contactmap/)服务器上得到残基接触信息文件；

4)设置参数：种群大小np＝50，算法的最大迭代代数g＝1000，交叉因子cr＝0.5，选择因子sr＝0.5，温度因子β＝2，接触信息因子βcon＝2，置迭代代数g＝0；

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

6)对种群中的每个构象个体ci，i＝{1,2，…,np}看作目标构象个体进行如下操作生成变异构象

6.1)在1到np范围内随机生成三个正整数n1,n2,n3,且n1≠n2≠n3≠i；

6.2)在构象cn1上随机选择一个loop区的片段替换cn3相同位置的所对应的片段，然后在cn2上随机选择一个不同于cn1上选则的loop区的片段替换cn3相同位置的所对应的片段，然后对cn3进行一次片段组装生成变异构象

7)对每个变异构象i＝{1,2，…,np}执行交叉操作生成交叉操作生成测试构象过程如下：

7.1)生成随机数rand1，其中rand1∈(0,1)；

7.2)若随机数rand1≤cr，则从目标构象中随机选择一个loop区的片段替换变异构象相同位置所对应的片段，否则变异构象保持不变；

8)对每个目标构象和测试构象进行如下操作；

8.1)用rosettascore3能量函数分别计算与的能量：和

8.2)若则构象替换否则产生一个随机数rand2，其中rand2∈(0,1)，若rand2<sr，则进行步骤8.3)；否则，进行步骤8.6)；

8.3)根据距离谱中的残基对分别计算和的残基间距离和然后根据公式(1)和(2)分别计算和的距离误差系数dtrial和dtarget，其中，t表示距离谱中残基对的数量，和分别表示和构象中第t对残基m,n对应cα原子之间的距离，表示距离谱中第t对残基m,n对应的距离谱的峰值距离；

8.4)若dtrial<dtarget，则构象替换构象接收次数count1加1，否则进行步骤8.5)；

8.5)计算目标构象和测试构象的距离误差系数差▽d＝dtrial-dtarget，按照概率以蒙特卡洛准则接受构象若接受，则count1加1，其中β为温度因子；

8.6)根据残基接触信息，由公式(3)分别计算和的接触能量contrial和contarget，其中，n表示残基接触信息文件中残基对的总数，k是残基接触信息文件中的第k对残基对，是第p个残基和第q个残基的接触概率，是第p个残基和第q个残基的cβ原子之间的距离；

8.7)contrial>contarget，则构象替换构象接收次数count2加1，否则进行步骤8.8)；

8.8)计算目标构象和测试构象的距离误差系数差▽con＝contarget-contrial，按照概率以蒙特卡洛准则接受构象若接受，则count2加1，其中βcon为接触信息因子；

9)g＝g+1，sr＝count1/(count1+count2)迭代运行步骤5)～8)，至g>g为止；

10)输出能量最低的构象为最终预测结果。

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

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