分子动力学模拟经验谈

1、分子动力学攻略此文为dddc_redsnow发表于biolover上的关于分子动力学的系列原创文章，相当经典与精彩，特此将系列文章整合，一起转载，望学习动力学的新手们共同学习，提高进步，在此特向dddc_redsnow本人表示感谢。动力学系列之一（gromacs，重发）在老何的鼓励下，发一下我的gromacs上手手册（我带人时用的，基本半天可以学会gromcas）#Processproteinfilesstepbystep#pdb2gmx-f2th_cap.pdb-o2th_cap.gro-p2th_cap.top-ignh-ternedit2th_cap.topeditconf-f2th_c

2、ap.gro-o2th_cap_box.gro-d1.5genbox-cp2th_cap_box.gro-cs-p2th_cap.top-o2th_cap_water.gromake_ndx-f2th_cap_water.gro-o2th_cap.ndxgenpr-f2th_cap_water.gro-n2th_cap.ndx-o2th_cap_All.itpgenpr-f2th_cap_water.gro-n2th_cap.ndx-o2th_cap_M.itpgenpr-f2th_cap_water.gro-n2th_cap.ndx-o2th_cap_C.itpneditFlavo.itpg

3、rompp-fem.mdp-c2th_cap_water.gro-p2th_cap.top-oprepare.tprgenion-sprepare.tpr-o2th_cap_water_ion.gro-np1-pq1#Minimizestepbystep#1.minimizationfixingwholeprotein#2.minimizationfixingmaincharinofprotein#3.minimizationfixingCaofprotein#4.minimizationwithoutfix#grompp-np4-fem.mdp-c2th_cap_water_ion.gro-

4、p2th_cap.top-ominimize_water.tprmpirun-np4mdrun-nice0-sminimize_water.tpr-ominimize_water.trr-cminimize_water.gro-eminimize_water.edr-gminimize_water.log&grompp-np4-fem.mdp-cminimize_water.gro-p2th_cap.top-ominimize_sidechain.tprmpirun-np4mdrun-nice0-sminimize_sidechain.tpr-ominimize_sidechain.trr-c

5、minimize_sidechain.gro-eminimize_sidechain.edr-gminimize_sidechain.log&grompp-np4-fem.mdp-cminimize_sidechain.gro-p2th_cap.top-ominimize_sidechain_ex.tprmpirun-np4mdrun-nice0-sminimize_sidechain_ex.tpr-ominimize_sidechain_ex.trr-cminimize_sidechain_ex.gro-eminimize_sidechain_ex.edrminimize_sidechain

6、_ex.log&grompp-np4-fem.mdp-cminimize_sidechain_ex.gro-p2th_cap.top-ominimize_all.tprmpirun-np4mdrun-nice0-sminimize_all.tpr-ominimize_all.trr-cminimize_all.gro-eminimize_allx.edr-gminimize_all.log&editconf-fminimize_all.gro-ominimize_all.pdb#Raisetemperaturestepbystep#1.raisefrom0Kto100Kfixingwholep

7、rotein#2.raisefrom100Kto200Kfixingwholeprotein#3.raisefrom200Kto300Kfixingwholeprotein#4.balancefixingmaincharinofprotein#5.balancefixingCaofprotein#grompp-np4-fheat.mdp-cminimize_all.gro-p2th_cap.top-otemperature100K.tprmpirun-np4mdrun-nice0-stemperature100K.tpr-otemperature100K.trr-ctemperature100

8、K.gro-etemperature100K.edr-gtemperature100K.log&grompp-np4-fheat.mdp-ctemperature100K.gro-p2th_cap.top-otemperature200K.tprmpirun-np4mdrun-nice0-stemperature200K.tpr-otemperature200K.trr-ctemperature200K.gro-etemperature200K.edr-gtemperature200K.log&grompp-np4-fheat.mdp-ctemperature200K.gro-p2th_cap

9、.top-otemperature300K.tprmpirun-np4mdrun-nice0-stemperature300K.tpr-otemperature300K.trr-ctemperature300K.gro-etemperature300K.edr-gtemperature300K.log&g_energy-ftemperature300K.edr-stemperature300K.tpr-otemperature300K.xvggrompp-np4-fheat.mdp-ctemperature300K.gro-p2th_cap.top-oT300K_M.tprmpirun-np4

10、mdrun-nice0-sT300K_M.tpr-oT300K_M.trr-cT300K_M.gro-eT300K_M.edr-gT300K_M.log&grompp-np4-fheat.mdp-cT300K_M.gro-p2th_cap.top-oT300K_Ca.tprmpirun-np4mdrun-nice0-sT300K_Ca.tpr-oT300K_Ca.trr-cT300K_Ca.gro-eT300K_Ca.edr-gT300K_Ca.log&grompp-np4-fheat.mdp-cT300K_Ca.gro-p2th_cap.top-oT300K_MD.tprmpirun-np4

