数据准备：输入数据格式

在做 edgeR 批次校正之前，有一件事比模型本身更关键——数据准备。我自己踩过最大的坑，就是序列都跑完了才发现样本名写错了一位，导致重新读入花了一整天。所以这一章咱们把输入数据格式、读取方法、完整性检查一次说透。

基因表达矩阵：长什么样？

edgeR 最核心的输入就是一个 整数计数矩阵，行是基因，列是样本，每个格子是比对到该基因的 reads 条数。说白了，它不能是浮点数，不能是 FPKM 或 TPM。为什么？

因为 edgeR 的统计模型假设 reads 数服从负二项分布，归一化后再转化的 logCPM 是中间产物，但原始输入必须是整数。我曾经有个项目拿 TPM 直接跑 GLM，结果 p 值全部偏小，幸亏早期检查出来了。

一般格式如下（CSV 或制表符分隔）：

gene_id    sample1    sample2    sample3    sample4
ENSG001    152        89         234        178
ENSG002    0          5          12         3
ENSG003    456        321        678        543
...

注意事项：

  • 行名唯一：不能有重复的基因 ID，否则合并时 edgeR 会报错。
  • 列为纯样本名：列名不要包含额外描述，比如“sample1_treated”最好另加分组表。
  • 值为非负整数：允许 0，但不能是负数或 NA。

分组信息表：另一张关键表

基因表达矩阵只管“每个样本每个基因有多少 reads”，而 实验设计信息 需要单独的一张表。这两张表通过 样本名 关联。分组表至少包含三列：

  1. sample_id：与表达矩阵的列名完全一致。
  2. group：实验分组，比如 "treated"、"control"。
  3. batch：批次信息，比如 "batch1"、"batch2"。

举个实际例子：

sample_id    group    batch
sample1      control  batch1
sample2      control  batch1
sample3      treated  batch2
sample4      treated  batch2

我建议分组表用 Excel 或 CSV 维护，因为生物学家经常要在里面加注释列（比如性别、测序日期），后续用 model.matrix() 提取就会很方便。

⚠️ 常见陷阱：样本名在两张表中一定要完全一致，大小写、下划线、连字符都不能差。我曾经因为一个样本名叫 "Sample-01" 而在表达矩阵里写成 "Sample_01" ，合并时丢了一整列数据。可以用 setdiff(colnames(counts), pheno$sample_id) 快速检查差异。

读取 counts 文件：用 R 代码实操

我个人习惯把 counts 文件存为 制表符分隔的 .txt，因为基因名里偶尔会带逗号（比如有些数据库用逗号做分隔），CSV 容易引发解析错误。下面是标准读取流程：

# 读取 counts 矩阵
counts <- read.delim("counts_matrix.txt", 
                     row.names = 1, 
                     check.names = FALSE)
# 读取分组信息
pheno <- read.delim("sample_info.txt", 
                    stringsAsFactors = FALSE)
# 确认样本名一致
all(colnames(counts) %in% pheno$sample_id)  # 应返回 TRUE

这里 row.names = 1 把第一列设置为行名，check.names = FALSE 防止 R 把样本名里的短横线变成点。别忘了，你写的代码要给别人读，保留原始名称能省去很多沟通成本。

检查数据完整性：五步自查

只看两眼就往下跑模型，那是勇士的做法。我建议按以下顺序逐一验证：

  1. 行列对应：表达矩阵的列数是否等于分组表的行数？
  2. 样本名排序：两个表的样本名顺序是否一致？用 all.equal(colnames(counts), pheno$sample_id) 检查。
  3. 整数检查：有没有出现小数？用 any(counts != round(counts)) 快速排查。
  4. 缺失值：用 any(is.na(counts)) 或 sum(is.na(counts)) 统计。
  5. 基因过滤：把所有行里 reads 总和 ≤ 10 的基因删掉（下文会详解）。
💡 一个小技巧：如果你用 DGEList() 构建对象，它会自动检查一部分格式问题，但不会检查样本名是否匹配。所以最好在读入后先跑一步 all() 校验。

我自己曾遇到过一个情况：counts 表里样本列顺序跟分组表不一样，但 %in% 返回 TRUE，后来用 match() 重新排序才暴露问题。所以保险代码应该是：

# 按分组表顺序重排列
counts <- counts[, pheno$sample_id]

知识结构速览

下面这张图帮你理清从原始计数到开始建模前的数据流：

edgeR 数据准备流程 原始计数矩阵 分组信息表 完整性检查 构建 DGEList GLM 模型 关键检查点： ① 样本名完全一致 ② 行名唯一 ③ 整数 ④ 无缺失值 ⑤ 过滤低表达基因前先不删列，保持样本顺序一致

常见格式对比

为了让你一眼看清不同输入格式的差异，我把常见搭配列出来：

格式组合优点缺点推荐场景
CSV 表达矩阵 + CSV 分组表Excel 可直接编辑，方便非编程人员基因名带逗号时解析易错小型项目（<10 样本）
TSV 表达矩阵 + CSV 分组表避免逗号冲突，R 中兼容性好不容易用 Excel 单元格预览中等规模（10-50 样本）
RData/Saved RDS 对象保留数据类型，无需重复解析可移植性差，别人需相同 R 版本团队内部重复使用
📋 我的推荐：生产环境用 制表符分隔的 .txt + CSV 分组表。两个文件都放在同一个项目目录下，读取代码写在一个 R 脚本顶部，方便复现。我每次还会在脚本开头加一行 sessionInfo() 记录环境。

数据准备中的避坑清单

最后分享几个真实踩过的坑，希望你能绕过去：

  • 样本名里的特殊字符：比如 "Sample-1 (replicate)" 里的括号和空格，有的平台导出 counts 时自动替换成点，但分组表里没变。解决方法：读取时用 make.names() 统一转换。
  • 忽略基因注释列：有些 counts 文件除了基因名还有染色体位置、长度等列，记得 row.names = 1 只取第一列作行名，其他列用重复。
  • 忘记过滤低表达基因：edgeR 官方文档说得很清楚——low counts 会影响离散度估计。我习惯用 filterByExpr() 自动过滤，但前提是分组表已经正确读入。
  • 用原始 counts 绘图：有次学生画了 counts 的热图，发现样本间差异巨大，其实是因为没做标准化。edgeR 的 cpm() 或 logCPM 才是可视化用的。
← 上一章 📖 返回目录 下一章 →