CN111815036B - 一种基于多残基接触图协同约束的蛋白质结构预测方法 - Google Patents

Abstract

一种基于多残基接触图协同约束的蛋白质结构预测方法，基于Rosetta的框架，首先利用Rosetta第一、第二阶段对种群初始化，然后通过对目标构象进行变异、交叉生成新的测试构象；其次，根据四个contact服务器预测得到的残基接触图，设计基于残基接触图的余弦相似度指标来辅助Rosetta能量函数score3对构象进行更新，从而引导算法采样得到能量更低且结构更紧凑的构象。本发明提供一种预测精度较高的基于多残基接触图协同约束的蛋白质结构预测方法。

Description

一种基于多残基接触图协同约束的蛋白质结构预测方法

技术领域

本发明涉及生物信息学、计算机应用领域，尤其涉及的是一种基于多残基接触图协同约束的蛋白质结构预测方法。

背景技术

蛋白质结构预测是结构生物信息学的主要研究内容，也是分子生物学中心法则尚未解决的一个重大基础性科学研究课题。2018年12月初，在墨西哥坎昆举行的全球蛋白质结构预测竞赛(CASP13)中，谷歌旗下DeepMind团队研发的AlphaFold获得总分第一名。AlphaFold使得蛋白质结构预测这一前沿性基础研究问题从科学的殿堂进入民众的视野，成为目前的一个“热议”方向，注定其会成为结构生物信息学发展历程中的一个重要里程碑；该工作同时也表明了计算机技术、信息技术与生命科学领域的深度交叉融合，将会有效地驱动和加速科学新发现。

蛋白质结构预测的重要性源于目前实验测定方法存在的局限性。X射线晶体衍射是目前测定蛋白质结构最有效的方法，所能达到的精度是其它方法所不能比拟的，主要缺点是蛋白质晶体难以培养且晶体结构测定的周期较长；多维核磁共振(NMR)方法可以直接测定蛋白质在溶液中的构象，但是对样品的需要量大、纯度要求高，目前只能测定小分子蛋白质。对于药物靶标—膜蛋白而言，现有实验测定技术极难获得其三维结构；

蛋白质只有折叠成特定的三维结构才能产生其特有的生物学功能。因此，了解蛋白质的三维结构(天然态结构)是了解蛋白质的生物功能的关键。蛋白质三维结构可以通过核磁共振和X光晶体衍射等实验方法得到，然而这些实验测定方法不仅耗时而且极其昂贵，对于某些不易结晶的蛋白质来说不适用。因此，根据的Anfinsen的热力学假说(具有最低能量的构象被认为是天然态结构)，很多计算算法被提出用于蛋白质结构预测。

在理论探索和应用需求的双重推动下，根据Anfinsen法则，利用计算机预测蛋白质结构技术在20世纪末得到了蓬勃发展。1994年美国马里兰大学的科学家Moult发起的CASP竞赛是一个世界性的蛋白质结构预测评比活动，客观地反映了当前蛋白质结构预测领域发展的最新技术水平，被誉为蛋白质结构预测的奥林匹克竞赛。该比赛旨在吸引计算机科学、生物物理学等不同领域的专家参与到蛋白质三维结构预测这一极具挑战性的生物信息学问题中来，共同评估发展现状和讨论未来的趋势。

通过计算技术进行蛋白结构预测通常用一个很复杂的能量函数来评估，其能量函数曲面通常具有成千上万的自由度和大量的局部最优解，构象搜索空间极大。为了进行构象空间搜索，从头预测方法通常首先根据基于知识的粗粒度能量模型获得构象空间的全局最小解，然后对其对应的构象进行精修，从而得到预测结构。因此，从头预测方法需要解决两个问题：1.建立合适的能量函数来评价构象的合理性；2.提出有效的构象空间搜索方法来搜索全局最优解。第一个因素本质上属于分子力学问题，主要是为了能够计算得到每个蛋白质结构对应的能量值。第二个因素本质上属于全局优化问题，通过选择一种合适的优化方法，对构象空间进行快速搜索，得到与某一全局最小能量对应的构象。

差分进化算法(DE)由于其结构简单，易于实现，鲁棒性强和收敛速度快等优点已被成功应用于蛋白结构预测。然而，随着氨基酸序列的增长，蛋白质分子体系自由度也增大，利用传统群体算法采样获得大规模蛋白质构象空间的全局最优解成为一项挑战性的工作；其次，粗粒度模型虽然减小了构象搜索空间，但是也导致了相互作用力之间的信息丢失，从而直接影响预测精度。

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

发明内容

