生物大分子模拟要点

1、1、 computational biology计算机生物学是利用计算的方法对复杂生命现象和过程进行研究和预测的科学。它是理论与数据分析、数学建模和计算机模拟在生物学、行为学、社会群体系统中的应用与研究。它是数学、生物、物理、化学和计算机科学等高度交叉的学科，并与人类的生活与健康紧密结合。2、the four choices, define a model for a computer simulation 自由度、力场、外部环境、运动方程Definition of a model for a computer simulation2.力场Interaclions orFor anv syst

2、eE four choices have to be madeDegrees of freedom: 1、 自由度 what are (Ite demenlary enlitiea or uparlide,fI ii terart ion =phYsiro-cheinkal-Ijiologirlkiiowltdye2tfEnpri'ntutrptessuifvalkcxlrrnjil mrcfq3、外部环境Boundary conditionsMethods to configurations along degrees of freedom *rgr equations of mol

3、iou144、运动方程第二分子模型的代表：量子力学、量子/分子力学、分子力学、联合原子模型、粗粒化模型分子力学，又叫力场方法(force field method),是基于经典牛顿力学方程的一种计算分子的平衡结构和能量的方法。与量子力学不同，它求解的是Newton方程，而不是薛定渭方程。(用量子力学处理分子，计算量特别大，因此量子力学方法适用于处理小分子体系。分子力学计算量较小，更适合于大分子体系，可考虑额外效应)简单分子力场：分子力场是分子力学的核心；简单作用模型；力场的可移植性成键/非键相互作用： 键长、键角、二面角、范德华作用、静电作用、交叉相互作用( 五种 类型，老师说一定会考)3、

4、Potential energy functions由于分子力学是经验的计算方法，不同的分子力学方法会采用不同的势能函数(PotentialEnergy Function, PEF)表达式，而且力场参数值也会不同。一般将分子的PEF分解成五部分：Fbonded尸 + f + Fbond -3H'etchangle-beudrotate-along-bond键伸缩能健弯曲能二面角扭转能Fnon-bonded=F+ F1 xm-der-汀 aa/sel ecfrosiaric范镌华作用能静电作用能静电相互作用：计算方法 Electrostatic Contributions-点电荷法：通过

5、经脸规则或者量化计算确定每个原子上的 部分电荷(partial charge) ,两个原子一之间的静电作用用 库仑公式来计算行V K-红 飒 仃 两种方金在处理有机小分子体系的时候震季相似，但是当用-偶极量针蠲融僦宙题翳%麹的偶极矩，通 过计算偶极-偶极相互作用来描述静电相互作用4、 polarizability effect将电介质放入电场，表面出现电荷。这种在外电场作用下电介质表面出现电荷的现象叫做电 介质的极化。带电的配基会是蛋白质部分受到极化。5、force field, two most popular force field分子力场并不计算电子相互作用，它是对分子结构的一种简化模型

6、。一个分子的能量可以近似看作构成分子的各个原子的空间坐标的函数，简单地讲就是分子的能量随分子构型的变化而变化，而描述这种分子能量和分子结构之间关系的就是分子力场函数。传统力场、第二代力场、通用力场比较不同程序计算得到的能量值无意义AMBER力场：由Kollman课题组开发的力场，是目前使用比较广泛的一种力场，适合处理 生物大分子。AMBER力场的势能函数形势较为简单，所需参数不多，计算量也比较小，这 是这个力场的一大特色，但也在一定程度上限制了这个力场的扩展性。本力场用谐振子模型计算键长伸缩能和键角弯转能，用傅立叶级数的形式来描述二面角扭转能，选用 Lennard-Jones势来模拟范德华力；

7、用库仑公式来描述静电相互作用。CHARMM 力场：由Karplus课题组开发，对小分子体系到溶剂化的大分子体系都有很好的 拟合。第三6、 local minimum能量优化，寻找研究体系能量极小的状态，从数学角度看是处理多维函数的优化问题。优化的方法，降低体系能量到最靠近的一个最小。窄而深的最小值的统计学大小可能比高能下宽的最小值要少。ABCcnformiitional parameter7、global minimum全局优化研究的是多变量非线性函数在某个约束区域上的全局最优解的特性和构造寻求全 局最优解的计算方法。由于很可能在一个全局优化问题里存在多个局部最优解，且它们不同于问题的全局最优

8、解，因此人们无法借助于经典的局部优化方法求解这些问题。8、steepest descents最速下降法又称为梯度法，是 1847年由著名数学家 Cauchy给出的。其优点是工作量少， 存储变量较少，初始点要求不高；缺点是收敛慢，效率不高，有时达不到最优解。最速下降 法从目标函数的负梯度方向一直前进，直到到达目标函数的最低点。梯度下降法的计算过程就是沿梯度下降的方向求解极小值。9、 conjugate gradients methods, compare共轲梯度法是求解特定线性系统的数值解的方法，其中那些矩阵为对称和正定。共轲梯度法是一个迭代方法，所以它适用于稀疏矩阵系统，因为这些系统对于像乔莱

9、斯基分解这样的直 接方法太大了。这种系统在数值求解偏微分方程时相当常见。共轲梯度法是介于最速下降法与牛顿法之间的一个方法，它仅需利用一阶导数信息，但克服了最速下降法收敛慢的缺点，又避免了牛顿法需要存储和计算Hesse矩阵并求逆的缺点， 共轲梯度法不仅是解决大型线性方程组最有用的方法之一，也是解大型非线性最优化最有效的算法之一。 共轲梯度法是一个典型的共轲方向法，它的每一个搜索方向是互相共轲的，而这些搜索方向d仅仅是负梯度方向与上一次迭代的搜索方向的组合，因此，存储量少，计算方便。10、two methods for global energy minimization模拟退火是受金属热加工技术

10、的启迪而发展起来的一种随机搜索算法。将一个优化问题比拟成一个金属物体，将优化问题的目标函数比拟成物体的能量，问题的解比拟成物体的状态，问题的最优解比拟成能量最低的状态，然后模拟金属物体的退火过程，从一个足够高的温度开始，逐渐降低温度，使物体分子从高能量状态缓慢的过渡到低能量状态，直至获得能量最小的理想状态为止，从而得到优化问题的全局最优解。给定初温片.随机产生初始状态$=闻，令k=（hRepeatRepeat产生新状态?=Geneteiif min1 ,exp|-（C（sj）-C（s）ltk J >= random0 J s=Sj； Until抽样稳定准则满足；退温fHi=iipdate

11、（4）并令人=女+1 ；Until算法终止准则满足； 输出算法搜索结果。新状态接受概率仅依赖于新状态和当前状态，并由温度加以控制。质量高；简单、通用、易实现。由于要求较高的初始温度、较慢的降温速率、较低的终止 温度，以及各温度下足够多次的抽样，因此优化过程较长。遗传算法是一种以自然选择和遗传理论为基础，将生物进化过程中适者生存规则与种群内部染色体的随机交换机制相结合的随机化搜索算法。L选择问题的一个编码;给定种群规模杂交概率凡和变异概率Pm的值，给定终止条件.2 .随机产生由N个字符串（长度相等）组成的种群令尸0。3 .对p磔"）中每个个体pop"）,计算其适应度/管/叩江

12、）4.If满足终止条件then选择最佳个体作为GA的结果，输出结果，goto 6Else第四相空间，时间平均，系统平均11、Ergodic Hypothesi s遍历性假设是一个分子系统的足够长时间的演化（保持能量守恒）可遍历（或无限接近）任何微观状态。（绝大多数情况下不成立）不是真实的，因为采样时间非常短。但在一些特殊性情况是真实的。12、Metropolis approachmetropolis是一种采样方法，一般用于获取某些拥有某些比较复杂的概率分布的样本。采样最基本的是随机数的生成，一般是生成具有均匀分布的随机数。13、General steps of a Monte Carlo si

13、mulation蒙特卡洛方法( Monte Carlo method ),也称统计模拟方法， 是20世纪40年代中期由 于科学技术的发展和电子计算机的发明，而被提出的一种以概率统计理论为指导的一类非常重要的数值计算方法。 是指使用随机数( 或更常见的伪随机数)来解决很多计算问题的方法。MC方法的起源于布丰投针实验。针对待求问题，根据物理现象本身的统计规律，或人为构造一合适的依赖随机变量的概率模型，使某些随机变量的统计量为待求问题的解，进行大统计量(N-8)的统计实验方法或计算机随机模拟方法。(大数定理，中心极限定理)(1)选取所要研究的系统并建立适当的模拟模型。(2)设定模拟区域的边界条件，选

14、取粒子间作用势模型。(3)设定系统所有粒子的初始位置和初始速度。(4)计算粒子间的相互作用力和势能，以及各个粒子的位置和速度。(5)待体系达到平衡，统计获得体系的宏观特性。大数定理与中心极限定理与模拟退火的区别：模拟退火：先升温，再降温，不关心过程MC模拟：温度不变，关心过程，温度等于生理温度随机数，各种随机数随机数必备的两个特点：独立性和均匀性伪随机数的统计检验：均匀性，在【0,1】内等长度子区间中随机数的数量是一样的1. Define the system 定义系统2. Choate the energy runction 选择能量晶黄3. Decide on how to make ra

15、ndom perturbations 淡定热何使随机犹动4. Run the simulation 运在仿真 Choose an initial state (random, or based on some experimental input) evaluate energy 选择一个初始状态 评估能量 Iteratively:反复迭代make a random perturbationevaluate the energy of the perturbed conformationdecide whether the perturbed conformution ur the origi

16、nal fonforni»tionis added to the sampleend iteration ufter certain steps& Computing jvcritges urinterest 计算平均的利矗 ug: structures, binding energies, free energy, Entropy.第五14、the typical time step in a MD simulation时间步长是根据分子振动或转动频率来确定的，太长的时间步长会造成分子间的激烈碰撞，体系数据溢出；太短的时间步长会降低模拟过程搜索相空间的能力，因此一般选取的时

17、间步长为体系各个自由度中最短运动周期的十分之一。设为分子中最小伸缩振动周期的十分之一，一般碳氢的伸缩振动周期最短，是 10fs。所以速度限制在 1fs。15、why , choosing such a time step一般来说模拟总时间步越多T单步长越大，模拟结呆就按 近实际情况，但是由于受计算机硬件的限制,不可能时一个足 幅长的时间过程进行模拟，同时单步长枝大，对体系演变反映 的跨度越长，计算误差也就越大；加果单时间步长选择太短， 小于体系分子位置的演变，也会造成误差的增加，所以模拟时 要整合所选取的求解分子运动方程的算法.合理选择步长以减 小误差.求解中的另一个重要问题是时间步法的选择：

18、步长h展同 实现的M口步数，决慝了相空间有多大一部分被抽样,看起来 越大越好.这科可以对电大的部分抽样，但是由h所决定的时 间标尺和实际体系发生变化的时间标尺并不一定完全相同，所 以时间步长的选取应根据研究对象的具体情况而定16、The Leap-frog method蛙跳算法是一种全新的启发式群体进化算法，具有高效的计算性能和优良的全局搜索能力。作为一种新型的仿生物学智能优化算法，结合了基于模因进化的模因演算法和基于群体行为的粒子群算法。该算法具有概念简单，调整的参数少，计算速度快，全局搜索寻优能力强， 易于实现的特点。蛙跳算法的思想是：在一片湿地中生活着一群青蛙。湿地内离散的分布着许多石头

19、，青蛙通过寻找不同的石头进行跳跃去找到食物较多的地方。每只青蛙个体之间通过文化的交流实现信息的交换。每只青蛙都具有自己的文化。每只青蛙的文化被定义为问题的一个解。湿地的整个青蛙群体被分为不同的子群体，每个子群体有着自己的文化，执行局部搜索策略。在子群体中的每个个体有着自己的文化，并且影响着其他个体，也受其他个体的影响，并随着子群体的进化而进化。当子群体进化到一定阶段以后，各个子群体之间再进行思想的交流（全局信息交换）实现子群体间的混合运算，一直到所设置的条件满足为止。边界条件，真空边界条件，周期边界条件，截断半径的边界条件17、Box types in PBC立方体、十二面体、六方柱、截断正八

20、面体，格子的选取必须足够大，考虑分子的旋转。正十二面体和截断正八面体用得比较多18、why ,electrostatic interactions are problematic in MD长程力：随距离的增加而缓慢减少。如静电力。这种静电想回作用会消退的。19、 the most efficient algorithm to compute electrostatic interactions so farEwald求和在k空间求和那一部分的级数会发散。通常为了避免这种问题，分子模拟的程 序在处理非电中性系统的Ewald求和是都会给系统加上一个均匀分布的背景电荷密度来中和这部分电荷。如果你的

21、体系净电荷不大（比如说只有-2 e）,那么这个背景电荷的影响是比较小的，反之，如果你的体系净电荷很大，那么 Ewald求和带来的误差就会很大。（使 用快速傅里叶变换）20、microcanonical, canonical, and isothermal-isobaric ensembles微正则系综（microcanonical ensemble）:系综里的每个体系具有相同的能量（通常每个体 系的粒子数和体积也是相同的）。正则系综 （canonical ensemble）:系综里的每个体系都可以和其他体系交换能量（每个体 系的粒子数和体积仍然是固定且相同的），但是系综里所有体系的能量总和是固

22、定的。系综内各体系有相同的温度等温等压系综（isothermal-isobaric ensemble）:正则系综的推广， 体系间可交换能量和体积， 但能量总和以及体积总和都是固定的。（系综内各体系有相同的粒子数。）正如它的名字，系综内各个体系有相同的温度和压强。21、how, couple the system at constant temperature and pressure in MD（原）NVT =一个明显的方式来改变系统的温度是调节速度比例。简单的速度扩展方法是通过因子入控制。NVT系综就是通过控制体积（密度）来控制体系的压力。压力可以保持在一个常数值通过简单地伸缩体积，即改变模

23、拟细胞的体积。模拟盒子的体积被因子入缩放。入=lf （Pith）TP22、Particle decompositionPD是一种简单的并行算法类型。PD算法比较容易实现，但是由于它需要各处理器不停交换粒子的位置信息，通信开销相当大，难以实现大规模 MD模拟。 同时它对硬件尤其是网络性能要求高，否则并行效率会大打折扣。在开始模拟，粒子被分配到处理器。粒子之间的力需要分配给处理器，这样力负载均衡。这种分解要求每个处理器知道系统中至少一半的粒子坐标，从而达到高N的处理器数量。因为这，PD没有很好的伸缩性。现在采用长程作用才是有用的，尤其 NMR距离或方向约束。 23、domain decomposi

24、tion.DD算法的原理是：模拟盒子被分解为一些相同大小的区域小盒子盒子的大小可以根据截断 半径得到，一般取边长等于或者大于截断半径。每个处理器分担一个盒子内粒子的作用力、位置和速度等，这样粒子可以在小盒子之间自由移动。这样的好处在于一个处理器只要考虑其自身分担区域与相邻区域的粒子间作用力，从而大大减少处理器之间的通信量。DD通常可以得到很高的并行效率，适用于大量复杂粒子的模拟。但由于处理器存储的是特定范围内的粒子， 在运行过程中难免会发生有些处理器存储粒子多，有些存储少的情况，这样会造成处理器负载不均。粒子数目少的处理器早早完成运算，然后等待负荷较重的处理器运算完才能交换信息，这会大大影响并

25、行效率。所以DD算法中负载平衡是关键，其编程相当复杂。24、The two major limitations of MD 详细的平衡并不容易。遍历性并不总是得到保证。第六25、RMSD and RMSFRMSD表示的是分子结构变化的程度。适合同一个蛋白的对比RMSD是计算在某一时刻的构象与目标构象所有原子偏差的加和，对原子数的平均。用来表是蛋白质结构之间差异的参数是两个结构之间原子位置的RMSD以一个结构中的原子与另外一个结构中对应原子为计算标的，因此，如果两个分子在座标系统中以不同的位置开始计算， 那么不管其结构是否相似， 这两者之间的 RMSD必定相当大。 也因为这样，我们为了要计算有意

26、义的RMSD ,两者的结构要尽可能的重叠。可以通过计算 RMSD来当作评估蛋白质结构的可信度 ：在模拟过程中，分子会不断的发生（fluctuation较小时）再进一步进变化，而对于我们而言，必须等到分子结构在稳定的状态下 行分析，这样才比较有意义。啊|口出)-ri(ta) A' r .L t=lRMSF值表示的是分子中各个原子运动的自由程度RMSF计算的是一个粒子在时间T内，其位置的偏移量之平方和随时间的平均后再开方。RMSF的计算可以写成积分的形式，时域上（0,T）区间的积分。RMSF和粒子的B-factor是相关的RM4226、Algorithms for cluster anal

