手机版

CADD计算 — 计算机辅助药物设计的分子模拟全管线实战

发布时间：2026-07-21 来源：科研学术网

小中大

字号：

CADD计算（Computer-Aided Drug Design）说的不是单一技术——是一条从靶点结构到候选化合物的完整管线。我接过一个案例：客户用CADD计算筛选了50万化合物库，挑出Top 100做实验验证，结果体外活性命中率不到5%——远低于CADD计算承诺的20-30%。问题不在CADD方法本身——在于管线中每个环节的精度损失会累积放大：对接的打分函数偏向某些化学型（scaffold bias）、MD模拟时间太短没看到配体从口袋中解离、结合自由能计算省略了熵贡献。这篇文章复盘CADD计算管线中每个环节的最优实践和常见陷阱。

一、靶点蛋白的准备——结构质量决定管线天花板

CADD计算的第一步是靶点蛋白三维结构的获取和准备。首选来源是PDB数据库中的X射线晶体结构（分辨率<2.5 Å）或高分辨率Cryo-EM结构（<3.5 Å）。同源建模（SWISS-MODEL或AlphaFold2预测结构）在CADD计算中只在无实验结构时作为备选——因为建模结构的主链偏差可能达到1-2 Å，侧链朝向的偏差更大——足以让对接程序把配体放到错误的亚口袋中。

蛋白准备工作流：用PDB2PQR或Maestro的Protein Preparation Wizard补全缺失的侧链和loop区（短缺失<5残基用Prime填补，长缺失需要同源建模修补）、确定质子化状态（H++或PROPKA在pH 7.4下预测可滴定残基的pKa）、添加氢原子并做约束能量最小化（重原子RMSD<0.3 Å，仅优化氢原子和添加的水分子，保持蛋白骨架接近晶体结构）。一个容易被忽略的细节：催化位点中His残基的互变异构体选择——对接前必须根据氢键网络手动确认His是ND1还是NE2质子化。

二、分子对接的打分函数与精度边界

CADD计算中的分子对接环节对接精度决定了虚拟筛选的命中率。AutoDock Vina是目前最常用的开源对接程序——速度快（约1秒/分子），打分函数的Pearson相关系数在PDBbind基准集上约0.55-0.65。但Vina的打分函数有明显的系统偏误：倾向于给高分子量（>500 Da）和高度疏水的分子打高分——因为疏水接触的面积贡献在Vina中被高估了约15-20%。

改进策略是在对接后加一轮基于力场的重打分（rescoring）：用AMBER或CHARMM力场对对接pose做短MD模拟（1 ns），从MD轨迹中用MM-GBSA或MM-PBSA计算结合自由能。这个两步策略将相关系数从Vina的~~0.6提升到MM-PBSA的~~0.75——但计算成本增加约100倍（从1秒/分子到约100秒/分子）。

对接的搜索算法精度取决于exhaustiveness参数。Vina默认exhaustiveness=8——对于筛选可以接受。对于lead optimization阶段的精细对接，exhaustiveness=32可以更充分地搜索配体的构象和位置空间——但CPU时间增加约4倍。对于含10+可旋转键的柔性配体，exhaustiveness至少设为24——否则搜索可能遗漏关键结合模式。

三、分子动力学验证——从静态对接到动态结合

CADD计算中最大的精度跳跃发生在从静态对接到MD验证这一步。静态对接告诉你”这个pose能不能放进口袋”——MD验证告诉你”这个pose在溶液中能不能稳定存在”。

做法是：把对接pose放入TIP3P水盒子中，加Na⁺/Cl⁻离子中和电荷并达到生理离子浓度0.15 M，做标准MD平衡（最小化→NVT 100 ps→NPT 500 ps），然后跑100-500 ns的生产模拟。核心分析指标是配体重原子的RMSD——如果RMSD在最后50%轨迹中稳定在<3 Å，说明对接pose在动力学上是合理的。如果RMSD持续漂移或超过5 Å——配体可能正在脱离口袋或迁移到错误的亚口袋。

一个典型的翻车案例：对接程序把配体放到了激酶的ATP结合口袋——打分-9.8 kcal/mol，看起来很完美。但MD跑到200 ns时配体RMSD突然跳到8 Å——配体从口袋中解离了。原因是对接的打分函数没有考虑结合位点中一个关键水分子的桥接作用——那个水分子在晶体结构中介导了配体-蛋白的氢键网络，但对接时被删掉了。重做MD时保留结合位点的结构水分子（用Dowser或WaterDock预测水分子位置），配体的RMSD稳定在2.1 Å，200 ns内未解离。

四、MM-PBSA结合自由能——精度与成本的平衡

CADD计算中MM-PBSA（或MM-GBSA）是计算结合自由能最实用的方法——比对接打分准，比FEP/TI自由能微扰快。MM-PBSA的典型精度：对于结构类似（congeneric）的配体系列，R²约0.5-0.7；对于跨化学型（diverse scaffolds）的比较，R²可能降到0.3-0.5。

MM-PBSA的实用技巧：不要对整个MD轨迹做MM-PBSA——只取RMSD平稳的后半段（如最后50 ns），每100 ps提取一帧（共500帧）。熵贡献（通过简正模分析或准简谐分析计算）的收敛极其缓慢——500帧也只能给出±3 kcal/mol的误差——如果不做熵校正直接比较焓值（ΔH_MM-PBSA），精度反而更稳定。

一个容易误判的场景：当两个配体的ΔG_bind（含熵）分别为-8.0和-6.5 kcal/mol时，它们的热力学排名可能是反的（如果排名第二的配体有更大的熵惩罚）。在CADD计算中对类似化学结构的配体做ranking时，可以先只看焓贡献ΔH，忽略熵——因为类似配体的熵差异通常<2 kcal/mol（刚性配体）或<5 kcal/mol（柔性配体）。

五、ADMET预测——别让好分子死在药代动力学上

CADD计算的最后一步是ADMET（吸收、分布、代谢、排泄、毒性）预测——筛出来一批体外活性好的分子，还要保证它们能口服吸收、不被肝脏快速代谢、没有心脏毒性（hERG抑制）。常用的预测工具有SwissADME（免费在线）和ADMET Predictor（商业软件）。

ADMET预测的核心指标：Lipinski五规则（分子量<500、logP<5、氢键供体<5、氢键受体<10）是口服药物开发的基础过滤——但它是必要条件不是充分条件。更关键的指标是Caco-2渗透性（预测肠道吸收，>10×10⁻⁶ cm/s为高渗透）和肝微粒体稳定性（t₁/₂>30 min为代谢稳定）。

一个经验：SwissADME对logP的预测在训练集覆盖的化学空间内（分子量200-500）较准（MAE约0.5 log单位），但对于含氟、硼、硅等非典型元素的分子——预测偏差可能达到1-2 log单位——因为在训练数据中这些元素出现较少。对这些特殊分子需要更高精度的量子化学方法（如COSMO-RS）。