GROMACS使用教程

GROMACS教程

一 Gromacs基本模拟流程 .......................................................................................... 3

1 下载pdb文件..................................................................................................... 3 2 用pdb2gmx 处理 pdb 文件 ............................................................................ 3 3 建立盒子............................................................................................................. 3 5 设置能量最小化................................................................................................. 4 6 用grompp程序进行文件处理 .......................................................................... 6 7 使用 genion 和tpr文件添加离子 ................................................................... 6 8 用fws_ion.pdb来产生能量最小化的输入文件 ............................................... 6 9 在后台运行能量最小化（在命令后加&） ..................................................... 7 二 设置位置限制性动力学模拟.................................................................................. 7 三 设置非限制性动力学模拟...................................................................................... 9

1 如何重启一个计算........................................................................................... 11 2 如何延长一个计算........................................................................................... 11 3 如何设置并行计算........................................................................................... 11 五 模拟结果分析........................................................................................................ 12

1 如何将特定帧的轨迹保存成*.pdb文件......................................................... 12 2 用ngmx观察轨迹文件（也可以用VMD观察轨迹文件) ........................... 12 3 比较常用的分析工具....................................................................................... 14

3.3 g_covar 计算斜方差 ............................................................................... 16 3.4 g_energy 能量数据作图，如压力、体积、密度等 ........................... 16 3.5 g_gyrate 测量回旋半径........................................................................ 17 3.6 g_rms 与 g_rmsdist 计算结构的RMSD 值 ................................... 17 3.7 g_rmsf 计算原子位置的根均方波动（ rmsf ） ................................. 18 3.8 do_dssp 计算模型的二级结构 .............................................................. 20 3.9 g_hbond 计算模拟过程中分子间的氢键的数目、距离或角度........ 21 3.10 g_saltbr 分析模拟中残基间的盐桥 ..................................................... 21

1

GROMACS 是一个使用经典分子动力学理论研究蛋白质动力学的高端的高效的工具。GROMACS是遵守GNU许可的免费软件，可以从以下站点下载：http://www.gromacs.org，并且可以在linux和 Windows上使用。

在本教程中，将研究一个从漏斗形蜘蛛的毒液中分离的毒素。我们将使用显性溶剂动力学的方法来进行研究。首先比较真空中和溶解的模型。我们将把毒素肽溶在水盒子里，紧接着用牛顿运动定律加以平衡。我们还将比较偿离子在显性溶剂动力学中的影响。

更全面的用法指导请参考官网的GROMACS 用户手册http://www.gromacs.org

注意：在本教程中，将要生成的gromacs（*.gro）结构文件，可以用VMD（下http://www.ks.uiuc.edu/Research/vmd/ ）查看。

2

一 Gromacs基本模拟流程

1 下载pdb文件