27、ysis聚类分析是一种分类的多元统计分析方法。按照个体或样品的特征将它们分类，使同一类别内的个体具有尽可能高的同质性，而类别之间则应具有尽可能高的异质性。从大量的构象选择一个小数量的“代表”做后续的分析。27、Normal mode analysis简正模分析qiMj最小化构象的势能海赛矩阵的计算海赛矩阵对角化28、elastic network model粗粒化弹性网络模型，生物分子是相当于一个三维弹性网络ENM仍和NMA 一样，但蛋白质模型大大简化。原子是由弹性的网络连接的。首先，没有对能源最小化必要，因为所有的弹性连接的距离被他们最小能量长度取代。第二，减少与NMA的标准对比。29、pr

28、incipal component analysis主组分分析主成分分析是设法将原来众多具有一定相关性(比如P个指标)，重新组合成一组新的互相无关的综合指标来代替原来的指标。主成分分析，是考察多个变量间相关性一种多元统计方法， 研究如何通过少数几个主成分来 揭示多个变量间的内部结构，即从原始变量中导出少数几个主成分， 使它们尽可能多地保留 原始变量的信息，且彼此间互不相关.通常数学上的处理就是将原来 P个指标作线性组合， 作为新的综合指标。计算集合坐标的方法 简正模分析(Normal mode analysis)CI =Idx. dXj 粗粒化弹性网络模型(Mastic network mod

29、el)Hessian矩阵的对角化: C = UKUT 主组分分析(Principal Component Analysis) G=”,心-卜加第七30、 Why it is generally difficult to compute free energy高能状态有一个非常低的概率和对积分微不足道的贡献，所以我们可以通过 MD或MC生成准确的估计。高能状态对积分做出重大贡献 ，所以结果F通过MD或MC聚合不准确。1( F( N、F = -k£IInQ - ksThiexp + ° " p(pI"I k演J J采样充分得到准确的自由能数值，自由能误差来源

30、主要是采样不够31、Free energy perturbation method自由能微扰给出了两种状态自由能的差异。自由能的差异是配分函数的比例，自由能微扰(FEP)是一种依靠统计学的计算自由能的方法。这种方法最早在1954年由Robert W. Zwanzig 提出来的。在这个方法中有个 Zwanzig 方程。自由能微扰的基本原理是从一个已知的体系出发，通过一系列微小的变化转到另一个体系，在每一个变化 步骤中做分子动力学模拟， 把每一步的体系的势能带入相应的公式中，就可以得到两步之间的自由能变化，把所有的自由能变化加在一起，就能得到两体系的自由能变化。 如果有两个配体分子 A和B,它们和

31、受体S 形成了复合物AS和B&为了求算这两个配体分子和受体之间结合自由能的差值，则需要通过热力学循环。先求算配体分子A和B之间的自由能变AG,然后再 求算出复合物AS和BS之间的自由能变 G,二者之差 就是我们需要得到的自由能变化 AG在一般的自由能微扰的计算中，求算的都是两个配体分子和受体之间结合能的差值。在有些应用中，也可 以采用这种方法来计算单个配体和受体之间的结合自由能。在这种情况下，突变的过程不是一个配体到另 外一个配体，而是一个配体到相应的溶剂分子 （一般为水分子）。微扰：跨度大的过程分 成几个小过程 “ a f af 100j die®。1=Jjre'W

32、i-c；o）o（r）032、potential of mean forcePMF是指自由能曲线.自由能表面沿着选择的坐标。计算出一个物理可实现的过程。描述一个体系的不同状态，需要有一个比较有区分度的变量，这个变量叫“反应坐标"，正确的选取反应坐标非常重要，需要对所研究体系的本质有比较深入的了解，比如蛋白质模拟中振动频率比较小的二面角，拉伸问题中的位移等等。体系的配分函数是 Q,他是对玻尔兹曼因子在所有自由度的上的积分，lnQ对应自由能，但这个量表示的是整个体系处在平衡态的热力学性质，我们可能更关心体系的不同状态的区另1J,所以在积分时加上一个delta函数：delta （z-x）其中

33、z是反应坐标，得到的就是平均力势PMF ,可以理解为反应坐标为 x时体系的自由能。33、reaction coordinates （选反应坐标）反应坐标是一个抽象的一维坐标代表进步沿反应途径。它通常是一个几何参数变化转换成一个或多个分子实体。在分子动力学模拟中，反应坐标叫做集体变量。34、MM/PBSA method在amber中，可以通过 mm/pbsa来计算这个相对自由能。MM/PBSA 全称是 MolecularMechanics/Poisson-Boltzmann Surface Area 。该方法成功运用于生物大分子体系，包括 DNA 的构象变化、蛋白-蛋白、蛋白-DNA、蛋白-小分

34、子的相互作用。从其名字便可看出，该方 法将结合自由能拆成了多项分别计算，如Molecular Mechanics代表动力学项。其具体拆法如下，受体R与配体L的结合自由能A Gb可写成如下形式， AGb=AEMM+ AGsolv-T A SAGb=AEMM+ AGsolv-T A S其中，A EMM是真空中受体-配体结合前后动力学能量的差值，动力学项又包含三项Eint、EedW和Eelec,其中int指键、键角及二面角，这些在分子动力学理论里都有详细介绍。相 比较而已，动力学项也是最好计算的，对体系跑完MD后该项便可从记录的能量文件中直接得到。 第八35、Two general strategi

35、es to enhance sampling efficiency on energy surface粗粒化模型与分子动力学模拟相结合提高采样效率的方法根据热浴耦合的原理提高在低频模空间的采样效率36、give an example for each strategy弹性网络模型(考过)S肽类似物的折叠/去折叠 第九37、multi-scale modeling多尺度模拟 体系大小 (system size) 原子模型粗粒化模型“连续介质模型Atomic model(CG model)(continuum model)38、the purpose for coarse-graining 减少自

36、由度频率快速运动 流畅的潜能的表面 长时间的步骤可以采取 微妙模拟的可能All-atom MD:time step of 1-2 fs, time frame sampled ns/psL?oarse-grained MD: time step of 20-50 fs, time frame sampled |is/ms39、the most important issue in coarse-graining 粗粒度的水平 我们需要多少描述细节 ？ 我们可以描述多少细节 ？ 由于大幅削减粒子的数量和交互 ，加快了时间步长。虽然失去细节但也要保证准确性。40、at least two kind

37、s of coarse-grained models of biomolecules联合原子模型在united-atom模型中，CH3 group看成一个粒子，仅仅是把C-H键和H-C-H角给近似掉了，主要因为这两种势能是属于硬自由度范围，对于模拟结果影响很小，往往是可以忽略的(所以模拟中键长常常固定往往是基于这个考虑)。格点模型(Lattice model)快速/详尽的一代的低分辨率的结构Q H疏水域基亲水的极 性蛾基41、Simulations methods used for structure calculation/refinementrestrained molecular dynamics (rMD)约束的分子动力学：喙MR匕口出腺的,口咽e+ 'dihed十 vdW coulombsimulated annealing (SA)模拟退火hybrid methods42、structure prediction and design(假设)：蛋白质分子天然构象为处于热力学最稳定，能量最低状态根据Anfinsen热动力学原理，蛋白质在细胞中应该处在它与环境的自由能最低态。这意味着可以根据物理、化学、生物学等知识来设计蛋白质的能量函数，据此寻找这种最低自由能所代表的结构。大多数蛋白质的天然活性状态接近热力学上的最稳定态。小分子设计：基于构效关系的小