首页 > 随笔 > 课前准备--分子动力学Gromacs参数的详细介绍与选择（一）

课前准备--分子动力学Gromacs参数的详细介绍与选择（一）

简书 2026-09-08 16:15 1 阅读 查看原文
分子动力学 Gromacs 计算化学

作者，Evil Genius

生化小课，CADD计算机辅助药物设计已经开启，对分子对接/分子动力学感兴趣的同学可以报名参加。

这一篇我们来汇总关于分子动力学Gromacs关键参数的理解与选择。

对于 GROMACS 分子动力学（MD），真正难的不是记住命令，而是理解 .mdp 参数为什么这么选、不同参数之间有什么关系，以及什么情况下应该改参数。

如果应用主要是 蛋白–配体分子对接 → GROMACS MD → 稳定性/结合模式/MM-PBSA分析，可以按照下面这套思路学习。

一、GROMACS MD的完整流程

典型的蛋白–配体 MD：

蛋白结构准备 → 配体参数化 → 力场 → 建立模拟盒 → 加水 → 加离子 → 能量最小化 → NVT → NPT → Production MD → 轨迹分析

GROMACS官方教程也是按照 建盒、溶剂化、能量最小化、NVT/NPT平衡、Production MD 的思路进行。

二、最重要的参数：力场

这是整个MD最核心的选择之一。

常见：

力场 特点 适合
AMBER 蛋白模拟非常常用 蛋白、蛋白–配体
CHARMM36m 蛋白/膜体系常用 蛋白、膜蛋白
OPLS-AA 小分子和蛋白均可 药物设计
GROMOS GROMACS经典体系 老体系、特定研究

对于蛋白–小分子体系，不要只选择“最常用”的力场，而要考虑配体参数是否能够和蛋白力场保持一致。

例如：

AMBER蛋白力场 + GAFF/GAFF2配体

或者：

CHARMM36m蛋白 + CGenFF配体

这一点非常重要。

三、Water Model：水模型

常见：

  • TIP3P

  • SPC

  • SPC/E

  • TIP4P

  • TIP4P-Ew

选择水模型时必须考虑力场兼容性。

例如：

CHARMM36 → 通常使用TIP3P体系
AMBER体系 → 根据具体力场/工作流选择兼容水模型

不要单独考虑“哪个水模型最好”。

应该考虑：

Force Field + Water Model + Ligand Parameters

三者是否匹配。

四、模拟盒 Box

常见：

gmx editconf

主要参数：

-d -bt

-d

表示蛋白距离盒子边界的最小距离。

例如：

-d 1.0

意味着蛋白距离盒边界至少约 1.0 nm。

常见：

1.0–1.2 nm

对于普通蛋白–配体MD，这是比较常见的选择。

模拟盒不能太小，否则蛋白可能与自己的周期性镜像发生不合理相互作用；但盒子太大又会增加水分子数量和计算成本。

Box类型

常见：

cubic dodecahedron octahedron

对于球状/近球状蛋白：

rhombic dodecahedron（菱形十二面体）

通常比较节省溶剂分子。GROMACS官方教程也展示了这种盒型，因为相同最小距离下其体积明显小于立方盒。


五、离子参数

通常：

gmx genion

主要考虑：

1. 电荷中和

例如：

Protein总电荷 = -5

那么需要加入：

5 Na+

2. 生理盐浓度

很多蛋白体系会设置：

0.15 M NaCl

例如模拟人体生理环境。

因此一般流程：

先中和 → 再加入NaCl至目标浓度

不要简单理解成“所有体系都必须0.15 M”。

六、能量最小化 EM

.mdp 中最常见：

integrator = steep

也就是：

Steepest Descent（最速下降法）

常见参数：

emtol emstep nsteps

例如：

integrator = steep emtol = 1000.0 nsteps = 50000

emtol

能量最小化停止标准。

例如：

emtol = 1000

表示当最大力下降到指定阈值附近时停止。

一般来说：

不是追求“最小化次数越多越好”，而是消除严重原子碰撞和不合理几何。

GROMACS官方教程强调，能量最小化主要用于在正式动力学前解除体系中的空间冲突和不合理几何。

七、时间步长 dt

这是非常重要的参数。

dt = 0.002

单位：

ps

所以：

0.002 ps = 2 fs

目前蛋白MD最常见：

2 fs

如果使用约束算法限制H键振动，例如：

constraints = h-bonds

通常可以使用2 fs。

如果不使用相应约束，则时间步长可能需要更小。

八、NVT参数

NVT：

Number + Volume + Temperature

即：

粒子数固定 + 体积固定 + 温度固定

主要任务：

让体系达到目标温度。

典型参数：

define = -DPOSRES integrator = md dt = 0.002 nsteps = 50000 tcoupl = V-rescale tc-grps = Protein_LIG Water_and_ions ref_t = 300 300 tau_t = 0.1 0.1

例如：

300 K

是蛋白体系非常常见的选择。

九、为什么NVT之后还需要NPT？

这是很多初学者最容易混淆的地方。

NVT

控制：

温度

NPT

控制：

温度 + 压力

NPT主要用于让：

  • 密度

  • 体积

  • 压力

达到合理平衡。

所以典型流程：

Energy Minimization ↓ NVT ↓ NPT ↓ Production MD

官方GROMACS教程也采用这种NVT → NPT → Production的标准流程。

十、Temperature Coupling

常见：

