基于蒙特卡洛模拟的蛋白质热力学分析方法技术

本发明专利技术涉及生物计算与蒙特卡洛模拟技术领域，具体公开一种基于蒙特卡洛模拟的蛋白质热力学分析方法，包括步骤A：确定蛋白质能量模型；步骤B：模拟及计算蛋白质系统态密度。采用本发明专利技术提供的方法，可以高效地分析和研究蛋白质折叠的整个热力学过程，进而对蛋白质折叠过程进行探索和研究。

【技术实现步骤摘要】
基于蒙特卡洛模拟的蛋白质热力学分析方法
本专利技术涉及生物计算与蒙特卡洛模拟
，具体公开一种基于蒙特卡洛模拟的蛋白质热力学分析方法。
技术介绍
蛋白质工程是现代生物技术发展的前沿领域，其根本目的是要将天然存在的蛋白质按照人们的设想进行改造，或者根据需要设计出具有某种特殊功能的非天然的新蛋白质。而这种改造和设计的重要基础之一是蛋白质折叠结构的预测。蛋白质折叠结构的形态在很大程度上决定了其可能具有的生物学功能。这就是说，蛋白质的结构与功能之间具有一致性。因此对蛋白质折叠结构的预测和研究在蛋白质工程，医药领域等都有着极其重要的意义。蒙特卡洛模拟是一种被广泛使用的随机搜索算法，通过随机改变蛋白质系统的粒子位置去搜索构象空间，获得系统的态密度，从而可以计算系统的各种热力学性质，进而研究蛋白质折叠的整个热动力学过程。但是，由于蛋白质能量空间复杂，能量值局部最低区间与全局最低区间存在能垒，传统的蒙特卡洛算法在低温下的模拟容易陷入能量局部最低区间，很难跳出来找到全局最优解。由此产生了很多新的蒙特卡洛模拟算法，如多正则系综（multicanonical）蒙特卡洛算法、模拟回火（simulatedtempering）算法、副本交换（ReplicaExchangeMethod）算法和Wang-Landau算法。Wang-Landau算法的基本思想是使用非波尔兹曼分布函数实现在能量区间自由行走的模拟仿真，从而更为广泛地搜索构型空间；同时还可自动获得蛋白质系统的态密度，从而可以计算宽广温度范围内的许多正则系统热力学量，进而研究蛋白质折叠的整个热动力学过程。Wang-Landau由于直观简洁的特点，其应用范围已经扩展到磁系统、液态晶体、液体、原子簇、自旋玻璃模型以及蛋白质折叠中。但Wang-Landau算法在计算精度和速度上还有待进一步的改进。
技术实现思路
本专利技术旨在克服现有技术在蛋白质热力学分析中计算精度和速度不足的缺陷，提供一种基于蒙特卡洛模拟的蛋白质热力学分析方法。本专利技术的技术方案为1、一种基于蒙特卡洛模拟的蛋白质热力学分析方法，其特征在于，包括步骤A：确定蛋白质能量模型；以及B：模拟及计算蛋白质系统态密度。进一步地，步骤A进一步包括：采用ECEPP能量力场模型以及角度坐标系，ECEPP能量力场的表达形式为：EECEPP=EC+ELJ+EHB+ETor其中是两电荷之间的库伦作用力，rij表示原子i和j之间的距离，是两原子之间的兰纳-琼斯作用力，是氢键作用力，ETor=∑lUl(1±cos(nlξl))是两面角旋转作用力，ξl是第l个两面角。进一步地，步骤A进一步包括：对所使用的蛋白质能量区间进行离散化处理，若取k个能量bin区间值，则需对[Emin,Emax]平均划分k个bin区间，用每个区间中间的一个能量值代表该能量区间值。进一步地，步骤B进一步包括：通过主从进程模式的MPI并行程序算法，模拟及计算蛋白质系统的态密度。进一步地，在所述主从进程模式的N个分进程中，分进程1为主进程，其余分进程均为子进程。进一步地，所述主进程包括以下步骤：S11：初始化蛋白质系统态密度函数的对数S(E)=lng(E)=0，直方图H(E)=0（Emin≤E≤Emax），s=1，修正因子df|E=(κΘ(E0-E)+1)lnf，其中Θ(E0-E)是Heaviside分段函数，κ、E0、f是相关于模型的参数；S12：t=1；S13：在主进程中，对原来的构型进行随机变动，产生新的构型，计算能量Enew，根据Metropolis准则确定新构型被接受的概率，t=t+1,s=s+1，修正因子为df|E=(κΘ(E0-E)+1)lnf；步骤S13循环tmax次；S14：所有进程间相互通信，主进程收集所有从进程的Stmp(E)和Htmp(E)并累加计算出全局的S(E)和H(E)，即全局的S(E)=S(E)+所有从进程的Stmp(E)，全局的H(E)=H(E)+所有从进程的Htmp(E)，将全局的S(E)和H(E)的广播给所有从进程，求得Hreal(E)=H(E)×[1+κΘ(E0-E)]，判断直方图平缓条件：（0<φ<1）若不满足则返回到S12继续迭代；若满足则至S15；S15：改变修正因子f，再返回执行S12继续迭代，直到满足条件其中ΔE=Emax-Emin为系统的能量范围；S16：一旦满足条件则对原来的构型进行随机变动，产生新的构型，计算能量Enew，根据Metropolis准则确定新构型被接受的概率，s=s+1，修正因子为S17：继续S16迭代，直到满足进程终止条件通过对所有进程的S(E)求平均得进而得到蛋白质系统的相对的态密度进一步地，所述从进程包括以下步骤：S21：初始化系统的态密度的对数S(E)=lng(E)=0，Stmp(E)=lngtmp(E)=0，直方图H(E)=0，Htmp(E)=0(Emin≤E≤Emax)，s=1，修正因子df|E=(κΘ(E0-E)+1)lnf，其中Θ(E0-E)是Heaviside分段函数，κ、E0、f是相关于模型的参数；S22：t=1；S23：在从进程中，对原来的构型进行随机变动，产生新的构型，计算能量Enew，根据Metropolis准则确定新构型被接受的概率，t=t+1,s=s+1，修正因子为df|E=(κΘ(E0-E)+1)lnf；步骤S23循环tmax次；S24：所有进程间相互通信，从进程将Stmp(E)和Htmp(E)发送给主进程，然后接收经主进程计算得出的全局的S(E)和H(E)更新原来的S(E)和H(E)，将Stmp(E)和Htmp(E)初始化为0，求得Hreal(E)=H(E)×[1+κΘ(E0-E)]，判断直方图平缓条件：其中0<φ<1若不满足则返回执行S22继续迭代；若满足则执行步骤S25；S25改变修正因子f，再返回执行S22继续迭代，直到满足条件其中ΔE=Emax-Emin为系统的能量范围；S26：一旦条件满足，则对原来的构型进行随机变动，产生新的构型，计算能量Enew，根据Metropolis准则确定新构型被接受的概率；S27：继续步骤S26迭代，直到满足进程终止条件其中进一步地，步骤S13中根据Metropolis准则确定新构型被接受的概率进一步包括：若接受新构型，则：S(Enew)=S(Enew)+df|E,H(Enew)=H(Enew)+1；否则：S(Eold)=S(Eold)+df|E,H(Eold)=H(Eold)+1。进一步地，步骤S23中根据Metropolis准则确定新构型被接受的概率进一步包括：...

