Module 03 — EGFR 配体数据集：化学多样性、活性分类与 pIC50 回归

本报告基于 TeachOpenCADD（Volkamer Lab，固定提交 4efdbd36723b5d07960142e15df393d7c6df2ba1）公开数据， 用本机计算完成。所有指标都由本会话实际运行得到；没有引用任何未执行的数字。 图 1–13 同时提供 PNG（阅读）、SVG 与 PDF。不是每个图元都是矢量：文字、坐标轴与曲线为矢量， 而包含数千点的稠密散点层已显式栅格化（否则文件无法使用），因此放大到极点时这些点层会成为位图。 每张图的每个子图都有独立的旁表 CSV。

1. 一句话结论

一个无需拟合参数的最近邻基线就已经很强，而且它的强弱在两个任务上并不一致。 （该基线仍会读取训练集的指纹与标签/数值，只是不拟合任何参数。）

• 该基线（下称 1-NN Tanimoto）的做法是：用指纹相似度选出一个最像的已知分子， 然后把它自己的标签 / 活性转移给查询分子 —— 并不是把相似度数值直接当作活性分数。
• 分类：这个无训练基线的 AUC = 0.817、AP = 0.799；多数模型高于它 （RandomForest 0.927、SVM 0.925、神经网络 0.886、sklearn MLP(5,3) 0.873）， 提升幅度约 0.06–0.11。
• 回归：方向相反 —— Ridge（MAE 0.789）与神经网络（0.749）都不如这个无训练基线（0.682）， 只有 RandomForest（0.633）明显更好。也就是说"上了机器学习就一定更好"在这份数据上不成立。
• 换成不共享骨架的测试集后：分类 AUC 0.884、回归 MAE 0.744。但要注意， 该划分的测试分子与训练集的最近邻相似度中位数仍有 0.66（实测）， 所以这不等于"全新化学型"下的性能。

因此本模块把它当作在已测化学空间附近做候选优先级排序的工具；它的预测仍然是模型活性预测， 但不足以支撑精确的倍数判断（见 §7）。

2. 数据与语义（先看清楚一行记录代表什么）

两个队列是不同的表：它们共享 3063 个 ChEMBL ID，但行数与构成不同。报告与图注中已分别标注， 不会再出现"同一份数据同时用于两个任务"的说法（图 1 标题已按实际队列范围改写）。

一个必须指出的数据问题：T022 表里 47 行的单位字符串不是 nM，但全部 3906 行的 pIC50 都在数值精度内与 9 − log10(原始 IC50 数值) 一致（最大绝对差约 1×10⁻⁶，来自源文件的小数位） ——也就是说，单位字符串在换算中被忽略了。例如 CHEMBL490266 记录为 IC50 = 1.8 /uM、 pIC50 = 8.745；若单位真是 µM，正确的值应是 5.745，相差 3 个对数单位。 我们保留冻结基准不动，另外用一个"仅 nM"的敏感性子队列单独复核（见 §6.3）。

活性阈值 6.3 对应 IC50 ≈ 501 nM——这是教学时选定的一个阈值，不是自然边界； 图 2B 显示 pIC50 本身是一条连续谱，并没有两个分开的峰。

3. 方法（关键选择与理由）

分子表示：RDKit Morgan 指纹（radius 2, 2048 bit）用于相似度、建模与适用域； RDKit path 指纹（maxPath=5）用于复现 T005 的 Butina 聚类协议。 两种表示给出的相似度分布差别很大（图 2A：平均 0.16 vs 0.27），所以"相似"这句话必须说清用的是哪种指纹。

两种划分，回答两个不同问题（注意两个任务用的划分参数并不相同）：

（两个任务都在训练集内部再切 15% 作早停验证：分类 3151/557 与 3275/578，回归 2323/411 与 2407/425。 因此骨架划分的测试集比随机划分小——分类 782 vs 927、回归 1074 vs 1172。）

