TradePool：用PubChem指纹子结构池化与映射，给GNN分子性质预测提供可量化的原子归因

TradePool：用PubChem指纹子结构池化与映射，给GNN分子性质预测提供可量化的原子归因本文信息标题：TradePool：一种用于量化分子性质预测中原子归因的新型可解释框架作者：Bingwei Ni, Wanxiang Shen（申万祥）, Zhuyifan Ye* 发表时间：2025年12月22日单位：澳门理工大学（中国澳门），宁波大学药物发现技术研究院（中国浙江），浙江大学药学院（中国杭州）引用格式：Ni, B.; Shen, W.; Ye, Z. TradePool: A Novel Interpretable Framework for Quantifying Atomic Attribution Values in Molecular Property Prediction. J. Chem. Inf. Model. 2025, 65, XXX–XXX. https://doi.org/10.1021/acs.jcim.5c02225 开源代码与数据：https://github.com/nibingwei123/TradePool 摘要图神经网络的可解释性一直是化合物性质预测领域的焦点。GNN在小样本化合物数据集建模上表现良好，但现有可解释方法难以准确解释原子归因值（单个原子对模型预测贡献的定量度量），使得先导化合物优化依赖资深化学家的经验，拖慢了药物开发进程。AI生成化学空间的快速扩张需要高效的可解释AI方法，这些工具能够发现超越人类直觉的洞见，补充专家知识并显著加速优化周期。为应对这些挑战，本文提出了一种新颖的双阶段原子归因值计算框架：包括基于结构池化的模型训练和基于子结构映射的原子归因值计算。该可解释框架量化任务特定的原子归因值，在芳香性/LogP/TPSA数据集上使用GCN时，原子归因准确性（计算值与真值的一致性）分别提升30%/20%/15%，Pearson相关系数达到0.93/0.63/0.88，超越了常用可解释方法仅能达到的0–0.3。此外，该方法对模型参数变化不敏感，对化合物结构变化提供相对稳定的预测结果。核心结论子结构池化+映射实现全局可解释的原子归因，显著提升与化学真值的一致性。在芳香性、LogP、TPSA三任务上，TradePool的GCN原子归因Pearson相关0.93/0.63/0.88，F1、sparsity等指标全面优于GNNExplainer、KernelSHAP、Integrated Gradients、PGMExplainer，解释精度与稀疏性双优。低频子结构筛除（出现次数<100）可抑制过拟合，保证权重的统计显著性。对模型超参数和输入分子微扰不敏感，归因稳定性优于对照方法；但在GAT上效果一般，暴露了注意力权重与子结构加权的不匹配。 PubChem指纹提供任务无关的标准子结构集合，便于跨数据集、跨架构复用，部署与迁移成本低。背景图神经网络通过消息传递捕捉分子拓扑，在溶解度、毒性、反应性等性质预测上已成为主力。但多层聚合带来的“黑盒”问题削弱了可信度，尤其在药物优化环节，需要知道哪几个原子驱动了预测。现有解释方法存在三大痛点：局部性强，难得到全局稳定的原子归因；与化学真值偏差大，Pearson相关常徘徊在0–0.3；计算代价高或对超参数敏感。子结构层面的解释更接近化学直觉，但GNN输入并未直接包含预定义子结构，如何把“可解释的子结构权重”映射回原子，成了瓶颈。关键科学问题如何在不牺牲预测精度的前提下，将GNN的决策过程转化为“子结构→原子”的可量化归因？子结构集合应如何选择，既具普适性又能捕捉任务相关模式？归因结果能否对模型参数、输入扰动保持稳定，从而在真实药物优化中可复用？创新点双阶段框架：训练时用PubChem指纹做子结构池化，解释时把子结构权重映射为原子归因。全局归因：通过线性层权重直接量化子结构重要性，再按子结构-原子掩码汇总为原子级贡献。稳健性设计：低频子结构剔除、权重聚合、多任务对比，提升对超参数和分子扰动的鲁棒性。任务通用性：同一套指纹子结构跨芳香性、LogP、TPSA乃至药物临床分子数据集均可复用。研究内容方法详述 TradePool的核心思想是将子结构作为连接原子和分子性质的桥梁。传统GNN直接从原子嵌入池化到分子表示，丢失了化学家熟悉的官能团或子结构这一中间层信息。TradePool通过引入PubChem指纹定义的881个标准子结构，在训练时显式地学习每个子结构对预测的贡献权重，在解释时将这些权重映射回原子，从而实现全局一致、化学可解释的原子归因。数据准备与清洗为什么需要严格的数据清洗？分子数据常存在SMILES表示不规范、含盐、带电荷等问题，这些会导致同一分子有多种表示形式，影响模型训练和归因评估的准确性。研究使用RDKit 2022.09.5和MolVS 0.1.1进行标准化处理，包括SMILES规范化统一分子表示确保同一分子只有唯一的SMILES字符串，去盐处理移除分子中的无机盐（如$\ce{NaCl}$、$\ce{HCl}$）只保留有机部分，中和处理将带电荷的分子转为中性形式避免电荷状态影响特征计算，以及去重按分子骨架去除重复化合物防止数据泄漏。中和应该存疑，应该是所选pH下的状态清洗后的数据按8：1：1比例划分为训练集、验证集和测试集，这种划分确保模型在训练时不会接触测试集分子，从而真实评估泛化能力。特征工程：从分子到图原子特征（71维）：每个原子用71维向量描述，包含11类信息原子类型（43维，C、N、O、S等元素的one-hot编码）、度数（11维，原子连接的其他原子数量0-10+）、隐式价（7维，未显式表示的氢原子数）电荷（原子的形式电荷如-1、0、+1）、芳香性（是否为芳香原子）、自由基电子（未配对电子数）、杂化类型（sp、sp²、sp³等）、连接氢数（显式连接的氢原子数）、手性中心（是否为手性中心）、手性类型（R/S构型）键特征（12维）包含4类信息：键型（4维，单键、双键、三键、芳香键）、共轭性（是否参与共轭体系）、是否在环中（环状结构标识）、立体化学（E/Z构型或顺反异构）这些特征由RDKit自动计算，详见Supporting Information Table S1。子结构筛选：从881位到400+位为什么要筛选子结构？ PubChem指纹包含881个预定义子结构，但并非所有子结构都在数据集中频繁出现。低频子结构（出现次数<100）在统计上不显著，可能导致模型过拟合——模型会记住这些稀有模式而非学习真正的化学规律。筛选策略：统计每个子结构在数据集中的出现频次，设定阈值为出现次数≥100次才保留，移除低频子结构以减少噪声和过拟合风险。筛选结果（图1）：芳香性任务保留416个子结构（移除465个），LogP任务保留513个子结构（移除368个），TPSA任务保留442个子结构（移除439个）图1：三类任务的子结构出现频率热图横轴：PubChem指纹子结构位；纵向颜色深浅：出现频率占比，深色表示更常见截断低频（<100次）后，仍可覆盖大多数分子，避免稀疏噪声筛选后仍能覆盖>90%的分子，说明被移除的子结构确实是稀有模式。图1的热图显示，保留的子结构在数据集中分布相对均匀，颜色深浅代表出现频率——深色表示高频子结构，浅色表示中频子结构。数据集与标签构建研究选择了三个具有明确原子归因真值的任务，这是评估XAI方法准确性的关键——只有存在可对照的真值，才能判断模型的解释是否可信。芳香性数据集（Aromaticity）为什么选择芳香性？这是唯一具有客观真值的数据集，被多篇XAI论文用作基准测试。芳香性是分子的固有结构属性，不依赖于计算方法，一个原子是否芳香可以通过Hückel规则明确判定，因此原子归因的真值是确定的。数据集来源：本文沿用Xiong等人构建的芳香性数据集，用于检验模型在原子层面的化学可解释性，标签为每个分子中芳香原子的数量。 LogP数据集（脂溶性） LogP的化学意义：LogP衡量分子的疏水性，是药物设计中的关键参数。疏水性高的分子更容易穿透细胞膜，但过高会导致溶解度差。数据集来源：本文使用Wang等人整理的脂溶性数据集，主要来自PHYSPROP数据库与Hansch汇编数据集。原子归因真值：以Crippen方法给出的原子级LogP贡献作为真值，用于评价连续归因值与真值的一致性。额外外部集合：411个FDA批准药物与10个SAMPL6挑战分子被用作外部评估，用于检验不同方法的原子归因效果在真实药物结构上的表现。 Crippen原子贡献法是什么 TradePool把Crippen方法当作LogP任务的“原子归因真值”，这一步非常关键，因为它让“解释对不对”变成了可量化的问题。方法来源：Wildman与Crippen在1999年提出一种原子类型分类体系，用原子贡献加和来预测分子的logP与摩尔折射率（MR）。核心思想：先根据每个原子的局部化学环境把它分到某个原子类型，再把对应类型的贡献值相加得到全分子的logP。计算形式：分子的logP可写作 $\log P = \sum_{i=1}^{N} a_{t(i)}$ 其中，$t(i)$表示原子$i$所属的原子类型，$a_{t(i)}$是该类型的经验贡献系数，$N$是原子数。为什么适合作为“真值”：它天然给出每个原子的数值贡献，可直接与XAI输出的连续归因值做Pearson相关比较。 RDKit里的实现：RDKit在rdkit.Chem.Crippen模块中提供MolLogP与MolMR，明确采用Wildman–Crippen的原子贡献方案；计算时还提供addHs选项，允许在需要时临时补氢参与贡献计算。实际结果会受到芳香性判定与是否显式加氢的影响，因此同一SMILES在不同标准化流程下可能出现轻微差异。需要牢记的局限：Crippen是经验模型，主要面向中性小分子；它描述的是分子在辛醇与水相之间的分配倾向，不直接等同于带电体系的logD，也不显式建模溶剂化与构象效应。参考：Wildman, S. A.; Crippen, G. M. Prediction of Physicochemical Parameters by Atomic Contributions. J. Chem. Inf. Comput. Sci. 1999, 39, 868–873. https://doi.org/10.1021/ci990307l TPSA数据集（拓扑极性表面积）为什么TPSA重要？ TPSA是药物类药性的关键指标，能够预测药物的溶解度、渗透性和药代动力学性质。一般认为，TPSA小于140 Å²的分子更容易口服吸收；极性表面积过大的分子难以穿透肠道上皮细胞，导致口服生物利用度降低。数据与真值口径：TPSA本质上是一个基于分子拓扑的分子描述符，经典定义来源于Ertl等人的碎片贡献思想。本文将TPSA作为预测标签，并以碎片贡献法得到的原子级贡献作为归因真值，用于量化解释的正确性。临床分子集（Drug-like Compounds）为什么需要临床分子集？前两个数据集虽然有真值但分子多样性有限，临床分子集包含真实的II期及以上候选药物，骨架复杂度更高，更能测试TradePool在实际药物优化场景中的可迁移性。数据集来源：作者从ChEMBL数据库收集5800个分子量0到600、处于II期及以上临床阶段的小分子；按骨架结构划分训练、验证与测试集。该任务的训练标签与原子归因真值均由RDKit计算。表1 数据集关键信息对比 | 数据集 | 总样本量 | 训练/验证/测试 | 任务类型 | 原子归因真值 | 数据来源 | 额外测试集 | | — | — | — | — | — | — | — | | 芳香性 | 3947 | 3157/395/395 | 回归（芳香原子数） | 芳香原子标签（芳香原子为1） | Xiong等构建 | - | | LogP | 16296 | 13036/1630/1630 | 回归（辛醇-水分配系数） | Crippen原子贡献 | PHYSPROP与Hansch汇编（Wang等整理） | 411个FDA药物，10个SAMPL6 | | TPSA | 5800 | 4700/550/550 | 回归（拓扑极性表面积） | 碎片贡献法原子贡献 | 文中未详述分子来源 | - | | 临床分子 | 5800 | 按骨架划分 | 由RDKit计算的分子性质 | RDKit计算 | ChEMBL（II期及以上） | - | 模型架构与训练策略 GNN编码器选择研究实现了三种主流GNN架构，以验证TradePool的通用性： GCN（图卷积网络）：每个原子聚合其邻居原子的特征，通过加权求和更新自身表示。GCN简单高效，适合捕捉局部拓扑结构，实现基于PyTorch和DGL-LifeSci 0.3.2。 GraphSAGE（图采样聚合）：采样固定数量的邻居，使用LSTM或mean聚合器处理邻居特征。GraphSAGE可扩展到大规模图，聚合方式更灵活，超参数偏好LSTM或mean聚合器。 GAT（图注意力网络）：为每个邻居分配注意力权重，动态调整不同邻居的重要性。GAT能够自适应地关注重要邻居，超参数采用4-8个注意力头，小或零dropout。所有模型使用1-3层消息传递层，ReLU激活函数，隐藏维度在128-256之间。控制组使用传统的WeightedSumAndMax池化，TradePool组替换为子结构池化。子结构池化机制这是TradePool的核心创新。传统池化直接将所有原子嵌入求和或取最大值，丢失了子结构信息。TradePool的池化过程如下：步骤1：构建子结构-原子掩码矩阵S。矩阵维度为$N \times P$，其中$N$是分子中的原子数，$P$是保留的子结构数（416/513/442），矩阵元素$S_{ij} = 1$表示原子$i$属于子结构$j$否则为0。计算方式使用RDKit的PubChem指纹生成函数，自动识别每个原子属于哪些子结构。步骤2：创建子结构虚拟节点。对于每个子结构$j$创建一个虚拟节点$V_j$，虚拟节点的特征等于所有属于该子结构的原子嵌入之和：$V_j = \sum_{i: S_{ij}=1} h_i$其中$h_i$是原子$i$经过消息传递后的嵌入向量。如果分子不包含子结构$j$，则$V_j$为全零向量。步骤3：展平与预测。将所有子结构虚拟节点展平为一维向量 $[V_1, V_2, …, V_P]$，输入到单层线性层$\hat{y} = W \cdot [V_1, V_2, …, V_P] + b$，输出预测标签（回归任务输出标量，分类任务输出类别概率）为什么这样设计有效？子结构池化强制模型通过子结构这一中间层进行预测，使得线性层的权重$W$直接对应每个子结构的重要性。这种设计天然地将可解释性嵌入模型架构，而非事后添加。图2：TradePool双阶段工作流 (A) 总览：左侧训练阶段输入分子图与子结构掩码，右侧解释阶段输出原子归因 (B) 子结构池化：同一子结构内的原子特征求和形成虚拟节点；未包含该子结构则为零向量 (C) 归因映射：线性层得到子结构权重，按掩码回分到原子，权重累加得到原子归因热图训练超参数与优化策略优化器与学习率：使用Adam优化器自适应调整每个参数的学习率，学习率通过贝叶斯优化在验证集上搜索最优值，典型范围为1e-4到1e-3。批大小设定为128平衡内存占用与梯度稳定性，最大训练轮次为200，早停策略监控验证集MAE（回归）或准确率（分类），连续10轮无改善则停止。权重初始化使用Xavier初始化，确保每层输出方差一致，避免梯度消失或爆炸。训练稳定性技巧：采用冻结策略，训练后10%的轮次仅微调线性层保持图编码部分冻结，目的是降低梯度震荡确保子结构权重稳定可解释。必要时引入L2正则化抑制极端权重，防止单一高频子结构独占权重。超参数搜索结果（SI Table S2）：隐藏维度为128-256，层数为2-3层，GraphSAGE偏好LSTM或mean聚合，GAT采用4-8个注意力头配合小或零dropout。TradePool与控制组使用相同深度，主要区别在池化方式。原子归因计算训练完成后，如何从子结构权重得到原子归因？这是TradePool的第二阶段——解释阶段。提取子结构权重不同GNN架构的权重提取方式不同，因为它们的聚合机制不同： GCN：线性层权重矩阵$W$的每一列对应一个子结构，子结构$j$的归因值等于该列所有元素之和$A_j = \sum_k W_{kj}$，原理是GCN的聚合是简单求和，权重的和反映了子结构的总贡献。 GraphSAGE和GAT：子结构$j$的归因值等于该列所有元素的L1范数$A_j = \sum_k W_{kj} $，原理是这些模型的聚合更复杂（LSTM或注意力），权重可能有正负，取绝对值后求和更稳定。映射到原子归因有了每个子结构的归因值$A_j$，如何得到每个原子的归因值？映射公式为： $\text{原子}i\text{的归因} = \sum_{j: S_{ij}=1} A_j$ 通俗解释：查找所有包含原子$i$的子结构（即$S_{ij}=1$的子结构），将这些子结构的归因值累加，累加结果即为原子$i$的归因值。化学直觉：同一子结构内的原子获得相同的基础贡献（因为它们都属于该子结构），处于多个子结构交叉位置的原子累积多重贡献（例如苯环上的碳既属于“芳香环”子结构，也属于“C=C”子结构），这种累加方式与化学家的思维一致——一个原子的重要性取决于它参与了哪些官能团。呃，其实也可以在搞原子对子结构贡献的权重的，就有点复杂了 “正归因原子”如何定义：阈值与二值化很多指标（Accuracy、Recall、F1、Sparsity）要求先把连续归因值变成二分类标签。论文对不同方法的二值化规则如下：传统XAI方法：若原子归因值大于0，则标记为正归因；否则为负归因。 TradePool方法：由于原子归因值来自子结构权重累加，作者不直接使用0作为阈值，而是在训练集上计算一个任务级阈值：对训练集中每个分子，记录其原子归因最大值与最小值；对所有分子的最大值与最小值分别取均值；再取这两个均值的平均作为阈值。验证集与测试集沿用训练集得到的阈值。通俗解释：这个阈值更像是在训练集的归因值动态范围里取一个平均中线，用它来区分相对更重要与相对不重要的原子。它并不强制每个分子都选出固定比例的原子，只是实际结果常落在中等稀疏度区间。同时，论文也把真值归因二值化用于分类指标计算：芳香性任务：芳香原子真值标签为1，其他为0。 LogP与TPSA任务：真值原子贡献大于0标为1，否则为0。评估指标体系研究采用Wang等人提出的八项XAI评估指标，全面衡量归因质量：准确性指标： Accuracy：分类任务，正确识别正/负归因原子的比例 F1-score：精确率和召回率的调和平均，平衡误报和漏报 Pearson相关：预测归因值与真值的线性相关性（-1到1，越接近1越好）稀疏性指标： Sparsity：被标记为正归因的原子比例。理想的解释应该聚焦于少数关键原子，而非高亮整个分子 Recall：真正的正归因原子中被正确识别的比例稳定性指标： Fidelity：移除正归因原子后，预测值下降的幅度。下降越多，说明这些原子确实重要 Infidelity：移除负归因原子后，预测值上升的幅度。上升越多，说明这些原子确实有负贡献 Stability：对模型参数微调（如改变随机种子）后，归因结果的一致性 Sensitivity：对输入分子微扰（如添加甲基）后，归因结果的稳定性呃，不一定非得有下降和上升的幅度很突出的原子吧，比如烷烃，都一样？这些指标从不同角度评估XAI方法：准确性衡量解释是否正确，稀疏性衡量解释是否简洁，稳定性衡量解释是否可靠。只有在所有维度都表现优秀，才能称为真正好的XAI方法。总结：TradePool双阶段框架 graph TB subgraph S1["阶段1 训练：子结构池化"] A1["分子图输入 原子71维+键12维"] --> B1["消息传递 GCN / GraphSAGE / GAT"] B1 --> C1["子结构池化 PubChem掩码求和虚拟节点"] C1 end subgraph S2["阶段2 解释：权重映射"] D1["线性层预测标签 同时产生子结构权重"] --> E1["选取子结构权重 不同模型取和或L1范数"] E1 --> F1["按掩码回分原子 归因=所有含该原子的子结构权重之和"] F1 --> G1["输出原子归因热图 稳定、可化学解释"] end S1 --> S2 图2详细展示了TradePool的核心工作流程。TradePool的创新在于将子结构作为可解释的中间层，连接原子级输入和分子级预测。训练阶段（图2A左侧）：分子图经过消息传递层后，不是直接进行全局池化，而是根据预先计算的PubChem指纹子结构掩码，为每个子结构创建虚拟节点。这一步骤（图2B）通过将属于同一子结构的所有原子嵌入求和来实现——如果某个分子不包含某个子结构，则对应的虚拟节点为零向量。这些子结构表示随后被展平并输入到线性层以预测最终标签。解释阶段（图2A右侧）：TradePool的优势得以充分体现。由于线性层的权重直接对应于每个子结构对预测的贡献，研究者可以提取这些权重作为子结构归因值。图2C展示了如何将子结构归因映射回原子：对于每个原子，查找所有包含该原子的子结构，将这些子结构的权重累加，即得到该原子的归因值。这种从子结构到原子的映射策略确保了归因的全局一致性——同一子结构内的原子获得相同的基础贡献，而处于多个子结构交叉位置的原子则累积多重贡献，这与化学直觉高度一致。结果与分析主任务预测性能：并未牺牲预测精度在GCN上，TradePool的原子归因Pearson相关：芳香性0.93，LogP 0.63，TPSA 0.88；常见解释方法多在0–0.30之间。 F1与sparsity均优于GNNExplainer、KernelSHAP、Integrated Gradients、PGMExplainer，说明归因更聚焦、冗余更少，解释“准而简”。 GAT上表现一般，源于注意力权重为标量，难与子结构权重对齐，提示池化假设需与注意力机制协同。预测精度与对照组相当或更优：TPSA任务GCN的MAE 1.157、RMSE 1.569（对照3.367、4.846），LogP任务GCN的MAE 0.299与对照0.296持平，显示可解释性增强未牺牲主任务性能（SI Tables S3–S5）。主任务预测性能是可解释性的基础。图3展示了TradePool在三个任务上的预测精度散点图，每行对应一个任务（芳香性、LogP、TPSA），每列对应一个GNN架构（GCN、GraphSAGE、GAT）。从图中可以看到，所有模型在三个任务上都取得了较高的$R^2$值，数据点紧密分布在对角线附近，表明预测值与真实值高度吻合。值得注意的是，TradePool引入的子结构池化机制并未牺牲预测精度——在大多数情况下，TradePool的$R^2$与使用传统WeightedSumAndMax池化的对照组相当，甚至在TPSA任务上表现更优。这证明了子结构池化不仅提升了可解释性，还通过结构化的中间表示增强了模型对任务相关模式的学习能力。图3：三任务的真值与预测散点行：芳香性、LogP、TPSA；列：GCN、GraphSAGE、GAT 颜色区分训练/验证/测试；对角越集中表明拟合越好，右上角图例标示$R^2$ 原子归因精度全面领先在保证预测性能的前提下，TradePool在原子归因质量上实现了显著突破。图4对比了TradePool与四种主流XAI方法（GNNExplainer、KernelSHAP、Integrated Gradients、PGMExplainer）在原子归因准确率和Pearson相关性上的表现。图4A显示，TradePool在芳香性、LogP和TPSA三个任务上，特别是在GCN架构上，原子归因准确率均达到最高，芳香性任务甚至接近0.9的准确率。图4B的Pearson相关更是揭示了TradePool的优势：在GCN上，TradePool的相关系数达到0.93（芳香性）、0.63（LogP）和0.88（TPSA），而其他方法大多徘徊在0–0.3之间，部分方法甚至出现负相关。这种量级的提升表明，TradePool生成的原子归因不仅在统计上与化学真值一致，而且能够准确捕捉到原子对分子性质的实际贡献。图4：原子归因总体表现对比 (A) 原子归因准确率：TradePool在三任务、特别是GCN上最高 (B) Pearson相关：TradePool显著领先，其他方法多在0–0.3之间图5想回答的问题是：如果一个方法把很多原子都判成关键原子，它当然容易拿到高召回，但这不一定是好解释。作者强调要同时看F1、召回与稀疏性，避免靠把整分子都高亮来“刷指标”。结论1：TradePool的召回不是靠过度归因换来的。论文指出，TradePool在芳香性与LogP任务的F1与召回都表现良好，更重要的是稀疏性维持在0.4–0.5，意味着大约只有40%–50%的原子被标为正归因，解释更聚焦、信息密度更高。结论2：KernelSHAP与PGMExplainer存在明显的过度归因倾向。它们在部分任务上召回较高，但对应稀疏性很低，说明方法倾向把接近90%的原子都判为正归因，从而抬高召回。作者认为稀疏性过高或过低都意味着解释存在缺陷：太低会导致解释冗余，太高又容易漏掉关键结构片段。图5：F1、召回与稀疏性对比 (A) F1-score：TradePool在LogP/TPSA上优势明显。 (B) 召回率：KernelSHAP与PGMExplainer高召回但伴随过多正归因。 (C) 稀疏性：TradePool保持0.4–0.5的稀疏度，解释更集中。为了直观展示不同XAI方法的归因质量，研究团队随机选取了三个任务（芳香性、LogP、TPSA）测试集中的分子，使用GCN模型生成原子归因热图并进行对比。图11中，绿色高亮表示被预测为正归因的原子，每一行对应一个不同的任务。从可视化结果可以清晰看到，TradePool的原子归因与化学真值高度吻合：在芳香性任务中，TradePool准确高亮了芳香环上的原子；在LogP任务中，疏水性基团（如芳环、烷基链）被正确识别；在TPSA任务中，含氧、含氮的极性原子得到强调。相比之下，GNNExplainer、KernelSHAP、Integrated Gradients和PGMExplainer等方法存在明显的误高亮问题——它们倾向于高亮更多的原子，包括一些与任务无关的位点，导致归因稀疏性降低、解释冗余增加。这种可视化案例验证了前面定量指标的结论：TradePool不仅在Pearson相关、F1等数值指标上优于对照方法，在实际化学结构解释的视觉一致性上也表现更佳，更符合化学家的直觉判断。图11：不同XAI方法在GCN模型上的原子归因可视化对比每行对应芳香性、LogP和TPSA三个任务之一。绿色高亮表示被预测为正归因的原子。 TradePool的高亮区域与真值最吻合，误高亮最少；其他方法存在明显的过度归因或归因偏差。子结构化学合理性提取各任务权重前十的SMARTS子结构（表3），与化学常识一致：芳香性任务40%含芳香键； LogP任务突出甲基、芳环、卤素； TPSA任务高频出现含氧、含氮片段O−H、N−O、S(=O)(=O)。子结构权重跨随机种子保持一致（SI Table S6），N=N、C=S、N−S、O(:C)(:C)等始终位列前十，支持归因的可重复性。在SMARTS里，: 表示芳香键（aromatic bond），:C 表示芳香碳原子（aromatic carbon） ~ 表示任意键（any bond）：不限定是单键、双键、三键还是芳香键，只要两原子之间“有键”就匹配。表3 三个任务权重最高的10个SMARTS子结构芳香性 LogP TPSA N#N C−I ≥1 O C=S ≥1 Br O(:C)(:C) N−S N#N N−H C(∼N)(:N) ≥1 Cl ≥1 P N−C:O:C ≥1 S O−H C(:N)(:N) C:C−N−C:C N#N ≥1 Cl C−Br ≥1 N O(:C)(:C) C(∼H)(∼H)(∼H) C#N N−N C(∼F)(∼F) S(=O)(=O) C#N ≥1 F N=O 稳定性与鲁棒性对模型参数微调或输入分子小幅扰动，TradePool的fidelity/infidelity方差最低，归因热图变化最小。归因稳定性的原因：子结构集合固定、权重全局学习、低频子结构滤除减少噪声。对指纹掩码随机置零10%或对分子增加单键旋转等扰动，原子归因排名的Spearman相关仍高于0.85，而对照方法掉到0.5以下，说明结构微扰下解释更稳。在超参数网格搜索（学习率、隐藏维度、层数）中，TradePool的Pearson相关标准差小于0.03，显著优于对照，超参敏感性低。 Fidelity/Infidelity统计（SI Tables S7–S8）：LogP任务TradePool的fidelity_mean=3.38、infidelity_mean=0.69，明显优于其他方法；TPSA任务TradePool保持正向fidelity 36.18，而KernelSHAP虽fidelity高但infidelity为负且方差大，说明TradePool稳定性更高。图6对应论文的稳定性实验：作者在芳香性测试集里对369个含苯环的分子做结构微扰，在苯环上添加1–2个甲基；不含苯环的分子不做修改。随后用同一个预训练模型分别对加甲基前后分子计算原子归因，并比较每个分子的原子归因准确率变化量。结论：TradePool对结构小改动更稳。论文报告所有方法的变化总体接近0，但TradePool的变化显著更小，说明当分子发生轻微修饰时，TradePool给出的关键原子集合更不容易漂移。对于药物优化而言，这意味着解释可以跨相邻类似物复用，降低“每做一次修饰就要重新理解解释”的成本。为什么这能叫稳定：芳香性任务的真值关键原子主要是芳香环原子，给苯环加甲基不会改变原来的芳香原子标签。理想的解释应继续高亮芳香环，而不是被新增甲基带跑偏。图6：小扰动下的稳定性 (A) 在分子上随机添加1–2个甲基的示例。 (B) 各方法扰动前后原子归因值变化，TradePool波动最小，说明对结构微扰不敏感。图7是图6的可视化证据：同一批分子在加甲基前后的归因热图对比。这里的绿色代表被方法预测为正归因的原子，也就是它认为的关键原子。在芳香性任务里，这些原子理想情况下应与真实的芳香原子位置一致。结论1：TradePool的高亮区域更贴近化学真值且更一致。加甲基前后，TradePool主要持续高亮芳香环原子，新增甲基不会导致模型把大量非芳香原子误判为关键。结论2：部分对照方法会把几乎所有原子都判为芳香。论文特别指出KernelSHAP与PGMExplainer会把所有原子都分类为芳香原子，这会造成“看似召回很高、实际毫无区分度”的解释，和图5中稀疏性异常的问题相呼应。把图5–7连起来读：TradePool不仅在指标上避免过度归因，还能在结构微扰后保持解释形状；而某些方法的高召回来自过度归因，导致热图失去化学可读性。图7：扰动前后的原子热图对比绿色高亮：被判定为正归因的原子；每两行对应同一分子扰动前后。 TradePool在扰动后保持高亮区域一致性，对照方法出现更多误高亮。图8：跨随机种子的敏感性。五个随机种子训练的GCN模型在芳香性测试集的原子归因准确率箱线图；TradePool方差最小，鲁棒性最佳。图9：不同随机种子下的原子热图每行对应一个随机种子训练的模型，绿色为正归因原子。 TradePool跨种子保持高亮模式一致，对照方法高亮位置漂移更大。 Fidelity的实验验证是通过移除原子来测试归因质量：将模型预测为正归因的原子张量置零后重新预测，若预测值显著下降，说明这些原子确实对预测有正贡献；反之，移除负归因原子应使预测上升。图10展示了在芳香性任务上，TradePool移除正归因原子后预测值下降最多，移除负归因原子后预测值上升最稳定，证明其归因方向与化学机制高度一致。相比之下，其他方法如GNNExplainer和Integrated Gradients移除原子后预测值变化较小或方向不一致，表明归因质量不佳。纵轴：模型对芳香性任务的输出（预测的芳香原子数量），数值越大表示预测的芳香原子越多。横轴（1–4）与图例一一对应：1为真值，2为完整分子预测，3为仅保留关键原子后的预测，4为移除关键原子后的预测。若一个方法真的找到了关键原子，那么黄色箱线图应明显低于绿色（拿走关键后预测下降），而红色应仍接近绿色（只看关键也能维持预测）。图例解释：蓝色为真值分布；绿色为完整分子预测；红色表示把非关键原子特征置零后的预测；黄色表示把关键原子特征置零后的预测。图10：移除正/负归因原子后的预测值变化在芳香性任务上，移除不同XAI方法预测的正/负归因原子后，模型的预测值变化。 TradePool移除正归因原子后预测值下降幅度最大（KernelSHAP也还行？），移除负归因原子后预测值上升幅度最小且稳定，验证其归因方向正确。其他方法移除原子后预测值变化较小或方向不一致，说明归因与模型决策机制不匹配。版面结构：每一列对应一种解释方法（Random、TradePool、GNNExplainer、KernelSHAP、IG、PGMExplainer），每一行对应一种GNN架构（GCN、GraphSAGE、GAT）。箱线图显示分布范围，箱体中线为中位数，三角形为均值。临床分子集迁移在5800个II期及以上候选药物上，TradePool在GCN与GAT上均优于传统基线，显示其对真实药物骨架的可迁移性。典型案例：含卤代芳环的口服候选物，TradePool高亮芳环与卤素原子，与疏水性主导的LogP真值一致；对照方法偏高亮杂原子，解释偏差较大，显示对实际药物骨架的解释可靠性。在LogP与TPSA任务中，去除正归因原子会导致TradePool预测下降幅度最大，去除负归因原子则上升最小（SI Figures S1–S2），进一步验证其归因方向符合化学机制。需要注意的是：正文与Supporting Information未给出该临床分子集的完整数值表，仅给出数据集构建方式与文字性结论描述。结果逻辑图（方法—结果—局限） graph TB subgraph S0["问题与设计"] Q1("难以获得稳定原子归因") --> Q2("采用指纹子结构做全局池化") end subgraph S1["实验管线"] D0("数据清洗与指纹筛选") --> D1("GCN/GraphSAGE/GAT训练") D1 --> D2("线性层子结构权重") D2 --> D3("掩码回分原子归因") end subgraph S2["核心结果"] R1("Pearson相关0.93/0.63/0.88") --> R2("F1与sparsity领先基线") R2 --> R3("扰动下归因稳定性最佳") R1 --> R4("权重前十子结构符合化学直觉") end subgraph S3["局限与改进"] L1("与GAT注意力不匹配") --> L2("计划加入可学习子结构生成") L3("小数据高复杂任务R^2偏低") --> L2 end Q2 --> D0 D3 --> R1 D3 --> R3 R3 --> L1 方法局限与改进方向对注意力模型支持不足：子结构等权假设与GAT的原生注意力冲突。数据集较小（<3000）或任务复杂时，$R^2$偏低，子结构权重难以学到任务相关性。未来计划：在训练中加入“生成-筛选”子结构模块（类似GAN），替换低权重子结构，提升任务相关性与多样性。化学与工程解读化学角度：子结构权重凸显芳香键、卤素、含氧氮片段，与芳香性、疏水性、极性表面积的主导因素一致，提升了模型的化学可信度。工程角度：使用标准指纹可避免任务特定规则，部署时只需计算指纹与权重矩阵，无需逐分子重新训练，适合大规模虚拟筛选。 Q&A Q1：为什么用PubChem指纹而不是ECFP或规则切分？ A1：PubChem指纹是公开字典，881位覆盖常见官能团，跨分子可比；数量适中，便于全局权重学习；规则切分在多数分子下碎片数<10，统计显著性不足。 Q2：子结构权重如何转成原子归因？ A2：训练后从线性层取每个子结构的权重（GCN取和，GraphSAGE/GAT取L1范数），再用子结构-原子掩码，将包含该原子的所有子结构权重相加，即为该原子的归因值。 Q3：为什么对参数和分子扰动更稳？ A3：归因依赖全局训练得到的固定子结构权重，而非逐样本优化；子结构数量大、权重聚合降低单一掩码变化带来的波动；低频子结构被剔除减少噪声。 Q4：数据清洗如何保证标签一致性？ A4：使用RDKit与MolVS标准化SMILES、去盐和中和，重复分子按骨架去重；标签计算遵循Crippen原子贡献或拓扑表面积分拆，保证训练与真值口径一致。关键结论与批判性总结潜在影响：为分子GNN提供全局、量化的原子归因路径，能直接指导先导优化与毒性定位，降低对专家经验的依赖。局限性：与注意力类模型存在机制不匹配；小数据、高复杂任务下权重难学；对子结构词表的覆盖度仍依赖预定义指纹。未来方向：引入可学习的子结构生成与淘汰机制；探索与GAT兼容的子结构加权方式；将方法拓展到蛋白-配体复合物、材料晶格等更大图结构。小编锐评：做可解释性分析的一种尝试了。我的体会是，如果更贴近人类语言，那还得是基团，但到底谁贡献多，会不会有相关，本身就是有点复杂的，case by case的解释是避免不了的。现在这样有解释已经不错了。做可解释性分析可以水这么多图，学到了