11、mdrun-nice0-sT300K_MD.tpr-oT300K_MD.trr-cT300K_MD.gro-eT300K_MD.edr-gT300K_MD.log&至此，全部搞定，可以跑了。PJ的东西放出来敏感，原创的教学东西也可以收精华吧动力学系列之二（amber基础篇）上次写了一个，结果一下在因为发贴的原因丢了，加上最近忙的要命，不过还是吃完饭干活前，静下心写一下amber,amber远比gromacs要复杂（不过还比不上charmm,估计第一次看charmm脚本的人一定觉得自己选错了行当，类fortran的语言加上自创的格式，唉，又是一段心酸的往事）。我把amber分成两部分讲，第一部分

12、讲基础，第二部分讲应用。这里先讲第一部分。关于整体概括.amber的文件格式有三种最为重要：top,crdandpdb。pdb可以生成topandcrd,topandcrd也可以转换成pdb，这些都是ascii码文件，可以手动编辑（这也是DNA+Protein时候有时不得不作的事情，巨大的工作量）。top文件中是拓扑文件，相当于gromacs的top，不过直接把lib的参数搞过来，不要再调用lib了。crd文件是坐标和速度文件，类似于gro文件，pdb不用我说了吧。amber中的概念：所有的单元从原子,小分子，残基，碱基等等一直到大分子都是unit,熟悉C+JAVA的朋友想一下object就理

13、解了，amber的xleap的操作就像一个外包体，利用object的独立性进行关联，组合或者分拆object,实现不同的功能，每一个unit都是独立单元。而标准的氨基酸，碱基就像class，你自己的蛋白就是object。每个class对应多个有自己名字的object。这些都是在xleap中实现的，xleap就是利用这种关系来导入，修改蛋白或者DNA/RNA。3.amber自带gaff力场（比gromacs那个网页算的小分子正规多了），适用于大部分有机小分子，当然修改力场文件可以让你获得更多的支持，尤其是特殊的原子比如说卤素，金属原子。电荷是采用的resp拟和，而不是原子的partialchar

14、ge,这个的优点见JACS的原始文献，不再赘述，主要是考虑极化。通过xleap搞定了大分子，小分子，水,溶剂离子，跑得时候就是sander了，当然8以后的pmemd是更好的选择，不过只能跑pme方式的动力学。轨迹文件amber做的不好，太大，没有压缩。断点续跑也不如gromacs，半个坐标的情况时有出现,只见原子叉叉满天飞啊，那叫一个壮观。分析工具首推ptraj（不是我推的，是他们，我个人喜欢carnal），加各种辅助工具可以实现gromacs的各种功能。优点：上次说过了，轨迹稳定，重复性好，可信度高。有些朋友非要跑出个地动山摇，惊涛骇浪的轨迹,不推荐您使amber,估计等到它跑到那一天，您也

15、毕业了。缺点：真的打算用amber吗，打算跑10ns吗？打算有什么特别有意思的想法吗？恭喜您，您可以考虑买硬盘了。amber的存储量大概是1ns-1.5GB,如果想发一个好一点的文章，阐述一个精细的机理，没有10条轨迹搞不定的，算算是多少硬盘.我的一个题平均120G轨迹。洋洋洒洒写了这么多，不过amber的真正强大之处还不在这里，它的强大在于它的应用广泛，有着各种的“人无我有，人有我精”的套件，下次我在应用篇里面详细讲一下ttsv发表于2007-9-2105:21动力学系列之三（amber应用篇）上海的天气真是一天一变啊，又是大雨滂沱了。反正出不去了，索性把amber的第二部分写完，大家周末也

16、有个消遣。上回说amber的基础应用谈到了他有很多应用的套件，这里我把最常用的以及比较容易发好文章的几个部分讲解以下：1。首先要说的就是用于药物设计平衡的sander/pmemd，一般来说，用amber作一下平衡是同源模建的重要可选步骤之一，尤其适用于结构不确定性大的loop和缺乏明确模板的地方，当然，如果您要是有macromodel，还是用那个，最近的文章发现申稿人很喜欢那个，可能是因为那个公司比较得到大家的认可吧2。MMPBSA,计算自由能的很好选择，KUNZ他们组就喜欢拿这个说事，作为虚筛的终结器。当然了方法的想法是很好的，不过还是有一些缺点的，比如对共轭体系，可极化严重的分子计算结果很

17、不准确。MMPBSA流行的时候只做一个这个旧可以发ProteinSci.现在基本不可能了，不过还是很好的一个计算自由能的方法，平均半个小时一个分子3。NMODE，有的朋友问要什么样的分析工具，NMODE绝对是一个可选的东东，具体的原理我不讲了，很简单，矩阵二次求导做频率，底频高幅运动就是它分析的对象，Maccoman他们组只做这个就在JMB上发了无数的文章（人家是开山鼻祖，我们做，呵呵，估计就是BJ也很难啊）4。Alan扫描：预测突变影响的工具，很粗糙，但是如果你计算资源不是很多的话，这个是推荐的方法。原理没什么可讲的，想用的去看manual，不到1页5。TMD,我最近做的东西，和什么SMD了