tcoupl = V-rescale

或者：

tcoupl = Nose-Hoover

V-rescale

优点：

  • 稳定

  • 常用

  • 适合平衡阶段

Nose-Hoover

更接近严格的canonical ensemble动力学。

实际蛋白MD中经常看到：

NVT：V-rescale

然后：

Production：Nose-Hoover

但具体方案需要和压力耦合、约束及研究目的整体匹配。

十一、Pressure Coupling

NPT最重要。

常见：

pcoupl = Parrinello-Rahman

常用于：

Production MD

而NPT平衡阶段也可能使用其他稳定的压力耦合方法。

对于普通蛋白溶液体系：

pcoupl = Parrinello-Rahman

是非常常见的选择。

目标：

ref_p = 1.0

即：

1 bar

十二、constraints

非常重要。

常见：

constraints = h-bonds

或者：

constraints = all-bonds

作用：

限制某些化学键的振动。

最常见的原因是：

固定高频H键振动 → 可以使用2 fs时间步长。

通常：

constraints = h-bonds

已经足够常见。

十三、长程静电：PME

蛋白–配体体系中非常重要。

常见：

coulombtype = PME

PME：

Particle Mesh Ewald

用于处理长程静电相互作用。

通常不要随便改成简单的cut-off。

典型：

coulombtype = PME rcoulomb = 1.0

具体cutoff要与所用力场推荐设置一致。

十四、van der Waals

常见：

vdwtype = cut-off rvdw = 1.0

但这里尤其要注意：

不同力场对非键相互作用的推荐参数可能不同。

所以不要简单复制网上的：

rvdw = 1.2

然后所有体系都使用。

应该：

优先遵循所选力场的官方推荐参数。

十五、Production MD最关键参数

真正的数据采集阶段：

integrator = md dt = 0.002

例如：

nsteps = 50000000

那么：

50000000 × 0.002 ps = 100000 ps = 100 ns

所以：

nsteps dt 模拟时间
500,000 0.002 ps 1 ns
5,000,000 0.002 ps 10 ns
50,000,000 0.002 ps 100 ns
100,000,000 0.002 ps 200 ns

十六、到底应该模拟多久？

这是科研中非常重要的问题。

不要简单认为：

100 ns = 标准答案

实际上应该根据体系判断。

初步研究

50–100 ns

可以观察：

  • RMSD

  • RMSF

  • 配体是否脱离

  • 结合模式变化

比较正式的蛋白–配体研究

100–200 ns

比较常见。

更深入研究

可以：

200–500 ns甚至更长

尤其是：

  • 构象变化

  • 蛋白变构

  • 配体解离

  • 多状态转换

十七、轨迹输出频率

例如：

nstxout-compressed = 5000

如果：

dt = 0.002 ps

那么：

5000 × 0.002 = 10 ps

也就是每：

10 ps

保存一次轨迹。

但需要注意：

轨迹保存太频繁 → 文件巨大。

GROMACS官方教程也特别提醒，输出频率需要根据最终分析需要决定；过于频繁会明显增加轨迹文件大小和运行负担。

十八、最常分析的MD指标

完成Production MD以后，不是看一条轨迹就结束了。

通常分析：

① RMSD

判断：

整体结构是否稳定


② RMSF

判断：

哪些氨基酸波动最大

特别适合观察：

  • 活性位点

  • Loop

  • N/C端

  • 配体结合区域

③ Radius of Gyration

即：

Rg

判断蛋白整体：

紧密程度 / 构象变化

④ Hydrogen Bonds

观察：

配体与蛋白之间氢键数量随时间变化

对于蛋白–配体体系非常重要。

⑤ SASA

Solvent Accessible Surface Area

观察：

蛋白表面暴露程度变化。

⑥ Ligand RMSD

非常重要。

普通RMSD主要看蛋白：

Protein RMSD

而配体还应该分析：

Ligand RMSD

这样才能判断：

配体是否一直稳定保持在结合口袋中。

十九、如果目的就是“分子对接 + MD”

建议先掌握这一套：

参数 推荐理解
Force field AMBER / CHARMM36m
Water 与力场匹配
Box dodecahedron
Box distance ~1.0–1.2 nm
Ion Neutralize + 必要时0.15 M
EM Steepest Descent
dt 0.002 ps
NVT 约300 K
NPT 约1 bar
Coulomb PME
Constraints h-bonds
Production 100–200 ns作为常见起点
Replicates 最好多个独立重复
RMSD Protein + Ligand
RMSF Residue level
Rg Protein
H-bond Protein–Ligand
SASA Protein / complex
MM/PBSA 结合自由能辅助评价

最重要的一点

GROMACS参数不是“背一套mdp文件”。

真正需要建立的是这条逻辑：

力场 → 水模型 → 配体参数 → cutoff/PME → 温度耦合 → 压力耦合 → 时间步长 → NVT/NPT → Production → 轨迹分析

尤其是蛋白–配体体系，最容易出问题的往往不是 md.mdp，而是：

① 配体参数化 → ② 电荷 → ③ 力场兼容性 → ④ 配体拓扑 → ⑤ 初始构象

这些问题如果前面处理错了，后面即使跑出 200 ns 的轨迹，结果也可能没有可靠的物理意义。蛋白–配体MD的标准工作流同样把配体参数化、拓扑、溶剂化、EM、NVT/NPT和Production作为连续环节处理。

生活很好，有你更好。

相关推荐