Machine Learning & AI · 2026-01-11

黑箱的透明化：通过原子敏感性分析实现可解释的pKa预测模型

“黑箱”的透明化：BCL-XpKa通过原子敏感性分析实现可解释的pKa预测模型 Title: Interpretable Deep-Learning pKa Prediction for Small Molecule Drugs via Atomic Sensitivity Analysis Authors: Joseph DeCorte,* Benjamin Brown, Rathmell Jeffrey, and Jens Meiler https://doi.org/10.1021/acs.jcim.4c01472 Cite This: J. Chem. Inf. Model. 2025, 65, 101-113 一、论文整体概览 1. 摘要、背景与科学问题摘要翻译机器学习（ML）模型如今在预测药物研发所必需的性质方面扮演着至关重要的角色，例如药物的对数尺度酸解离常数（pKa）。尽管近期在架构上取得了进展，但由于缺乏基准真实数据，这些模型在面对新化合物时常常泛化能力不佳。此外，这些模型也缺乏可解释性。为此，通过精心设计的分子嵌入，可以通过观察模型对输入分子进行原子扰动后的响应，来获取化学结构中的原子级分辨率信息。在此，我们提出了BCL-XpKa，一个基于深度神经网络（DNN）的多任务分类器，用于pKa预测，它通过Mol2D描述符来编码局部原子环境。BCL-XpKa为每个分子输出一个离散分布，该分布存储了pKa预测值以及模型对该分子的不确定性。BCL-XpKa能很好地泛化到新的小分子上，其性能与现代ML pKa预测器相当，在泛化任务中优于多个模型，并能准确模拟常见分子修饰对分子可电离性的影响。然后，我们通过原子敏感性分析（ASA）利用BCL-XpKa的精细描述符集和以分布为中心的输出，该分析无需重新训练模型即可将分子的预测pKa值分解为其各自的原子贡献。ASA揭示了BCL-XpKa已经隐式地学习到了关于分子亚结构的高分辨率信息。我们进一步通过在93.2%的复杂小分子酸和87.8%的碱中识别电离位点，展示了ASA在为蛋白质-配体对接准备结构方面的效用。最后，我们应用带有BCL-XpKa的ASA方法，识别并优化了一款最近发表的KRAS降解PROTAC的物理化学缺陷。背景在计算辅助药物研发领域，准确预测化合物在体内的行为（如生物利用度、溶解度等）对于节约研发时间和成本至关重要。其中，分子的酸解离常数（pKa）是一个决定其在生理pH下电离状态的关键物理化学性质，深刻影响着药物的吸收、分布、代谢、排泄和毒性（ADMET）。传统的预测方法中，量子力学（QM）计算能够提供与实验相当的精度，但其巨大的计算成本使其无法应用于药物发现早期阶段对数以亿计化合物的虚拟高通量筛选（vHTS）。因此，机器学习（ML）方法，特别是定量结构-活性/性质关系（QSAR/QSPR）模型，因其极高的预测速度而成为主流。这些模型通过分子指纹或图神经网络（GNNs）等方式将化学结构转化为数学表示，并学习结构与性质之间的关系。本文解决的科学问题总结尽管ML方法取得了巨大成功，但仍面临两大核心挑战，这也是本文着力解决的科学问题：性能与泛化问题：现有的ML模型大多依赖于数量有限的高质量实验数据进行训练，这常常导致模型在面对训练集中未见过的、新颖的化学骨架时泛化能力差，容易过拟合。可解释性问题：大多数先进的ML模型（尤其是深度学习模型）如同一个“黑箱”，我们很难理解模型是基于分子的哪些具体结构特征做出某一特定预测的。这种可解释性的缺乏阻碍了我们对模型预测结果的信任，也使得我们难以从模型的“智慧”中获得化学洞见来指导后续的药物设计。本文旨在通过创新的模型架构（BCL-XpKa）和新颖的可解释性分析方法（ASA）来同时应对这两个挑战。 mindmap root(可解释性pKa预测分析思路) )为可解释性服务的模型架构( ::icon(fa fa-cogs) **多任务分类（MTC）架构** **核心思路** 预测pKa落入离散区间的概率分布 **关键产出** 输出**概率分布**而非单个值 为ASA提供可比较的分布 **附加价值** 分布的标准差可作为**模型不确定度** **局部原子描述符（Mol2D）** **核心思路** 仅编码原子及其一阶邻居 使模型对单原子扰动更敏感 **关键特性** **可逆性**：描述符可直接映射回化学亚结构 是实现归因分析的基础 )原子敏感性分析（ASA）( ::icon(fa fa-atom) **核心方法：基于扰动的归因** **扰动方式** 将目标杂原子替换为同构的碳原子 **差异量化** 用**KL散度**衡量扰动前后 模型输出的pKa概率分布差异 **分数计算** 通过指数函数放大KL散度 得到最终ASA分数 **应用一：识别关键电离位点** **做法** 寻找分子中ASA分数最高的原子 **结论** 高分原子大概率是主要质子化或去质子化位点 **价值** 快速标注分子质子态 用于对接或MD模拟前的结构准备 **应用二：洞察模型学习到的化学知识** **做法** 分析特定亚结构在不同化学环境下的ASA分数变化 **化学原理验证** **诱导效应**：邻近吸电子或给电子基团 会相应降低或提高碱性氮的ASA分数 **上下文理解** 模型能区分局部环境相同但整体有别的基团 （如咪唑 vs. 吲哚） **研究启发** 验证模型是否学到真实的化学规则 增加对“黑箱”预测的信任 **应用三：指导先导化合物优化** **完整工作流** 1. **定位缺陷**：用ASA找到导致不良pKa的原子 2. **提出方案**：对高分原子进行生物电子等排替换 3. **快速验证**：用模型预测新分子的pKa 4. **结构确认**：用对接等方法确认活性 2. BCL-XpKa 模型简介作者首先构建了一个名为 BCL-XpKa 的pKa预测模型，其核心是一个多层感知机（MLP）。该模型的设计巧妙，集成了几个关键特性：分类而非回归（Multitask Classification, MTC）：不同于传统模型直接预测一个连续的pKa值，BCL-XpKa将pKa范围划分为多个离散的“桶”（bins），并预测分子的pKa值落入每个“桶”的概率。最终的pKa值是这个概率分布的期望值。这种做法的好处是：可以直接从输出分布的标准差中读出模型对预测的不确定度。通过识别模型在哪些分子上表现出“高不确定性”或“高置信度但高错误率”，可以指导训练数据的优化。在性能上与回归模型相当，甚至略优。双模型架构：为了处理既有酸性基团又有碱性基团的复杂分子，作者分别训练了 BCL-XpKaAcid 和 BCL-XpKaBase 两个模型，用于分别预测一个分子中酸性最强和碱性最强的pKa值。图1：BCL-XpKa的架构评估（A） BCL-XpKa 使用独立的模型来预测分子的酸性和碱性pKa值。它使用 Mol2D 局部原子环境描述符来嵌入分子，然后使用一个多层感知机（MLP）来对pKa值所属的1-pKa单位区间进行分类。区间边缘交替包含和不包含端点。极值区间（pKa≤0, pKa>12）在其无界的一侧是开放的。（B）用于pKa预测的多任务分类误差随“桶”尺寸的增加而变化。小的“桶”允许更高的精度，但每个桶的数据更少；而大的“桶”精度较低，但每个桶的数据更多。（C） BCL-XpKa与使用相同分子描述符和训练集训练的最佳性能回归架构在两个外部测试集上的性能对比。（D） “留下一类” （leave-class-out, LCO）方法，其中一种分子亚结构被从模型训练中移除，并在之后用作结构新颖的测试集。（E）模型误差由LCO亚结构和描述符类型决定。（F） LCO亚结构的误差与包含该亚结构的TS-Acid或TS-Base分子数量的关系。 3. Mol2D 描述符为何对 ASA 至关重要？ BCL-XpKa模型选择使用 Mol2D 描述符而非更复杂的GNN，这是实现原子级别可解释性（ASA）的基石。参考其原始论文 BCL::Mol2D—a robust atom environment descriptor for QSAR modeling and lead optimization，Mol2D 的核心优势在于其设计上的简洁性与可逆性。核心定义：Mol2D的核心是原子环境（Atom Environment, AE）。一个 AE 是以某个原子为中心，包含其周围一定化学键距离内的原子及其成键信息。BCL-XpKa使用的是 height=1 的AE，这意味着它只考虑中心原子和与它直接相连的邻居原子。与传统指纹的关键区别：计数而非存在与否：传统指纹（如 Molprint2D）通常是二进制的，只记录某种AE是否“存在”。而 BCL::Mol2D 是一个计数向量，它记录了分子中每种特定AE出现的次数。这提供了更丰富的信息，例如可以区分五元环和六元环。细粒度的原子类型：Mol2D 不仅考虑元素类型，还考虑了原子的杂化状态/轨道构型（’Atom type’ 编码），这使得它能够区分同样是氮原子，但在不同化学环境下的细微差别。通用 AE 库：BCL::Mol2D 的描述符向量的每一个维度都对应一个从大型化合物库（超过90万个类药分子）中预先构建好的“通用AE库”中的特定AE。这意味着描述符的索引是固定的，任何分子都可以被映射到这个统一的向量空间中。可逆性（Reversibility）——实现ASA的关键：这是 Mol2D 最重要的特性。由于描述符向量的每个索引都唯一地、固定地对应着一个具体的化学亚结构（即一个AE），我们可以从描述符向量反推回它所代表的化学结构。这种清晰的“描述符-结构”对应关系，使得当我们扰动一个原子时，我们能精确知道是哪些维度的描述符发生了变化。这为衡量模型对特定原子变化的敏感度提供了直接、无歧义的途径。相比之下，许多复杂的GNN模型其内部表示（节点嵌入）是经过多轮信息传递后高度抽象化的向量，难以直接映射回具体的、独立的原子或化学键贡献，从而使原子级别的归因分析变得非常困难。 4. BCL-XpKa 模型性能与表现总结 BCL-XpKa模型尽管采用了相对简单的多层感知机（MLP）架构，但在多个基准测试中展现了极具竞争力的性能。与主流预测器的性能对比：在多个外部标准测试集（如Novartis、SAMPL6-8）上，BCL-XpKa的平均绝对误差（MAE）与包括ChemAxon、QupKake以及基于GNN的MolGpKa和Uni-pKa等在内的多种现代pKa预测器不相上下。例如，在Novartis-Acid测试集上，其MAE为0.79。优秀的泛化能力：模型的核心优势在于其对新化学骨架的泛化能力。在“留下一类”（Leave-Class-Out, LCO）的交叉验证中，模型需要预测从未在训练集中见过的、特定化学亚结构分子的pKa。结果显示，使用Mol2D描述符的BCL-XpKa显著优于使用传统MACCS和Morgan指纹的同等模型，平均MAE分别为1.1（BCL-XpKa）、1.46（MACCS）和1.20（MFP2）。准确捕捉化学趋势：模型不仅能预测绝对pKa值，更能准确地再现微小化学修饰所引起的pKa变化趋势。例如，在包含71对仅有细微结构差异的分子测试中，BCL-XpKa能够正确预测pKa变化方向的比例高达81.7%。这对于指导药物化学中的先导化合物优化尤为重要。数据策略的有效性：该研究还表明，尽管模型主要在预测数据（来自ChEMBL）上进行训练，但其性能全面优于仅使用少量实验数据训练的同等模型（BCL-MLP-MTC-EO），证实了在大规模预测数据基础上进行训练策略的有效性。二、原子敏感性分析（ASA）方法细节与应用这部分是该研究的核心。原子敏感性分析（Atomic Sensitivity Analysis, ASA）的设计初衷是：在不重新训练模型的情况下，将模型对整个分子的pKa预测值“分解”到每个原子上，从而理解哪个原子或基团对最终的预测贡献最大。 1. ASA的核心原理 ASA的核心思想是“基于扰动的敏感性分析”。它通过系统性地、有物理意义地扰动分子中的每一个原子，并观察模型预测结果的变化剧烈程度，来判断该原子对原始预测的重要性。如果对某个原子的微小改动导致了模型预测结果的巨大变化，那么这个原子就被认为是“敏感的”或“重要的”。 2. ASA的具体实施步骤 graph TD subgraph "ASA 核心流程（针对单个原子）" A["**第1步：父本预测** 将原始分子输入BCL-XpKa 获得pKa概率分布 P_parent"] --> B; B["**第2步：原子扰动** 将分子中的杂原子 a 替换为同构的碳原子 生成'扰动分子'"] --> C; C["**第3步：扰动预测** 将'扰动分子'输入模型 获得新的pKa概率分布 P_perturbed"] --> D; D["**第4步：量化差异** 计算两个分布的差异 使用Kullback-Leibler (KL)散度 D_KL(P_perturbed || P_parent)"] --> E; E["**第5步：计算ASA分数** 通过经验公式放大信号 ASA(a) = exp[S * D_KL] - 1 得到原子 a 的敏感性分数"] end style A fill:#e3f2fd,stroke:#1e88e5,stroke-width:2px style B fill:#fff3e0,stroke:#ef6c00,stroke-width:2px style C fill:#e3f2fd,stroke:#1e88e5,stroke-width:2px style D fill:#e8f5e9,stroke:#2e7d32,stroke-width:2px style E fill:#f3e5f5,stroke:#6a1b9a,stroke-width:2px 第1步：获取父本分子的预测分布：将原始的、未经修改的“父本分子”（parent molecule）输入到BCL-XpKa模型中，获取模型输出的pKa概率分布 Pparent。第2步：对单个原子进行扰动：遍历分子中的每一个杂原子（非碳、氢原子），将其替换为一个保持价态和杂化状态正确的碳原子。第3-4步：获取扰动分布并量化差异：将这个新的“扰动分子”输入到同一个BCL-XpKa模型中，获取其pKa概率分布 Pperturbed，并使用Kullback-Leibler（KL）散度来衡量 Pparent 和 Pperturbed 这两个概率分布的差异。第5步：计算最终的ASA分数：原始的KL散度值需要经过一步经验性的去噪和放大，才能得到最终的ASA分数。其计算公式为： $\text{ASA}(\text{atom } a) = e^{[S \cdot D_{KL}(P_{\text{perturbed}} || P_{\text{parent}})]} - 1$ 这个公式通过指数函数非线性地放大差异，使得影响显著的原子的分数远高于影响微弱的原子。 3. ASA的分析思路与应用识别关键功能位点（Ionization Site Identification）做法是直接找出分子中ASA分数最高的原子。这通常就是模型认为的、决定该分子pKa值的主要电离/质子化位点。论文在Novartis测试集上进行了验证，该测试集中的酸性分子平均有2.93个潜在电离亚结构，碱性分子平均有2.61个。ASA方法在识别最酸性原子时达到了灵敏度96.6%和特异性82.9%，在识别最碱性亚结构时也表现出色。这种方法的直接应用是在药物研发早期，可以快速、批量地为化合物库中的分子标注质子化状态，为后续的对接、MD模拟等步骤提供更准确的输入结构。图3：用于分子分析的原子敏感性（A） ASA协议示意图。“扰动pKa分布”和“父本pKa分布”指的是由BCL-XpKa输出的离散分布。（B）一个由BCL-XpKaAcid评分的酸的ASA分数示例。在这里，磺酰胺的氮原子被正确地选择为比其他潜在的酸性亚结构更具酸性。（C）一个由BCL-XpKaBase评分的碱的ASA分数示例。（D）一个碱，其中酰胺的氧原子在存在更具碱性的氮原子的情况下，主导了ASA分数。这种情况在61个含酰胺的碱性化合物中出现了4次。（E）用于BCL-XpKaBase分解的阳性（蓝色）和阴性（红色）对照亚结构的ASA分数。（F）通过添加一个胺基来调节吡啶氮的ASA分数，示例显示在x轴下方。ns = 不显著。*** = p < 0.001。洞察模型的“化学知识”与学习机制（Probing Model Learning）做法是比较同一官能团在不同化学环境下的ASA分数，或比较局部环境相同但整体化学性质迥异的基团。论文发现，邻近的吸电子基团（EWG）会显著降低吡啶氮的ASA分数（即降低其碱性贡献），而给电子基团（EDG）则会提高其分数（见图4F）。例如，在图2D的分子系列中，将哌啶（pKa 11.2，预测10.45）芳构化为吡啶（pKa 5.20，预测5.45），其碱性显著下降，模型准确地捕捉了这一趋势。这证明了模型虽然只学习了局部原子环境，但隐式地捕捉到了上下文依赖的化学规则。这种分析可以用来验证模型是否学到了正确的化学知识，而不是仅仅记住了某些表观特征，从而增加我们对模型预测的信任度。图4：亚结构的原子敏感性分析（A-D）常见亚结构在作为分子的主要电离位点时与存在更主要电离位点时的ASA分数小提琴图。（E）常见亚结构在作为主要电离位点时的ASA分数小提琴图。（F）相邻的吸电子基团（EWG）和给电子基团（EDG）对吡啶氮ASA分数的影响，通过ASA分数的变化（ΔASA）来衡量。（G）分子对称性对ASA分数的“掩蔽效应”。指导先导化合物的理性优化（Lead Optimization）这是一个非常实用的应用场景，也是可以借鉴的完整工作流。论文以一个已知的KRAS降解剂PROTAC（P-1, PDB: 8QU8）为例，该分子与靶蛋白形成的复合物中，其连接臂上的一个叔胺与KRAS的Q62残基形成了关键的盐桥相互作用（图5A-B）。流程：问题定位：ASA分析显示，这个叔胺氮原子具有最高的ASA分数（12.1），确认了它就是导致PROTAC在生理条件下可能质子化的主要碱性位点（图5C）。而这种质子化状态不利于细胞膜通透性。提出优化方案：针对这个高分子的叔胺，进行生物电子等排替换，例如将其替换为酰胺，设计出候选分子P-2（图5D）。快速虚拟验证：BCL-XpKa模型预测P-2的pKa值显著降低至3.23，成功消除了碱性。结构验证：通过对接模拟发现，新的P-2分子依然能够与KRAS的Q62残基形成一个关键的氢键，保持了必要的结合模式（图5F-H）。结论：这个流程展示了如何利用ASA精确定位分子的物化性质缺陷来源，并指导进行高效、理性的化学修饰，从而在保持活性的前提下优化类药性。图5：用于药物设计的原子敏感性分析（A）泛KRAS降解PROTAC P-1与VHL和KRAS形成的三元复合物的晶体结构（PDB: 8QU8）。（B） PROTAC P-1，其pKa由BCL-XpKaBase计算为6.51。（C） P-1连接臂中氮原子的ASA分数。（D）提出的P-1连接臂生物电子等排替换修饰及其由BCL-XpKa预测的pKa值。（E） P-1和P-2对接到8QU8中VHL-KRAS蛋白-蛋白相互作用界面的三元复合物模型全局视图。（F-H） 8QU8晶体结构和P-2酰胺修饰的结合位点视图，显示它们支持相似的PROTAC构象，并保留了与KRAS Q62的氢键。三、 ASA的局限性与未来方向作者在论文的讨论部分明确提到了当前框架的一些局限性，这对于我们借鉴和改进该方法至关重要。仅限于原子级别：ASA旨在识别单个原子对预测的贡献，但不能直接输出“官能团级别”的贡献。一个原子的影响往往与它所在的整个官能团或药效团紧密相关，而ASA目前无法直接解耦这种集体效应。对非直接电离原子的影响处理不完美：一些本身不电离但能显著影响pKa的原子（例如，通过强诱导效应或共振效应）偶尔会得到异常高的分数，从而干扰对真正电离位点的判断。例如，在一个含有酰胺的碱性分子中，有少数情况是酰胺的氧原子（本身不质子化）得到了最高的ASA分数，这可能是因为它被扰动后对分子整体电子云的改变过大，从而“掩盖”了真正的质子化位点。扰动方式单一：将杂原子替换为碳是一种有意义但简化的扰动方式。对于某些复杂的化学环境，这种替换可能无法完全反映该原子在真实化学修饰中的作用。未来方向：指导数据增强：通过ASA识别出模型预测不佳或不确定的化学结构类型，可以指导性地扩充训练集，从而提升模型的性能和泛化能力。整合到药物发现工作流：作者致力于将ASA整合到更大的药物发现工作流中，例如用于超大规模虚拟筛选（vHTS）的分子库预处理，以确保正确的分子质子化状态，提高筛选的命中率。拓展到其他性质预测：论文提出，未来可以将ASA的思想应用于ADMET（吸收、分布、代谢、排泄、毒性）等更复杂性质的预测模型中，以理解和优化这些关键的药物属性。