1OMB.pdb (http://www.rcsb.org/pdb/)

2 用pdb2gmx 处理 pdb 文件

pdb2gmx –ignh –ff G43a1 –f 1OMB.pdb –o fws.pdb –p fws.top –water spce

pdb2gmx 此命令将pdb文件转换成gromacs文件并产生拓扑文件。

-ignh 因为本pdb文件是由 NMR产生的，含有氢原子，因此用-ignh选项

忽略文件中的氢原子。

-ff 指定力场（G43a1是Gromos96力场，一个通用原子力场）。 -f 读入pdb文件，

-o 指定一个新产生的pdb文件（也可以是其它多种类型文件）的文件名。 -p 指定新产生的拓扑文件名。拓扑文件包含了所有力场参数（基于一开

始选择的力场），因此非常重要。

-water 来指定水模型研究表明SPC/E 水模型在水盒子模拟中表现最好。用

SPC/E 水模型研究长程静电相互作用较好。

#注：对于下面将要用到的任何命令，都可以使用“-h”查看该命令的使用方法，比如，对于命令pdb2gmx 可以使用： pdb2gmx –h

3 建立盒子

editconf -bt cubic –f fws.pdb –o fws.pdb –d 0.9

用上面的命令建立了一个简单的立方体盒子.

-d 决定了盒子的尺寸，即盒子边缘距离分子边缘 0.9nm （9?）。理论上在绝大多数系统中，-d 都不能小于0.85nm。

注：editconf 也可以用来进行gromacs文件（*.gro）和pdb 文件（*.pdb）的相互转化。

例如：editconf –f file.gro –o file.pdb 则将file.gro 转换为 file.pdb

现在就可以用产生的文件进行真空模拟了。真空模拟就是先能量最小化，然后进行动态模拟。

4 在盒子中放入溶剂

3

genbox –cp fws.pdb –cs spc216.gro –o fws_b4em.pdb –p fws.top

genbox命令在editconf产生的盒子基础上生成水盒子。上面的命令行指定了SPC水盒子。

genbox命令可以在给定尺寸的盒子中加入正确数目的水分子。

5 设置能量最小化

em.mdp文件：Gromacs用*.mdp 文件指定所有计算的参数。

它用最速下降法消除原子位置碰撞。编辑文件，将 nsteps 变成400。如果最小化不能收敛，就用nsteps=500 再做一次。（最小化在400步内一般是能收敛的，但不同的平台可能结果会不一样。）要重做的话，必须重新运行grompp（注意：预处理器的位置在你的机器上可能不同，用which命令来定位，即 which cpp）

em.mdp文件内容：

4

title –标题随便取（最长64 个字，简单点好） cpp –指定预处理器的位置

define –传递给预处理器的一些定义。–DFLEXIBLE 告诉 grompp将spc水

模型而非刚性 SPC包含进拓扑文件，以便用最陡下降法进一步最小化能量。

constraints –设置模型约束 integrator – steep，告诉 gompp 使用最速下降法进行能量最小化。cg则代表使

用共轭梯度法。

dt –能量最小化用不用。只在动力学模拟中用（如 md）。 nsteps –在能量最小化中，指定最大运行步数。

nstlist –更新邻居列表的频率。 nstlist = 10表示每10 步更新一次。 rlist –短程邻居列表的阈值。

coulombtype –告诉gromacs如何计算静电。PME为particle mesh ewald 法（参

见 Gromacs用户手册）

rcoulomb –指定库仑力阈值

vdwtype –告诉Gromacs如何计算范德华作用(cut-off, Shift 等） rvdw –指定 LJ 或Buckingham势能距离阈值 EM Stuff

emtol –最大的力如果小于此值则能量最小化收敛（结束）（单位kJ mol–1 nm–1）

emstep –初始步长（nm）

5

3.2 特性研究 g_confrms

要比较最后结构和初始PDB 文件的差异，用g_confrms ( 用g_confrms –h 查看详细信息 )。此程序计算两个结构的最小二乘拟合。

g_confrms –f1 1OMB.pdb –f2 md.gro –o fit.pdb

你将被提示选择一个组（两次都选（组4 ））。程序将报告RMSD 值，并产生一个输出文件（ fit.pdb ）。输出文件中包含两个位置重叠的结构。

3.3 g_covar 计算斜方差

也可用于从动态轨迹计算平均结构。如计算1ns动态模拟的后200ps的平均结构：

g_covar –f traj.xtc –s topol.tpr –b 801 –e 1000 –av traj_avg.pdb

警告- 平均结构往往较粗糙，需进一步执行能量最小化。

3.4 g_energy 能量数据作图，如压力、体积、密度等

g_energy –f md.edr –o fws_pe.xvg

首先要选择输出（*.xvg）的数据。输出文件是一个电子数据表，可以用Xmgr或Grace打

开。它是一个文本文件，在进行一些小的改动后可以用 Microsoft Excel打开。 如用上面的命令，你将看到如下结果（你的可能不同）：

如计算势能，输入“Potential”，回车 再按一次回车

我们得到一个平均势能和RMSD 的总结（单位kJ/mol ）

16

输入如下命令用Grace打开 *.xvg文件：

xmgrace -nxy fws_pe.xvg

可以到以下地址下载Grace http://plasma-gate.weizmann.ac.il/Grace/ .Grace 只能在linux和unix 上运行。如果没有Grace或Xmgr，可以作为空格分隔文件导入 MS Excel 。

3.5 g_gyrate 测量回旋半径

这个指标用于度量结构的紧密度。此程序计算某（些）原子质量与分子重心的关系。

g_gyrate –f md.trr –s md.tpr –o fws_gyrate.xvg

3.6 g_rms 与 g_rmsdist 计算结构的RMSD 值

用g_rms 计算动态模拟过程中的结构与初始结构的结构偏差。（-dt 10 选项告诉程序每10 帧计算一次）

g_rms –s *.tpr –f *.xtc –dt 10

计算与NMR结构的rmsd值用如下命令：

g_rms –s em.tpr –f md.trr –o fws_rmsd.xvg

选择4 组（ Backbone）计算最小二乘拟合。程序生成一个 rmsd 随时间变化的图（ rmsd.xvg ）。以空格分隔文件的形式导入 Excel。

选择一个 RMSD 图上平衡的范围（用g_rms计算）。例如：

17

上面的实例是一个1ns模拟（你的结果可能不同）。这个模拟需要延长到完全平衡。

3.7 g_rmsf 计算原子位置的根均方波动（ rmsf ）

与g_covar相似，此程序也可以计算平均结构。例如，计算一个 2ns（2000ps ）模拟的后500ps的平均结构，用如下命令：

g_rmsf –f traj.xtc –s topol.tpr –b 1501 –e 2000 –o traj_rmsf.xvg –ox traj_avg.pdb

上面的例子我们用 200-500ps范围计算平均结构是因为我们看到这段比较稳定（与原始结构比）。命令如下：

g_rmsf –s md.tpr –f md.trr –b 200 –e 500 –ox fws_avg.pdb

提示时选择组1“Protein”。

推荐的真空能量最小化的em.mdp文件。先用最陡下降法，再用共轭梯度法。警告！需要用pdb2gmx 重新生成拓扑文件，尤其是当你选择特定组（而非整体系统）计算平均结构时。

18

程序g_rmsf 也可以用来计算温度因子。计算的温度因子可以和 X光晶体结构的温度因子比较。

g_rmsf –s md.tpr –f md.xtc –o rmsf.xvg –oq fws_bfac.pdb

19

仍选择“Backbone” 组。

3.8 do_dssp 计算模型的二级结构

前提是你必须在电脑中（/usr/local/bin）安装了dssp 程序( http://swift.cmbi.ru.nl/gv/dssp/ )。

do_dssp –s md.tpr –f md.trr –o fws_ss.xpm

选择计算组1（Protein ）。用xpm2ps将xpm 文件转成 eps格式。然后用ImageMagick 转化程序将eps文件转成 png文件或其他格式文件。

xpm2ps –f fws_ss.xpm –o fws_ss.eps convert fws_ss.eps fws_ss.png

残基数在y 轴，时间（ps ）在x轴。看下面的NMR结构：

从上面的 dssp 图上，我们看到3个红色区域代表3 个beta 片层（下图中的黄色部分）。中间的较短区是最不稳定的。下边的图是用pymol 程序( http://pymol.sourceforge.net/ )做的。

20

3.9 g_hbond 计算模拟过程中分子间的氢键的数目、距离或角度

g_hbond –f md.trr –s md.tpr –num fws_hnum.xvg

Gromacs 4.0的缺省值为： r≤30.35nm a≤30o

用-r 和-a 标签设定其它限值。默认地，g_hbond 计算共体受体距离（r DA）。你也可以用-da no 标签改变这种计算，设置 da为no 则计算rHA距离。

3.10 g_saltbr 分析模拟中残基间的盐桥

程序会输出一系列*.xvg文件。给出-/- ，+/-（盐键）和+/+残基间的距离。

g_saltbr –f md.trr –s md.tpr

21

3.9 g_hbond 计算模拟过程中分子间的氢键的数目、距离或角度

g_hbond –f md.trr –s md.tpr –num fws_hnum.xvg

Gromacs 4.0的缺省值为： r≤30.35nm a≤30o

用-r 和-a 标签设定其它限值。默认地，g_hbond 计算共体受体距离（r DA）。你也可以用-da no 标签改变这种计算，设置 da为no 则计算rHA距离。

3.10 g_saltbr 分析模拟中残基间的盐桥

程序会输出一系列*.xvg文件。给出-/- ，+/-（盐键）和+/+残基间的距离。

g_saltbr –f md.trr –s md.tpr

21

本文来源：https://www.bwwdw.com/article/7p1w.html