RNA‑seq 表达数据与机器学习概述

做 RNA‑seq 数据分析的人，迟早会遇到“特征太多、样本太少”这件事。几年前我第一次拿到转录组数据，心想不就是个基因表达矩阵嘛，扔进分类器不就好了？结果模型全在过拟合，验证集的 AUC 只有 0.52，跟抛硬币差不多。后来才意识到，RNA‑seq 数据的结构特点，决定了它的机器学习套路和图像、文本完全不同。

这一章我们先捋清楚：RNA‑seq 数据长什么样、它能解决哪些生物问题、以及常见的陷阱在哪里。你不需要一下子记住所有细节，但得对这个基本面有个直觉。

RNA‑seq 数据的 DNA（数据特性）

一个标准的 RNA‑seq 表达矩阵，行是基因（通常 2 万～6 万个），列是样本（几十到几百个）。你想想看，这就是典型的“高维小样本”场景。再加上基因表达值往往有大量零值、批间效应、测序深度差异……说白了，就是一个噪声大、分布偏、稀疏度高的矩阵。

关键数字一览：

特性 典型值范围 对机器学习的影响
特征数（基因） 20,000 – 60,000 极易过拟合，需要强正则化或特征选择
样本数 50 – 2,000 样本稀少，交叉验证稳定性要求高
零值比例 20% – 70% 稀疏矩阵，部分模型（如 GBDT）效率下降
动态范围 0 – 10⁶（CPM） 需要对数变换或归一化
批次效应 几乎必然存在 需做 ComBat / limma 校正，否则模型学的是批次

我个人的习惯是，拿到数据后先做三步检查：

  1. 看表达值的分布（直方图），有没有严重的长尾。
  2. 检查样本的聚类情况（PCA 或 UMAP），如果不同批次完全分开，就得先校正。
  3. 统计每个基因的表达中位数，去掉表达量一直趋近于零的基因。

这三步走完，数据才算“可以用”。

三个常见应用场景

1. 分类：识别疾病亚型或组织类型

最经典的任务——拿 RNA‑seq 表达谱判断样本是肿瘤还是正常、是基底型还是 luminal 型。这里模型输出是离散类别。我记得有一次用公开的 TCGA 数据做乳腺癌亚型分类，特征有 2.5 万个基因，我直接扔给随机森林，训练集准确率 99%，测试集直接掉到 65%。为什么？因为没做特征选择，模型记住了几十个基因的噪音。后来改成先选出前 1000 个差异基因，再训练，测试集准确率升到 88%。所以，分类任务的瓶颈不是算法，而是特征工程。

2. 生存分析：预测患者生存时间或风险

生存分析用的不是普通回归，而是 Cox 比例风险模型或者它的深度变体（DeepSurv、Cox‑NNet）。输入是基因表达，输出是风险评分。这里有个坑：Cox 模型本身不能处理高维特征，必须配合 L1/L2 正则化（Lasso‑Cox）。我项目中遇到过一个泛癌生存预测的数据，直接用全基因做 Cox 回归，结果 90% 的系数都是零——Lasso 自己帮你选了不到 100 个关键基因。嗯，这说明基因表达量的预测信号其实是稀疏的。

# R 示例：用 glmnet 构建 Lasso‑Cox 模型
library(glmnet)
cv_fit <- cv.glmnet(x = expr_matrix, y = Surv(time, status), family = "cox")
plot(cv_fit)
coef(cv_fit, s = "lambda.min")  # 提取非零系数基因

3. 聚类：挖掘新的分子亚群

当没有标签时，用无监督聚类（k‑means、层次聚类、或非负矩阵分解 NMF）发现表达模式相似的样本组。RNA‑seq 聚类最大的陷阱是维度灾难——在高维空间中，欧式距离会失效。我建议先用 PCA 或 t‑SNE 降到 30～50 维再做聚类，否则聚类结果几乎等于随机。另外，聚类数的确定可以用轮廓系数或者 consensus clustering。

我的小技巧： 在做聚类之前，先用 PCA 做一个Elbow图（解释方差 vs 主成分数），选一个刚好让曲线“拐弯”的维度。再降维到这个维度上跑聚类，结果比直接全基因跑稳定得多。

机器学习处理 RNA‑seq 的核心流程

下面的图把整个过程串起来了。你可以把它想象成数据从“原始计数”到“业务结论”的管道。

RNA‑seq 机器学习全流程 原始 FASTQ 表达矩阵 (TPM/CPM) 质控+批次校正 (ComBat, limma) 特征选择 (差异基因/方差筛选) 模型训练 (Lasso, RF, Cox) 评估与应用 (AUC, C-index, 聚类) 💡 关键洞察 • 特征选择是 RNA‑seq 性能提升的第一杠杆。 • 批次效应不消除，模型学到的可能不是生物学差异。 • 模型选择不重要，数据质量和特征质量才重要。

一个快速起步的代码框架

不管后面用什么模型，数据预处理的骨架都差不多。下面是一个极简的 Python 流程，加载 TPM 矩阵并做初步特征筛选：

import pandas as pd
import numpy as np
from sklearn.preprocessing import StandardScaler
from sklearn.feature_selection import VarianceThreshold

# 加载表达矩阵 (row=基因, col=样本)
expr = pd.read_csv("expression_tpm.csv", index_col=0)

# 1. 去除低表达基因 (至少 20% 样本中表达 > 1)
mask = (expr > 1).sum(axis=1) / expr.shape[1] > 0.2
expr = expr.loc[mask]

# 2. 对数变换 (log2(TPM+1))
expr_log = np.log2(expr + 1)

# 3. 过滤低方差特征 (防止全是常数)
selector = VarianceThreshold(threshold=0.1)
expr_filt = selector.fit_transform(expr_log.T)  # 转置为样本×特征

print(f"基因数从 {len(expr)} 减少到 {expr_filt.shape[1]}")
⚠️ 我曾经踩过的坑： 早年我直接用原始 counts 做 PCA，结果第一个主成分反映的是测序深度，而不是生物学差异。后来养成了习惯：先用 CPM/TPM 归一化，再做 log2 变换，最后看是否需要标准化。另外，特征选择一定不能在训练集和测试集合并后整体做——不然数据泄露，测试集评估就废了。

小结

RNA‑seq 表达数据的特点是“胖而杂”，但机器学习在很多场景下仍然很实用。这一章我们覆盖了三个主要应用（分类、生存分析、聚类），也提到了批处理、特征选择、高维稀疏性这些实际中绕不开的问题。后面的章节我们会深入每一个超参数（惩罚系数、树深度、学习率……），一步步看它们怎么影响模型表现。

嗯，先消化这些基础。接下来第二章，我们会拿一个真实的 TCGA 乳腺癌表达数据集，从头搭建第一个分类基线模型。

← 上一章 📖 返回目录 下一章 →