⚠️ 两个划分不能相互混用比较。分类与回归的 seed、测试比例、队列和样本数都不同， 因此不能把"随机划分 vs 骨架划分"的分数变化归因于某一个因素（例如模型容量）。 两种划分下训练样本量、测试成员、化学空间覆盖同时改变 —— 这是一个多因素同时变化的比较， 只能描述"在这种划分下观察到什么"，不能识别原因。

例如 sklearn MLP(5,3) 在骨架划分下 AUC 下降 0.105、而 RandomForest 只下降 0.043； 这看起来像容量差异，但本模块没有做"固定划分、只改容量"的受控对照， 因此不能据此断言"小容量导致更依赖相似分子"。这是一个尚未验证的假设。

⚠️ 骨架划分只保证训练集与测试集不共享骨架，不保证跨簇没有近邻。 实测骨架划分的测试分子，其与训练集的最近邻相似度中位数仍有 0.66，最高到 1.0。 因此本文一律用"实测的最大训练相似度"来解释结果，而不把它称为"完全陌生的分子"。

模型：Dummy（最频繁类）与 1-NN Tanimoto（无需拟合参数）作为对照。 1-NN 的机制要说清：它不是把相似度数值当成活性分数，而是用相似度挑出最像的一个训练分子， 再把那个分子的标签（分类）或 pIC50（回归）转移给查询分子。它没有拟合任何参数。 其余为 RandomForest、SVM-RBF、sklearn MLP(5,3)（T007 示例），以及一个 PyTorch 神经网络基线 （分类 2048-512-128-1 + Dropout + 早停；回归 2048-64-32-1，即 T022 的架构，Adam + MSE）。

评价纪律：所有预处理/早停只使用训练集内部切出的验证集；测试集只用于最终打分，不用来调参。 报告 AP 时明确写作 average precision，不称其为 PR 曲线积分。 SVM 的 sigmoid(decision_function) 只是分数，不是校准概率（旁表 score_kind 列有记录）。

4. 化学空间覆盖（图 1–5）

• 活性跨越约 10 个对数单位，且是连续谱（图 1A、2B）。
• 骨架高度集中：4635 个分子只有 1475 个 Murcko 骨架（31.8%）；最大的一条—— S1 苯胺基喹唑啉（图 4 左上，c1ccc(Nc2ncnc3ccccc23)cc1）——单条就占 504 个分子（10.9%）。 这正是 EGFR 抑制剂领域最经典的铰链结合骨架。
• 骨架大小 ≠ 活性高低：S2 只是一个苯环（94 个分子，平均 pIC50 仅 4.19，2% 活性）， 而 S4/S5（51/47 个分子）活性比例高达 98%/96%、平均 pIC50 8.5/8.1。图 4 用 RDKit 固定键长绘制， 保证六个骨架同一尺度可比。
• Butina 聚类（d=0.20）：1128 个簇，592 个单例（52.5%），最大簇 145 个分子（图 3A/B）。
• 近邻普遍存在：5.6% 的分子与另一个分子拥有完全相同的 Morgan 指纹；中位最近邻相似度 0.809。
• 相似不等于活性相同（图 2B，本模块最重要的教学点）。这里的统计口径要写清楚： 每个分子只与它自己那一个最近邻配对，因此这是一组逐分子的近邻关系， 不是唯一分子对的计数，也不是整批分子的预测错误率。 在最近邻相似度 ≥ 0.7 的那 3887 个分子中，有 808 个（20.8%） 与自己的最近邻相差 ≥ 1 个对数单位。 作为参照，在全部 4635 个分子上，这个比例是 1053 个（22.7%）。 也就是说：对"最近邻已经很像"的分子，直接照搬最近邻的活性会有约五分之一的分子偏差 ≥ 10 倍。

5. 留出预测（图 6–11）

5.1 活性分类

读法（注意这是实测，不是预设）：换成不共享骨架的测试集后，所有可训练模型在 AUC 上都下降。 逐模型变化如下（只统计可训练模型，不含未拟合的 Dummy 与 1-NN 基线）：