18、，FMD了。一大堆兄弟姐妹啊，不过也快做烂了，可以模拟蛋白大幅运动的过渡态6。TI。目前公认的计算自由能最好方法，缺点是相当的耗时间，一个师兄做n年,结果发现不适合他的体系，韶华逝水啊，没有时间的朋友千万不要碰上面的这些方法都是最常用的，也是文章中暴光频率最高，精通2到3个，一个好体系，综合起来，可以冲击JMB或者JACS。呵呵，外面雨终于小了，希望大家看了能有帮助，会一个技术不难，难的是把这个技术用的恰到好处，下次有时间聊聊怎么综合使用各种工具分析问题动力学系列之四（能量分析）明天就是5.1，可是手里还是一大堆的活，实验也是诡异的要命，时光真是快啊。前面对两个基本软件介绍了一下，本来想讲讲c

19、harmm和macromodel，后来想想，没有什么意思，还是讲一下分析方法对大家更适用，如果大家真的想听这两个软件，后面我再补充。动力学分析方法多种多样，但是能量永远是不变的话题，很多问题的说明都是以能量为基础的，这里说的能量包括各种力场能量，溶剂化能，自由能等。我主要对自由能讲解一些如何进行计算。在amber中，可以通过mmpbsa来计算这个相对自由能，在免费的软件gromacs中，能量的计算（不知道有心的朋友是否自己加和过他提供的能量项，是不是缺点儿什么啊相当的粗糙，虽然提供了计算自由能的工具，但是源码中却是空，不知道到哪个版本才能实现。amber是比较好用，但是也是一个收费软件，特别是

20、想发文章的朋友要注意了。如果是不想买amber的朋友，一个中间路线就是用gromacs跑轨迹，用autodock算能量（自由能），虽然不如mmpbsa原理上准确，但是忽略部分熵效应也在变动不大的体系也还不错。下面这个脚本是我自己写的，可以实现这一功能，希望能对大家有个帮助：#!/bin/sh#$1issetfordifferentlistname#$2issetfortheligandname#preparethe$2.bndinadvancewithcommand:#/root/autodock/dist305/bin/pdbtoatmlig.pdbq|root/autodock/dist3

21、05/bin/atmtobndlig.bndforsnapshotincat$1dorec=protein_mn_hoh1.$snapshotlig=$2.$snapshotecho$rec#assigntheatomicsolvationparameterstoPDBQ-formattedversionofyourmacromolecule,root/bin/addsol$rec.pdbq$rec.pdbqsroot/bin/autotors-A$2.bnd$lig.pdbq$lig.out.pdbq$lig.out.$rec.root/autodock/dist305/bin/autodo

22、ck3-p$lig.out.$rec.dpf-l$lig.out.$rec.epdb.log-c$1_Bindingroot/bin/Docked.awk$lig.out.$rec.epdb.log$1_Dockedroot/bin/Inter_score.awk$lig.out.$rec.epdb.log$1_Inter_Energyroot/bin/Intra_score.awk$lig.out.$rec.epdb.log$1_Intra_Energy#rm-f*.map*.fld*.xyz*.err*.pdbqs*.glg*.log*.com*.gpf*.dpfBinding.awk#!

23、root/bin/gawk-fBEGINname=ARGV1split（name,array,.）if（$0/EstimatedFreeEnergyofBinding/）printarray1$9$10Docked.awk#!root/bin/gawk-fBEGINname=ARGV1split（name,array,.）if（$0/FinalDockedEnergy/）printarray1$7$8动力学系列之五（外一篇,1-4楼加在一起）上次发了系列四的能量计算后，不少朋友pm我，既然可以用多个大分子构象计算能量曲线，反过来是否可以用一个大分子对小分子数据库做类似dock,gold那种虚拟

24、筛选那？答案是可以的，不过这个要通过脚本来实现，虽然也有一个程序，但是那个程序是我写给ddgrid这个工程的，这半年暂时不能放出来，年底给大家共享（主要原因是怕阿）。言归正传，不过这个脚本也可实现一样的功能，就是时间上慢一些，空间上占用多一点，以前在别的地方放过早期的版本，这次放一个完全版，而且包含了我们常用的参数设置。虽然和动力学这个系列关系不是非常紧密，但是在进行多个阳性物对比动力学的时候结合上一个脚本可以把大量工作简单化。废话少说，放东西先：主程序,由于论坛屏蔽原因，我用了替代字符，而且好像总是砍掉我的字符啊不得不分开发，请大家见谅&请转换为EOF后使用/准备文件#!/bin/shreceptor=receptor.mol2rec=$receptor%.mol2mol2topdbqs$receptorforligandincatmol2.listdolig=$ligand%.mol2deftors$ligand&/选