【技术保护点】
一种基于蒙特卡洛模拟的蛋白质热力学分析方法，其特征在于，包括步骤A：确定蛋白质能量模型；以及B：模拟及计算蛋白质系统态密度。

【技术特征摘要】
1.一种基于蒙特卡洛模拟的蛋白质热力学分析方法，其特征在于，包括步骤A：确定蛋白质能量模型；以及B：模拟及计算蛋白质系统态密度；步骤B进一步包括：通过主从进程模式的MPI并行程序算法，模拟及计算蛋白质系统的态密度；在所述主从进程模式的N个分进程中，分进程1为主进程，其余分进程均为从进程；所述主进程包括步骤：S11：初始化蛋白质系统态密度函数的对数S(E)＝lng(E)＝0，直方图H(E)＝0，Emin≤E≤Emax，s＝1，修正因子df|E＝(κΘ(E0-E)+1)lnf，其中Θ(E0-E)是Heaviside分段函数，κ、E0、f是相关于模型的参数；S12：t＝1；S13：在主进程中，对原来的构型进行随机变动，产生新的构型，计算能量Enew，根据Metropolis准则确定新构型被接受的概率，t＝t+1,s＝s+1，修正因子为df|E＝(κΘ(E0-E)+1)lnf；步骤S13循环tmax次；S14：所有进程间相互通信，主进程收集所有从进程的Stmp(E)和Htmp(E)并累加计算出全局的S(E)和H(E)，即全局的S(E)＝S(E)+所有从进程的Stmp(E)，全局的H(E)＝H(E)+所有从进程的Htmp(E)，将全局的S(E)和H(E)广播给所有从进程，求得Hreal(E)＝H(E)×[1+κΘ(E0-E)]，判断直方图平缓条件：若不满足则返回到S12继续迭代；若满足则至S15；S15：改变修正因子f，再返回执行S12继续迭代，直到满足条件其中ΔE＝Emax-Emin为系统的能量范围；S16：一旦满足条件则对原来的构型进行随机变动，产生新的构型，计算能量Enew，根据Metropolis准则确定新构型被接受的概率，s＝s+1，修正因子为S17：继续S16迭代，直到满足进程终止条件通过对所有进程的S(E)求平均得进而得到蛋白质系统的相对的态密度2.如权利要求1所述的蛋白质热力学分析方法，其中，步骤A进一步包括：采用ECEPP能量力场模型以及角度坐标系，ECEPP能量力场的表达形式为：EECEPP＝EC+ELJ+EHB+ETor其中是两电荷之间的库伦作用力，rij表示原子i和j之间的距离，是两原子之间的兰纳-琼斯作用力，是氢键作用力，ETor＝∑lUl(1±cos(nlξl))是两面角旋转作用力，ξl是第l个两面角。3.如权利要求1所述的蛋白质热力学分析方法，其中，步骤A进一步包括：对所使用的蛋白质能量区间进行离散化处理，若取k个能量bin区间值，则需对[Emin,Emax]平均划分k个bin区间，用每个区间中间的一个能量值代表该区间能量值。4.如权利要求1所述的蛋白质热力学分析方法，其中，所述从进程包括步骤：S21：初始化系统的态密度的对数...

【专利技术属性】
技术研发人员：彭丰斌，魏彦杰，张慧玲，
申请(专利权)人：深圳先进技术研究院，
类型：发明
国别省市：广东;44