即可训练模型的 ΔAUC 介于 −0.105 ~ −0.031（均值 −0.055）。 注意降幅最小的是 PyTorch 神经网络（−0.031），不是 RandomForest 或 SVM； Dummy 的 0.000 是不变量的必然结果，不纳入"可训练模型"的区间。 但如 §3 所述，两个划分的训练样本量、测试成员与化学空间覆盖同时改变， 所以这一比较不能识别下降的原因；把"小模型掉得多"解释成容量效应， 在本模块中没有受控对照支持，只是一个待验证的假设。

图 7 的纵轴刻意设为完整 0–1 并画出 0.5 机会线，让初学者能看见 Dummy 基线就在 0.5—— 若像通常那样把轴从 0.5 起画，这根最重要的参考柱会完全不可见。

5.2 pIC50 回归

单位：MAE/RMSE 均为 pIC50 对数单位；简单基线是"预测训练集均值"。

两个必须如实说的观察：

1. 无训练的 1-NN 基线（0.682）优于神经网络（0.749），只有随机森林明显胜过它（0.633）。 在这份高度同源的类似物数据上，"查最相似的已知分子"是非常强的对手， 训练一个神经网络的边际收益必须先跨过这个基线才谈得上有意义。
2. MACCS（166 bit）比 Morgan（2048 bit）差（MAE 0.714 vs 0.633）。 表示的选择在这里是可测量的一阶效应，不是细节。

6. 误差与适用范围（图 8、9、12、13）

6.1 误差与"离训练集多远"的关系（是趋势，不是单调曲线）

• 分类（图 9A）：把测试分子按其到训练集的最大 Tanimoto 相似度分箱后，平均 |P(active) − 真实标签| 在低相似度一端更高。
• 回归（图 12A）：总体是"越相似误差越低"的趋势，但并不是单调的。实测分箱值：

可以看到 0.3–0.4 那一箱反而是局部高点（骨架划分 1.590），而最低相似度箱只有 1–3 个分子， 平均值不稳定。因此"低相似度一定更差"这句话本模块不支持； 能支持的只是"高相似度端误差更低，两端之间关系不单调，且低相似度分箱样本很少"。 另外，中等 R² / 中等相关性并不等于模型没有区分信息——见 §7 对 Spearman 的说明。 - 图 9B / 12B 说明骨架划分的测试分子整体更远离训练集 （分类中位数 0.661 vs 0.797；回归 0.675 vs 0.771）。

6.2 "斜率小于 1"到底意味着什么（一次重要纠正）

初学者常把"拟合斜率 < 1"直接读成"预测被压缩到均值附近"。这两个是不同的事：

斜率 = Pearson r × sd(预测)/sd(实测)

骨架划分的斜率 0.69 主要来自相关性弱（r = 0.685），而不是预测被压缩 （sd(pred) 1.511 vs sd(measured) 1.500）。它的平均偏差 mean residual 为 −0.14 对数单位 （略偏低）。图 12C 把这三个量分别画出来，读者可以自己判断属于哪一种情形。

注意图 12C 横轴上 1.0 在三行的含义并不相同：对拟合斜率是"预测与实测一比一变化"， 对 Pearson r 是"完美线性相关"（线性相关，不是秩相关）， 对 sd 比是"预测与实测一样宽"。图中只保留诊断点与参考线，长解释见 figures/CAPTIONS.md 与本节。

因此本模块不再使用"predictions collapse toward the mean / why predictions are conservative" 这类笼统说法，改为陈述具体量。

6.3 一个混淆的比较（如实标注）

"仅 nM"敏感性子队列的 RandomForest MAE 为 0.593，低于冻结基准的 0.633。 但这不能归因于"单位质量更好"：去掉 47 行非 nM 记录的同时，测试集成员也变了 （样本数、划分成员都不同），两个因素同时变化。这是一个混杂比较，只能说明 "排除单位可疑行后性能没有变差"，不能说明原因。见 results/regression_unit_sensitivity.json。

6.4 应用到无标签分子（图 13）

