摘要：本文用Flare GIST分析了PDE10A的apo结合位点，用计算的水合位点与水合自由能解释了苗头化物与先导化合物优化过程中的SAR，回溯性地证明了Flare GIST在基于结构中指导分子设计的作用。根据GIST分析结果，选择感兴趣的修饰位置，用Hit Expander自动、快速地生成衍生物，实现对高能水进行置换的设计，其中包括文献报道活性提高30-300倍的设计，还包含尚未见文献报道的、另人感兴趣的设计。Flare GIST与Hit Expander联合使用是强大且可靠的基于结构的药物设计工具。

作者：肖高铿
日期：2024-05-17

前言

在从苗头到先导的优化过程中或先导化合物的优化通常会用口袋填充策略。当填充的空腔较大时，我们可以用基团替换往配体上添加新的片段（三个或三个以上原子）；当填充的空腔较小时，可以每次往配体上添加一个原子进行细小的修饰。这种小修饰得当的话，可带来结合亲和力的巨大提升。比如，Merck公司在PDE10A抑制剂MK-8189的发现过程中，如图1所示，在苗头化合物9系列的噻唑环上引入甲基或甲氧基，带了活性高达30-300倍的提升¹。

图1. PDE10A抑制剂9系列化合物及其K_i

在先导化合物优化时，如图2所示，在R₁的吡啶环上引入甲基或甲氧基，也带来10倍的活性提升²。最终，化合物18（MK-8189）作为候选药物进入临床研究²。

图2. PDE10A抑制剂13-18的化学结构式及其K_i

本文的目的是以上述的9系列苗头化合物以及先导化合物的优化为例，演示如何使用Flare^TM的水分析工具GIST^3,4对PDE10A的apo结合位点进行分析，以利用水合位点信息建立先导化合物优化策略；然后利用Flare^TM的Hit Expander⁵对苗头化合物与先导化合物进行快速修饰、生成新的结构；最后选择满足优化策略的化合物进入下一轮测试。

结果

分析PDE10A的结合位点

图3. 化合物9h

当9系列6位取代基为丙基时的化合物9h（图3）与PDE10A的共晶结构已经解释，PDB代码为5C29。分析结合位点，如图4所示，可以发现除了9h重原子所在位置之外，整个空腔几乎都被水分子占据。其中HOH1027、927所处的子口袋为先导化合物13-18的R₁的结合位点，HOH1015、909、907所在的子口袋为先导化合物13-18的R₂的结合位点，也是9系列6位噻唑取代基的结合位点。了解这些水所在位置apo结构的水合自由能以及其它热力学性质对于基于结构的设计具有非常重要的作用。

图4. 化合物9h与PDE10A的共晶结构PDB 5C29。

PDE10A结合位点apo结构的GIST分析

采用Flare^TM GIST⁴在默认参数下分析PDB 5C29结合位点的apo结构，我们只对水密度显著高（\(Density\ge\)4）且水合自由能显著高或低（\(ΔG\ge0.5或\le-0.5\)）的水合位点感兴趣，也就是local \(
ΔG\ge0.5或\le-0.5\)的水合位点，结果如图4所示。

图4. 化合物9h与PDE10A共晶结构PDB 5C29结合位点apo结构的GIST分析结果。其中黄色化合物：9h；分子表面：PDE10A；等值面：GIST计算的水密度与热力学性质

观察水的密度图，可以发现，结合口袋里，几乎全部的水都被water density=4的等值图覆盖，这说明GIST预测的水合位点与实验值非常一致。

图5. 化合物18与PDE10A-9h共晶结构apo结合位点GIST的Water Density=4等值图。绿色分子：化合物18；飘带：PDE10A（PDB 5C29）；灰色网格状：GIST计算的density=4等值面图

将化合物18（MK-8189）与PDE10A的共晶结构PDB 8DI4按C_α原子叠合到9h与PDE10的共晶结构PDB 5C29，可以获得化合物18与PDE10A apo结合位点水的关系。如图5所示，可以发现化合物13-18系列的R₁基团（见图2）与水929、1027重合，这意味R₁的芳香基团与甲基将这两个水置换；同样，我们也清晰的看到R₂基团（见图2）与水1015重合，甲基置换了水1015。

