GAFF立场在Gromacs中的实现

GAFF立场在GMX软件中的实现步骤

1 小分子的GAFF参数可以用Antechamber来构建，电荷目前有2个比较合理的处理小分子，AM1-bcc和RESP电荷。AM1-bcc可以用antechamber直接得到，也可以用 chimera来实现。具体可以结合自己会的来做即可。Resp需借助第三方程序来处理得到（高斯+amebr）。拿到小分子GAFF参数后，需要用脚本来进行转化成GMX识别的参数文件，因为GMX和amber 格式的参数文件很不一样，这样就导致程序识别不了GAFF的参数文件。转化可以用amb2gmx.pl和acpype.py脚本.这样就可以得到小分子的gro和top文件。gro文件GMX可以识别，但top文件GMX还是不能识别，需要自己手动改动。需要对GMX的立场和参数文件有一定的熟悉度。最后得 gro与itp文件。实现命令如下：

a.产生小分子的高斯输入文件：（假设小分子名字为lig.mol2)

antechamber -i lig.mol2 -fi mol2 -o https://www.360docs.net/doc/4e428942.html, -fogcrt -pf y -gm "%mem=4096MB" -gn "%nproc=8" -nc 1 -gk "#HF/6-31G* SCF=tight Test Pop=MK iop(6/33=2) iop(6/42=6)"

b.用高斯计算小分子的静电势：

nohup g09 https://www.360docs.net/doc/4e428942.html, > lig.log &

c.把高斯的输出结果做resp电荷拟合：

antechamebr -ilig.out -fi gout -o lig.mol2 -fo mol2 -c resp -rn LIG

d.小分子参数的检查与坐标改动（因为高斯计算结果不会改动小分子的构想，但会改变小分子的移动，即相对运动。消除相对移动只需要把坐标改成原始坐标即可，然后再产生GAFF缺失的参数文件）

parmchk -i lig.mol2 -f mol2 -o lig.frcmod

lig.frcmod即是小分子的缺失的参数文件

e.生成小分子的amber参数文件（lig.prmtop与lig.inpcrd）

先写个leap.in脚本文件，内容如下：

source leaprc.ff99SB

source leaprc.gaff

loadamberparamslig.frcmod

lig=loadmol2 lig.mol2

check lig

saveamberparmlig.prmtoplig.inpcrd

savepdblig lig.pdb

quit

然后运行tleap命令即可：

tleap -f leap.in

这样就拿到了小分子的参数文件（lig.prmtop与lig.inpcrd）

f.调用脚本程序（amb2gmx.pl和acpype.py），尽量用acpype.py,比较稳定，需要安装python2.7版本环境，才能运行。命令如下

python acpype.py -p lig.prmtop -x lig.inpcrd -d

这样就产生不太合法的lig.gro与lig.top文件，gro文件GMX可以识别，但是top文件GMX 识别不了，这样需要对GMX的参数和立场文件非常熟悉才行，不然不知道该如何处理才能生成itp文件，其实itp文件与top文件

大概一样，但是top文件比较广泛，可以包含itp文件，小分子需要包含在蛋白质里面，所以必须把lig.top文件处理成lig.itp文件。具体做法就是除去表头，改动原子类型，再出去后面的附加信息即可。这样就拿

到了lig.gro与lig.itp2个GMX必须的输入文件。

2、拿到小分子参数之后（lig.gro与lig.itp），蛋白文件（例如为rec.pdb）下面步骤为把蛋白和小分子叠合到一起，这样就可以做分子动力学来研究小分子与蛋白之间作用。

a.处理蛋白，加水模型，忽略蛋白自带的氢原子，得到gmx的受体文件

pdb2gmx -f rec.pdb -o rec_gmx.pdb -water tip3p -ignh

b.得到小分子的pdb文件，以便合并

editconf -f lig.gro -o lig.pdb

c.合并蛋白与小分子

grep -h ATOM rec_gmx.pdb lig.pdb > rec_lig.pdb

d.把小分子参数加到蛋白里面

在蛋白的top文件加上lig.itp文件即可，还要在系统说明中加上小分子

3、加水盒子、离子(盒子大小、类型、离子类型根据研究体系而定)

