引言
在保证安全性和生物学功能的前提下，寻找高亲和力的配体是小分子药物开发的目标。因此，精确预测结合自由能一直是计算机辅助药物设计中的一个重要方向，而使用建立在严密物理基础上的蒙特卡罗或分子动力学模拟方法一直被认为是解决此问题最严谨的方法。然而，高昂的计算成本及采样算法和力场局限性仍阻碍着这类方法的广泛应用。过去十年间，硬件及算法的不断发展，越来越多药企将相对自由能计算工具（如薛定谔的FEP+）纳入其药物研发流程。
最近，来自德国默克集团的研究者们展示了他们从2016年起的一项对FEP+的大型前瞻性研究工作。在这项研究中，作者们旨在将FEP+应用到所有合适的内部在研药物发现项目中，以达到三个目的：（1）对计算工具进行前瞻性盲评，（2）评估药物研发中普遍存在的时间限制、资源限制和信息限制对使用该方法效用的影响，（3）对此前2015年以来FEP+加入的新功能进行基准测试。从2016年到2019年，作者前瞻性地将FEP+应用于12个靶点和23个化合物系列，进行了超过35,000个独立扰动计算。最终，作者获得了6,000多个化学实体的有效预测，合成和测试了400多个预测得到的新颖分子，产生了大量的前瞻性数据，提供了该方法在典型小分子药物发现的实际流程中准确性的详细评估。
项目中自由能计算工作流程
在过去三年中，作者建立了一个用于在项目中部署自由能计算的工作流程（图1）。首先，作者通过收集可用的蛋白结构数据和实验测得的结合亲和力，对给定的靶点和感兴趣的化合物系列使用FEP的一般可行性进行的评估。在这一阶段，作者通常需要至少一个高分辨率的共晶结构，且其配体应属于靶点化合物类似物。这一严格的要求源于作者在三个项目中的失败经验，作者尝试在没有X射线衍射结构的情况下使用同源模建模型，结果均未获得令人满意的结果。虽然作者使用同源模建模型未能成功，但当作者使用蛋白晶体结构和对接技术获得复合物结构时，两个项目中一成功一失败。与实验结果不一致的那个项目中，后来获得共晶结构表明蛋白结合位点有相当大的柔性。在第一个项目中，作者无法与实验数据取得良好的一致性。后来获得的共晶体结构在结合位点显示出相当大的蛋白质柔性，尽管预测的结合姿势与晶体结构相对相似。而成功的项目中，后续的共晶结构与预测时的一致性很强。
图1 德国默克的自由能计算流程。
图片来源：Journal of Chemical Information and Modeling
一旦获得足够的结构数据，作者将收集具有实验结合亲和力的同类配体的数据集（至少10个配体，最好20个），以及有关生化和生物物理测定的所有可用信息。注意，关于配体数据集大小的建议是一个“经验法则”，通常在早期项目中其实很难获得更大的数据集。如果配体数据集足够大，则可以根据被修饰的分子的R-基团进行拆分，因为不同位点可能准确性不同。然后根据这些数据集进行回顾性自由能计算评估，将预测值与实验值进行比较。作者将这些回顾性的计算实验称为验证研究（System Validation）。在这一阶段，通常会评估不同的输入结构和体系设置，以便找到最佳参数用于以后的前瞻性计算。实际流程中，由于时间限制这一阶段通常只能评估3种可能的模型。
针对出现的较大的异常值（| ΔGpred – ΔGexp | > 2 kcal / mol），需要进行详细的分析以便了解的其产生原因。如果获得的RMSE pw < 1.3 kcal / mol且存在较大的异常值能够充分解释的话，通常可以认为验证研究成功了。尽管FEP所需的确切准确度取决于具体应用场景，如通常认为2 kcal / mol的准确度可用于对大型化合物库打分，但作者发现验证研究中使用更严格的阈值（RMSE pw < 1.3 kcal / mol）可以更大概率的保证在前瞻性预测中获得足够准确度（RMSE < 2 kcal / mol）。
理想情况下，如果数据集的动态范围合适则FEP预测也应产生良好的排名。但实际上，经常只有动态范围有限的数据集可用。在这种情况下，作者以“试用”方式进入了前瞻性预测阶段。作者预测了所有正在合成的分子，并在合成以后评估这些前瞻性预测的准确性。最后，作者再决定是否在生产模式（Production）下将FEP应用于该项目。
在成功完成验证阶段之后，FEP项目进入生产模式，对复合物构象进行前瞻性计算。这些新化合物必须与验证阶段的化合物足够相似。对于新的骨架和新的晶体结构信息，必须进行新的验证研究。作者在整个项目中密切监视前瞻性预测的准确性，并跟踪已合成了哪些化合物。所有数据存在数据库中，使用自动化的工作流程。根据作者的经验，这种对合成分子的持续监控以及及时更新的前瞻性预测准确性对于在项目团队中建立最初的信任以及随后监测衍生化合物何时超出模型的适用范围是必不可少的。
内部药物研发项目中的
FEP可行性和验证性研究
在三年的过程中，作者评估了28个靶点上FEP的可行性（图2 A）。作者对18个靶点和44个化合物系列进行了验证研究，并对14个靶点和25个化合物系列进行了前瞻性计算。大部分无法进行验证研究的靶点的主要原因是缺乏相关的结构数据（7个靶点）。总体而言，一旦有足够的结构和结合亲和力数据可用于进行验证研究，作者观察到FEP验证研究的失败率相对较低。图2 B显示了18个靶点在验证研究中的准确性。总的来说，作者获得了针对14个靶点和21个化合物系列的高精度（RMSE pw <1 kcal / mol）和可接受精度（RMSE pw <1.3 kcal / mol）的预测。在该计划的早期阶段，作者判断验证研究成功的标准较宽，因此了RMSE pw大于1.3 kcal / mol的一些系列也进入了生产模式。后续发现由于验证研究中准确率较低总会导致前瞻性预测准确度更低。预测准确率不仅在不同靶点之间而且在同一靶点蛋白的不同化合物系列之间都有差异。此外，在在研项目中前瞻性地使用FEP时，作者经常会面临各种各样可能影响方法准确性的挑战。图2 C显示了对这些挑战的定性评估。几乎所有项目都有至少一个方面可能影响在其中应用自由能计算。毫无疑问，现实生活中的药物发现项目并不是理想的案例方案。
在作者尝试进行验证研究的项目中，作者遇到的最常见挑战是至少一部分配体的结合模式存在不确定性，并且由于怀疑的蛋白质构象变化而导致蛋白质结构存在不确定性（分别为66％和44％） 。在六个项目中，作者发现实验数据的来源会影响作验证研究是否成功的判断。在一种情况下，作者最初将预测的亲和力与功能分析的结果进行了比较，发现偏差较大。但是，在将相同的预测亲和力与SPR数据进行比较时，作者发现了很好的一致性，因此决定将系列推向生产模式。在四个项目中，作者发现小分子力场的参数可能无法准确描述相互作用。在其中两个项目中，使用更高版本的OPLS3e力场提高了验证研究中使用的配体组的准确性。在一个项目中，力场的变化与部分电荷有关。在另一个项目中，该化合物具有取代脂肪环。而最近版本的Force Field Builder中对该环中扭转电势进行重新参数化，提高了预测的准确性。然而对于三个出现较大离群值的项目，作者使用较新版本的力场进行重新计算也并未改善结果。
图2 内部项目的FEP可行性、验证结果和挑战。（A）对28个靶点进行FEP评估的结果。（B）使用FEP +进行验证研究的结果。（C）在所有项目中都可能遇到预测准确性方面的挑战。
数据来源：Journal of Chemical Information and Modeling
内部项目的前瞻性FEP+预测结果
对于12个靶点的19个化合物系列，作者获得了包含至少五个数据点的前瞻性预测数据集（图3）。与验证研究相比，有几个普遍的原因可以解释前瞻性中较大的误差。首先，前瞻性数据集大于原始验证集。样本量越大，RMSE的估计就越可靠。其次，在整个项目过程中，新设计的化合物往往与验证研究中使用的化合物不太相似，与晶体结构中配体也不太相似。这导致在配体的结合方式和质子化状态方面有更高的不确定性。第三，作者发现作者内部化合物库中的多种分子仍然是小分子力场面临的挑战。对于几乎每个新的化学系列，即使使用了具有非常高的扭转电势的新型OPLS3e力场，作者也必须重新设定一些扭转电势。第四，使用前瞻性自由能预测时，作者倾向于关注极端预测（例如排名最高的化合物）。可以通过使用选择偏差校正来减轻关注极端预测值产生的偏差。然而，作者发现这对最大离群值的影响不大。最后，许多作者非常感兴趣的化学结构修饰都对该方法产生了固有的挑战（例如，从芳环系统到脂肪链、电荷变化、使用柔性连接结构的引入新官能团等）。尤其是最后一种对构象采样产生了巨大的挑战，但实际上它经常出现在早期的苗头化合物优化和片段优化中。
图3 来自12个靶点的19个化学系列的前瞻性FEP +结果。数值上下角标为90%置信区间。
数据来源：Journal of Chemical Information and Modeling
将从FEP+获得的排名与从Glide对接和Prime MMGB-SA打分中获得的排名进行比较时，作者发现整体上FEP+的性能优于这些常规的基于结构的药物设计方法。在四种情况下（靶点1 /系列4、靶点4 /系列1、靶点5 /系列3和靶点6 /系列1），Prime MM-GBSA似乎产生了更好的效果，比基于FEP+的排名要好。甚至在其中两种情况下（靶点1 /系列4和靶点4 /系列1），Glide打分也好于FEP+。但是由于数据集都非常小（<10个配体），因此置信区间非常大，所以很难就不同方法的相对性能得出最终结论。尽管如此，这可能仍暗示着有机会在某些情况下使用更简单的评分方法，而像FEP +这样的计算量大的方法可以用在其他方法无法对配体进行准确排名的时候。作者还将前瞻性内部数据集上的FEP+性能与通过简单的描述符（例如分子量和log P）进行排名进行比较，FEP+也优于这些“无效模型”。
新基准数据集的构建与评估
基于作者在无法公开的内部项目中自由能计算的丰富经验，作者决定构建一个新的基准数据集，由8个具有挑战性的、最近发布的、具有成药潜力的相关靶点及小分子组成，总共包含264个配体。此基准数据集中的蛋白质靶点和配体化学空间可以代表作者内部项目的情况，但配体与内部数据集没有重叠。总体而言，此数据集很好地说明了作者在内部项目中面临的许多挑战。与先前发布的基准数据集相比，配体结构变化包括分子总电荷和电荷分布变化、及开环和母核跃迁（图4）。同时，配体组还显示出结构多样性略有增加。总体而言，作者在这些数据集实现了良好的相关性（图5）。鉴于此数据集中包含的挑战，这是一项了不起的成就。分析不同类型的结构变化时，作者发现涉及净电荷或电荷位置变化或分子核心/骨架变化的转化时会显示出较低的准确性。FEP+软件已经认为这些变化从本质上讲更加难以计算并使用了特殊的采样设置。尽管如此，出现这些变化时，计算结果仍显示较大的偏差。有趣的是，对于其余类型的结构变化，作者未发现误差与变化程度（即改变的重原子数）的强烈相关性。
图4 新基准集中不同类型的转换的示例。（A）在Eg5中添加柔性链。（B）HIF-2α中的闭环转化。（C）SHP-2中带电荷的胺的移动。
数据来源：Journal of Chemical Information and Modeling
图5 新基准集的FEP +结果。
数据来源：Journal of Chemical Information and Modeling
有趣且令人震惊的是，根据作者此前得到的验证实验成功标准，作者只能将基准中的8个化合物系列中的一个判为成功。作者发现将λ从5ns提高至20ns可以降低RMSEpw，从而可以使成功的系列扩大为3个，但这对预测亲和力和实验亲和力之间的相关性没有影响。为了研究FEP+的最新版本是否在新基准数据集上会有更好的性能，作者使用Schrödinger 2020-1版重新计算了PFKFB3和SYK数据集。对于PFKFB3，作者获得了略微改善的性能，对于SYK其准确率较低。接下来，作者将更详细地讨论c-Met和SYK 的FEP+结果。
c-Met数据集上的准确率适中。对于c-Met，FEP+预测重现了从氨基甲酸酯结构变为各种芳杂环的构效关系（图6）。它成功地预测出嘧啶比两种噻唑、咪唑、恶二唑和哒嗪活性更好。但是当从吡啶变为嘧啶时，它未能成功再现活性的增强。作者最初假设吡啶化合物的质子化情况的存在会负面影响与铰链区结合，因此吡啶和嘧啶化合物质子化状态的不同可能是活性差异的原因。但是使用Jaguar计算得到吡啶化合物的pKa为4.5，这使得质子化不太可能发生。对于吡啶化合物的预测值与实验值的偏差的另一种可能的解释可能是在结合位点缺乏杂环的旋转异构体的采样。由于其对称性，嘧啶化合物可以两种旋转状态结合，通过N原子接触蛋白质的铰链区。相反，吡啶只能以一种构象结合以建立这种相互作用。因此，如果在模拟过程中仅取样一种构象，则会低估嘧啶的结合力。确实，嘧啶和吡啶化合物的轨迹表明，在复杂的模拟中旋转受到阻碍，而在溶剂中两种状态均被采样。
数据来源：Journal of Chemical Information and Modeling
对于SYK基准案例，作者发现准确度较低，这导致预测的亲和力与实验亲和力之间的相关性较差。一组离群值与化合物CHEMBL3265015相关。作者观察到，涉及该化合物的六个扰动计算中有四个显示绝对误差> 1 kcal / mol（一个示例在图7 C中显示）。化合物CHEMBL3265015在邻位的苯环上带有甲基，这可能会影响环的旋转。与c-Met情况相似，在复合物中沿该扭转角的采样不完整。从化合物CHEMBL3265009到CHEMBL3265003转化时，作者还注意到一个大的离群值（图7 C）。在这里，分子是通过两个延伸到溶剂中的芳环生长的。同样，根据实验，从CHEMBL3264999到CHEMBL3265003的预测相对亲和力被高估了2 kcal / mol以上。对于每个λ窗口以20 ns采样进行计算时，观察到了相同的结果。但是，FEP+却准确预测了从CHEMBL3265003到大小相似的分子CHEMBL3265004的扰动。作者在其他几个项目中也观察到了类似的高估。这种对结合袋入口处生长和置换溶剂的基团的高估可能是由于水模型高估了水的流动性或低估了水与蛋白质的相互作用。
图7 SYK基准案例的FEP+结果。（A）预测亲和力与实验亲和力相关性较差。（B）| ΔΔGpred- ΔΔGexp |误差直方图。（C）误差> 1.5 kcal / mol的示例。
数据来源：Journal of Chemical Information and Modeling
FEP+对项目的影响和实操中的挑战
尽管获得了总体上令人鼓舞的前瞻性FEP +结果，并且与简单的SBDD方法相比具有明显的优势，但在整个计划中，作者已经清楚地发现，预测的有效性和获得有意义的排名所需的预测准确性在很大程度上取决于所研究的化合物系列和实验测得的亲和力分布范围。例如，对于靶点3 /系列1和靶点3 /系列3，尽管RMSE大于1.5 kcal / mol，但作者发现FEP +具有良好的排名。这有助于确定化合物的合成优先顺序，并且被项目中的化学家认为是有价值的。另一方面，在靶点1 /系列1-3中，尽管就RMSE而言具有相似的准确性，但作者没有获得任何预测性排名。根据作者的经验，在探索的化学空间对应的活性分布范围较大的项目中，可以容许相对较低的预测准确性。一旦优化达到了化学空间中的亲和力“峡谷”，此时小的结构变化不会对亲和力产生很大的影响，那么此时应用自由能计算对项目的价值就非常有限（就像靶点1 /系列1-3和靶点5 /系列1）。有趣的是，作者发现这种情况（亲和力峡谷）经常出现在（后期）先导物优化项目中。在这一计划开始之初，作者认为先导物优化应该是FEP的主要应用场景，因为通常在此阶段进行的小规模化学修饰最适合使用该方法。然而，与作者最初的预期相反，当要探索的化学空间更广阔并且活性仍然是主要的优化目标时，FEP表现更好，例如苗头-先导（hit-to-lead）优化，片段（fragment）优化。此外，在这些情况下，这些化合物的合成通常是“木桶中的最短板”，化学家面临的合成挑战限制了可以进行实验测试的分子数量。在这种情况下，对候选分子进行排序以专注于最有前景的分子被认为是非常有价值的，因此此时使用这些计算量大的方法并不会被认为没有必要。
强有力的沟通对于成功实施项目中的自由能计算也至关重要。作者亲身体会到了在项目团队中充分了解该方法的功能和局限性是多么重要。这有助于选择计算适用范围内的化合物，以FEP+进行排序。最初，作者在项目中使用FEP+时总会遇到药物设计者想要计算的分子不在计算适用范围内；在计划实施的后期阶段，作者专注于对定制库进行计算排名，这些库在设计之初就将适用范围考虑在内。这种方法的缺点是，在某些情况下此类库的计算结果会由项目团队用于设计新分子。但是这种方式下，由于未合成FEP预测的确切分子，因此很难评估该方法的影响。为避免此问题，作者则在之后使用FEP+对此类分子进行计算，以便能够通过自动工作流程评估预测的质量（请参见上文）。总的来说，作者发现当结果可以被解释或合理化时，化学家对FEP预测的接受程度更高，例如，分析相互作用或配体柔性。因此，作者建议在进行FEP预测的同时提供此类分析，尤其要关注那些被预测为最佳和最差的化合物。
总之，要在项目中充分利用自由能计算，必须仔细权衡，在预测准确性、适用范围、关键优化目标、所关注化合物系列的合成可及性和项目时间之间取得平衡。作者发现，筛选大型定制库是提供增加FEP附加值的有效方法。对于这些库，作者筛选至少50–100个设计的新化合物，一般筛选数量为实际可合成数量的5-10倍。
结论
自由能计算在制药行业中越来越常用，并且已成为计算化学家工具箱中的一项强大功能。在这里，作者描述了为在项目中使用自由能计算而建立的一般工作流程，报告了在多个内部药物发现项目中前瞻性使用FEP+计算的宝贵数据。作者同时提供了一个新的基准数据集，该基准数据集的可供其他研究者使用，可能推动相关方法的进一步开发。除了预测的准确性外，作者还确定了影响该方法在项目中的效果的多个重要的实操因素。作者期待在不久的将来FEP+将作为一种专家工具，通过大规模计算支持相关药物研发项目。
基准数据集获取地址
www.github.com/MCompChem/fep-benchmark
参考文献
Schindler, Christina EM, et al. Large-scale assessment of binding free energy calculations in active drug discovery projects. Journal of Chemical Information and Modeling (2020). (Article ASAP)。DOI: 10.1021/acs.jcim.0c00900