☰
药物相互作用预测模型复现:从摩根指纹到图神经网络的深度学习实践
2026/10/27 7:39:41 网站建设 项目流程

简介:一份聚焦药物相互作用预测的深度学习项目资源,面向生物信息、医药数据挖掘方向的研究者与Python开发者。压缩包以decagon1-master项目为主体,包含模型定义、数据处理、训练评估等完整模块,适合希望复现基于图神经网络与深度学习方法的药物相互作用预测流程的学习者。包内共19个文件,以Python脚本为主,配合Jupyter Notebook用于探索性分析和测试,另含架构示意图与依赖说明,整体仅621KB,轻量便于下载。已有192人学习。内容覆盖多类型深度学习模型、分子数据预处理、特征工程、模型训练优化及评估指标等关键环节,同时包含模型解释与部署相关思路,可帮助快速搭建实验环境、理解药物分子表征与相互作用建模的核心步骤。

1. 基于深度学习的药物相互作用预测.zip:解开压缩包前先搞清楚这几件事

如果同事丢给你一个基于深度学习的药物相互作用预测.zip,你的第一反应多半不是兴奋,而是先确认里面到底有没有能跑起来的东西。药物相互作用(DDI)预测是临床用药里的老问题:两种药一起吃,会不会出现毒性叠加或疗效抵消?过去靠数据库枚举和药师手工查,现在更常见的做法是把药物分子变成向量或图,交给深度学习模型自动学交互模式。这类zip一般装的是数据预处理脚本、模型代码、训练好的checkpoint和一份很简的readme,适合想把公开数据跑出基线、再往药企场景迁移的工程师或研究生。开头这步别急着训练,先把数据链路看明白。

2. 数据与特征工程:把 DrugBank 的原始表变成能训练的药物对样本

2.1 DrugBank 之外还有哪些数据源,DDI 阳性数据怎么挑

DDI 预测最常见的建模方式是把任务定义成药物对二分类:给定 drug A 和 drug B,判断它们是否存在相互作用。正样本好办,凡是数据库里明确写了"存在相互作用"的药物对都算;负样本反而是麻烦事,因为"没有记录"不等于"没有相互作用"。先看数据源怎么选。

我一般会优先用 DrugBank 起步,它免费版就能导出药物条目 CSV 和相互作用条目,里面每个相互作用还带一段 description,比如"增加QT间期延长的风险",这为后面做细粒度分类留了退路。TWOSIDES 是另一个常见选择,它来自 FDA 不良事件报告系统,用药物对和副作用共现频次做统计推断,偏挖掘性质,条目很噪。ChEMBL 质量高,但 DDI 不是它的主力赛道,assay 数据整理起来更费劲。下表是三条路线的基本情况:

数据源原始形式DDI 标注颗粒度适合的用途
DrugBank药物条目 + 相互作用对药物对级 + 文本描述二分类、多分类、规则分析
TWOSIDES药物对 + 副作用共现药物对级挖掘新候选、补充负样本
ChEMBLassay 活性数据条目级,需自行判定特征校验、外部验证集

整理 DRUGBANK 免费版时,通常拿到的是两个文件:drugs.csv 里包含 drugbank_id、name、smiles(有的版本叫 SMILES 列),interactions.csv 里是 drug1_id、drug2_id、description。有的 zip 里已经帮你把这两个 CSV 合并成可以直接训练的样本表,但合并方式未必严谨,常见做法是只保留两列 drugbank_id 都存在于药物表的相互作用对,防止下游解析时遇到空指纹。

2.2 用 RDKit 把 SMILES 转成 Morgan 指纹,参数怎么定

药物分子必须变成模型能吃的张量。主流路线有三条:基于 Morgan 指纹的向量、基于 SMILES 序列的 Transformer embedding、基于分子图的节点特征。做 DDI 二分类,我一般选指纹起步,因为它快、可复现、不依赖额外的语料模型。RDKit 的 GetMorganFingerprintAsBitVect 把分子编码成固定长度的二进制向量,官方叫 ECFP 系指纹,radius=2 等价于 ECFP4 的感受野。

