简介面向生物信息学研究者、深度学习初学者及相关领域开发者这份42页PDF完整演示了如何以PyTorch为工具、以Transformer为核心构建基因表达谱分类模型。文档从实际痛点出发既涵盖高维数据清洗、标准化、特征选择等预处理环节也深入讲解输入层设计、位置编码、多头注意力机制、前馈神经网络、全连接层与输出层的实现思路并给出超参数调整、正则化、早停等优化策略。书中还设置了癌症亚型分类、疾病诊断辅助、药物反应预测三个案例结合与传统方法的对比实验让读者能直接借鉴完整建模流程。资源仅有1个PDF文件容量为2.15MB支持目录章节跳转与阅读器大纲定位结构条理清晰方便按章节查阅。目前已有98人浏览学习适合希望快速掌握基因表达谱分类模型构建方法并了解PyTorch与Transformer实践应用的读者。1. 基因表达谱分类为什么需要 PyTorch 与 Transformer基因表达谱分类的任务并不新鲜但它的数据形态一直和常规表格分类不太一样特征数量通常在数万甚至十万量级而样本量往往只有几十到几百。传统机器学习遇到这种场景第一反应是做特征选择或降维例如用方差过滤后接 SVM 或随机森林这在小样本下确实稳定。但问题在于基因表达谱里真正有价值的信号是基因与基因之间的共表达关系单个基因的读数本身很难独立支撑分类决策而手工构造共表达特征又依赖太多先验知识。Transformer 这类自注意力架构最初为序列建模设计但它真正吸引人之处在于能自动学习任意位置之间的关系。把每个基因看作序列中的一个 token表达量看作 token 的特征自注意力天然就能建模基因之间的依赖不需要像 CNN 那样依赖局部窗口也不需要像 GNN 那样显式构造基因调控网络。PyTorch 在这里起到的作用是提供灵活的自动微分和成熟的 Transformer 层实现让模型设计可以快速迭代而不必从零实现反向传播和优化器。本文面向的读者是已经跑得动 PyTorch 基础框架、但还没把 Transformer 用到基因数据上的算法工程师或生信工程师。文章会先讲如何把表达谱转换成 Transformer 的输入再给出一套可运行的最小模型代码随后讨论训练时的批次效应、类别不平衡等实际问题最后落在一个很实用的技巧通过注意力权重做基因层面的可解释性分析。2. 把基因表达谱变成 Transformer 能吃的序列输入2.1 基因表达谱的数据特点决定了 token 化方式基因表达谱通常是一个样本 × 基因的矩阵值代表转录组测序如 RNA-seq得到的表达丰度。RNA-seq 数据常见的问题包括值域跨度极大从 0 到数万、存在大量零值、样本间总读数差异明显。如果直接把原始 count 或 TPM 喂给模型梯度会被极端值主导注意力权重也会偏移到表达量大的基因上而不是关系上有意义的基因上。我见过不少新手在这里踩坑认为 Transformer 对输入分布不敏感直接丢原始值进去。实际上自注意力里的矩阵乘法对尺度非常敏感不同基因之间的读数差异可以达到 3 到 4 个数量级这会直接破坏 QK 点积的稳定性和梯度回传。因此 token 化之前标准化优先级最高。常见的预处理管线是# 假设表达矩阵 expression_matrix.csv行是样本列是基因 python preprocess.py \ --input expression_matrix.csv \ --log1p \ --standardize per_gene \ --save_npy expression_processed.npy也可以用 Python 直接写import numpy as np from sklearn.preprocessing import StandardScaler # 读入表达矩阵形状 (n_samples, n_genes) X np.load(expression_raw.npy) # 第一步log1p 压缩动态范围处理 RNA-seq 的偏态分布 X np.log1p(X) # 第二步按基因标准化每个基因在样本维度上均值为 0方差为 1 scaler StandardScaler() X scaler.fit_transform(X) np.save(expression_processed.npy, X)这段代码里 log1p 的作用是把 count 值的乘法差异变成加法差异标准化则是让每个基因处于同一量纲避免高表达基因主导注意力得分。需要注意StandardScaler 必须在训练集上拟合再用同一组均值和方差去变换验证集和测试集。整个预处理管线的效果往往比后面换模型结构影响更大。2.2 基因排序问题顺序无关数据如何做位置编码文本序列的 token 有天然顺序所以 Transformer 需要位置编码。基因表达谱里的基因顺序则完全由芯片或探针排列决定没有生物学意义。这时如果你直接用标准正弦位置编码等于人为强加了一个不存在的顺序关系模型会去学习这个虚假顺序干扰真实的基因共表达信号。我一般会做两件事。第一把基因按染色体坐标或功能注释重新排序这样相邻位置的基因在生物学上具有潜在关联位置编码就带有先验信息。第二使用可学习的位置嵌入learnable positional embedding让模型自己决定基因顺序中什么是重要的。PyTorch 里的nn.Parameter可以很方便地实现这一点import torch.nn as nn class GenePositionalEmbedding(nn.Module): def __init__(self, n_genes: int, d_model: int): super().__init__() # 可学习的位置嵌入而不是固定的三角函数 self.position_embeddings nn.Parameter( torch.randn(1, n_genes, d_model) * 0.02 ) def forward(self, x): # x: (batch, n_genes, d_model) return x self.position_embeddings这里的n_genes是输入基因数d_model是特征维度。位置嵌入和基因嵌入直接相加前者负责告诉模型当前基因在序列中的位置上下文后者负责编码基因本身的表达特征。相比正弦位置编码可学习版本在大规模数据上更灵活但需要保证训练数据量足够否则容易过拟合到训练集的顺序上。2.3 基因数量过大时的维度压缩策略人类基因组中常见的表达谱芯片覆盖约 20000 个基因RNA-seq 则更多。直接把这 2 万个基因全部作为 token 输入 Transformer计算复杂度是 n²显存和训练时间都会很吃紧。实际操作中需要用特征选择做必要的剪枝。我常用的过滤策略有三个层次策略操作适用场景方差过滤保留方差最大的前 n 个基因快速粗筛丢掉低信息量基因差异表达过滤用 limma 或 DESeq2 筛选组间差异显著的基因监督分类任务利用标签信息通路先验筛选按 KEGG / GO 通路注释只保留通路内基因有明确的生物学假设时比如我最近在处理一个乳腺癌数据集时原始矩阵有 21000 个基因先用方差过滤保留 top 5000再用差异表达分析筛到 2000最终输入维度从batch × 21000降到了batch × 2000训练时间减少了约 8 倍验证集准确率反而小幅提升因为去掉了大量噪声基因。3. 用 PyTorch 构建基因表达谱分类的 Transformer 模型3.1 模型骨架输入投影、编码器、分类头基因表达谱最后一个维度是标准化后的表达值而 Transformer 期望的输入是维度为d_model的向量序列。所以第一步是用一个线性层把单个基因的表达值映射到高维空间类似于 NLP 中的 embedding 层。然后把基因序列送入 TransformerEncoder最后取序列的全局表示做分类。下面是一个可以直接跑通的最小模型代码使用 PyTorch 自带的TransformerEncoderLayer来避免自己实现多头注意力import torch import torch.nn as nn class GeneTransformerClassifier(nn.Module): def __init__( self, n_genes: int, d_model: int 256, nhead: int 8, num_layers: int 4, dim_feedforward: int 512, num_classes: int 2, dropout: float 0.1 ): super().__init__() self.n_genes n_genes self.d_model d_model # 输入投影将每个基因的表达值映射到 d_model 维 self.input_proj nn.Linear(1, d_model) # 可学习位置嵌入 self.pos_embedding nn.Parameter( torch.randn(1, n_genes, d_model) * 0.02 ) # 分类 token保存在序列最前面用于聚合全局信息 self.cls_token nn.Parameter( torch.randn(1, 1, d_model) * 0.02 ) # Transformer 编码器层 encoder_layer nn.TransformerEncoderLayer( d_modeld_model, nheadnhead, dim_feedforwarddim_feedforward, dropoutdropout, activationgelu, batch_firstTrue ) self.encoder nn.TransformerEncoder(encoder_layer, num_layersnum_layers) # 分类头 self.classifier nn.Sequential( nn.LayerNorm(d_model), nn.Linear(d_model, num_classes) ) def forward(self, x): # x: (batch, n_genes)标准化后的表达值 # 1. 每个基因独立投影到 d_model 维度 x x.unsqueeze(-1) # (batch, n_genes, 1) x self.input_proj(x) # (batch, n_genes, d_model) # 2. 添加位置嵌入 x x self.pos_embedding # 3. 拼接分类 token cls_tokens self.cls_token.expand(x.size(0), -1, -1) x torch.cat([cls_tokens, x], dim1) # (batch, n_genes1, d_model) # 4. Transformer 编码 x self.encoder(x) # 5. 取分类 token 对应位置的输出进行分类 cls_output x[:, 0, :] logits self.classifier(cls_output) return logits逻辑说明input_proj是一个nn.Linear(1, d_model)把每个基因的表达值映射成向量表示。这样做的合理性在于单个基因的表达量本身有意义但更重要的是它与其他基因的组合模式高维映射为后续注意力计算提供了更丰富的表示空间。cls_token是借鉴 BERT 的做法在序列最前面拼接一个可学习的向量经过多层编码后这个位置聚合了整个序列的信息。分类时只需取第一个位置的输出避免了设计池化策略的问题。nn.TransformerEncoderLayer内部已经包含多头自注意力、前馈网络、残差连接和层归一化参数nhead控制注意力头数num_layers控制编码器堆叠层数dim_feedforward控制中间层宽度。3.2 数据集与 DataLoader不要让数据加载成为训练瓶颈基因表达谱数据量不大但格式往往比较混乱。有的是 CSV有的是 HDF5有的存在临床信息需要按样本 ID 对齐。PyTorch 的Dataset抽象非常适合处理这类问题。下面给出一个带标签的简单数据集实现import torch from torch.utils.data import Dataset, DataLoader class GeneExpressionDataset(Dataset): def __init__(self, X, y): self.X torch.tensor(X, dtypetorch.float32) self.y torch.tensor(y, dtypetorch.long) def __len__(self): return len(self.y) def __getitem__(self, idx): return self.X[idx], self.y[idx] # 使用示例 dataset GeneExpressionDataset(X_train, y_train) dataloader DataLoader( dataset, batch_size32, shuffleTrue, num_workers2, drop_lastTrue )参数说明batch_size受显存限制基因数多时需要调小drop_lastTrue保证最后一个不完整 batch 不会进入训练避免 BatchNorm 统计量偏差num_workers在 Windows 上建议设为 0否则可能报错Linux 下可以适当调大加速数据读取。3.3 训练循环与优化器选择Transformer 的训练对优化器比较敏感。AdamW 是常见选择需要设置 weight decay 来抑制过拟合。学习率不宜过大我一般从 1e-4 起步配合 warmup 策略让模型在前几个 epoch 内逐步适应梯度尺度。import torch.optim as optim from torch.optim.lr_scheduler import OneCycleLR model GeneTransformerClassifier( n_genes2000, d_model128, nhead4, num_layers3, dim_feedforward256, num_classes2 ) criterion nn.CrossEntropyLoss() optimizer optim.AdamW(model.parameters(), lr1e-4, weight_decay1e-2) # 一个周期学习率调度warmup 阶段线性上升之后余弦下降 scheduler OneCycleLR( optimizer, max_lr1e-4, steps_per_epochlen(dataloader), epochs50 ) for epoch in range(50): model.train() total_loss 0.0 for X_batch, y_batch in dataloader: optimizer.zero_grad() logits model(X_batch) loss criterion(logits, y_batch) loss.backward() # 梯度裁剪防止深层 Transformer 梯度爆炸 torch.nn.utils.clip_grad_norm_(model.parameters(), max_norm1.0) optimizer.step() scheduler.step() total_loss loss.item() print(fEpoch {epoch1}/50, Loss: {total_loss/len(dataloader):.4f})这里的关键参数weight_decay1e-2对小数据集很重要基因表达谱样本量少不加正则化很容易在训练集上过拟合clip_grad_norm_设置梯度裁剪阈值我在训练基因数据时经常碰到某个 batch 的梯度异常偏大裁剪后训练稳定很多OneCycleLR的max_lr通常可以比普通学习率略高因为 warmup 阶段提供了缓冲。4. 训练实战批次效应、样本不均衡与超参数调优4.1 批次效应基因数据里最容易被忽视的坑基因表达数据来自不同的测序批次、不同实验室、甚至不同平台这种技术差异带来的系统性偏移叫批次效应。如果训练集和验证集来自不同批次模型学到的可能不是疾病信号而是批次差异。这一点被忽略时验证集准确率会虚高换个新样本直接就崩。处理批次效应的常见做法有两种一种是数据层面用 ComBat 或 limma removeBatchEffect 回归掉批次协变量另一种是模型层面把批次信息作为额外的 token 输入让注意力机制学会忽略批次相关的模式。数据层面处理代码示例使用 Python 的 combat 库from combat.pycombat import pycombat # data: 基因表达矩阵行是基因列是样本 # batch: 每个样本所属批次标签 data_corrected pycombat(data, batch) # 用校正后的数据替代原始数据在模型层面把批次 ID 映射成 embedding 拼到输入中class BatchAwareGeneTransformer(GeneTransformerClassifier): def __init__(self, n_batches: int, **kwargs): super().__init__(**kwargs) self.batch_embedding nn.Embedding(n_batches, self.d_model) def forward(self, x, batch_ids): x x.unsqueeze(-1) x self.input_proj(x) x x self.pos_embedding # 把批次信息加到每个基因的表示上 batch_vec self.batch_embedding(batch_ids).unsqueeze(1) x x batch_vec cls_tokens self.cls_token.expand(x.size(0), -1, -1) x torch.cat([cls_tokens, x], dim1) x self.encoder(x) return self.classifier(x[:, 0, :])两者的取舍ComBat 简单有效适用于批次信息已知且各批次样本量足够的情况batch embedding 则让模型端到端学习如何处理批次相关特征在批次信息不完整或需要保留更多原始信息时更灵活。4.2 类别不平衡准确率不可靠改用加权损失与 AUPR 评估基因表达谱分类经常遇到类别不平衡。比如某种亚型占总样本的 90%模型只要全预测成多数类就有 90% 准确率看起来不错但没有实际价值。我用两个手段应对损失函数层面做类别加权评估指标层面用 AUPR 而不是 accuracy。# 统计各类别样本数 from collections import Counter counts Counter(y_train) total sum(counts.values()) # 权重与频率成反比 weights torch.tensor([total / counts[i] for i in range(num_classes)], dtypetorch.float32) criterion nn.CrossEntropyLoss(weightweights)CrossEntropyLoss的weight参数会在计算损失时给少数类更大的惩罚让梯度更新更偏向少数类。不过要注意权重不能过于极端否则多数类完全学不好。一般权重范围控制在 1 到 5 倍之间具体可以看验证集上的 AUPR 调整。评估阶段用以下脚本from sklearn.metrics import average_precision_score, roc_auc_score model.eval() with torch.no_grad(): logits model(X_val) probs torch.softmax(logits, dim1)[:, 1].numpy() aupr average_precision_score(y_val, probs) auroc roc_auc_score(y_val, probs) print(fAUPR: {aupr:.4f}, AUROC: {auroc:.4f})AUPR 对少数类的性能更敏感它的含义是查准率-查全率曲线下的面积不会因为多数类占比高而被虚高。样本量只有几十个时AUPR 比 AUROC 更容易暴露模型的真实泛化能力。4.3 超参数选择从小模型出发用验证集做决策Transformer 在基因数据上的表现对超参数并不特别敏感但有几个参数有明确的调整方向。以下是我常用的初始配置基于多个数据集的经验总结超参数推荐初始值调整方向d_model64 ~ 128基因数多时适当增大但超过 256 后收益不明显nhead4 ~ 8d_model 必须能被 nhead 整除头数过多在小数据集上易过拟合num_layers2 ~ 4样本量少时不要超过 4 层深层会放大噪声dropout0.1 ~ 0.3过拟合时调大欠拟合时调小学习率1e-4 ~ 3e-4配合 warmup 使用观察损失曲线下降速度一个实用的判断方法训练集 loss 和验证集 loss 的差距。差距过大说明过拟合优先调大 dropout 或 weight decay两个 loss 都下不去说明模型容量或数据预处理有问题这时候换模型不如回去检查标准化是否合理。我通常会固定一个 3 层、d_model64 的基线模型跑通流程后再逐项调参。基因数据样本量小调参空间有限过度搜索反而容易在验证集上过拟合。5. 注意力权重提取关键基因验证模型学到的生物学信号5.1 从多头注意力中提取基因对关系Transformer 在基因数据上的一个独特优势是训练完成后注意力权重直接给出了基因对之间的关系强度这些权重可以被反向映射回基因 ID。对于分类任务可以统计与 cls token 注意力得分最高的基因列表这些基因就是模型认为最有助于分类判别的基因。attention_weights [] # 存储每个图层的注意力权重 def get_attention_hook(module, input, output): # 在 forward 过程中捕获注意力权重 attention_weights.append(output[1].detach().cpu()) # 注册 hook 到第一层 encoder model.encoder.layers[0].self_attn.register_forward_hook(get_attention_hook) model.eval() with torch.no_grad(): logits model(X_test_batch) # attention_weights[0] 的形状: (batch, heads, seq_len, seq_len) # seq_len n_genes 1位置 0 是 cls token cls_attention attention_weights[0][:, :, 0, 1:].mean(dim1).mean(dim0) # 形状: (n_genes,)表示每个基因与 cls token 的平均注意力得分 gene_ranking torch.argsort(cls_attention, descendingTrue) top_k gene_ranking[:20].numpy()这段代码通过register_forward_hook捕获第一层自注意力模块的权重矩阵然后取出所有样本、所有注意力头中 cls token 对每个基因的平均注意力。得分最高的基因可以被看作模型的判别依据。需要注意的是output[1]是nn.MultiheadAttention返回的注意力概率矩阵但只有在模型内部调用时才会返回注册 hook 时如果模型封装过需要确认返回值格式。得到 top 基因列表后用注释数据库做富集分析# 保存基因列表 with open(top_genes.txt, w) as f: for idx in top_k[:50]: f.write(gene_symbols[idx] \n) # 之后可以用 gseapy.enrichr 做 GO / KEGG 富集分析如果模型学到的 top 基因在已知通路中显著富集比如免疫相关通路或细胞周期通路说明模型不仅过拟合数据还学到了有生物学意义的模式。这一步在高影响因子审稿流程中几乎是必需内容。5.2 与基线模型的对比验证Transformer 并不总是比传统模型更好尤其在样本量极少的情况下。做对比实验时不能只和 SVM、Random Forest 比还要和 MLP 比这样才能确认收益来自注意力机制而非单纯的非线性模型容量。from sklearn.linear_model import LogisticRegression from sklearn.metrics import average_precision_score, roc_auc_score # Logistic RegressionL2 正则C 需要调参 lr LogisticRegression(C0.1, max_iter2000) lr.fit(X_train, y_train) lr_probs lr.predict_proba(X_val)[:, 1] # 与 Transformer 的预测结果做对比 print(fLR AUPR: {average_precision_score(y_val, lr_probs):.4f} fAUROC: {roc_auc_score(y_val, lr_probs):.4f})对比实验有一个重要细节基因筛选或标准化时如果使用全量数据含验证集来拟合 scaler 或筛选基因会导致信息泄漏让验证指标偏乐观。正确的做法是只在训练集上做特征选择和 scaler 拟合验证集和测试集用相同参数变换。常见做法是把这部分封装成一个 Pipeline避免手工操作时遗忘。5.3 模型鲁棒性的一个快速检查方法最后一个实用技巧用基因置换测试评估模型是否真的依赖基因关系。把验证集中每个样本的基因顺序随机打乱然后输入模型。如果模型的预测概率发生显著变化说明模型依赖于基因间的相互作用而不是单个基因的表达量如果概率几乎不变说明模型可能在走捷径只看了高表达基因的绝对值。X_val_shuffled np.random.rand(*X_val.shape) * 0 # 伪代码示意 for i in range(len(X_val)): X_val_shuffled[i] np.random.permutation(X_val[i]) model.eval() with torch.no_grad(): orig_probs torch.softmax(model(torch.tensor(X_val)), dim1)[:, 1].numpy() shuffled_probs torch.softmax(model(torch.tensor(X_val_shuffled)), dim1)[:, 1].numpy() # 计算预测相关性的变化 from scipy.stats import spearmanr rho, _ spearmanr(orig_probs, shuffled_probs) print(fSpearman correlation: {rho:.3f}) # rho 越低说明置换后预测变化越大模型越依赖基因间关系当 rho 大于 0.9 时模型很可能只是在利用单基因的判别能力Transformer 的注意力机制没有真正学到基因关系。这时需要回去检查位置编码和 cls token 的设计或者增大d_model看能否捕捉更复杂的模式。这个检查做一次不会耽误超过几分钟但能省下之后分析模型失败原因的半天时间。本文还有配套的精品资源点击获取