Mendelevium
【笔记整理|2023-09】Amber和GROMACS分子动力学模拟实用指南
本文总结了在使用Amber、GROMACS和NAMD进行分子动力学模拟时的实用技巧、常见问题和最佳实践。
AmberTools使用经验
版本更新和兼容性
AmberTools22改进
AmberTools22解决了早期版本的Python兼容性问题
参数生成工具改进
parmchk2 vs parmchk：
- parmchk2（Amber14引入）比parmchk更优秀
- parmchk2对所有子结构进行搜索打分，比较所有参数后选择最适合的参数
- parmchk只检查某几个子结构的参数文件来获取缺失参数
# 使用parmchk2生成缺失参数
parmchk2 -i ligand.mol2 -f mol2 -o ligand.frcmod
AmberTools更新管理
# 更新AmberTools
./update_amber --update
# 检查可用的bug修复
# 参考：[Amber Bug修复页面](https://ambermd.org/BugFixes.php)：https://ambermd.org/BugFixes.php
小分子参数化
antechamber使用
# 从Gaussian输出文件生成mol2文件
antechamber -i bay.log -fi gout -o bay.mol2 -fo mol2
# acpype工具依赖关系问题
# acpype依赖于AmberTools但Amber不包含acpype
# 通过conda安装会获取另一个ambertools版本
# 解决方案：在base环境中使用pip安装
pip install acpype
GROMACS使用技巧
性能优化
GPU使用限制
GROMACS大部分体系用多GPU，和单GPU比很难获得有效的提升
- GROMACS 4.6.x后支持CPU/GPU混合模式
- 短程非键相互作用在GPU上计算，长程和键相互作用在CPU上计算
- 通过调整短程相互作用截断距离来优化GPU/CPU负载平衡
建议GROMACS版本选择
# 对于PLUMED用户，建议使用GROMACS 2022.5而非2023版本
gq says use gmx 2022.5 instead of 2023 for plumed
常见操作命令
基础模拟运行
# 能量最小化
gmx mdrun -deffnm em_tpr
# 自由能计算脚本示例
bash gmx_fep_re_sep_conti.sh WT-M132-re quick 3 2>error.log
力场和膜体系
CHARMM36力场移植
- CHARMM36 GROMACS移植讨论：https://gromacs.bioexcel.eu/t/newest-charmm36-port-for-gromacs/868/9
- 注意力场兼容性和参数一致性问题
膜体系模拟设置
推荐设置来避免生物分子跑出盒子：
# 在mdp文件中设置
comm-grps = protein
comm-mode = angular
这样可以持续消除蛋白质的平动和转动。
膜体系构建最佳实践
构建工具对比
PACKMOL的局限性
虽然也可以用Packmol构建蛋白质、核酸浸在溶剂环境中的体系，但是这样做明显不如用动力学程序自带的专用工具好，因为：
- Packmol产生的水的密度偏低
- 水的分布特征和实际体相水相差较大
- NPT模拟后盒子变形、收缩得厉害
- 可能出现溶质与其镜像最近距离太近的问题
推荐构建方法
使用MD程序专用的溶剂化工具：
# GROMACS推荐使用gmx solvate
# 使用事先NPT平衡好的溶剂盒子（如spc216.gro）
# 通过平移复制来填充真空区，溶剂分布更理想
Amber膜体系构建
可用工具和力场
构建Amber膜体系的工具选择：
- AMBAT：Amber自带工具
- CHARMM-GUI：图形界面，支持多种力场
- DABBLE：第三方工具
- PACKMOL-Memgen：最新推荐工具
LIPID21力场：
LIPID21 is the latest and recommended lipid force field.
力场兼容性
Stockholm lipids (SLipids)：
Parameters are available for saturated and unsaturated PC, PS, PE, PG, SM lipids and cholesterol.
They are supposed to work with AMBER99SB/AMBER99SB-ILDN/AMBER03/GAFF FF for proteins and small molecules.
在CHARMM-GUI中使用Amber力场
回答”setup a lipid bilayer full of popc in Amber force field with charmm-gui”的问题：
在Force Field Options步骤中可以选择Amber力场，这是在任何构建模块的最后一步（通常是输入生成步骤）。
磷脂分子理解
sn-2位置含义
sn-2 hydrocarbon in phospholipid指磷脂分子甘油骨架上第二个碳原子所连接的脂肪酸链。
sn来自stereochemical numbering（立体化学编号），用于区分甘油分子的三个碳原子位置。
高级功能和技巧
牵引和约束
GROMACS Pull Code
使用pull code在配体和脂质双分子层质心之间添加距离约束：
# 在mdp文件中设置pull参数
pull = yes
pull_ngroups = 2
pull_group1_name = ligand
pull_group2_name = membrane_com
pull_coord1_type = distance
pull_coord1_geometry = distance
PLUMED集成
# PLUMED使用与GROMACS相同的内部单位
PLUMED internal units: the same as gromacs
# 在PLUMED中添加约束的示例
RESTRAINT ARG=d1 KAPPA=1000 AT=2.0
力场开发和修改
GROMACS力场扩展性问题
rtp文件其实并不难写，和rtf的复杂度几乎相同，扩展参数的复杂度和prm也基本相同。
问题是gmx建模的可扩展性极差，频繁更改力场文件令人难以接受，所以也没人开发自动转化为rtp等格式、自动加入gmx格式力场的程序。
解决方案
- 对非聚合物体系，暂且忍受现有限制
- 对特殊聚合物，往往需要用VMD/tleap建模再转换
- 对偶尔使用的residue，手动添加到GROMACS力场中
常见错误和解决方案
编译和安装问题
Boost库依赖
# 检查Boost版本和组件
Found Boost: /path/to/anaconda3/envs/AMBER22/lib/cmake/Boost-1.78.0/BoostConfig.cmake
(found version "1.78.0") found components: thread system program_options iostreams regex timer chrono filesystem graph
构建工具链问题
# cgenff工具编译
pyinstaller -F cgenff_charmm2gmx_py3_nx2.py
文件格式和拓扑问题
GROMACS vs Amber拓扑差异
只有GROMACS在.top文件中可能有moleculetype（Amber/NAMD：列出所有原子），
所以从其他程序转换的拓扑只能列出所有原子，使得复杂约束生成非常困难！
sed脚本处理拓扑
# 在topol.top中添加包含文件
sed -i "/\#endif/a\#include \"LIG.itp\"" topol.top
sed -i "/\#endif/a\n\#include \"LIG.itp\"" topol.top
资源和参考
官方教程
- Amber基础教程4b：https://ambermd.org/tutorials/basic/tutorial4b/
- Amber膜体系教程：https://ambermd.org/tutorials/MembraneSystems.php
- Amber高级教程16：https://ambermd.org/tutorials/advanced/tutorial16/
- Amber高级教程38：https://ambermd.org/tutorials/advanced/tutorial38/index.php
第三方资源
- AMBER antechamber指南：https://emleddin.github.io/comp-chem-website/AMBERguide-antechamber.html
- PACKMOL用户指南：https://m3g.github.io/packmol/userguide.shtml
- GROMACS伞型采样教程：https://group.miletic.net/en/tutorials/gromacs/5-umbrella/
社区讨论
总结
选择合适的MD程序和工具组合是成功进行分子模拟的关键：
- Amber: 适用于生物分子体系，参数化工具成熟
- GROMACS: 高性能，适合大规模并行计算
- NAMD: 灵活的参数控制，适合复杂体系
建议根据具体研究需求和计算资源选择最合适的工具组合。
本文基于2023年9-12月技术讨论记录整理，涵盖实际模拟中遇到的问题和解决方案