editconf -f rec_lig.pdb -o rec_lig_box.gro -c -d 1.2 -bt cubic

genbox -cprec_lig_box.gro -cs spc216.gro -o com_wat.gro -p topol.top

grompp -f ions.mdp -c com_wat.gro -p topol.top -o ions.tpr

genion -s ions.tpr -o model_wat.gro -p topol.top -pname NA -np 3

这样就加好了水盒子和离子来中和体系，特殊体系，需要特殊处理，还可以调节PH和生理盐浓度等等（自己慢慢琢磨~）

4、下面就做连续的能量优化，来放松模拟的体系，因为蛋白加氢原子和加水分子都是按照特定的排序加上的，可能会产生不合理的碰撞等等

第一步优化：优化水分子和加入的离子，限制蛋白和小分子（限制文件需要自己产生，因为蛋白的限制文件有，但小分子的限制文件没有，可以用用make_ndx+genrest这个命令组合

来处理，这里面就不详细说了）

grompp -f min1.mdp -c model_wat.gro -p topol.top -o min1.tpr

mdrun -deffnm min1(程序会产生默认的名字)，但模拟的体系都比较大，所以需要并行处理，如下用mpirun 24个核来做并行计算

mpirun -np 24 mdrun -deffnm min1

第二步优化：优化蛋白侧链，限制蛋白骨架和小分子（限制文件需要自己产生，因为蛋白的限制文件有，但小分子的限制文件没有，可以用用make_ndx+genrest这个命令组合来处理，这里面就不详细说了）

grompp -f min2.mdp -c min1.gro -p topol.top -o min2.tpr

mdrun -deffnm min2(程序会产生默认的名字)，但模拟的体系都比较大，所以需要并行处理，如下用mpirun 24个核来做并行计算

mpirun -np 24 mdrun -deffnm min2

第三步优化：全局优化

grompp-f min3.mdp -c min2.gro -p topol.top–omin3.tpr

mdrun -deffnm min3(程序会产生默认的名字)，但模拟的体系都比较大，所以需要并行处理，如下用mpirun 24个核来做并行计算

mpirun -np 24 mdrun -deffnm min3

5、加热体系（300K为例），NVT过程，加热过程需要限制蛋白和小分子的过度摆动，需要做限制，限制可以利用优化过程的参数文件

grompp -f nvt.mdp -c min3.gro -p topol.top -o nvt.tpr

mdrun -deffnmnvt (程序会产生默认的名字)，但模拟的体系都比较大，所以需要并行处理，如下用mpirun 24个核来做并行计算

mpirun -np 24 mdrun -deffnmnvt

6 平衡体系，npt过程

grompp -f npt.mdp -c nvt.gro -p topol.top -o npt.tpr

mdrun -deffnmnpt (程序会产生默认的名字)，但模拟的体系都比较大，所以需要并行处理，如下用mpirun 24个核来做并行计算

mpirun -np 24 mdrun -deffnmnpt

6 prod MD，npt过程(也可以是nvt过程，根据研究需要来做，但npt比nvt要稳定)

grompp -f prod.mdp -c npt.gro -p topol.top -o prod.tpr

mdrun -deffnm prod (程序会产生默认的名字)，但模拟的体系都比较大，所以需要并行处理，如下用mpirun 24个核来做并行计算

mpirun -np 24 mdrun -deffnm prod

7 动力学结果分析（。。。）分析过程是动力学的精髓部分，不同人的出来不同的结果，要抓住模拟的假设，根据目的去分析想要的数据，来支持自己的假设。没有假设就去做动力学，基本不会得到什么好的结果。关于这一部分，你需要多看文献才能知道为什么要做分子模拟。以上是如何利用GAFF立场在GMX软件里实现分子模拟，之所以这个比较重要，因为GAFF 立场处理有机小高分子有很大的优势，但在GMX没有GAFF立场，所以需要利用各种软件组合来实现。这些以后做多了，你就会了！后面

我会教你如何实现膜蛋白的模拟，慢慢来，还要做实验。

https://www.360docs.net/doc/4e428942.html,/s/blog_3d102e0601014xld.html