第3章 数据准备：输入数据格式要求，count矩阵构建，样本信息表的准备
讲DESeq2之前，我得先强调一件事：数据准备比跑差异分析本身更容易翻车。我接过不少"为什么跑出来全是NA"的问题，一查十有八九是输入格式出了岔子。说白了，DESeq2的要求不复杂，但稍微不注意，就给你报一堆看不懂的错。
本章咱们就把输入数据这件事掰开揉碎。主要说三块：
- count矩阵到底长啥样？
- sampleTable怎么建才能让DESeq2认出来？
- 常见报错和我的避坑经验。
3.1 count矩阵：整数、无缺失、行为基因
DESeq2接收的是原始整数计数（raw counts），不是FPKM/TPM，也不是log后的值。为什么？因为它内部要用负二项分布建模，小数或者归一化后的数据会破坏统计假设。
3.1.1 基本格式
我通常把count矩阵存成一个data.frame或matrix，行是基因（或转录本），列是样本。比如这样：
| Gene | Sample1 | Sample2 | Sample3 |
|---|---|---|---|
| GeneA | 120 | 95 | 134 |
| GeneB | 0 | 5 | 2 |
| GeneC | 456 | 389 | 501 |
- 所有值必须是非负整数（0可以，但不要有NA）
- 行名必须是基因ID（Ensembl ID、symbol等，但不能有重复）
- 列名必须是样本ID，与sampleTable中的样本名完全一致
核心：DESeq2只认整数。小数点、负数、NA进来，它会直接罢工。
3.1.2 从哪里获取count矩阵？
常用的上游流程：
- featureCounts：输出一个表，包含Geneid和每个样本的count列。直接read.table读进来，取需要的列即可。
- HTSeq-count：每个样本一个文件，需要合并。我习惯用
merge或cbind手动拼。 - Salmon / Kallisto / RSEM：输出的是估计count或TPM，需要先取"estimated counts"列，然后合并。
举个例子，从featureCounts输出构建count矩阵：
raw_counts <- read.table("featureCounts_results.txt",
header=TRUE, row.names=1)
# 一般第1列是Geneid，后面列是样本名
# 取出需要的样本列（假设是第7到第9列）
count_matrix <- raw_counts[, 7:9]
colnames(count_matrix) <- c("Control1","Control2","Treat1")注意：很多工具生成的count矩阵最后会有几列统计信息（如Chr、Start、End、Strand），记得去掉。
3.1.3 常见坑：非整数与缺失值
我曾经接手一个项目，跑DESeq2时报错“每个基因至少有一个样本大于0”，结果查了半天发现有一列是TPM（带小数）。还有个坑是Salmon输出的estimated counts有时候是浮点数，虽然大部分接近整数，但DESeq2还是会报warning。我的做法是：
# 四舍五入取整
count_matrix <- round(count_matrix)
# 确保没有NA
count_matrix[is.na(count_matrix)] <- 0注意：不要用round处理salmon的估计count？我一般不推荐，因为DESeq2官方说可以用tximport直接生成整数count。后面第5章会细讲tximport。如果你手头只有浮点count，可以试试DESeqDataSetFromMatrix时设置countData = round(countData)，但最好用专门的工具。
3.2 样本信息表（sampleTable）：设计矩阵的“说明书”
DESeq2需要你提供一个data.frame，告诉它每个样本属于什么条件、批次等。这就是sampleTable。我习惯存成一个.csv文件，读进来保证列名规范。
3.2.1 基本结构
至少有两列：rownames（样本ID，与count矩阵列名一致）和因子列（如condition、batch）。
| sample | condition | batch |
|---|---|---|
| Control1 | Control | 1 |
| Control2 | Control | 1 |
| Treat1 | Treat | 2 |
| Treat2 | Treat | 2 |
- 样本ID必须唯一，且与count矩阵的列名一一对应（连拼写、大小写都不能错）
- 因子列尽量写成字符或因子类型，不要用数字1/2代表条件，除非你清楚自己在做什么
- batch等协变量列可以是非数值（如"batch1","batch2"），也可以是数值（但要作为因子处理）
3.2.2 构建示例
假设你有6个样本，3个对照、3个处理，来自两个批次：
sampleTable <- data.frame(
sample = c("C1","C2","C3","T1","T2","T3"),
condition = c("Control","Control","Control","Treat","Treat","Treat"),
batch = c("batchA","batchA","batchA","batchB","batchB","batchB")
)
rownames(sampleTable) <- sampleTable$sample我自己写代码时，习惯把rownames(sampleTable)设为样本名，方便后面DESeqDataSetFromMatrix直接按行名匹配。
3.2.3 与count矩阵对齐
运行DESeq2前，务必保证count矩阵的列顺序与sampleTable的行顺序一致。DESeq2内部会按rownames匹配，但如果有样本在count里不在sampleTable里，会报错。
我常用的检查：
all(colnames(count_matrix) %in% rownames(sampleTable))
all(rownames(sampleTable) %in% colnames(count_matrix))小技巧：直接用match重排count矩阵的列，让它们按sampleTable的顺序排好：count_matrix <- count_matrix[, rownames(sampleTable)]。
3.3 一张图理清数据准备流程
下面这张SVG图概括了从原始测序数据到DESeq2输入的完整准备步骤。你可以照着这个流程自查。
3.4 常见报错及我的避坑指南
3.4.1 "Error in DESeqDataSet(se, design) : every gene contains at least one zero"
这个错其实不是真的说所有基因都是零——而是提示你的count矩阵里有缺失值（NA）或者非整数。DESeq2内部要求所有值非负整数且没有缺失。我曾经遇到一个case，是因为HTSeq输出某个样本的计数文件里有一行"no_feature"，被我读进来变成了NA。只要用is.na排查一下就能找到。
处理办法：
- 检查count矩阵是否有NA
- 检查是否有小数（用
any(count_matrix %% 1 != 0)） - 如果确实有基因在所有样本都是0，不用担心，DESeq2会自动过滤掉低表达基因。但报这个错通常不是全零的问题。
3.4.2 "Error: object 'condition' not found"
这个我碰到过好几次：设计公式里写了~ condition，但sampleTable里没有叫"condition"的列。有时候列名是"Condition"（大写C），或"group"。DESeq2是大小写敏感的。
建议：写设计公式之前，先colnames(sampleTable)看一眼。
3.4.3 "Error in checkForExperimentalErrors(object) : the model matrix is not full rank"
这是多重共线性问题。常见于：你把批次和条件完全混淆了（比如所有对照在batch1，处理都在batch2），设计公式写成~ batch + condition时，DESeq2无法区分。解决方法要么去掉batch项，要么重新设计实验让批次和条件交叉。
重要：批次校正不是万能的。如果你的实验设计本身就有混淆（batch和condition完全一致），任何统计方法都拆不开。设计实验时就要考虑交叉或部分重复。
3.4.4 count矩阵行名重复
有时候合并多个样本的计数时，基因名重复（比如不同转录本有相同gene symbol）。DESeq2会报错“duplicate rownames”。我习惯合并前先检查一下：any(duplicated(rownames(count_matrix)))。如果有重复，可以用aggregate按基因名求和，或者用make.names(unique=TRUE)加后缀区分。但我个人更推荐前者——同一个gene symbol的计数应该合并。
3.5 小结
数据准备阶段，花10分钟检查count矩阵和sampleTable，可以省掉后面几小时debug时间。我每次跑新项目都会跑下面这几行代码，算是我的个人习惯：
# 基础检查
stopifnot(is.matrix(count_matrix) || is.data.frame(count_matrix))
stopifnot(all(count_matrix == floor(count_matrix), na.rm=TRUE)) # 整数
stopifnot(!any(is.na(count_matrix))) # 无缺失
stopifnot(all(colnames(count_matrix) %in% rownames(sampleTable))) # 列名匹配
stopifnot(!any(duplicated(rownames(count_matrix)))) # 无重复基因只要这几行过了，DESeq2的输入数据基本不会出幺蛾子了。