蛋白模拟分子模拟

我做过一个蛋白-配体复合物的显式溶剂分子动力学模拟项目，客户想确认他的靶点蛋白在水溶液里能不能保持天然构象，以及配体结合后哪些区域发生柔性变化。我这个项目让我把蛋白 MD 的平衡判断标准重新讲了一遍。

左侧图是蛋白与水盒子的示意图。蛋白骨架用管状表示，α-螺旋和 β-折叠清晰可见，周围浅蓝色小点是水分子。我做全原子显式溶剂 MD，把蛋白放在一个长方体水盒子里，水盒子边界距离蛋白表面至少 10 Å，避免蛋白周期性镜像相互作用。水模型我选的是 TIP3P，力场用 AMBER ff19SB，配体参数用 GAFF2 加 RESP 电荷。

右侧图是 50 ns 模拟过程中的 RMSD 和回转半径 Rg。蓝色实线是主链 RMSD，前 20 ns 快速上升，在 20-25 ns 达到平台，稳定在约 3.2 Å；绿色虚线是 Rg，迅速达到 2.4 Å 后基本平稳。这个图告诉我：体系在 25 ns 左右已经达到热力学平衡，后面的轨迹可以用来做构象采样和统计分析。如果 RMSD 一直漂移，说明力场或初始结构有问题；如果 Rg 突然下降，可能是蛋白发生塌缩或非特异性聚集。

模拟流程上，我先用 PDB 文件准备蛋白结构，加氢、补全缺失 loop、质子化状态根据 pH 7 确定；然后能量最小化 5000 步，NVT 升温到 300 K，NPT 平衡 1 ns，最后生产运行 50 ns。时间步长 2 fs，氢键用 SHAKE/SETTLE 约束。温度控制用 Langevin 阻尼，压力控制用 Parrinello-Rahman。

踩坑经验：第一个是起始结构。如果初始 PDB 分辨率低或 loop 缺失，RMSD 会假性偏高。我用了 MODELLER 补 loop，并用 Ramachandran 图检查合理性。第二个是配体电荷：很多小分子用半经验方法算电荷会偏差很大，我这次用 RESP 拟合静电势，和实验结合能更接近。第三个是模拟时长：50 ns 对球蛋白来说通常够看局部柔性，但大的构象变化需要数百 ns 到 μs，我这个项目结束时我也给客户说明了这个局限。

我的经验是，蛋白 MD 的平衡不能只看 RMSD 一个指标，要同时看 Rg、二级结构含量、氢键数目、能量漂移。如果这几个指标都平稳，才能放心做后续分析。另外，显式溶剂虽然贵，但比隐式溶剂更能准确描述水介导的氢键网络和离子效应。

更多 GROMACS/AMBER 蛋白分子动力学模拟的实战技巧，我汇总在 [蛋白模拟分子模拟](https://www.keyanxueshu.com/category/md/) 栏目里。[蛋白模拟分子模拟](https://www.keyanxueshu.com/)

判断蛋白 MD 是否平衡，除了 RMSD 和 Rg，我通常还会看三个辅助指标。第一是二级结构含量随时间变化，用 DSSP 分析，如果 α-螺旋和 β-折叠比例在平衡段波动小于 5%，说明结构稳定。第二是蛋白-配体氢键数目，如果配体和结合位点的关键氢键持续存在，结合模式可靠。第三是溶剂可及表面积 SASA，如果 SASA 突然下降，可能是蛋白塌缩或聚集。这三张图和 RMSD/Rg 一起放在补充材料里，审稿人基本不会再质疑平衡问题。另外，50 ns 对大多数酶和受体足够，但如果是 intrinsically disordered protein，建议至少 500 ns 以上才能采样到多态平衡。