下面这段代码演示从药物表读出 SMILES 并批量转成指纹,是所有后续模型共用的前置步骤:

import csv from pathlib import Path from rdkit import Chem from rdkit.Chem import AllChem import numpy as np data_dir = Path("./data") drug_file = data_dir / "drugs.csv" inter_file = data_dir / "interactions.csv" feature_map = {} # drugbank_id -> 1024 维指纹数组 valid_ids = [] with open(drug_file, "r", encoding="utf-8") as fp: reader = csv.DictReader(fp) for row in reader: smi = (row.get("smiles") or "").strip() if not smi: continue mol = Chem.MolFromSmiles(smi) if mol is None: continue bit_vec = AllChem.GetMorganFingerprintAsBitVect(mol, radius=2, nBits=1024) feature_map[row["drugbank_id"]] = np.array(bit_vec, dtype=np.float32) valid_ids.append(row["drugbank_id"]) print(f"valid drugs: {len(feature_map)}")

逻辑说明很直接:先按行读 CSV,丢掉没有 SMILES 的条目,再用 MolFromSmiles 解析,解析失败说明这行的分子表示有问题,干脆跳过;最后把 RDKit 的 BitVect 转成 numpy 数组。这里有几个参数值得你亲手调:radius 控制特征感受野,radius=2 是 ECFP4,radius=3 是 ECFP6,交互预测任务里 ECFP4 通常更稳定,因为过大半径会把官能团信息混在一起;nBits 对应向量长度,1024 是常见起点,小数据集 512 也够,你要是拿几千对样本硬上 2048,意义不大还浪费显存。

从 BitVect 到 numpy 数组这一步很多人忽略:直接np.array(bit_vec)拿到的其实是对象的引用,必须像我这样显式转换,否则后面进 PyTorch 会报类型错误。还有一点,指纹是稀疏的 0/1 向量,不需要额外归一化,标准化反而会扭曲二值语义,Dropout 和 BatchNorm 比归一化更管用。

2.3 负样本构造:随机配对与“未知即负”的取舍

特征就绪后,最大的坑是负样本。现实世界中绝大多数药物对从不被检验,所以“数据库里没记录”不能直接当成安全标签。但入门阶段、包括很多已发表论文,实际都在用“未知即负”的假设:把不属于已知正对的药物对当作负样本。这个假设有代价,评估出来的 AUPR 会略微虚高,但工程上完全可接受,因为你后续还要用真实病例数据做二次验证。

我建议负正比从 5:1 到 10:1 之间选。比例太小模型见不到足够负样本,比例太大正样本被稀释,AUPR 反而难看。负采样一定要固定随机种子,否则每次跑出来的指标都不一样,调参时就会怀疑人生。下面是一个可复用的负采样脚本:

import numpy as np known_positive = set() with open(inter_file, "r", encoding="utf-8") as fp: reader = csv.DictReader(fp) for row in reader: a, b = row["drug1_id"], row["drug2_id"] known_positive.add(tuple(sorted([a, b]))) rng = np.random.default_rng(42) # 固定种子,保证可复现 neg_pairs = [] pos_count = len(known_positive) neg_target = pos_count * 10 while len(neg_pairs) < neg_target: a, b = rng.choice(valid_ids, size=2, replace=False) key = tuple(sorted([a, b])) if key in known_positive or a == b: continue neg_pairs.append((a, b, 0))

while 循环里做了两件事:避免自配对、避免把正样本重复放进负样本集合。replace=False保证一次随机不会抽到同一个药物,这是最容易被忽略的细节——同一个药和自己配对在化学上没意义,放进训练集只会教模型记住“同一药物就判负”。生成完负样本后,把所有正负对合并,再依照 feature_map 把两边的指纹拼成一个 2048 维向量,标签是 1/0,就能直接喂给下一章的模型了。

负采样还有更严格的变体:要求负样本的药物对在已知的相似药物结构上也不与任何正样本冲突,这需要额外计算,入门阶段不必做。先把简单的跑通,再看指标是否需要修正。

3. 模型选型:先用 MLP 跑基线,再用图神经网络吃结构信息

