6.1 建模工作流与经典方法

6.1 建模工作流与经典方法

本节摘要：QSAR 建模的骨架是一条四段流水线：带标签数据、结构描述符、学习算法、严格验证。本节沿方法阶梯从汉施方程讲到三维场方法再到机器学习，并交付一个 EGFR 类靶点的最小可运行案例——重点不在算法新旧，而在每个环节的工程判断与指标的真实量感。

荐书引擎开造。第3章备好了数据，第4章备好了特征，本节把它们装配成预测机器；装配顺序本身就是知识体系的见证——漏掉任何一段上游，下游都会以过拟合的方式讨债。

从结构到数字到预测的四段流水线

第一段，数据。活性标签的口径决定问题形态：连续 pIC50 走回归，"大于某阈值为活性分子"走分类。数据量决定方法层级——几十个分子只配汉施方程，几百个可以上机器学习，几十万才轮到深度学习展示身手。

第二段，特征。描述符矩阵（4.1、4.3 节的产出）是机器学习路线的口粮；指纹展开或指纹计数向量是相似度友好的替代；三维描述符与场变量留给构象可信的场景。特征工程（4.3 节）在这里兑现价值。

第三段，算法。见下文方法阶梯。

第四段，验证。本节只做最小动作（留出测试集），完整纪律留给 6.2 节——但请现在就接受一个观念：没有测试集的模型分数都是自说自话。

方法阶梯从低到高有四层。汉施方程（1964 年奠基）把活性写成疏水、电子、立体三类参数的线性组合：疏水项 π 反映分子穿越脂水界面的倾向，电子项 σ 捕捉取代基对反应中心的极化，立体项刻画基团体积——三项参数对应药代、键合能与空间适配三种直觉，是"可解释建模"的百年祖师。它样本效率极高，几十个分子就能拟合，但表达力止于线性。偏最小二乘与经典二维 QSAR扩展到上百描述符，靠潜变量消化共线性，上世纪九十年代的主力。三维场方法（CoMFA 与后续变体）换了一个视角：把一组对齐的分子放进公共网格，逐格计算探针感受到的立体与静电场能量，场值经偏最小二乘映射到活性——产出著名的等势图（绿色区域放正电基团有利、红色区域相反），第一次把"定量预测"与"下一轮怎么改结构"连成一句话。它的软肋在对齐：分子构象与叠合方式稍变，结果随风摇摆。机器学习阶梯：随机森林、梯度提升、支持向量机吞下全部描述符，非线性交互全权交给模型；再往上就是第8章的图神经网络——直接在分子图上端到端学习，连描述符都可以不要。

图 6-1：QSAR 建模的四段流水线与两层反馈

案例：一个靶点的最小可运行 QSAR

背景。用 3.1 节清洗好的 EGFR 类数据集（约四千个独特分子，pIC50 从四到十）训练一个基线模型，目标是给后续虚拟筛选提供排序引擎，并搞清"这个模型可信到什么程度"。

操作。描述符、骨架划分、随机森林三步走：

import numpy as np
import pandas as pd
from rdkit import Chem
from rdkit.Chem import rdMolDescriptors as rd
from rdkit.Chem.Scaffolds import MurckoScaffold
from sklearn.ensemble import RandomForestRegressor
from sklearn.metrics import mean_squared_error, r2_score

# 1) 描述符矩阵（取若干二维描述符）
def feats(smi):
    m = Chem.MolFromSmiles(smi)
    return [rd.CalcExactMolWt(m), rd.CalcTPSA(m),
            rd.CalcCrippenDescriptors(m)[0],
            rd.CalcNumRotatableBonds(m), rd.CalcNumRings(m)]

X = np.array([feats(s) for s in data["smiles"]])
y = data["pic50"].values

# 2) 按 Murcko 骨架划分：同骨架分子整体进训练或测试
scafs = pd.Series(
    [MurckoScaffold.MurckoScaffoldSmiles(smiles=s) for s in data["smiles"]])
train_mask = ~scafs.isin(test_scaffolds).values   # 测试骨架清单预先指定

model = RandomForestRegressor(n_estimators=300, random_state=0)
model.fit(X[train_mask], y[train_mask])
pred = model.predict(X[~train_mask])

print("测试集 RMSE：%.2f" % mean_squared_error(y[~train_mask], pred) ** 0.5)
print("测试集 R2：%.2f" % r2_score(y[~train_mask], pred))
​

结果。典型量感：随机划分时测试 RMSE 约零点七，骨架划分后涨到零点九上下，R2 从零点六掉到零点四——差距本身就是信息（6.2 节细说）。对四千样本、五列描述符的基线而言，这份成绩算诚实合格。

解读。RMSE 零点九个 pIC50 单位的物理含义要翻译给项目同事：pIC50 差一个单位，活性差十倍；也就是说这个模型的典型预测误差约十倍上下，排序比绝对值可靠（排序任务对系统性偏移不敏感）。这个量感解释了 QSAR 在项目里的真实角色——不是替你定剂量，而是替你排优先级：从五万分候选里把最可能强的百分之几浮出水面，精准度交给实验。另一个实践要点：特征只用了五列，不是省事，是基线纪律——先用最少的可信特征立基线，再加特征看指标是否真的变好，变不好就退回，防的是 4.3 节说的过拟合通道。

变式。沿阶梯各升一级的取舍：换成梯度提升通常能把骨架划分 RMSE 再压零点几，代价是需要更认真的调参与更严的交叉验证；样本量不足时（几百以内）回退偏最小二乘或贝叶斯线性回归，防过拟合反而更稳；等对齐构象可信赖时，上 CoMFA 换取等势图的可解释产出——给药物化学同事的沟通价值有时超过精度本身。分类形态（活性/非活性阈值化）在数据跨度大、极端值多时更稳健，指标也要换成富集度（6.2 节）。

本节要点回顾

• 四段流水线：数据、特征、算法、验证，缺一段下游还债；验证意识要在训练前就位。
• 方法阶梯按数据量选层：几十个分子用汉施，几百上千用机器学习，几十万才轮到深度学习。
• 三维场方法贡献等势图：预测与"怎么改结构"连成一句话，软肋在构象与对齐。
• RMSE 要翻译成倍数：零点九个 pIC50 单位约等于十倍误差——模型排优先级，不定剂量。
• 下一站：造好的引擎要过审讯室——6.2 节的三件套让模型自己交代能信到哪、不能信到哪。