T022 附带的 60 个无标签分子（test.csv）没有真实标签，因而无法打分，只能给出预测： 预测 pIC50 范围 4.33–8.40，其中 26/60（43%） 与训练集的最大相似度 < 0.4， 即接近一半落在模型经验较少的区域。图 13 用随机森林 200 棵树之间的散布作为 离散度提示（这不是校准过的置信区间）。该关系按实测描述，不预设方向： 本次实测为负相关——Spearman(相似度, 树间散布) = −0.660（n = 60，p < 0.001）， 即相似度越高、树间散布越小；分箱中位数为 0.4 以下 1.208（n=26）、0.6 以上 0.646（n=28）。 这只是这 60 个分子上的观察，样本量小且没有标签可验证，不能外推为一般规律。

7. 这些模型能预测什么，对陌生分子有哪些局限

能做（本次有量化支持）： - 在已测化学空间附近给分子做候选优先级排序：随机划分下分类 ROC-AUC 0.927 / AP 0.937， 骨架划分下 0.884 / 0.900，可作为排序信号。 但本模块没有做独立的筛选富集验证（没有诱饵集、没有前瞻验证、没有富集因子计算）， 因此"用富集方式挑候选是有效的"这句话尚未被本次数据检验，不能声称。 - 判断"某个新分子像不像已知活性物"：无需拟合参数、只用训练集指纹挑最近邻再转移其标签/数值的 1-NN 基线，就已提供分类 AUC 0.817 / 回归 MAE 0.682 的信号。 - 对性质接近已知类似物的分子给出量级估计：随机划分 MAE ≈ 0.63 对数单位（约 4 倍 IC50 量级误差）； 骨架划分下 MAE ≈ 0.74。

不能做 / 局限（按本次实测量化）： 1. 不能声称在新化学型上保持精度。 不共享骨架时 MAE 升到 0.744、R² 从 0.66 降到 0.55， AP 从 0.937 降到 0.900。而且这个骨架划分仍保留了中位 0.66 的最近邻相似度； 真正全新化学型（相似度 < 0.4）的表现本次没有可评分的标签，属于尚未验证的范围。 2. 不足以支撑精确的倍数判断。 骨架划分 R² = 0.356、Spearman = 0.707、回归斜率 0.690。 但要注意：R² 与相关性中等，不等于模型没有区分能力 —— Spearman ≈ 0.71 表示 预测仍能对分子强弱做出有信息的排序；不能支持的是把预测 pIC50 = 8 与 6 直接读成 "实验活性差 100 倍"。它仍然是模型的活性预测，只是精度不足以做定量倍数推断。 3. 不能跨靶点 / 跨实验条件外推。 全部数据来自 ChEMBL EGFR IC50 测定， 不同实验的条件、构体与突变背景并不完全可比；47 行非 nM 单位的内部矛盾也提醒需人工核对。 4. 分类的"活性"定义是一个选择。 pIC50 ≥ 6.3 是教学阈值（IC50 ≈ 501 nM）； 换阈值会同时改变类别比例与全部指标。 5. 不能假定"低相似度一定更差"。 图 12A 是趋势而非单调曲线（§6.1）， 且最低相似度分箱只有 1–3 个分子。 6. 不能假定模型越复杂越好，也不能假定越简单越好。 神经网络在分类上接近但未超过 RandomForest；在回归上不如 1-NN 基线（§5.2）。 7. 相似度不能替代实验。 图 2B 显示，即便最近邻相似度 ≥ 0.7，仍有 20.8% 的分子 与其最近邻相差 ≥ 1 个对数单位。

下一步最有信息量的实验：合成或采购一批与现有训练集相似度 < 0.4 的分子并实测活性， 才能量化"新化学型"下的真实误差，并顺带做一次带诱饵集的富集验证。 在此之前，本模块的所有数值只在已测化学空间附近有效。

8. 复现方式

环境：本机 omicos-env Python 3.11（RDKit 2025.09.6、scikit-learn 1.7.2、torch 2.14.0+cpu）。 脚本位于 code/，全部接受 --input/--out 参数，不绑定任何个人路径。

