3. 蛋白配体体系建模与结构准备

结合自由能计算的第一步不是跑模拟，而是把实验结构“翻译”成力场能识别的拓扑。一个可计算体系需要同时满足四个条件：蛋白链完整、质子化态明确、配体构象合理、力场参数齐备。本节以PDB 1IEP（Abl激酶–伊马替尼复合物）为例，走完从原始结构到GROMACS可读体系的完整链路。

一、结构获取与清洗

从RCSB下载1IEP.pdb后，先做一次“体检”：检查分辨率（1.75 Å）、缺失残基、非标准残基、水分子与金属离子。常用命令：

grep -E "^HETATM|^HETNAM" 1IEP.pdb | sort -u
pdb4amber -i 1IEP.pdb -o 1iep_clean.pdb --nohyd --dry

清洗原则：保留与配体直接配位或催化必需的水和金属（如激酶中的Mg²⁺），删除远离结合口袋的结晶水；删除多余链、交替构象（altLoc）只保留占有率最高的一个；对组氨酸、半胱氨酸等可电离残基暂不处理，留待质子化态分配阶段统一决策。

二、缺失残基补全

1IEP的A链在激活环附近存在若干缺失残基。补全策略取决于缺失位置：

• 远离结合口袋的柔性环：可用Modeller或SWISS-MODEL同源建模补全，或直接截断为短连接，避免引入虚假构象。
• 靠近口袋的缺失：必须补全，否则口袋形状失真，直接影响ΔG。推荐用同源结构（如2HYY）作为模板，Modeller建模后取DOPE评分最低的模型。
• 末端缺失：若不影响口袋，可保留为带电末端；若影响，需加ACE/NME封端。

补全后务必检查：主链二面角是否落在Ramachandran允许区、侧链是否与配体冲突、缺失区是否引入不合理电荷。可用MolProbity或PROCHECK快速评估。

三、质子化态分配

质子化态错误是自由能计算中最隐蔽的系统误差来源。以1IEP为例，需重点检查：

• His：A链His361位于催化位点，通常为HIP（双质子化）或HIE；用H++或PDB2PQR在pH 7.4下预测。
• Cys：Abl激酶无关键半胱氨酸配位，但需确认无游离Cys被误质子化。
• Asp/Glu：生理pH下默认去质子化，但埋藏在疏水环境中的Asp需单独评估。
• 配体：伊马替尼的哌嗪环N在pH 7.4下为质子化态（pKa≈8.1），必须带+1电荷。

推荐流程：pdb2pqr --ff=AMBER --with-ph=7.4 1iep_fixed.pdb 1iep_pqr.pdb，再人工核对每个可滴定残基的预测结果，尤其是口袋内残基。

四、配体构象生成与力场参数化

配体参数化是整条链路中最容易出错的一环。以伊马替尼为例，标准流程如下：

1. 从PDB中提取配体坐标，加氢：obabel -ipdb lig.pdb -omol lig.mol -h。
2. 用Open Babel或RDKit生成多个构象，选取与晶体结构RMSD最小的构象作为起始。
3. 计算RESP电荷：先用HF/6-31G*优化，再用antechamber -i lig.mol2 -fi mol2 -o lig_bcc.mol2 -fo mol2 -c bcc -nc 1生成AM1-BCC电荷。
4. 用parmchk2检查缺失参数，生成lig.frcmod。
5. 在tleap中组装：source leaprc.gaff2; LIG = loadmol2 lig_bcc.mol2; saveamberparm LIG lig.prmtop lig.inpcrd。

关键检查点：配体总电荷是否与质子化态一致（伊马替尼为+1）、原子类型是否被GAFF2正确识别、键长键角是否合理。若使用GAFF，需确认芳香环、酰胺键等关键基团的参数覆盖。

五、体系组装与拓扑生成

蛋白与配体分别参数化后，需合并为统一拓扑。推荐用tleap或pdb2gmx+gmx editconf组合：

# 合并蛋白与配体
cat protein.pdb lig.pdb > complex.pdb
# 生成拓扑
gmx pdb2gmx -f complex.pdb -o complex.gro -ff amber99sb-ildn -water tip3p
# 添加配体参数
# 在topol.top中include lig.itp，并在[ molecules ]段添加LIG

随后定义盒子（建议配体最小距离1.2 nm）、溶剂化（TIP3P）、加离子（0.15 M NaCl中和至电中性）。最终检查：体系总电荷为0、无原子重叠、配体未逃出口袋。

六、常见陷阱与验证

完成组装后，不要急于跑模拟。先做三项快速验证：

• 能量最小化前检查：用gmx grompp生成tpr，观察是否有警告（如原子重叠、电荷异常）。
• 配体位置验证：计算配体与晶体结构的RMSD，应小于0.5 Å；若大于1 Å，说明构象生成或对接环节有问题。
• 质子化态复核：用gmx editconf查看体系总电荷，确认与预期一致（伊马替尼+1，蛋白净电荷需计算）。

至此，一个可用于后续自由能计算的体系已经搭建完成：蛋白链完整、质子化态明确、配体参数齐备、溶剂与离子环境合理。下一阶段将在此基础上进行能量最小化与平衡，为生产模拟做准备。