Machine Learning & AI · 2025-10-08

解构血脑屏障渗透性：一个可解释的多模态深度集成框架

解构血脑屏障渗透性：一个可解释的多模态深度集成框架一、论文整体概览 1. 论文基本信息标题：Interpretable Multimodal Deep Ensemble Framework Dissecting Blood–brain Barrier Permeability with Molecular Features 中文译名：使用分子特征解构血脑屏障渗透性的可解释多模态深度集成框架期刊：The Journal of Physical Chemistry Letters DOI: 10.1021/acs.jpclett.5c01077 发表年份：2025 Citation: J. Phys. Chem. Lett. 2025, 16, 5806-5819 2. 摘要、背景与科学问题摘要翻译血脑屏障渗透性（BBBP）预测在药物发现过程中扮演着关键角色，特别是对于靶向中枢神经系统（CNS）的化合物。尽管机器学习（ML）已显著推动了BBBP的预测，但目前仍迫切需要能够揭示调控BBB渗透性的物理化学原理的可解释性ML模型。在本研究中，我们提出了一个多模态ML框架，该框架整合了分子指纹（Morgan, MACCS, RDK）和图像特征以改进BBBP预测。分类任务（BBB可渗透 vs. 不可渗透）通过一个结合了多个基础分类器的堆叠集成模型来解决。在可比较的评估设置下，所提出的框架与近期的方法相比，展示了有竞争力的预测稳定性、泛化能力和特征可解释性。除了预测性能，我们的框架还结合了主成分分析（PCA）和沙普利加性解释（SHAP）分析，以突显对预测有贡献的关键指纹特征。回归任务（logBB值预测）则通过一个多输入深度学习框架来解决，该框架包含一个用于处理指纹的Transformer编码器，一个用于提取图像特征的卷积神经网络（CNN），以及一个用于增强特征交互的多头注意力融合机制。从多模态特征中提取的注意力图（Attention maps）揭示了分子表示内部的令牌（token）级关系。这项工作提供了一个可解释的框架，用于以增强的透明度和机理洞察力来建模BBBP，并为未来结合透明描述符和物理信息特征的研究奠定了基础。背景与科学问题血脑屏障（BBB）是保护中枢神经系统（CNS）的关键生理屏障，但它也成为CNS药物研发的巨大障碍。准确预测一个候选药物能否穿透BBB，是其成药性的决定性因素之一。传统的实验方法成本高昂且耗时，因此开发快速、准确且可靠的计算模型至关重要。近年来，机器学习（ML）和深度学习（DL）模型在BBBP预测上取得了很高的准确率。然而，这些高性能的模型往往像一个“黑箱”，研究人员难以理解其做出特定预测的具体依据。这种可解释性的缺失不仅阻碍了我们对模型预测的信任，更重要的是，我们无法从模型学到的知识中提炼出清晰的、指导性的化学规则来辅助新药的理性设计。因此，本文的核心科学问题是：如何在保证高预测精度的前提下，构建一个透明、可解释的BBBP预测框架，从而不仅“知其然”（预测结果），更能“知其所以然”（揭示分子结构与BBB渗透性之间的构效关系）？ mindmap root(可解释性分析思路与实践) )特征集质量初评( ::icon(fa fa-flask) **PCA降维可视化** **核心思路** 在建模之前快速评估特征集的质量与判别能力 **关键发现** MACCS指纹比Morgan指纹 能更有效地分离BBB正负样本 **研究启发** 这是筛选有效分子表示方法的一种重要且高效的前置步骤 )模型归因分析( ::icon(fa fa-search-plus) **SHAP值分析** **核心思路** 定量计算每个“分子亚结构”对最终预测的贡献度 **高贡献度（高SHAP值）的关键亚结构** MACCS_43 极性官能团（氢键供体/受体、磺酸盐） MACCS_39 亚硫酸（酯） MACCS_37 氨基酰胺（如脲结构） **化学原理验证** **极性表面积（PSA）原理** 极性基团（如MACCS_43）增加PSA从而**降低**穿透脂质血脑屏障的能力 （表现为高的负SHAP值） **卤化效应** MACCS_46（溴代基团）的负贡献可能源于分子量增加或代谢不稳定 **上下文依赖性** **结论** 亚结构的最终效果受到 分子整体拓扑和周围化学环境的共同调节 **具体表现** 同一亚结构（如MACCS_38）在不同分子中可产生相反的SHAP贡献 **研究启发** 为药物化学家提供**可操作的优化线索** 指导基于构效关系的**理性药物设计** )深度模型内部机制探索( ::icon(fa fa-project-diagram) **注意力图可视化** **核心思路** 揭示Transformer等深度模型在预测时“关注”的特征区域 **两种分析模式** **内部结构注意力** 分析指纹序列内部各部分的重要性 **跨模态注意力** 分析“指纹特征”与“图像空间区域”的关联 **关键发现** 模型注意力会**从模糊逐渐聚焦** 到化学上有意义的区域 （例如`C=C`和`C-O-C`官能团） **研究启发** 验证深度模型是否抓住了**正确的物理化学特征** 为理解复杂模型的**内部工作机制**提供直观窗口 3. 模型框架总结为解决上述问题，作者提出了一个多模态深度集成框架，其核心是融合不同来源的分子信息来提供更丰富的表征。多模态特征输入：模型不依赖于单一的分子表示，而是同时整合了多种信息。其中分子指纹（Morgan、MACCS 和 RDK）和2D分子图像特征是两大核心输入模态。指纹与模型的使用方式：论文中的框架分别针对分类和回归两个任务设计了不同的模型。在评估时，Morgan、MACCS和RDK这三种指纹是分开独立使用的，即用每一种指纹分别训练和评估模型，以比较不同分子表示方法的效果。它们并未融合成一个单一的巨大特征向量。双任务模型架构：图2中展示的(a)和(b)是针对两个不同任务的两种独立模型。分类模型（图2a，BBB+ vs. BBB-）：该模型仅使用分子指纹作为输入。它采用一个堆叠集成模型（Stacking Ensemble Model）。该模型将多个基础分类器（如逻辑回归、随机森林、XGBoost等）的预测结果作为元特征（meta-features），再由一个最终的分类器进行综合决策，以提高模型的稳定性和泛化能力。回归模型（图2b，logBB值预测）：该模型采用了多模态输入，即同时使用分子指纹和2D分子图像。它是一个更复杂的多输入深度学习网络，使用Transformer处理序列化的指纹特征，用CNN处理图像特征，最后通过多头注意力机制（Multi-Head Attention）将这两种不同模态的特征进行深度融合。图像特征的价值：论文明确提出，通过多模态融合来丰富分子表示是其核心策略之一。在回归模型中，作者专门设计了CNN和注意力模块来处理和融合图像特征。结论部分也强调，多样化分子模态的融合（即指纹+图像）结合透明的归因技术，能够提供更准确和有意义的预测。PCA分析（图8 e,f）显示，在与MACCS和RDK指纹融合后，特征空间的解释方差有所提升或保持高位，这表明图像特征确实为模型提供了有价值的互补信息，特别是在通过跨模态注意力分析揭示两种特征的关联时，其价值更为凸显。二、核心可解释性方法与发现本文的亮点在于系统性地应用了多种前沿的可解释性技术来剖析其模型，从不同维度揭示了BBBP的分子层面的驱动因素。 1. 特征空间分析：PCA降维可视化在进行复杂的模型解释之前，作者首先使用了主成分分析（PCA）这一经典的无监督降维方法，来直观地评估不同分子指纹对BBB+/BBB-两类分子的区分能力。做法与发现：将所有分子的Morgan指纹和MACCS指纹分别通过PCA降到二维空间进行可视化。图6(a)显示，在使用Morgan指纹时，BBB+（红色）和BBB-（蓝色）两类分子的数据点严重重叠，难以区分，且前两个主成分仅能解释总方差的极小部分（PC1: 1.36%, PC2: 1.16%）。这表明Morgan指纹生成的特征向量虽然信息量大，但可能过于稀疏或其线性组合难以捕捉到类别间的清晰界限。相比之下，图6(b)显示，在使用MACCS指纹时，两类分子的数据点形成了相对清晰可分的簇，且前两个主成分解释了更多的方差（PC1: 11.31%, PC2: 8.0%）。这说明MACCS指纹定义的166个预设化学亚结构，能够更有效地捕捉与BBB渗透性相关的结构差异。应用与价值：PCA分析虽然简单，但它是在建模之前快速评估特征集质量和判别能力的有效手段。通过这种方法，作者在早期就得出结论：MACCS指纹在这种二元分类任务中比Morgan指纹更具信息量，这为后续选择MACCS作为主要特征进行SHAP分析提供了依据。 2. SHAP分析：量化分子指纹的贡献 SHAP（Shapley Additive Explanations）是一种源于合作博弈论的模型解释方法，它可以为单个样本的预测结果计算出每个输入特征的贡献值（SHAP值）。一个正的SHAP值表示该特征将预测推向正类（如BBB+），负值则推向负类（如BBB-）。做法与发现：作者对表现最好的分类模型（基于MACCS指纹的堆叠模型）进行了SHAP分析。图6(c)的蜂群图（Beeswarm plot）直观地展示了所有测试样本中，对模型影响最大的前几个MACCS指紋特征。关键特征识别：分析发现，MACCS_43（通常代表富含氢键供体/受体和磺酸盐的极性官能团）、MACCS_39（O-S(=O)O）、MACCS_37（N-C(=O)N）和MACCS_38（N-C(-C)-N）等特征具有最高的平均SHAP值，表明它们对模型的预测有决定性影响。化学意义的验证：图6(d)展示了包含这些关键亚结构的具体分子示例。例如，MACCS_43 在一个分子中贡献了+0.130的正SHAP值，而在另一个分子中贡献了-0.169的负SHAP值。这与化学直觉相符：极性基团通常会增加分子的极性表面积（PSA），从而降低其穿透富含脂质的血脑屏障的能力（对应负的SHAP值，预测为BBB-）。MACCS_46（代表溴代亚结构）的SHAP值为负，这可能是因为过度卤化会增加分子量或引入代谢不稳定性，从而整体上降低了渗透性。上下文依赖性：SHAP分析还揭示了亚结构贡献的上下文依赖性。例如，MACCS_38（对称的二胺结构）在某些分子中贡献为正（+0.108），而在另一些分子中为负（-0.057），这表明其最终效果受到分子整体拓扑结构和周围化学环境的调节。应用与价值：SHAP提供了一种强大的、定量的手段，可以将抽象的模型预测归因于具体的、化学家可以理解的分子亚结构。这使得模型的决策过程不再是“黑箱”，而是可以被验证和理解的。论文指出，这种由SHAP引导的分析为以结构功能关系为基础的CNS靶向药物设计提供了可操作的见解。虽然本文未直接展示用SHAP结果指导模型调优或新实验，但其揭示的关键特征无疑可以用于指导特征工程（例如，构建仅包含最重要特征的简化模型）或提出需要通过实验验证的化学假说（例如，系统性地修饰MACCS_43代表的基团来验证其对渗透性的影响）。图6：特征分析与SHAP可解释性（a）使用Morgan特征的BBB+/BBB-样本的前两个主成分（PC）的PCA得分图。（b）使用MACCS特征的PCA得分图。（c）展示了使用MACCS指纹的模型中，样本级别的SHAP值分布的蜂群图。（d）包含关键亚结构的代表性分子的结构可视化，并标注了其对应的SHAP值（pos代表正贡献，neg代表负贡献）。 3. 注意力机制可视化：揭示模型内部焦点对于用于logBB值预测的多输入深度学习模型，作者利用其核心组件——注意力机制（Attention Mechanism）——来探索模型在进行预测时，其“注意力”集中在哪些特征上。做法与发现：内部结构注意力（Intra-Structure Attention）：图11(a)展示了Transformer编码器内部的注意力热图。这张图揭示了模型在处理一个分子的指纹序列时，不同指纹“令牌”（tokens）之间的相互依赖关系。跨模态注意力（Cross-Modal Attention）：图11(b)展示了从训练的第1个周期到第50个周期，分子指纹特征和CNN提取的图像特征之间的跨模态注意力图的演化。可以清晰地看到，随着训练的进行，模型逐渐学会将指纹中的特定信息（符号/化学特征）与图像中的特定空间区域（视觉特征）对应起来。模型的注意力从最初的模糊、分散状态，逐渐锐化并聚焦于化学上有意义的区域。具体案例分析：以一个BBB+的分子divinylether为例，模型的高度注意力权重区域，无论是内部结构注意力还是跨模态注意力，都准确地对应于其分子结构中的关键官能团，如烯烃（C=C）和醚键（C-O-C）。应用与价值：注意力可视化为理解深度学习模型（特别是基于Transformer的模型）的内部工作机制提供了一个直观的窗口。它能告诉我们模型在做决策时“正在看哪里”，从而验证模型是否抓住了正确的物理化学特征，而不是依赖于数据中的某些伪影或噪声。这种方法为以一种更具机理性的方式理解BBB渗透性提供了支持。图11：注意力权重的可视化（a）代表指纹内部结构的注意力热图。（B）从第1个训练周期到第50个训练周期，结合图像特征的跨模态注意力热图的演化过程。三、本文的局限性与未来展望作者在论文的结论部分坦诚地指出了当前工作的局限性，并对未来研究方向进行了展望。需要更先进的特征选择技术：尽管当前框架表现良好，但未来可以引入更先进的特征选择方法（如LASSO、SHAP-RFE）来进一步优化输入特征，可能会提升模型性能和可解释性。需要更广泛的外部验证：目前模型的验证主要基于B3DB数据集。为了证明其更广泛的适用性，未来需要在更多、更多样化的外部数据集上进行验证。计算预测与实验验证的鸿沟：模型最终需要与真实的实验结果相结合。未来的工作需要整合实验分析，以弥合计算预测与药理学现实之间的差距。展望：作者希望这个可解释的框架能够为药物发现早期阶段CNS活性化合物的设计和优先级排序做出贡献，并为未来融合更多透明化描述符和物理信息特征的研究铺平道路。

Machine Learning & AI · 2025-10-08

Mendelevium

Contact

Interpretability