GROMACS 蛋白质结构整理与缺失坐标修复

蛋白模拟的第一步是确认结构完整、原子身份和化学连接。pdb2gmx可以按力场补氢并生成拓扑，但不能直接为任意缺失环区预测可靠结构。本篇使用小蛋白1CRN与含血红素的1A6N，整理PDB后进入二硫键模拟和血红素建模。

1 下载并查看真实结果#

下载蛋白质、二硫键与血红素完整实例包，包含8组真实短模拟的输入、生产轨迹、TPR、EDR和绘图脚本。

unzip gromacs-protein-disulfide-heme-wave232*.zip -d protein-example
cd protein-example
sha256sum -c MANIFEST.sha256
bash view.sh
ls analysis/figures

view.sh只处理参考数据，不提交计算，也不需要下载力场。它使用集群现有的MDAnalysis2.9.0、NumPy、Matplotlib环境；可用 ANALYSIS_PYTHON 指定自己的Python。环境路径中带有qvasp名称，但这里不调用qvasp或VASP。8个实例是两种体系、两个种子、两类GPU的组合；硬件重复不是新增独立种子。

2 先看PDB里有什么#

module load ase/3.29.0 gromacs/2026.4-p100-threadmpi
head -20 inputs/1CRN.pdb
rg '^SSBOND' inputs/1CRN.pdb
rg '^HETATM' inputs/1A6N.pdb | head
ls prepared-structures
cat prepared-structures/PREPARATION.json

包里保留了RCSB的1CRN结构和1A6N结构，来源地址与哈希见 UPSTREAM-PDB.json。1CRN有46个残基和三对二硫键。这里使用的1A6N解析构建体有151个蛋白残基，保留其原编号，His93与血红素Fe相连；不要把其他肌红蛋白结构的残基数量和辅因子编号照搬过来。

结构整理脚本会做以下工作：

• 为同一残基一致地选择A替代构象，并保留没有替代标记的公共骨架原子，避免拼出A、B混合侧链。
• 保留1CRN的二硫键信息；保留1A6N的血红素重原子，准备后续对应拓扑。
• 本例移除晶体水和硫酸根，然后重新加水和离子。实际研究若依赖结构水或离子，应单独判断其作用并保留合适拓扑。
• 统一链标识便于程序识别，同时保留残基编号；血红素与蛋白通过同一个分子拓扑连接。

3 缺失坐标修复到底做了什么#

为了给出可检查的操作例子，脚本从完整1CRN中主动删去20–22号残基，以及Pro5的侧链原子，再对坐标做旋转和平移。随后用剩下的43个共同Cα原子拟合完整模板，只补回缺失原子，保留原有可见原子的坐标。

本例43个Cα拟合RMSD为 0.0003544 Å，完整重建坐标相对已知答案的RMSD为 0.0000782 Å。数值很小是因为模板来自同一个完整结构，这不是未知环区预测的精度证明。

重新做一次结构整理与已知缺失控制，可在独立子目录运行，避免覆盖已经生成的文件：

mkdir repair-demo
cp -r inputs repair-demo/
cp prepare_structures.py repair-demo/
cd repair-demo
python3 prepare_structures.py
ls prepared-structures
cat prepared-structures/PREPARATION.json
cd ..

脚本输出 crambin-known-gap-control.pdb、crambin-repaired.pdb 和 myoglobin-heme.pdb。这份修复脚本针对附件中的固定案例；要换体系，需先准备对应完整模板和残基映射。序列缺失、未知环、无法对应的侧链不能靠复制其他残基坐标解决，应先完成同源/结构建模并检查局部几何，再最小化和平衡。

4 补氢、拓扑与非标准组分#

普通标准蛋白可先这样生成AMBER99SB-ILDN/TIP3P拓扑：

gmx pdb2gmx -f prepared-structures/crambin-repaired.pdb   -o processed.gro -p topol.top -ff amber99sb-ildn -water tip3p -ignh
head -12 processed.gro
rg 'CYX|SSBOND' topol.top

-ignh丢弃输入中的氢并按所选力场重新补氢；它不会自动修复缺失的主链。生成时检查端基、组氨酸质子化和二硫键。pH、金属配位、催化位点和修饰残基可能需要具体判断。

小分子、金属中心和辅因子必须有匹配的参数来源，不能因为PDB里有原子坐标就认为已有力场。这里血红素使用官方CHARMM36端口的预定义参数和原始PHEM连接补丁，没有调用CGenFF程序为未知配体自动生成参数。下一步按二硫键完整流程加水、加离子并提交计算。