图6. 关键水合位点的热力学性质。绿色分子：化合物18；飘带：PDE10A（PDB 5C29）；红色等值图：local unhappy water，ΔG\(\gt\)1.0kcal/mol；绿色等值图：local happy water，ΔG\(\lt\)-1.0kcal/mol。

图6呈现了上述化合物13-18系列的R₁与R₂取代基甲基化水合位点的自由能。R₁甲基置换了水1027，该甲基所在位置的水合自由能ΔG\(\gt\)1.50kcal/mol，这意味在R₁的吡啶4位引入甲基，可以提高化合物的结合自由能，这与化合物13-18系列的SAR一致（图2）：引入一个甲基或甲氧基，K_i降低一个数量级（～10倍）。

图7. 结合水1015的自由能以及9系列化合物的结合模式。飘带：PDE10A（PDB 5C29）；红色等值图：local unhappy water，ΔG\(\gt\)1.0kcal/mol；绿色等值图：local happy water，ΔG\(\lt\)-1.0kcal/mol。

化合物18的R₂取代基末端甲基（见图2）置换HOH1015的SAR可以从9系列化合物6-位取代的SAR（图1）获得。将9s与PDE10A共晶结构PDB 5C2A按C_α原子叠合到9h与PDE10的共晶结构PDB 5C29，然后在结合位点编辑获得9z、9aa、9y的结合模式，结果如图7所示。可以发现，在9系列的噻唑环上引入置换高能unhappy的HOH1015基团之后活性提高了30-300多倍。比如，9a引入甲基（9y）与甲氧基（9bb）后活性分别提高了38与100倍；9z在噻唑环上引入甲基（9s)后，活性分别提高了333倍。

综上所述，用GIST分析apo结合位点可以为基于结构的设计提供非常丰富、明确的水合位点及其热力学信息，可以定性地解释化合物18（MK-8189）的R₁与R₂的SAR，用来指导新化合物的设计。

用GIST定量计算高能水置换的自由能变化

在上一小节中，我们用GIST计算的水合位点上的ΔG等值图来定性说明哪些水是高能的“unhappy water”，哪些是低能的“happy water”，并据此来评估哪些水被配体置换后可以获得结合自由能的增加。为了了解一个配体原子占据水分子位置的去溶剂化自由能，需要求解被占据格点上ΔG的累加值。

图8. 用GIST来计算当配体原子占据结合水分子的位置时发生的去溶剂化自由能。蓝色球：配体的原子𝑖；红色网格：“unhappy water“网格k。

当配体原子占据了结合位点里水分子的位置时，会发生对应的蛋白去溶剂化自由能变化，配体获得对应的结合自由能增益，表示为ΔG_watdisp。为了评估受体在结合过程中的去溶剂化自由能，首先识别被配体原子置换的体素（Voxel），如公式1所示：

$$
V = \left\{v_{i} | v_{i} \in ligand\right\} \cdots(1)
$$

图8解释了计算方法，在配体结合到蛋白口袋的过程中，如果一个体素包含在一个配体原子的范德华半径内，则认为该体素被配体置换。累加这些体素的能量值，并将总和乘以体素的体积，得到以kcal/mol为单位的值，如公式2所示：

$$
ΔG_{watdisp} = \alpha \times vol \times \sum_{v_{i} \in V}G_{GIST}(v_{i})\cdots(2)
$$

其中α是缩放因子。在Balius等人⁶的研究中，探索了多种α值，当α=-0.5时，比起没有使用基于GIST的蛋白去溶剂化的打分方法，大部分情况下虚拟筛选的性能得到显著地提高。其中，在GIST计算时格点为0.5 Å（见方法部分，space=0.5 Å），因此这里的vol = 0.125 ų。而且，在计算的时候，对体素v_i上的G_GIST值设置了上限截断值+3kcal/mol/ų。

