蛋白质结构和网格的蒙特卡罗模拟

Posted

技术标签:

【中文标题】蛋白质结构和网格的蒙特卡罗模拟【英文标题】:monte carlo simulation of protein structure and grid 【发布时间】:2013-09-24 08:21:55 【问题描述】:

我正在处理蛋白质结构的蒙特卡罗模拟脚本。在蒙特卡洛脚本编写之前,我从未做过。我将大规模扩展该计划。 根据蛋白质 xyz 坐标,我必须定义盒子大小。这个盒子将被分成一个大小为 0.5 A 的网格。根据距离和角度标准,我必须根据玻尔兹曼概率分布来分配点。

我的程序应该通过 0.5 A 的网格在每个方向上移动并生成随机点并检查距离和角度的情况。如果满足条件放置点,则根据玻尔兹曼概率分布丢弃该点。

这是我生成随机点的代码

from __future__ import division    
import math as mean    
from numpy import *   
import numpy as np   
from string import *    
from random import *    

def euDist(cd1, cd2):# calculate distance
    d2 = ((cd1[0]-cd2[0])**2 + (cd1[1]-cd2[1])**2 + (cd1[2]-cd2[2])**2)
    d1 = d2 ** 0.5
    return round(d1, 2)

def euvector(c2,c1):# generate vector
    x_vec = (c2[0] - c1[0])
    y_vec = (c2[1] - c1[1])
    z_vec = (c2[2] - c1[2])
    return (x_vec, y_vec, z_vec)


 for arang in range(1000):  # generate random point
        arang = arang + 1
        x,y,z = uacoord
        #print x,y,z

        x1,y1,z1 = (uniform(x-3.5,x+3.5), uniform(y-3.5,y+3.5), uniform(z-3.5,z+5))
        pacord = [x1,y1,z1]                 # random point coordinates
        print pacord

我完全震惊于从蛋白质结构的 xyz 坐标生成盒子大小以及如何定义大小为 0.5 的网格。如何检查框中的每个点。 任何帮助都将是可观的。

【问题讨论】:

那么,你有一堆 3D 点,你想生成一个包含所有点的框? 顺便说一句,我真诚地希望您知道自己在做什么:蛋白质折叠通常具有多个最小值,而幼稚的 MC 模拟往往会卡住。 定义盒子和网格大小我会生成点,然后会做进一步的计算。不知道box和grid怎么定义??????我的循环将如何从每个方向移动?????? 出于好奇,您为什么不使用当前免费且支持 GPU 的 MC 程序之一?我使用的那些会以您希望的任何方式随机化起点。 你看过github.com/ndexter/Aeolotopic-Monte-Carlo-Simulation/blob/… 吗?虽然它是用 C 语言编写的,但我相信你会对实现有所了解,并且它会帮助你编写 python 代码。 【参考方案1】:

喜欢您的主题、问题和方法。我不知道它会在这里停留多久。

矩形空间中的 3D 网格在有限元分析和所有其他分析物理问题的技术中很常见。看看网格是如何生成的。

通常有两部分:空间中 (x,y,z) 点的网格和连接它们的框。具有两个元素的简单 2D 网格如下所示:

D               E               F
o (1, 0) ------ o (1, 1) ------ o (1, 2)
+               +               +
+               +               +
+               +               +
o (0, 0) -------+ (0, 1) -------+ (0, 2)
A               B               C

这个网格中有六个点:

A (0, 0)
B (0, 1)
C (0, 2)
D (1, 0)
E (1, 1)
F (1, 2)

还有两个盒子:

1 - (A, B, E, D)
2 - (B, C, F, E)

因此,一种可能的方法是遍历所有框,检查质心的位置,并进行相应的调整。

我将从您的代码中外部化网格定义并从文件中读取它。这样您就可以使用相同的代码处理不同的网格。

如果我理解正确,网格将在空间中保持固定。你试图随机化蛋白质本身的运动。这就像一个流体问题,一种欧拉方法:您正在跟踪材料在固定控制体积上的运动。除了它是一种蛋白质而不是一种液体。

因此,您还将对空间中的蛋白质的初始条件有一个单独的定义。您计划及时增加它,以查看蛋白质如何根据您的规则折叠。

我接近了吗?

【讨论】:

感谢达菲莫。我不是在看蛋白质的运动,而是试图分配一堆点,它们靠近蛋白质并与蛋白质原子形成相互作用。我是通过随机点生成来完成的,但是每次我得到新结果时都会使用那个程序。所以我想到了这种方法。【参考方案2】:

您是否考虑过使用 PyRosetta?它更容易使用,因为您需要的许多功能已经内置。您可以在 PyMol 中实时可视化您的输出。我在 PyRosetta 中写了一个类似的脚本,相当容易编写和修改,它完成了它应该做的事情。

【讨论】:

我知道 PyRosetta。我正在做一个大项目并开发一个程序。所以我不能使用这样的可用程序。如果可能的话,您能否向我提供您的脚本或告诉我如何生成框和定义网格。谢谢。【参考方案3】:

我的代码是为蛋白质折叠应用程序编写的,但总体思路是相同的。它从某个温度开始并逐步降低温度,如果新位置/结构的能量得分低于前一个,则蛋白质(或在您的情况下为“点”)随机移动,如果不是姿势将根据 Metropolis 分布进行评估。你必须定义你的 scorefunction(),一个定义随机起始位置的函数和一个将你的点从其原始位置移动的移动器。下面的代码是根据我原来的脚本修改的,只是给大家一个大概的思路。

kT_lower=0.1
kT_upper=100
ktemp=kT_upper

max_iterations=15000

i=-1

#create random start point
pose=create_pose()

#evaluate start point
starting_score=scorefunction(pose)

while i<max_iterations:
    i=i+1
    new_pose=random_move(pose)
    if scorefunction(new_pose)<scorefunction(pose):
        pose=new_pose
    else:
        #linear decrease of kT
        ktemp=kT_upper-i*(kT_upper-kT_lower)/max_iterations
        #exponentatial decrease of kT
        #ktemp=math.exp(float(i)/float(max_iterations)*float(-5))*kT_upper+kT_lower

        try:
            p=math.exp(DeltaE/ktemp)
        except OverflowError:
            p=-1

        if random.random()<p:
            pose=new_pose
            print str(i)+'; accept new pose, metropolis'
        else:
            print str(i)+'; reject new pose!'

【讨论】:

看起来和听起来很像使用模拟退火进行优化。 它是模拟退火,或者至少我是如何理解模拟退火的概念的。 嗨 Ashafix。我又被困住了。您能否提供您的电子邮件地址,以便我可以在没有垃圾邮件的情况下获得您的帮助。我有同样的查询。我得到了这个想法,但发现实施起来有困难。 您好 Awanit,我很乐意为您提供帮助。只需检查我的个人资料并通过其他个人资料页面向我发送消息。

以上是关于蛋白质结构和网格的蒙特卡罗模拟的主要内容,如果未能解决你的问题,请参考以下文章

数学模型:4. 蒙特卡罗模拟

使用蒙特卡罗模拟多线程计算 Pi

蒙特卡罗方法

蒙特卡罗(洛)模拟

蒙特卡罗算法(Monte Carlo method)

蒙特卡罗(Monte Carlo) 模拟