为了克服现有的蛋白质结构预测方法采样效率低、预测精度较低的不足，本发明基于Rosetta，引入多张残基接触图来指导构象空间优化，提出一种采用效率高、预测精度高的基于多残基接触图协同约束的蛋白质结构预测方法。

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

一种基于多残基接触图协同约束的蛋白质结构预测方法，所述方法包括以下步骤：

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

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

3)根据目标蛋白序列，利用RaptorX-Contact服务器(http://raptorx.uchicago.edu/ContactMap/)、ResPRE服务器(https://zhanglab.ccmb.med.umich.edu/ResPRE/)、DeepMetaPSICOV服务器(http://bioinf.cs.ucl.ac.uk/psipred/)、NeBcon服务器(https://zhanglab.ccmb.med.umich.edu/NeBcon/)得到四张contactmap，分别为contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon；

4)设置参数：种群大小NP，算法的迭代次数G，交叉因子CR，置迭代代数g＝0；

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

6)将种群中的构象个体C_i,i∈{1,2,3,…,NP}看作目标构象个体

执行以下操作生成变异构象

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

6.2)在构象C_n1位置上随机选取一个9片段替换构象C_n3的相同位置所对应的片段，再从构象C_n2位置上随机的选取一个与构象C_n1选取位置不相同的9片段替换偶像C_n3的相同位置所对应的片段，然后用对构象C_n3进行3片段组装生成变异构象个体

7)对变异构象

i∈{1,2,3,…,NP}执行交叉操作生成测试构象

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

7.2)若随机数rand1≤CR，则从目标构象

中随机选择一个3片段替换到变异构象

的相应位置,否则变异构象

不变；

8)对每个目标构象

和测试构象

进行选择操作；

8.1)用Rosetta score3能量函数分别计算

和

的能量：

和

8.2)若

则构象

被拒绝，否则继续执行步骤8.3)；

8.3)首先，把

和

的L×L二维contactmap转化为长度为L×L的一维向量

和

把contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon转化为4个长度为L×L的一维向量

和

其中L为蛋白质序列的长度；然后分别计算

和

与

余弦相似度并求和得到

和

计算方法如下：

8.4)若

则构象

替换构象

且转到步骤9)；

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

10)输出结果。

本发明的技术构思为：基于Rosetta的框架，首先利用Rosetta第一、第二阶段对种群初始化，然后通过对目标构象进行变异、交叉生成新的测试构象；其次，根据四个contact服务器预测得到的残基接触图，设计基于残基接触图的余弦相似度指标来辅助Rosetta能量函数score3对构象进行更新，从而引导算法采样得到能量更低且结构更紧凑的构象。本发明提供一种基于多残基接触图协同约束的蛋白质结构预测方法。

本发明的有益效果为：首先，通过结合不同服务器预测的残基接触图信息，缓解了单一残基接触图所存在的召回率、精确度不足的问题；其次，设计了基于残基接触图的余弦相似度指标辅助Rosetta能量函数score3对构象进行更新，从而引导算法采样得到能量更低且结构更紧凑的构象。

附图说明

图1是四张预测的残基接触图。

图2是基于多残基接触图协同约束的蛋白质结构预测方法对蛋白质1TEN采样得到的构象分布图。

图3是基于多残基接触图协同约束的蛋白质结构预测方法对1TEN蛋白结构预测得到的三维结构。

具体实施方式

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

参照图1～图3，一种基于多残基接触图协同约束的蛋白质结构预测方法，所述方法包括以下步骤：

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

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

3)根据目标蛋白序列，利用RaptorX-Contact服务器(http://raptorx.uchicago.edu/ContactMap/)、ResPRE服务器(https://zhanglab.ccmb.med.umich.edu/ResPRE/)、DeepMetaPSICOV服务器(http://bioinf.cs.ucl.ac.uk/psipred/)、NeBcon服务器(https://zhanglab.ccmb.med.umich.edu/NeBcon/)得到四张contactmap，分别为contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon；

4)设置参数：种群大小NP，算法的迭代次数G，交叉因子CR，置迭代代数g＝0；

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

6)将种群中的构象个体C_i,i∈{1,2,3,…,NP}看作目标构象个体

执行以下操作生成变异构象

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

6.2)在构象C_n1位置上随机选取一个9片段替换构象C_n3的相同位置所对应的片段，再从构象C_n2位置上随机的选取一个与构象C_n1选取位置不相同的9片段替换偶像C_n3的相同位置所对应的片段，然后用对构象C_n3进行3片段组装生成变异构象个体

7)对变异构象

i∈{1,2,3,…,NP}执行交叉操作生成测试构象

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

7.2)若随机数rand1≤CR，则从目标构象

中随机选择一个3片段替换到变异构象

的相应位置,否则变异构象

不变；

8)对每个目标构象

和测试构象