cd <本模块目录>
PY=/home/phys/src/omicos-full/omicos-env/.venv/bin/python3
PACK=/mnt/d/Omicos-dev/data/cadd-tutorial-inputs
T007=$PACK/T007_compound_activity_machine_learning/data/EGFR_compounds_lipinski.csv
T022=$PACK/T022_ligand_based_screening_neural_network/data/CHEMBL25_activities_EGFR.csv
T022T=$PACK/T022_ligand_based_screening_neural_network/data/test.csv

# 1) 化学多样性（只产出紧凑结果，不写大矩阵）
$PY code/run_diversity.py --input $T007 --out .

# 2) 活性分类 --torch-threads 1 见 §9 的线程实测
$PY code/run_classification.py --input $T007 --out . --torch-threads 1 --jobs 2

# 3) pIC50 回归
$PY code/run_regression.py --input $T022 --unlabeled $T022T --out . --torch-threads 1 --jobs 2

# 4) 补齐主运行未持久化的交付物（确定性重放，见脚本 docstring）
$PY code/export_replay_artifacts.py --t007 $T007 --t022 $T022 --out .

# 5) 从已存结果重画全部图（不重新计算，不重新训练）
$PY code/redraw_figures.py --t007 $T007 --t022 $T022 --out .
#    也可只重画一部分：  --only div,clf,reg

# 6) 【仅在缺失交付物需要恢复时运行】补齐分类骨架划分的权重与学习曲线
#    只写这两个文件，不动已记录的指标/预测；脚本会打印重算与记录值的差
$PY code/recover_classification_scaffold_nn.py --input $T007 --out . --torch-threads 1

重绘入口：code/redraw_figures.py 只读 results/ 与两个输入 CSV。 指纹、Tanimoto、聚类、UMAP 与模型训练都不会被触发；改动只影响样式。 code/compact_similarity.py 是已废弃的开发产物——它存在单侧近邻与只对首块抽样两个缺陷， 不是推荐入口；当前交付的 similarity_sample.npz / nearest_neighbour_pairs.csv 由 run_diversity.py 直接生成（对称矩阵求最近邻、对每种表示在全量上抽样），未经过该脚本。

需要区分两件事： 1. 完整运行上表的 1–4 步时，在全新目录里要得到全部 15 张图，第 1 步必须带 --umap（图 5 需要该嵌入）； 2. 重绘已有结果不需要任何重算。

验证范围如实说明：本次只对重绘入口做了独立进程的 fresh-output 验证 （复制 results/ 到全新空目录 → 15×{PNG,SVG,PDF}，且未创建 models/）。 完整流水线本身未做过同等的 fresh-output 验证，它是在本模块目录内按上述顺序执行的。 recover_classification_scaffold_nn.py 是开发期补齐脚本，不属正常全新训练后的必跑步骤。

目录

module-03/
├── REPORT.md                 本报告
├── README.md                 运行说明与交付清单
├── DEVELOPMENT-NOTES.md      本轮开发过程中的作图/运行指导与修正记录
├── code/                     common.py + 6 个可执行脚本（见上，含恢复脚本）
├── inputs/input_manifest.json  输入路径 + SHA256 + 行数 + 角色
├── results/                  全部数值表（紧凑格式，见下）
├── figures/                  15 张图 × PNG/SVG/PDF + 38 个子图旁表 CSV
└── models/                   已训练神经网络的最佳权重（.pt）

交付体积说明：早期版本一度写出 similarity_distributions.csv， 包含 21 478 590 行两两相似度、约 571 MB。现已替换为语义等价的紧凑产物： similarity_distribution_bins.csv（100 分箱频数表，13 KB）、 similarity_sample.npz（每种表示 20 万条可复现随机样本，记录样本量与种子）、 nearest_neighbour_pairs.csv（每分子一行最近邻）。三个文件的样本量、分箱边界与来源 都记录在 similarity_pairwise_stats.json 中。

