MCPB.py制备金属离子力场

用MCPB.py的方法生成锌离子力场整体思路：按照tutorial（ambermd.org/tutorials/advanced/tutorial20/mcpbpy.htm）去做。用到的软件：Amber(AmberTools15orHigher)及Gaussian09制作人：月牙泉tutorial具体操作过程中遇到的问题及解决方法1、准备PDB、mol2文件制作ZN.pdb:新建文件，从原pdb文件中复制出ZN的pdb文件信息，保存成pdb格式。运行：antechamber-fipdb-fomol2-iZN.pdb-oZN_pre.mol2-atamber-pfy手动修改生成的ZN_pre.mol2文件中的的ZN离子类型为ZN（一定要注意大小写），并修改电荷为2.0，然后另存为ZN.mol2。2、加氢（网址：biophysics.cs.vt.edu）,保存top，crd文件。运行：ambpdb-p0.15_80_10_pH7.2_3spu.top-c0.15_80_10_pH7.2_3spu.crd3spu_Hpp.pdb修改3spu_Hpp.pdb中与ZN形成配位键的氨基酸类型，保存为3spu_Hpp_fixed.pdb。通过可视化软件进行修改，一般用chimera进行手动修改。修改方法：一般主要只将His修改成相应的HIE,HID类型。（后文将详细讲解）3、合成新的PDB文件运行：cat3spu_Hpp_fixed.pdbZN.pdb3spu_H.pdbKeyError:”HIE-HD1”这个错误所提示的是：HIE在δ位有H；注意：HIE和HID的区别，HIE是在ε位有H，δ位没有H；HID是在δ位有H，ε位没有H；HIP在δ、ε位都有H。解决方法：通过VMD、Chimera等确定His的类别，然后修改pdb文件，并删除多余的H。分别是HIS、HIE、HID特别注意利用加H网站加H后，要特别注意用chimera查看与zn形成配位键的His的C、N互换的情况，如果不注意那么就会使后面的结果无法进行。在对3q6x进行操作时，加氢之后只修改了HIS的类型，然后继续进行操作，到了MCPB.py-i3q6x.in-s3b这一步时就出现如下的错误：***Checktheresiduenamesprovided...Yougavetheresiduenames:['HE1','HD1','AP1','HE2','CS1','HE3','ZN1','ZN2','Z71']Databasehadthem:['91-HIE','93-HID','95-ASP','179-CYS','221-HIE','242-ZN','243-ZN','244-ZZ7']Traceback(mostrecentcalllast):File/usr/software/Amber/amber14/bin/MCPB.py,line594,inmodulepremol2fs,mcresname,1,chgfix_resids,g0x,lgchg)File/usr/software/Amber/amber14/lib/python2.7/site-packages/mcpb/resp_fitting.py,line529,inresp_fittingraisepymsmtError('Error:Thenumberoftheresiduenamesgivenis'pymsmtexp.pymsmtError:Error:Thenumberoftheresiduenamesgivenismismatchthedatabase!错误的原因及解决方法出现氨基酸数目不匹配的情况的原因是：我们在输入文件3spu.in中规定了cut_off=2.8但是在生成的3spu_Hpp.pdb中，因为120位His的C、N互换使得原来的ε-N偏离到了C的位置上，那么此时ZN与120位上的ε-N的距离就大于2.8埃了，那么就不会将120位的His算为ZN的sidechain里面了，但是可能是在MCPB.py-i3q6x.in-s2s这一步时程序会计算进去，导致了这样的错误。解决方法：通过chimera查看原来的pdb文件与生成的3spu_Hpp.pdb中互换的部分，然后手动修改pdb文件，进行保存即可。一定要在修改完pdb文件后再通过chimera或者时vmd检查确认，预防眼花修改错误了。4、生成Guassian准备文件MCPB.py-i3spu.in-s1a（注意：输入命令时不要随便加减空格键哦）生成的是侧链相关的文件，主要有**_sidechain_opt.com,**_sidechain_fc.com,**_large_mk.com等相关文件。3spu.in文件#original_pdb3spu_fixed_h.pdbgroup_name3spucut_off2.8ion_ids33523353（这里是重新合成的3spu_H.pdb文件中ZN的原子序号）ion_mol2filesZN.mol2（注意文件名一定要相同，文件名在unix系统中要求比较苛刻）gau_versiong09（有时候这个会报错的，所以加后面的就好）或者是software_versiong09(g09必须加上)#（#必须加上，如果不加会报错）报错及解决方法当运行：MCPB.py-i3q6x.in-s1a时会报错Traceback(mostrecentcalllast):File/usr/software/Amber/amber14/bin/MCPB.py,line556,inmodulecutoff,watermodel,2,largeopt,sqmopt,scchg,lgchg)File/usr/software/Amber/amber14/lib/python2.7/site-packages/mcpb/gene_model_files.py,line1718,ingene_model_filesbdedatms,libdict,autoattyp)File/usr/software/Amber/amber14/lib/python2.7/site-packages/mcpb/gene_model_files.py,line1297,inbuild_standard_modelattype=libdict[resname+'-'+atname][0]KeyError:'ZN-ZN'解决方法：按照所报的错误可知是zn的原因，其实很简单，修改一下ZN.mol2中的zn的原子的大小写就可以了,统一改为ZN，还有把DU改为ZN。（当时为了解决这个错误，真是愁死个人，找了N多个网页也没有找到，最后好像是看的程序代码才半懂非懂的修改大小写，这才更解决的，当时可花费了不少不少的时间呀）5、运行G09：g093spu_sidechain_opt.com3spu_sidechain_opt.logg093spu_sidechain_fc.com3spu_sidechain_fc.logformchk3spu_sidechain_opt.chk3spu_sidechain_opt.fchkg093spu_large_mk.com3spu_large_mk.log运行g093spu_sidechain_opt.com3spu_sidechain_opt.log出错一般出的错误是会在3spu_sidechain_opt.log中显示出来。检查这个运行是否正常运行结束的方法有三个：1、看时间，这个运行时间大约在30个小时左右（体系有124个原子，调用了30个线程大约是30个小时，如果你忘记修改3spu_sidechain_opt.com中的%NProcShared=30的话，那你就可以放清明假了，因为生成的com文件中默认的%NProcShared=2，所以你就可以想象这得算多久了，虽然这个不是倍数级别的累加，但是也要运行个四五天的样子，如果后面的两个com文件中都没有修改，那么恭喜你，你可以乐乐呵呵的放个五一假期了），如果提前结束，那么你就要查看log文件中最后的几行信息了。2、运行命令echo$?如果返回值是0，那么就是正常结束，表明没有什么错误（这个命令在Unix中是通用的，专门查看命令是否正常结束哦），如果非0，那么你就必须查看log文件啦。3、命令一结束，直接查看log文件，看最后几行有没有令人讨厌的Error！！如果有,恭喜你，又中奖了，那么一般会是：Errorininternalcoordinatesystem.ErrorterminationviaLnk1ein***/g09/l103......(***表示的是g09的路径)（有时候l103也可能是l202）错误原因及处理方法如果在3spu_sidechain_opt.log文件最后出现：Errorininternalcoordinatesystem.ErrorterminationviaLnk1ein***/g09/l103......说明你的体系是不合理或者说是有错误的。那么检查的方法就是用gaussview查看3spu_sidechain_opt.com文件，这时候你会看到，其中的成键状态是比较乱的，主要是单双键成键不规则，还有环的成键不规范，这时候你就需要修改过来了，然后保存成.com文件。同样报l202错误时也是这个原因，一定要仔细看看成键的情况哦。保存com文件后的注意事项因为gaussview目前实验室里只有window版本，所以在保存转移到centos系统里之后，还需要修改以下地方：1、在com文件的第一行输入$RunGauss（注意这个是单独成为一行，成为新的第一行）2、修改%Chk=后面的路径，一般修改为%Chk=3q6x_sidechain_opt.chk，因为我们习惯上在运行终端命令时，会进入到该文件所在的文件夹下进行操作，所以可以把=后面的路径直接去掉，直接跟文件就好。这里所说的修改%Chk=后面的路径指的是在window下的保存时自动生成的路径，所以在centos下自然找不到这个路径，所以就直接删除跟上文件名就OK啦。3、当然是修改前面提到过的%NProcShared=，我习惯上设置为30，因为咱们小超级计算机的线程数是32，留下两个差遣干其它活，如果你设置为28也没什么问题，时间上也不会相差太多，因为g09的效率不是与线程数成线性关系的，所以28与30也差不了太多。6、生成最后的文件MCPB.py-i3spu.in-s2s生成3spu.frcmod文件MCPB.py-i3spu.in-s3b生成相应的mol2文件MCPB.py-i3spu.in-s4生成3spu_mcpbpy.pdb、3spu_tleap.in文件tleap-s-f3spu_tleap.in3spu_tleap.out生成用于tleap的prmtop、inpcrd文件运行到现在，出现的问题就比较少了如果不放心，最好检查以下这个模型的构建效果。(我在这里没有出现什么错误，所以就把网页上的东西稍加整理放进来了)1、通过vmd进行检查。运行：vmd-parm73spu_solv.prmtop-rst73spu_solv.inpcrd选择Graphics---Representations输入“sameresidueaswithin3ofresnameZN1”然后在DrawingMethod中选择Licorice然后点击Apply,查看我们生成的ZN是否有四个键相连，如果是，那么就要成功了。2、运行cpptraj-p3spu_solv.prmtop如果没有error那么就可以开开心心地进行MD啦但是如果出现了error：AtomXXXwasassignedalowermolecule#thanpreviousatom........那么我们就用cpptraj进行修补原子顺序