3.1 药物对二分类的最小模型:特征拼接 + MLP

拿到拼接后的 2048 维特征向量,最常见、也最值得先做的模型就是 MLP。MLP 的好处是黑匣子嫌疑小、训练快、能马上给出一个可对照的基线分数。很多论文里所谓“深度学习 DDI 预测”惊讶地提升几个点,其实就是因为基线 MLP 没有认真调。所以你自己做实验时,先把 MLP 当基准,后面模型能超过它再谈创新。

一个够用的 MLP 结构通常是三层,第一层 256 维,第二层 128 维,输出单节点。中间要夹 BatchNorm 和 Dropout:指纹输入是高维稀疏二值,BatchNorm 能防止某几个 bit 的取值偏移主导梯度;Dropout 是防止模型直接把训练集里的药物对死记下来。PyTorch 实现长这样:

import torch import torch.nn as nn class MLPDDI(nn.Module): def __init__(self, in_dim: int = 2048, hidden_dim: int = 256): super().__init__() self.net = nn.Sequential( nn.Linear(in_dim, hidden_dim), nn.BatchNorm1d(hidden_dim), nn.ReLU(), nn.Dropout(0.3), nn.Linear(hidden_dim, hidden_dim // 2), nn.ReLU(), nn.Dropout(0.2), nn.Linear(hidden_dim // 2, 1) ) def forward(self, x): return self.net(x)

forward 直接返回未过 sigmoid 的 logits,激活留给损失函数处理。这里有几个参数值得解释:hidden_dim 256 是我在小样本 DDI 任务里的习惯,药物对样本量通常只有几千到几万,过大的 hidden_dim 会让模型快速过拟合;Dropout 0.3 是第一层后的常规值,第二层降到 0.2 是防止过度正则化把信号也压没了。BatchNorm 放在 Linear 之后、激活之前,别放反,否则激活后的分布会破坏 BN 的统计假设。

3.2 升级到图网络:把药物集合建模成图,预测缺失边

MLP 的局限在于只看两个药物的指纹,完全不利用“已知相互作用网络”的结构。例如阿司匹林和 A 药相互作用,A 药和 B 药相互作用,这些信息对预测阿司匹林-B 药有没有帮助?MLP 说没有,图神经网络说有。

常见做法是把所有药物当作节点,已知相互作用当作边,整个数据库变成一张大图,DDI 预测就是在这张图上做“缺失边预测”。节点特征仍然用 Morgan 指纹,或者进一步用原子级特征生成。下面的代码片段用 PyTorch Geometric 定义了一个两层的 GCN 链路预测器:

import torch import torch.nn as nn from torch_geometric.nn import GCNConv class GCNLinkPredictor(nn.Module): def __init__(self, in_dim: int, hidden_dim: int = 128, n_layers: int = 2): super().__init__() self.convs = nn.ModuleList() self.convs.append(GCNConv(in_dim, hidden_dim)) for _ in range(n_layers - 1): self.convs.append(GCNConv(hidden_dim, hidden_dim)) self.decoder = nn.Sequential( nn.Linear(hidden_dim * 2, 128), nn.ReLU(), nn.Dropout(0.2), nn.Linear(128, 1) ) def forward(self, x, edge_index, edge_label_index): for conv in self.convs: x = torch.relu(conv(x, edge_index)) src, dst = edge_label_index[0], edge_label_index[1] return self.decoder(torch.cat([x[src], x[dst]], dim=-1))

forward 分两段:先用 GCN 在图结构上做消息传递,得到每个药物的新向量;再针对候选边(训练时是正负样本对)把两端节点向量拼起来,过一个小 MLP 打分。注意 GCNConv 的 edge_index 是训练用的边集合,而 edge_label_index 是你想预测的候选对,两者必须分开。GCN 层数这里只放两层,原因很实际:药物图直径通常不大,两层已经能传播到二阶邻居,再多层会发生过平滑,所有节点向量趋向一致,预测对之间的差异就消失了。

3.3 输入输出形状和 Batch 组织:训练时最容易写错的维度

MLP 和 GCN 的输入差异看起来只是“向量拼起来”和“图结构”,但真正写训练循环时报错几乎都出在维度上。MLP 的输入形状是 (batch, 2048),输出 (batch, 1),损失函数需要 y 是 (batch,) 的一维标签,所以你在用 BCEWithLogitsLoss 时一定要做一次logits.flatten(),否则 PyTorch 会默默广播,指标还算得出来,只是结果完全没意义。

GCN 的输入形状则完全不同:x 是 (drug_num, in_dim),edge_index 是 (2, edge_num),前向传播出来后再用x[src]做索引,得到的是 (candidate_edge_num, hidden_dim),拼接后进 decoder。最容易翻车的点是 edge_label_index 的构建:训练集的正样本边必须从图中移除,否则 GCN 在消息传递时已经“看见”了要预测的边,测试指标会虚高一截。这个叫 label leakage,做链路预测时是头号错误,比数据划分泄漏还隐蔽。

Batch 组织上,MLP 用普通 DataLoader 就行;GCN 不能直接对整个图做 batch,常见做法是每次采样一个子图,或者用 NeighborSampler 对每个候选边采样固定阶邻居。入门阶段用全图 + 小 batch 候选边也能跑,只是数据量大后显存吃紧,这点放到避坑章节细说。

4. 训练、验证与超参数:用 AUPR 盯早停,别被 Accuracy 骗了

4.1 按药物分桶做数据划分,避免训练和测试互相泄露

DDI 数据集划分比普通分类任务讲究得多。如果你直接用 sklearn 的 train_test_split 随机划分配对样本,同一药物可能同时出现在训练集和测试集的不同药物对里,模型相当于在考自己背过的药,测试 AUPR 会虚高到不真实。正确做法是按药物 ID 分组,让训练集和测试集没有任何共同药物。

实现上不需要依赖复杂库,手动按药物划分就行:

all_drugs = list(feature_map.keys()) rng.shuffle(all_drugs) split = int(len(all_drugs) * 0.8) train_drugs = set(all_drugs[:split]) test_drugs = set(all_drugs[split:]) train_pairs = [] test_pairs = [] for a, b, label in all_pairs: if a in train_drugs and b in train_drugs: train_pairs.append((a, b, label)) elif a in test_drugs and b in test_drugs: test_pairs.append((a, b, label)) # 跨分组的边直接丢弃,保留纯净评估

逻辑说明:只保留药物完全落在同一集合里的样本对,跨组样本直接丢弃。代价是可用样本变少,尤其是测试集可能只剩 20% 药物内部的配对,但保住了评估可信度。更严格的做法是用分子骨架(Bemis-Murcko scaffold)分桶,防止结构相似的药物在训练和测试里重复出现,这在发表论文时是加分项,工程上先按 ID 分组就够。

4.2 BCEWithLogitsLoss + AUPR 的组合,类别权重怎么给

正负样本比 1:10 时,交叉熵损失会偏向多数类。直接套用公式负样本数 / 正样本数会得到权重 10,但实践里我通常会把 pos_weight 设在 3 到 6 之间,因为过高的权重会让模型把大量负样本误判成正样本,反而拉低 AUPR。PyTorch 里用 BCEWithLogitsLoss 时把 pos_weight 传进去就行,它内部会在正样本项上乘这个权重。

评估指标一定要用 AUPR,而不是 Accuracy 或 ROC-AUC。原因在于 DDI 是典型的“正样本稀少但更有价值”的任务:Accuracy 会被 90% 的负样本推高到 95% 以上,看着光鲜,实际上正样本可能一个都没抓住。AUPR 在正样本占比低时对召回率更敏感,能直接反映模型找阳性交互的能力。模型训练期间每个 epoch 都算验证集 AUPR,用早停保住最佳状态。

4.3 学习率、batch size、早停 patience 的一组可行起点

DDI 类任务的数据量不大,调参规律其实和《动手深度学习》这类书里讲的原则一致:先固定模型结构,再拿小数据集摸学习率量级。我常用的起点参数如下表,照着跑通常不会离谱:

超参数建议起点调整方向
学习率AdamW, lr=3e-4AUPR 前 3 个 epoch 不涨就降到 1e-4
batch size256显存不够先降到 128,别动 lr
最大 epoch80早停触发可提前终止
早停 patience10验证 AUPR 连续 10 轮不涨即停
dropout0.3 / 0.2过拟合明显时加到 0.5
pos_weight5高召回优先则调到 8~10

早停时要保存验证集指标最好的模型权重,而不是最后一轮的权重。这个细节救过我很多次:DDI 模型到第 60 轮左右验证集 AUPR 开始下降,训练集继续下降,如果你傻傻保存最后一个 checkpoint,拿到的就是过拟合版本。写训练循环时要在每个 epoch 结束后比较当前 AUPR 和 best_aupr,能更新就 torch.save 一份临时文件。调参过程里你很快会发现,同一个 lr 换一个随机种子可能就让 AUPR 掉两点,这种“玄学”波动不必过度焦虑,属于小样本任务的正常方差。

训练循环的核心代码大致是这样,注意 logits 和 y 的形状对齐:

import torch import torch.nn as nn from sklearn.metrics import average_precision_score def train_one_epoch(model, loader, optimizer, criterion, device): model.train() total = 0.0 for x, y in loader: x, y = x.to(device), y.to(device) optimizer.zero_grad() logits = model(x).flatten() loss = criterion(logits, y) loss.backward() optimizer.step() total += loss.item() * x.size(0) return total / len(loader.dataset) def evaluate_aupr(model, loader, device): model.eval() all_y, all_pred = [], [] with torch.no_grad(): for x, y in loader: x = x.to(device) logits = model(x).flatten() preds = torch.sigmoid(logits) all_y.append(y.numpy()) all_pred.append(preds.cpu().numpy()) return average_precision_score(np.concatenate(all_y), np.concatenate(all_pred))

evaluate 函数里用torch.no_grad()是必须的,否则评估阶段会消耗大量显存甚至 OOM。average_precision_score 是 sklearn 对 AUPR 的实现,直接用预测概率而不是阈值后的类别,这个区别很关键:如果你先按 0.5 阈值转成 0/1 再算 AP,指标会失真。

5. 复现这个 zip 的 5 个常见坑与排查:从解压乱码到数据泄露

5.1 坑一:zip 解压后路径带中文或空格,RDKit 与 torch 一起罢工

现象:代码在原机器上跑得好,换到你这边解压后,pd.read_csv("./data/药物 相互作用.csv")报 UnicodeDecodeError,或者 torch.load 找不到文件。原因分两半:Windows 下中文路径经过 zip 解压可能出现编码错位,路径里的空格在某些旧版 CUDA 扩展里会被误解析。解决:解压后第一件事就是把目录重命名成纯英文小写,例如ddi_project,所有子目录不要有空格;代码里统一用pathlib.Path拼接路径,不要手工写字符串加+。这个坑不解决,后面每一步都会在 IO 层反复翻车。

5.2 坑二:SMILES 解析返回 None,数量还不少

现象: feature_map 最终只有 800 个药物,但源表有 1000 行,你还没注意到。原因:DrugBank 免费版里有一部分词条是混合物、聚合体或带特殊标记的 SMILES,RDKit 的 MolFromSmiles 会返回 None,还有一些条目 SMILES 列本身就是空的。解决:解析循环里加一个失败计数器,把失败的 drugbank_id 和原始 SMILES 打印到一个日志文件里。如果失败比例超过 5%,说明数据版本或清洗方式有问题;如果是零星几条,直接过滤掉影响不大。千万不要静默吞掉异常,否则后面预测时发现缺了某个药,回头查数据又得从头跑。

5.3 坑三:随机划分数据集,测试指标虚高导致“模型作弊”

现象:MLP 测试 AUPR 高达 0.98,合理预期应该是 0.6 到 0.8。原因:直接用 train_test_split 按行切分样本,同一药物出现在训练对和测试对里,模型见过这个药物的向量,测试时相当于是开卷考试。解决:回到 4.1 节按药物分桶的划分方式。怎么快速判断当前数据集是否泄漏?统计一下训练集和测试集药物集合的交集大小,如果交集不是空集,就有泄漏。更隐蔽的泄漏是分子骨架相似性导致的结构泄漏,指标上表现为“模型好像学得太完美了”,这时候可以用 scaffold 划分重新评估。

5.4 坑四:类别不平衡让 Accuracy 虚高,AUPR 却崩盘

现象:训练日志里 Accuracy 升到 97%,但 AUPR 只有 0.2,模型几乎没识别出正样本。原因:负样本占 90% 以上,模型学到的策略是全判负就能拿到高准确率,交叉熵损失被海量负样本主导。解决:训练阶段用 pos_weight,把正样本的梯度权重提上来;评估阶段抛弃 Accuracy,只信 AUPR。排查时可以打印验证集正负样本的预测分数分布,如果正样本的分数普遍低于 0.5,说明损失权重还不够,把 pos_weight 往上调,或者改用 Focal Loss 对难分样本加权。

5.5 坑五:全图 GCN 训练显存 OOM,batch 一调小又报错

现象:GCN 前向传播时 CUDA out of memory,把 batch size 调小后,报错变成“index out of range”。原因:全图 GCN 在每层都保留整张图的中间特征,药物图只要有几千个节点和几万条边,中间张量就会暴涨。batch size 减小只影响候选边数量,不影响图规模,所以看起来没变化。解决:两招一起用。第一招是邻域采样,PyTorch Geometric 的 NeighborSampler 或 GraphSAINT 可以在每个 batch 只保留候选边周围 k 跳的子图,显存开销断崖式下降;第二招是训练和推理阶段都加torch.no_grad(),推理时不计算梯度,显存需求减半。要是权重文件加载时还报 pickle 或 weights_only 相关错误,多半是 PyTorch 版本升级后的兼容问题,用torch.load(path, map_location="cpu", weights_only=False)加载旧 checkpoint,再重新保存一次。

6. 从二分类到可解释的细粒度交互:用 Tanimoto 近邻做冷启动与类型预测

二分类只是起点。DrugBank 相互作用条目里那段 description 才是真正的业务价值,例如“增加华法林出血风险”和“降低降压药疗效”,这两者对应的干预措施完全不同。进阶做法是把单标签二分类升级成多标签分类:每个相互作用类型对应一个输出头,模型预测该药物对是否属于这个类型区间。损失函数换成多标签 BCE,评估指标改用 macro-F1,正负样本不平衡的问题会更突出,pos_weight 也要按类别分别统计。类型少时两三个标签还好,标签超过二十个时就要考虑相关性建模,比如“出血风险”和“凝血异常”天然互相关联。

第二个值得做的是冷启动场景。新药上市时没有任何相互作用记录,模型怎么给它一个初始风险提示?我一般会用 Tanimoto 相似度近邻投票:计算新药与已知药物的分子指纹相似度,取最相似的 K 个药物,看它们历史上被标注过的相互作用标签,按比例投票给出建议。下面是一个最小实现:

from rdkit import Chem, DataStructs from rdkit.Chem import AllChem def cold_start_prediction(smi: str, drug_mols, interaction_labels, k: int = 5): fp_ref = AllChem.GetMorganFingerprintAsBitVect(Chem.MolFromSmiles(smi), radius=2, nBits=1024) sims = [] for idx, (drug_id, fp_drug) in enumerate(drug_mols): sim = DataStructs.TanimotoSimilarity(fp_ref, fp_drug) sims.append((sim, idx)) top_k = sorted(sims, reverse=True)[:k] label_sum = np.zeros(interaction_labels.shape[1]) for sim, idx in top_k: label_sum += interaction_labels[idx] return label_sum / k

k 取 5 到 10 之间效果比较稳,k 太小会被单个相似药物带偏,k 太大则混入结构差异明显的噪声。这份近邻分数适合做成带警告的参考值,而不是自动决策。做这个课题时我最大的教训是:不要急着上复杂模型,先把数据拆分和负样本做对,否则后面所有创新都是空中楼阁,我也曾因为数据泄漏白高兴了整整一周。先把这些踩坑经验带进复现流程,你的指标和结论都会靠谱得多,希望帮到你。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询