$$
G_{GIST} = \left \{
\begin{aligned}
& 0\ &if\ |G_{GIST}| \le 0.5 \\
& +3\ &if\ G_{GIST} \ge +3 \\
& G_{GIST} &\ otherwise
\end{aligned}
\right.\cdots(3)
$$

在这里，使用α=-1，体素v_i上的G_GIST上限截断值=+3kcal/mol/ų来计算HOH1015的ΔG_watdisp，以及9z与9s的ΔG_watdisp差值，结果如表1所示。

表1. ΔG_watdisp的计算值、及其变化与实验值的比较

单位：kcal/mol；*:水合自由能，与置换自由能符号相反。

结果表明，距离HOH1015水分子1.4 Å范围之内，包含了40个体素：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45

x          y           z          distance   dG          dG_trunced
--------   ---------   --------   --------   ---------   ----------
10.19521    13.17006    41.7274    1.39791    5.19446     3.00000
10.19521    13.17006    42.2274    1.25321    2.33735     2.33735
10.19521    13.17006    42.7274    1.29882    0.03409     0.03409
10.19521    13.67006    41.7274    1.25108    3.21279     3.00000
10.19521    13.67006    42.2274    1.08702    0.83338     0.83338
10.19521    13.67006    42.7274    1.13930    0.01786     0.01786
10.19521    14.17006    41.7274    1.29471    0.10853     0.10853
10.19521    14.17006    42.2274    1.13695    0.01672     0.01672
10.69521    12.67006    42.2274    1.28051    0.79780     0.79780
10.69521    13.17006    41.7274    1.06506    20.36070    3.00000
10.69521    13.17006    42.2274    0.86646    7.79934     3.00000
10.69521    13.17006    42.7274    0.93121    0.15240     0.15240
10.69521    13.67006    41.2274    1.27633    2.22786     2.22786
10.69521    13.67006    41.7274    0.86337    17.14060    3.00000
10.69521    13.67006    42.2274    0.60151    5.86556     3.00000
10.69521    13.67006    42.7274    0.69153    0.16267     0.16267
10.69521    14.17006    41.2274    1.31912    0.19725     0.19725
10.69521    14.17006    41.7274    0.92546    1.19197     1.19197
10.69521    14.17006    42.2274    0.68766    0.22516     0.22516
11.19521    12.67006    41.7274    1.30518    1.86438     1.86438
11.19521    12.67006    42.2274    1.14887    0.28206     0.28206
11.19521    13.17006    41.2274    1.30314    1.33371     1.33371
11.19521    13.17006    41.7274    0.90253    6.82475     3.00000
11.19521    13.17006    42.2274    0.65648    2.68514     2.68514
11.19521    13.17006    42.7274    0.73984    0.07187     0.07187
11.19521    13.67006    41.2274    1.14421    1.21156     1.21156
11.19521    13.67006    41.7274    0.65240    4.87310     3.00000
11.19521    13.67006    42.2274    0.20500    1.47881     1.47881
11.19521    13.67006    42.7274    0.39803    0.03380     0.03380
11.19521    14.17006    41.2274    1.19176    0.13567     0.13567
11.19521    14.17006    41.7274    0.73259    0.24281     0.24281
11.19521    14.17006    42.2274    0.39126    0.08861     0.08861
11.19521    14.17006    42.7274    0.51912    0.01662     0.01662
11.69521    12.67006    41.7274    1.37248    0.06875     0.06875
11.69521    13.17006    41.2274    1.37054    0.10683     0.10683
11.69521    13.17006    41.7274    0.99738    0.17787     0.17787
11.69521    13.17006    42.2274    0.78178    0.04549     0.04549
11.69521    13.17006    42.7274    0.85298    0.02345     0.02345
11.69521    13.67006    41.2274    1.22042    0.09250     0.09250
11.69521    13.67006    41.7274    0.77835    0.06222     0.06222
Total unique voxels: 40
dG : 11.200 kcal/mol
dG_trunced: 5.291 kcal/mol

其中18体素的下限截断值高于0.5kcal/mol，它们的G_GIST范围在0-20kcal/mol/ų：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23

x           y           z          distance    dG            dG_trunced
--------   ---------   --------   ---------   -----------   -------------------
10.19521    13.17006    41.7274    1.39791     5.19446       3.00000
10.19521    13.17006    42.2274    1.25321     2.33735       2.33735
10.19521    13.67006    41.7274    1.25108     3.21279       3.00000
10.19521    13.67006    42.2274    1.08702     0.83338       0.83338
10.69521    12.67006    42.2274    1.28051     0.79780       0.79780
10.69521    13.17006    41.7274    1.06506     20.36070      3.00000
10.69521    13.17006    42.2274    0.86646     7.79934       3.00000
10.69521    13.67006    41.2274    1.27633     2.22786       2.22786
10.69521    13.67006    41.7274    0.86337     17.14060      3.00000
10.69521    13.67006    42.2274    0.60151     5.86556       3.00000
10.69521    14.17006    41.7274    0.92546     1.19197       1.19197
11.19521    12.67006    41.7274    1.30518     1.86438       1.86438
11.19521    13.17006    41.2274    1.30314     1.33371       1.33371
11.19521    13.17006    41.7274    0.90253     6.82475       3.00000
11.19521    13.17006    42.2274    0.65648     2.68514       2.68514
11.19521    13.67006    41.2274    1.14421     1.21156       1.21156
11.19521    13.67006    41.7274    0.65240     4.87310       3.00000
11.19521    13.67006    42.2274    0.20500     1.47881       1.47881
Total unique voxels: 18
dG : 10.904 kcal/mol
dG_trunced: 4.995 kcal/mol

用方程（2）加合法计算这些体素的自由能贡献为10.9kcal/mol，按方程（3）取上限截断值后为5.00 kcal/mol。可以预期配体原子对HOH1015的置换可带来-5.0kcal/mol的结合自由能增加。

计算的9s、9z的ΔG_watdisp分别为-41.599、-38.233kcal/mol，它们的差值ΔΔG_watdisp=3.37kcal/mol，也就是9s的甲基在置换HOH1015后比9z多获得的结合自由能，这与它们的实验结合自由能差值3.44kcal/mol几乎一致！注意到，9s因甲基带来的水置换自由能增值（-3.44kcal/mol）未达到预期对HOH1015置换带来的自由能增值（-5.0kcal/mol），这与此前的研究结论^10-12一致，即：水分子在结合位点里的自由能与置换该水分子的配体亲合力没有直接的相关性。

分子对接研究

用分子对接打分来预测化合物的结合亲和力是常用的方法之一。GNINA是流行的分子对接软件AutoDock Vina的深度学习打分版本^7-9。用GNINA将化合物16、18在PDB 8DI4的结合口袋里进行score_only模式的对接计算，计算时将蛋白里面的水全部保留，结果如表2所示。

表2. 化合物16、18的分子对接GNINA打分结果

注：AutoDock Vina的打分项为加权之前的值；实验值pKd与ΔG_binding是根据K_i计算转化而来。

16比18少了一个置换HOH1027的甲基，16与18的主要差异体现在与蛋白的形状互补与疏水相互作用上：16的gauss项与疏水相互作用项比18的略小，因此16比18预测的结合自由能（Vina score）略高，分别为-9.22、-9.74 kcal/mol，但并没有显著的差异。相似的，GNINA深度学习打分CNNaffinity预测16、18的结合亲和力分别为8.31与8.42（pK_d），也没有显著差异。这说明，无论是经典的Vina打分、还是深度学习亲和力打分都不能解释实验观察到的活性差异。

用GNINA将化合物9z、9s在PDB 5C2A的结合口袋里进行score_only模式的对接计算，计算时将蛋白里面的水全部保留，结果如表3所示。

表3. 化合物9z、9s的分子对接GNINA打分结果

注：AutoDock Vina的打分项为加权之前的值；实验值pKd与ΔG_binding是根据K_i计算转化而来。

9z比9s少了一个置换HOH1015的甲基，在Vina Score里9z与9s的差异也仅在与蛋白的形状互补与疏水上：9z的gauss项与疏水相互作用项比9s的略小，因此9z比9s预测的结合自由能（Vina score）略高，分别为-7.87、-8.62 kcal/mol，预测的相对活性0.84kcal/mol，远远低于实验的相对值3.44kcal/mol。GNINA深度学习打分CNNaffinity预测9z、9s的结合亲和力分别为7.79、8.28（pK_d），相对值为0.49pKd，远低于实验相对值2.51（pKd）。虽然GNINA的CNNaffinity准确地预测了化合物9s的结合亲和力，但是不能明确地、显著地与化合物9z区分开。无论是经典的Vina打分、还是深度学习亲和力打分都不能确切、有信心地解释实验观察到R₁与R₂置换水HOH1027、HOH1015带来的活性差异。

两种基于结构方法的比较：分子对接打分与水合位点信息对基于结构设计指导作用的比较

首先，在业务流程上两者有着显著的区别。分子对接需要同时利用蛋白与配体的信息，先要有设计出来的化合物，然后进行分子对接打分，最后比较打分值、进行优先级排序。而水合位点信息仅需蛋白结构信息，计算水合位点及其热力学性质，然后以水为中心，产生分子设计的想法与建议，需要构思、设计新化合物该是什么样子。从这一点上看，我更倾向于用虚拟筛选与基于结构设计来区分两者：虚拟筛选方法需要有现成的化合物、然后打分、进行优先级排序；而以水为中心的设计，不直接得到化合物，而是根据水合位点及其热力学的信息去构思、设计新化合物。

在这个例子中，显然分子对接打分对9z/9s，16/18这两个分子对没有足够的区分度，不能信心充足地排序、进行优先级排序；而根据GIST给出的水合位点与自由能，利用在R₁、R₂末端两处足够大且显著的水合位点自由能，可以充满信心地用来指导新化合物的优化设计，并与实验结果匹配良好。

用Hit Expander进行自动配体生成

从上述GIST分析结果可知，在R₁与R₂取代基末端的芳香环附近有高能的水可供置换，仅需引入一个单原子或双原子小基团可以实现水分子替换。Hit Expander⁵是Flare™一个组件，它通过通过添加小官能团或进行芳香单原子替换（例如，芳香碳原子转化为芳香氮原子）快速地生成衍生物库，随后可以使用Flare的其他组件进行分析。接下来我们开始用Hit Expander⁵自动生成潜在的水分子置换设计。

如图9所示，我们选中化物9aa 6位取代基末端噻唑环上靠近HOH1015的一个碳原子，然后点击Hit Expander，选择添加甲基、甲氧基、氯、腈基以及三氟甲基。

图9. 用Hit Expander对9aa末端噻唑选中的原子进行自动修饰。飘带：PDE10A（PDB 5C29）；红色等值图：local unhappy water，ΔG\(\gt\)1.0kcal/mol；绿色等值图：local happy water，ΔG\(\lt\)-1.0kcal/mol。

结果如图10所示，Hit Expander⁵生成了5个新化合物，它自动进行了几何优化、叠合到起始化合物9aa上，我们仅需要考察结果即可。我们可以看到，新生成的5个化合物都很好地实现了对HOH1015的置换。不仅包含了Shipe等人¹报道的甲基与甲氧基取代设计，还包含了令人感兴趣但未见文献报道的设计。

图10. 起始化合物9aa（青色）以及Hit Expander自动修饰的5个结果。飘带：PDE10A（PDB 5C29）；红色等值图：local unhappy water，ΔG\(\gt\)1.0kcal/mol；绿色等值图：local happy water，ΔG\(\lt\)-1.0kcal/mol。

计算9aa以及5个Hit Expander衍生物ΔG_watdisp，结果如表4所示。与起始化合物9aa相比，5个衍生物的ΔG_watdisp都获得了-3kcal/mol以上的增加（ΔΔG_watdisp），这表明这些衍生有望具有更好的结合亲合力。

表4. 化合物9aa以及5个衍生物的ΔG_watdisp

Unit：kcal/mol

至此，Hit Expander显示出与GIST水分析绝佳搭配的潜力，与GIST水分析结果组合可提供快速、自动地分子设计。新设计的分子可以进一步用进行互补性打分、FEP等进行进一步的分析、过滤。

方法

蛋白与配体结构准备

将共晶结构下载到Flare³中，并用Protein Prep工具进行结构准备以添加氢原子、优化氢键、消除原子冲突并给蛋白结构与配体分配最佳质子化状态，任何截短的蛋白质链被封端作为蛋白质准备的一部分。

GIST水分析

用Flare V8³对“湿”结合位点Apo结构（不包含配体与金属离子）进行GIST分析，共晶结构的A链水包含在GIST分析里面，具体的GIST分析条件如下：

• Calculation method: Normal
• Ligand: None
• Grid spacing: 0.5 Å
• Grid Definition：Ligand
• Chains: A Chain, A Water
• Simulation length: 20ns
• Solvent Model: explicit TIP4Pew Water

在这个计算中，在分子动力学模拟之前将共晶结构的水链包含在GIST分析里。

ΔG_{watdisp计算}

根据GIST计算的Local happy water的dx文件，通过pyflare的python脚本计算。

结论

在本文中，用Flare GIST分析了PDE10A的apo结合位点，用水合位点与水合自由能解释了苗头化物与先导化合物优化过程中的SAR，回溯性地证明了Flare GIST在基于结构中指导分子设计的作用。还根据GIST分析结果，选择感兴趣的修饰位置，用Hit Expander自动、快速地生成衍生物，实现对高能水进行置换的设计，其中包括文献报道活性提高30-300倍的设计，还包含尚未见文献报道的、另人感兴趣的设计。Flare GIST与Hit Expander联合使用是强大且可靠的基于结构的药物设计工具。

文献

1. Shipe WD, Sharik SS, Barrow JC, McGaughey GB, Theberge CR, Uslaner JM, Yan Y, Renger JJ, Smith SM, Coleman PJ, Cox CD. Discovery and Optimization of a Series of Pyrimidine-Based Phosphodiesterase 10A (PDE10A) Inhibitors through Fragment Screening, Structure-Based Design, and Parallel Synthesis. J Med Chem. 2015 Oct 8;58(19):7888-94. doi: 10.1021/acs.jmedchem.5b00983.
2. Layton ME, Kern JC, Hartingh TJ, Shipe WD, Raheem I, Kandebo M, Hayes RP, Huszar S, Eddins D, Ma B, Fuerst J, Wollenberg GK, Li J, Fritzen J, McGaughey GB, Uslaner JM, Smith SM, Coleman PJ, Cox CD. Discovery of MK-8189, a Highly Potent and Selective PDE10A Inhibitor for the Treatment of Schizophrenia. J Med Chem. 2023 Jan 26;66(2):1157-1171. doi: 10.1021/acs.jmedchem.2c01521.
3. Flare V8. http://www.cresset-group.com/software/flare
4. Flare^TM GIST. https://www.cresset-group.com/software/flare-3drism/#GIST
5. Flare^TM Hit Expander. https://www.cresset-group.com/software/flare-hitexpander
6. Balius, T.E. et al. (2017) “Testing inhomogeneous solvation theory in structure-based ligand discovery,” Proceedings of the National Academy of Sciences, 114(33), pp. E6839–E6846. Available at: https://doi.org/10.1073/pnas.1703287114.
7. Trott, O. and Olson, A.J. AutoDock Vina: Improving the speed and accuracy of docking with a new scoring function, efficient optimization, and multithreading. Journal of Computational Chemistry, 2009, 31(2), 455–461. Available at: https://doi.org/10.1002/jcc.21334.
8. Sunseri, J., & Koes, D. R. Virtual Screening with Gnina 1.0. Molecules, 2021, 26(23), 7369. https://doi.org/10.3390/molecules26237369
9. GNINA 1.0, https://github.com/gnina
10. Haider, K. and Huggins, D.J. (2013) “Combining solvent thermodynamic profiles with functionality maps of the Hsp90 binding site to predict the displacement of water molecules,” Journal of Chemical Information and Modeling, 53(10). Available at: https://doi.org/10.1021/ci4003409.
11. Abel, R. et al. (2008) “Role of the active-site solvent in the thermodynamics of factor Xa ligand binding,” Journal of the American Chemical Society, 130(9), pp. 2817–2831. Available at: https://doi.org/10.1021/ja0771033.
12. Nguyen, C.N. et al. (2014) “Thermodynamics of Water in an Enzyme Active Site: Grid-Based Hydration Analysis of Coagulation Factor Xa,” Journal of Chemical Theory and Computation, 10(7), pp. 2769–2780. Available at: https://doi.org/10.1021/ct401110x.

获取试用，请联系我们

想要在您自己的项目上尝试Flare GIST与Hit Expander，您可申请Flare试用版并了解如何将GIST与Hit Expander方法组合应用于其他基于结构设计的项目。