文章转自http://articles.biolover.com/

对动力学软件的一些认识

最近看到坛子上动力学问的人多起来了，我做过一点这方面的工作，谈一点自己的感受，对大家可能有些帮助吧
1.什么样的体系适合做动力学
不是所有的体系都适合做动力学，有的体系跑了分析了快一年才发现原来并不适合做动力学（组里有过这样的先例），轨迹扔了可惜，不扔是垃圾，还占空间，被人骂（组里硬盘空间太少，才几个T）。动力学最重要的是选择体系，这也就是说什么样的体系适合做动力学（或者说适合发文章的）.
以下几类都是比较好作的：
a. ligand induced conformation
b. open-closed mechanism
c. binding mechanism
d. large-scale mechnism
不好做的有（但是作出来会很好）
a. membrane protein
b. ion-channel protein
c. Protein-DNA/RNA interaction
千万不要去尝试（会让你做完欲哭无泪的）
a. homoloy-based MD(especillay <30%)，几乎文章做一个拒一个，除非你自己组里有人给作实验验证
b. independent stable protein，分析不出什么东西，小东西比方说水啊什么的，这个都已经过时了，不要扣细节了
c. Super large-scale MD, 这个我不用说了吧，估计我们作不了的呵呵，别的地方也很难做
d. 不知道自己为什么去作MD，先跑上再说的那种（最危险的一个），没有人指点你会迷失在分析的海洋中。
2.选择什么样的软件做动力学
动力学的软件多如牛毛，选择哪个去做哪？我只谈主流的自己摸过的几个（排名不分先后），Gromacs, amber, charmm, gromose和macromodel。NAMD我没有摸过，不敢评价。
a. Gromacs大家问的最多，我先说。一个字！快。其他的优点还有分析傻瓜，简单易学（我自己带人的话，一个上午之内保证学会，大家不要板砖啊，确实是缺乏内涵，傻瓜的很），可以pull（而且可以想怎么pull就怎么pull，还可以多个pull）。我说说缺点：结果可靠性差，重复性难，高档次分析工具欠缺，分析工具中bug多（仔细看就会知道，结果一看就知道不对），并行通信量大.....。这个软件只适合没有经济能力购买商业软件和需要长时间动力学(微秒）的同学使用
b. amber和charmm（偷懒了，放在一起），两者差不多，但是前者商业化程度更好，做得比较容易上手；后者还是老样子，脚本复杂，分支狂多，这也和karplus学生众多有关(孩子多了就是不好管啊)。优点都是结果相对准确，可信度高，具有高档次分析工具。缺点是上手难，跑的慢，理论基础要好，狂占空间（我都被管理员说了无数次了，唉）
c. gromose，拥有和CPMD的无缝接口是最大的优点。什么，不知道CPMD是什么，算我没说，跳过这个软件吧
d. macromodel，商业化最好，拥有能力最强，号称动力学黄金标准的MD软件。一个字，贵。但是人家薛定额做的就是好，结果被各大公司使用频率最高。
e. 其他，还有很多动力学软件，NAMD是比较出名的一个，但是我没有摸过，没有什么特色主要是。更小的大家如果感兴趣，一句话：小玩贻情，大作伤身。
3.应该做什么样的分析
rmsd，rmsf, t, e,dis, angle这些俗套就算了，最重要的是你对你的体系精通，各种实验结果了然于胸，机制脱口欲出。最好掌握一些发展前沿的分析手段，一两句话说不清楚，这个我下次有时间再开个贴说说。
洋洋洒洒写了一中午，后来写困了，可能前言不搭后语了，大家不要拍我啊，仅供大家参考
ps:老何，这个东西能不能加个精啊，我好有动力继续讲下去：）

动力学系统之二－amber基础篇

上次写了一个，结果一下在因为发贴的原因丢了，加上最近忙的要命，不过还是吃完饭干活前，静下心写一下amber，amber远比gromacs要复杂(不过还比不上charmm，估计第一次看charmm脚本的人一定觉得自己选错了行当，类fortran的语言加上自创的格式，唉，又是一段心酸的往事...)。我把amber分成两部分讲，第一部分讲基础，第二部分讲应用。这里先讲第一部分。
1. 关于整体概括. amber的文件格式有三种最为重要：top, crd and pdb。pdb可以生成top and crd, top and crd也可以转换成pdb，这些都是ascii码文件，可以手动编辑(这也是DNA+Protein时候有时不得不作的事情，巨大的工作量)。top文件中是拓扑文件，相当于gromacs的top，不过直接把lib的参数搞过来，不要再调用lib了。crd文件是坐标和速度文件，类似于gro文件，pdb不用我说了吧。
2. amber中的概念: 所有的单元从原子，小分子，残基，碱基等等一直到大分子都是unit，熟悉C++,JAVA的朋友想一下object就理解了，amber的xleap的操作就像一个外包体，利用object的独立性进行关联，组合或者分拆object, 实现不同的功能，每一个unit都是独立单元。而标准的氨基酸，碱基就像class，你自己的蛋白就是object。每个class对应多个有自己名字的object。这些都是在xleap中实现的，xleap就是利用这种关系来导入，修改蛋白或者DNA/RNA。
3. amber自带gaff力场（比gromacs那个网页算的小分子正规多了），适用于大部分有机小分子，当然修改力场文件可以让你获得更多的支持，尤其是特殊的原子比如说卤素，金属原子。电荷是采用的resp拟和，而不是原子的partial charge，这个的优点见JACS的原始文献，不再赘述，主要是考虑极化。
4. 通过xleap搞定了大分子，小分子，水，溶剂离子，跑得时候就是sander了，当然8以后的pmemd是更好的选择，不过只能跑pme方式的动力学。
5. 轨迹文件amber做的不好，太大，没有压缩。断点续跑也不如gromacs，半个坐标的情况时有出现，只见原子叉叉满天飞啊，那叫一个壮观。分析工具首推ptraj（不是我推的，是他们，我个人喜欢carnal），加各种辅助工具可以实现gromacs的各种功能。
6. 优点： 上次说过了，轨迹稳定，重复性好，可信度高。有些朋友非要跑出个地动山摇，惊涛骇浪的轨迹，不推荐您使amber，估计等到它跑到那一天，您也毕业了。缺点：真的打算用amber吗，打算跑10ns吗？打算有什么特别有意思的想法吗？恭喜您，您可以考虑买硬盘了。amber的存储量大概是1ns->1.5GB, 如果想发一个好一点的文章，阐述一个精细的机理，没有10条轨迹搞不定的，算算是多少硬盘....我的一个题平均120G轨迹。
洋洋洒洒写了这么多，不过amber的真正强大之处还不在这里，它的强大在于它的应用广泛，有着各种的“人无我有，人有我精”的套件，下次我在应用篇里面详细讲一下

动力学系列之三－amber应用篇

上海的天气真是一天一变啊，又是大雨滂沱了。反正出不去了，索性把amber的第二部分写完，大家周末也有个消遣。
上回说amber的基础应用谈到了他有很多应用的套件，这里我把最常用的以及比较容易发好文章的几个部分讲解以下：
1。首先要说的就是用于药物设计平衡的sander/pmemd，一般来说，用amber作一下平衡是同源模建的重要可选步骤之一，尤其适用于结构不确定性大的loop和缺乏明确模板的地方，当然，如果您要是有macromodel，还是用那个，最近的文章发现申稿人很喜欢那个，可能是因为那个公司比较得到大家的认可吧
2。MMPBSA, 计算自由能的很好选择，KUNZ他们组就喜欢拿这个说事，作为虚筛的终结器。当然了方法的想法是很好的，不过还是有一些缺点的，比如对共轭体系，可极化严重的分子计算结果很不准确。MMPBSA流行的时候只做一个这个旧可以发Protein Sci. 现在基本不可能了，不过还是很好的一个计算自由能的方法，平均半个小时一个分子
3。NMODE，有的朋友问要什么样的分析工具，NMODE绝对是一个可选的东东，具体的原理我不讲了，很简单，矩阵二次求导做频率，底频高幅运动就是它分析的对象，Maccoman他们组只做这个就在JMB上发了无数的文章（人家是开山鼻祖，我们做，呵呵，估计就是BJ也很难啊）
4。Alan扫描：预测突变影响的工具，很粗糙，但是如果你计算资源不是很多的话，这个是推荐的方法。原理没什么可讲的，想用的去看manual，不到1页
5。TMD，我最近做的东西，和什么SMD了，FMD了。。。。一大堆兄弟姐妹啊，不过也快做烂了，可以模拟蛋白大幅运动的过渡态
6。TI。目前公认的计算自由能最好方法，缺点是相当的耗时间，一个师兄做n年，结果发现不适合他的体系，韶华逝水啊，没有时间的朋友千万不要碰
上面的这些方法都是最常用的，也是文章中暴光频率最高，精通2到3个，一个好体系，综合起来，可以冲击JMB或者JACS。呵呵，外面雨终于小了，希望大家看了能有帮助，会一个技术不难，难的是把这个技术用的恰到好处，下次有时间聊聊怎么综合使用各种工具分析问题

动力学系列之四－能量分析

明天就是5.1，可是手里还是一大堆的活，实验也是诡异的要命，时光真是快啊。前面对两个基本软件介绍了一下，本来想讲讲charmm和macromodel，后来想想，没有什么意思，还是讲一下分析方法对大家更适用，如果大家真的想听这两个软件，后面我再补充。动力学分析方法多种多样，但是能量永远是不变的话题，很多问题的说明都是以能量为基础的，这里说的能量包括各种力场能量，溶剂化能，自由能等。我主要对自由能讲解一些如何进行计算。在amber中，可以通过mmpbsa来计算这个相对自由能，在免费的软件gromacs中，能量的计算（不知道有心的朋友是否自己加和过他提供的能量项，是不是缺点儿什么啊^v^）相当的粗糙，虽然提供了计算自由能的工具，但是源码中却是空，不知道到哪个版本才能实现。amber是比较好用，但是也是一个收费软件，特别是想发文章的朋友要注意了。如果是不想买amber的朋友，一个中间路线就是用gromacs跑轨迹，用autodock算能量（自由能），虽然不如mmpbsa原理上准确，但是忽略部分熵效应也在变动不大的体系也还不错。下面这个脚本是我自己写的，可以实现这一功能，希望能对大家有个帮助：
#!/bin/sh #'$1' is set for different list name #'$2' is set for the ligand name #prepare the $2.bnd in advance with command: #/~root/autodock/dist305/bin/pdbtoatm lig.pdbq| ~root/autodock/dist305/bin/atmtobnd >lig.bnd for snapshot in `cat $1` do rec=protein_mn_hoh1.${snapshot} lig=$2.${snapshot} echo ${rec} #assign the atomic solvation parameters to PDBQ-formatted version of your macromolecule, ~root/bin/addsol ${rec}.pdbq ${rec}.pdbqs ~root//bin/autotors -A $2.bnd ${lig}.pdbq ${lig}.out.pdbq > ${lig}.out.${rec}.epdb.com ~root/autodock/dist305/bin/autodock3 -p ${lig}.out.${rec}.dpf -l ${lig}.out.${rec}.epdb.log -c > $1_Binding ~root//bin/Docked.awk ${lig}.out.${rec}.epdb.log >> $1_Docked ~root//bin/Inter_score.awk ${lig}.out.${rec}.epdb.log >> $1_Inter_Energy ~root//bin/Intra_score.awk ${lig}.out.${rec}.epdb.log >> $1_Intra_Energy #rm -f *.map *.fld *.xyz *.err *.pdbqs *.glg *.log *.com *.gpf *.dpf Binding.awk #!~root/bin/gawk -f BEGIN{ name=ARGV[1] split(name,array,".") } { if($0 ~ /Estimated Free Energy of Binding/) {print array[1]" "$9" " $10} } Docked.awk #!~root/bin/gawk -f BEGIN{ name=ARGV[1] split(name,array,".") } { if($0 ~ /Final Docked Energy/) {print array[1]" "$7" " $8} }

ps：这段代码作者原文就是未经排版的...

动力学系列之五－虚拟筛选

上次发了系列四的能量计算后，不少朋友pm我，既然可以用多个大分子构象计算能量曲线，反过来是否可以用一个大分子对小分子数据库做类似dock, gold那种虚拟筛选那？答案是可以的，不过这个要通过脚本来实现，虽然也有一个程序，但是那个程序是我写给ddgrid这个工程的，这半年暂时不能放出来，年底给大家共享（主要原因是怕阿）。言归正传，不过这个脚本也可实现一样的功能，就是时间上慢一些，空间上占用多一点，以前在别的地方放过早期的版本，这次放一个完全版，而且包含了我们常用的参数设置。虽然和动力学这个系列关系不是非常紧密，但是在进行多个阳性物对比动力学的时候结合上一个脚本可以把大量工作简单化。废话少说，放东西先：
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
主程序, 由于论坛屏蔽原因，我用了替代字符，而且好像总是砍掉我的字符啊
不得不分开发，请大家见谅
"&&"请转换为"EOF"后使用
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
//准备文件
#!/bin/sh
receptor=receptor.mol2
rec=${receptor%.mol2}
mol2topdbqs $receptor
for ligand in `cat mol2.list`
do
lig=${ligand%.mol2}
deftors $ligand <<&&
//选项
1

c

c

&&
#if(-z ${lig}.pdbq) then
# deftors $ligand <<&&
# a
# 1
# c
# c
# &&
#endif
#vi ${lig}.pdbq << &&

#:%s/ Br/ b /g
#:%s/ Cl/ c /g
#wq
#&&

//continue
//第一步
modhalogen ${lig%.mol2}.pdbq > temp_pdbq
mv temp_pdbq ${lig%.mol2}.pdbq
mkgpf3 ${lig%.mol2}.pdbq ${rec}.pdbqs
mkdpf3 ${lig%.mol2}.pdbq ${rec}.pdbqs
vi ${rec}.gpf <<&&
:%s/npts 60 60 60/npts 80 80 80/g
:wq
&&
//第二步
vi ${lig}.${rec}.dpf <<EOF
:%s/tstep 2.0/tstep 0.2/g
:%s/qstep 50.0/qstep 5.0/g
:%s/dstep 50.0/dstep 5.0/g
:%s/ga_num_evals 250000/ga_num_evals 1500000/g
:wq
&&
//第三步 autogrid3 -p ${rec}.gpf -l ${rec}.glg echo epdb ${lig}.pdbq > ${lig}.${rec}.epdb.com echo stop >> ${lig}.${rec}.epdb.com autodock3 -p ${lig}.${rec}.dpf -l ${lig}.${rec}.epdb.log -c < ${lig}.${rec}.epdb.com #autodock3 -p ${lig}.${rec}.dpf -l ${lig}.${rec}.epdb.log -c log.list vi log.list <<&& :%s/\.\///g :wq &&
//循环总结
for result in `cat log.list`
do
score_list.awk $result >> score_list
done
rm -f log.list
sort -r +1n score_list >sort_list
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
辅助awk
----------------------------------------------------------------------------------------------------------------
#!～root/bin/gawk -f
BEGIN{
name=ARGV[1]
split(name,array,".")
}
{
if($0 ~ /Estimated Free Energy of Binding/)
{print array[1]" "$9" " $10}
}

linux的朋友可以直接把这个脚本放到/sh/bin目录下，然后设置autodock的环境变量，直接可以执行了。
顺便说一句，这个脚本不能改变金属的选择，而且注释的地方是根据不同需要可以使用的。好累啊，这么分开发帖子，以后一定写在文件里面上传
明天要出去玩了，回来把前一时间做的膜蛋白的一些流程和辅助工具和大家共享

动力学系列之六－膜蛋白的模拟

今天是假期的最后一天（不过好像大家都是8号上班了，我们特殊，不过后面是惨痛的10天连续工作），整理一下思路，写一些膜蛋白的模拟，膜蛋白的工作量大，需要进行的准备工作多，做一个膜蛋白相当于两个球蛋白的工作量。我先介绍一下膜蛋白的一些基本概念：
1。膜蛋白－顾名思义，就是存在于细胞膜和细胞核（线粒体膜也有，不过不属于主流，暂时忽略之）双层脂质膜上的大型蛋白质。膜蛋白从拓扑结构可以分为几种，酶（3500个左右），7次跨膜的G偶联受体（2000个左右），离子通道（1000个左右），核的激素受体（150个左右）
2。生物膜膜蛋白可分为外周膜蛋白和内在膜蛋白，后者约占整个膜蛋白的70%～80％，它们部分或全部嵌入膜内，有的则是跨膜分布，如受体、离子通道、离子泵、膜孔、运载体（transporter）以及各种膜酶等等。象水溶性蛋白质一样，要深入了解膜蛋白的功能必须解析它们的三维结构。第一个水溶性蛋白质——肌红蛋白的三维结构的解析是由英国Kendrew于1957年用X射线衍射法完成的，从而使他获得了诺贝尔奖。随着相关技术的日益改进与发展，迄今蛋白质解析出具有原子分辨率的三维结构已达13 000个左右，而且正以2000个/年的速率递增。但是其中内在膜蛋白只有26个左右（根据2000年6月的统计），换言之，仅占所有已经解析出三维结构的蛋白质的0.2%。因此，用MD理解膜蛋白的作用机制也成为2000年后的一个热点。
3。当获得了一个膜蛋白后，该怎么下手那（废话，当然是越快越好，世界上n多个组都瞄着新的膜蛋白那，JMB和PNAS的常客啊），准备过程我简要介绍一下，细节如果大家感兴趣，我今后再详细说一下，不过确实的繁琐啊，我花了两个月时间搞这个东东：
a. 修补残基（这个都知道怎么做了，不说了）
b. 选择合适的膜，膜的种类也是很多的，所以大家不要挑花眼啊(网上有，自己down，我也可以提供一个，呵呵，别人的，不过是我们修改了一点不合理结构的，获得方法, pm我)
c. 把蛋白正确的放到膜里面，去掉重叠部分，打通蛋白内部的通道（最关键的一步，决定了是否可以正确的模拟和结果的可靠性）
d. 膜两侧包水加box
e. MD
第三步是做膜蛋白的精髓，不同的组使用的方法不同，我自己写了一个简单的Interactive的程序，供组里的人使用，这次是第一次发布到外面（dos版，linux或者sgi的，可以pm我），可以非常简单的实现这一步，不要做繁琐的处理，google上可以搜到钾通道的处理方法，8页english, 相当的麻烦，gromacs给了一个命令，但是膜不合适（他给的膜也就适合做demo，缺乏合理性）。当然，即使搞定了上述5步，强大的cpu和足够的硬盘是基本的保障，否则2006年交的作业，2007年可以才能进行分析。
明天EMBO就要开始了，牛人我知道的来了两个，都是做的最前沿的东东，一个是free enregy landscape, 另一个是structrual genomics。去听听，有什么感想后面和大家共享

动力学系列之七－同源模建


最近去听了几场embo, 觉得自己还是有很多地方需要提高啊，而且今年的embo的特色很突出，可能代表了当今的潮流吧，等embo结束，我写一个自己的感受，大家分享。动力学里面核心是生物大分子，当然有NMR和X-RAY最好，如果没有，就需要同源模建技术来构建蛋白。同源模建技术在90年代风靡一时，简单的模建在一些高档次的文章中独撑一面，2000年以来，随着技术的成熟和应用的广泛，这项技术逐渐成为一个基本工具出现在各种文献中。同源模建是一个根据模板蛋白将一级序列转成3D结构的技术总称，包括了threading和homology两种，当让后来发展的fugue啊，geneatlas等等都是以次为基础。在这两种中，又以后者重要（因为简单而且准确性高），在很多软件中都提供了相应的功能，做的最好的商业软件是INSIGHTii和modeller（这个非常傻瓜，但是结果的可控制性不高）, 免费的有spdbvierwer和一些在线服务器.。如果买了这个软件的朋友一定首选insightII做，因为审稿人很偏爱这个；没有买的朋友最好的方法是找买了的单位合作（D版发文章会给你的导师和单位带来很多麻烦），如果找不到，只好用spdbviewer了，但是这样的结构最好不要发在纯计算的杂志上，被拒的危险很大。说了这么多废话，我这次主要介绍INSIGHTii的用法和一些使用细节：
1. insightII
最强大的同源模建软件（前提是你看了所有的manual，理解了各种算法和参数的意义），可以进行各种模建和验证预测。主要分以下几个步骤：
a. input sequence，alignment（pariwise and multi alignment，主要是对类似结构要很好了解，manual的排列原则写的非常好，建议大家仔细看）
b. refine （我个人使用的笔记）
1.
Check and correct potentials, partial charges, and hydrogens[url=] before submitting a Discover/CHARMm job[/url].
① Builder from the Module pulldown
② For hydrogens, use the Hydrogen command in the Modify pulldown
③ Potential atom types and partial charges, use the commands in the Forcefield pulldown.

2.
Check and correct steric[url=] overlaps[/url]
① Bump command from the Measure pulldown
② Relieve them by moving atoms with the Torsion command in the Transform pulldown.
3.
Find unacceptable steric[url=] overlaps caused by one to three interactions or short bond lengths.[/url]
① Using the Geometry command in the Modify pulldown in the Buildermodule.
4.
End Repair
5.
Splice Repair

6.
Energy Minimization


① Relax command provides for the selective minimization of sections of the model protein.

7.
Molecular Dynamics


① Explore command in the Refine pulldown can accomplish this task.
同源模建技术给我们带来了更多选择和探索新领域的机会，但是有一点请大家注意，如果你的蛋白同源性不高，波动性过大，关键位点差异显著的话，要谨慎的选择同源模建，因为这样构造出来的蛋白很可能误导你后面的工作，因为往往同源模建在我们这里都是一个课题的第一步，决定了后面这个课题是否成功的关键。

经典好贴！！

建议开个版面叫“老生常谈”，这种原创性的超级宝贝都收到里面去，瓦咔咔。


赞LZ！

大牛猛啊。。。。

Perfect!

建议开个版面叫“老生常谈”，这种原创性的超级宝贝都收到里面去，瓦咔0咔。
.

我太同意这句话了哈哈

这是药物所的人发的，确实是好东西，不过要有点道行才看的懂啊！

确实很基础的帖子
帮助很大
谢了！