进行选择操作；

8.1)用Rosetta score3能量函数分别计算

和

的能量：

和

8.2)若

则构象

被拒绝，否则继续执行步骤8.3)；

8.3)首先，把

和

的L×L二维contactmap转化为长度为L×L的一维向量

和

把contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon转化为4个长度为L×L的一维向量

和

其中L为蛋白质序列的长度；然后分别计算

和

与

余弦相似度并求和得到

和

计算方法如下：

8.4)若

则构象

替换构象

且转到步骤9)；

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

10)输出结果。

以序列长度为87的蛋白质1TEN为实施例，一种基于多残基接触图协同约束的蛋白质结构预测方法，其中包含以下步骤：

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

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

3)根据目标蛋白序列，利用RaptorX-Contact服务器(http://raptorx.uchicago.edu/ContactMap/)、ResPRE服务器(https://zhanglab.ccmb.med.umich.edu/ResPRE/)、DeepMetaPSICOV服务器(http://bioinf.cs.ucl.ac.uk/psipred/)、NeBcon服务器(https://zhanglab.ccmb.med.umich.edu/NeBcon/)得到四张contactmap，分别为contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon；

4)设置参数：种群大小NP＝100，算法的迭代次数G＝300，交叉因子CR＝0.5，置迭代代数g＝0；

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

6)将种群中的构象个体C_i,i∈{1,2,3,…,NP}看作目标构象个体

执行以下操作生成变异构象

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

6.2)在构象C_n1位置上随机选取一个9片段替换构象C_n3的相同位置所对应的片段，再从构象C_n2位置上随机的选取一个与构象C_n1选取位置不相同的9片段替换偶像C_n3的相同位置所对应的片段，然后用对构象C_n3进行3片段组装生成变异构象个体

7)对变异构象

i∈{1,2,3,…,NP}执行交叉操作生成测试构象

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

7.2)若随机数rand1≤CR，则从目标构象

中随机选择一个3片段替换到变异构象

的相应位置,否则变异构象

不变；

8)对每个目标构象

和测试构象

进行选择操作；

8.1)用Rosetta score3能量函数分别计算

和

的能量：

和

8.2)若

则构象

被拒绝，否则继续执行步骤8.3)；

8.3)首先，把

和

的L×L二维contactmap转化为长度为L×L的一维向量

和

把contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon转化为4个长度为L×L的一维向量

和

其中L为蛋白质序列的长度；然后分别计算

和

与

余弦相似度并求和得到

和

计算方法如下：

8.4)若

则构象

替换构象

且转到步骤9)；

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

10)输出结果。

以序列长度为87的蛋白质1TEN为实施例，运用以上方法得到了该蛋白质的近天然态构象，运行300代所得到的结构与天然态结构之间的平均均方根偏差为2.86，最小均方根偏差为2.01，预测得到的三维结构如图3所示。

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

Claims (1)

1.一种基于多残基接触图协同约束的蛋白质结构预测方法，其特征在于：所述方法包括以下步骤：

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

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

3)根据目标蛋白序列，利用RaptorX-Contact服务器、ResPRE服务器、DeepMetaPSICOV服务器、NeBcon服务器得到四张contactmap，分别为contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon；

4)设置参数：种群大小NP，算法的迭代次数G，交叉因子CR，置迭代代数g＝0；

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

6)将种群中的构象个体C_i,i∈{1,2,3,…,NP}看作目标构象个体

对每个

执行以下操作生成变异构象

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

6.2)在构象C_n1位置上随机选取一个9片段替换构象C_n3的相同位置所对应的片段，再从构象C_n2位置上随机的选取一个与构象C_n1选取位置不相同的9片段替换构象C_n3的相同位置所对应的片段，然后对构象C_n3进行3片段组装生成变异构象个体

7)对变异构象

i∈{1,2,3,…,NP}执行交叉操作生成测试构象

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

7.2)若随机数rand1≤CR，则从目标构象

中随机选择一个3片段替换到变异构象

的相应位置,否则变异构象

不变；

8)对每个目标构象

和测试构象

进行选择操作；

8.1)用Rosetta score3能量函数分别计算

和

的能量：

和

8.2)若

则构象

被拒绝，否则继续执行步骤8.3)；

8.3)首先，把

和

的L×L二维contactmap分别转化为长度为L×L的一维向量

和

把contactRaptorX、contactResPRE、contactDeepMetaPSICOV和contactNeBcon转化为4个长度为L×L的一维向量

和

其中L为蛋白质序列的长度；然后分别计算

和

与

余弦相似度并求和得到

和

计算方法如下：

8.4)若

则构象

替换构象

且转到步骤9)；

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

10)输出结果。

Images