已保存的 ID 与权重

• results/classification_split_ids.csv / results/regression_split_ids.csv： 每个划分的 train / validation / test 分子 ID（含 pIC50 与标签）。
• results/*_predictions.csv：每个已报告模型的逐分子测试预测 （分类 6 个模型 × 2 划分；回归 6 个模型，含 MACCS RF）。
• results/clf_mlp_learning_curves.csv / reg_mlp_learning_curves.csv：学习曲线。
• models/clf_mlp_random_seed22.pt、models/clf_mlp_scaffold_seed22.pt（后者由 §8 的恢复脚本补齐）、 models/reg_mlp_random_seed42.pt 与 models/reg_mlp_scaffold_seed42.pt：神经网络最佳权重。

缺失交付物的恢复（已补齐）

分类任务骨架划分的神经网络权重与学习曲线在首次运行时未被持久化。 按"可以为了补齐丢失的必要交付物做有记录的恢复"这一原则，我们用一个专用脚本 code/recover_classification_scaffold_nn.py 只补这两个缺失文件：

• 只恢复缺失部分：models/clf_mlp_scaffold_seed22.pt 与 results/clf_mlp_learning_curves.csv 中该划分的行。
• 保留已有结果：随机划分的权重与曲线、全部已记录指标、全部逐分子预测，均未改动。 交付的骨架划分指标与预测仍是首次记录的那一份。
• 协议完全相同：同一队列、同一 Murcko 骨架划分（seed 22）、同一内层验证划分 （train 内部 15%，seed 22）、同一架构（2048-512-128-1）与同一早停规则。
• 恢复结果与记录值一致（这是可核查的证据，不是"看起来差不多"）：

782 个测试分子的预测全部落在 5×10⁻³ 容差内（最大绝对差 2.4×10⁻⁷，相关系数 1.0000）。 完整记录见 results/classification_scaffold_nn_recovery.json。

一次性说明：第一次恢复尝试曾因训练集行序与 scaffold_split 返回顺序不一致而给出偏离的结果； 该问题在脚本中被修正（直接使用 scaffold_split 返回的训练/测试帧）。 上面表格是修正后的结果，脚本注释中保留了该注意事项。 除这两个文件外，其他模型一律未重跑。

也因此，本模块不再遗留可恢复的交付缺口。

9. 运行事实记录（含一次被中断的运行）

• 分类：torch.set_num_threads(2)，sklearn n_jobs=2（分类脚本默认值；首次执行使用默认线程池）。
• 回归：首次尝试使用 PyTorch 默认线程池（16 线程），运行约 26 分钟后被 开发复核控制器（development review controller）以 SIGINT 中断，未产出任何结果表 （该中断不是最终用户的手动操作）。请注意：由于命令通过管道 | tail 执行， shell 返回码为 0，但这不代表成功——判定依据是结果文件是否存在。 中断时仅留下了随机划分的权重文件。
• 重试：同一输入、同一目标、同一划分、同一 seed（42）、同一架构与超参，仅把进程内 Torch 线程数显式设为 1、 sklearn n_jobs=2，4 分 45 秒完成全部回归。这与开发复核控制器提供的吞吐探针一致 （batch 32、2048→64→32→1：1 线程 ≈ 0.017 s / 8 step，16 线程 ≈ 3.265 s / 8 step）—— 小 batch 下过大的线程池是净开销。线程配置与中断/重试事实已写入 results/regression_run_info.json。

10. 数据来源与致谢

• TeachOpenCADD（Volkamer Lab 及各章作者/贡献者），固定提交 4efdbd36723b5d07960142e15df393d7c6df2ba1，CC BY 4.0。 相关章节：T002 / T004 / T005（表示、相似度、聚类）、T007（活性分类）、T022（pIC50 回归）。
• ChEMBL 活性数据（含 CC BY-SA 3.0 条款），PubChem 查询结果按各节来源使用。
• 本报告为教学目的重新组织任务、执行计算并撰写中文解读；许可不替代底层数据条款。