3、数据清洗与过滤：去除低丰度OTU、低深度样本、质控步骤

数据清洗，说白了就是“去伪存真”。

你从测序公司拿到的 OTU 丰度表，里面混杂着大量“噪音”。低丰度的 OTU、深度不足的样本、还有各种技术误差，都会严重影响后续分析。我刚开始做项目时，就吃过这个亏——拿着未清洗的数据直接做 beta 多样性，结果样本聚类完全混乱，折腾了两天才发现是几个低深度样本在捣鬼。

所以这一章，我们重点解决三个问题：哪些 OTU 该删、哪些样本该去、哪些步骤必须做。

3.1 为什么要清洗数据？

原始 OTU 表本质上是“稀疏矩阵”。绝大多数 OTU 只在少数样本中出现，且丰度极低。如果把所有 OTU 都拿去做分析，会带来三个麻烦：

• 统计噪音：低丰度 OTU 的计数波动大，容易产生假阳性差异
• 计算负担：几万个 OTU 会让距离矩阵计算慢到怀疑人生
• 生物学意义模糊：低丰度 OTU 往往是测序误差或罕见污染，不代表真实物种

此外，样本深度（总 reads 数）差异过大时，直接比较丰度会产生偏差。你想想看，一个样本测了 5 万条 reads，另一个只测了 8 千条，后者会“看起来”物种更少，但这可能是采样不足导致的假象。

3.2 去除低丰度 OTU

这一步我会根据两个维度来过滤：总丰度和出现率。

基于总丰度的过滤：删除所有样本中总计数低于某个阈值的 OTU。我建议先从 5 或 10 开始尝试。如果数据量较大（比如超过 100 个样本），可以把阈值设到 50。嗯，这里要注意，阈值太小过滤不干净，太大又会误删稀有但真实的物种。

基于出现率的过滤：只保留在至少 X% 的样本中出现的 OTU。我个人习惯用 10% 作为默认值。举个例子，如果有 30 个样本，那么 OTU 至少要在 3 个样本中出现才保留。

下面是具体的 R 代码：

# 读取 OTU 丰度表（行为 OTU，列为样本）
otu_table <- read.csv("otu_table.csv", row.names = 1)

# 过滤前查看维度
dim(otu_table)  # 比如输出：[1] 15000    30

# 基于总丰度过滤：删除总计数小于 10 的 OTU
otu_table_filt <- otu_table[rowSums(otu_table) >= 10, ]

# 基于出现率过滤：保留在至少 10% 样本中出现的 OTU
min_samples <- ceiling(0.1 * ncol(otu_table))
otu_table_filt <- otu_table_filt[rowSums(otu_table_filt > 0) >= min_samples, ]

# 过滤后查看维度
dim(otu_table_filt)

我在项目中遇到过一种情况：过滤后 OTU 数量从 1.5 万骤降到 3000 左右。当时还有点慌，觉得删太多了，但后续分析反而更稳定、更可重复。说白了，剩下的才是“靠谱”的物种。

3.3 去除低深度样本

样本深度太低，表征的微生物群落往往不完整。深度不足 10,000 的样本，我会直接标记为“可疑”。如果这批样本整体深度都很低，比如中位数才 5,000，那我就会把阈值降到 3,000——但要记住，这种低深度数据后续分析的可靠性会打折扣。

R 中操作很简单：

# 计算每个样本的总 reads 数
sample_depths <- colSums(otu_table_filt)

# 查看深度分布
summary(sample_depths)

# 设定阈值，去除深度低于 10,000 的样本
keep_samples <- which(sample_depths >= 10000)
otu_table_clean <- otu_table_filt[, keep_samples]

# 同时更新元数据
metadata <- metadata[keep_samples, ]

3.4 质控步骤

除了去除低丰度和低深度，还有一些更精细的质控步骤。我个人习惯按顺序做这几件事：

步骤	说明	R 包/函数
1. 去除污染 OTU	删除线粒体、叶绿体等宿主 DNA 污染的 OTU	phyloseq 的 subset_taxa()
2. 去除低置信 OTU	基于 blast 或人工检查删除可能为嵌合体的 OTU	vsearch / chimera.uchime
3. 去除批次效应	检查不同测序批次或实验批次间的系统性偏差	vegan 的 adonis2() 检测
4. 检查空白对照	空白样本中高丰度的 OTU 往往是实验污染，需去除	手动筛选

去除污染 OTU 的代码示例：

# 假设有个 taxonomy 表，包含各 OTU 的分类信息
# 使用 phyloseq 对象操作
library(phyloseq)

# 去除线粒体（需根据实际分类信息调整筛选逻辑）
ps <- subset_taxa(ps, !grepl("Mitochondria", Rank5))
ps <- subset_taxa(ps, !grepl("Chloroplast", Rank5))

3.5 数据归一化与总结

清洗完成后，数据还需要归一化才能进行后续的 alpha/beta 多样性分析。我常用的方法有三种：

1. 抽平（Rarefaction）：将所有样本抽平到同一深度，简单直接，但会丢弃部分数据
2. TSS（Total Sum Scaling）：将每个样本的 counts 转为相对丰度，保留所有信息
3. CSS（Cumulative Sum Scaling）：metagenomeSeq 包的方法，适合高动态范围的数据

抽平操作：

# 先确定最小的样本深度
min_depth <- min(colSums(otu_table_clean))

# 使用 vegan 包抽平
library(vegan)
otu_table_rarefied <- rrarefy(t(otu_table_clean), sample = min_depth)

到这里，你的数据应该已经“干净”了。我在每次清洗后都会做两件事：一是检查样本间的总 reads 数分布是否更均匀了，二是随机挑几个 OTU 查看它的原始计数和过滤后计数，确认没有误删关键物种。

总结：数据清洗不是机械化的流程，而是需要结合实验设计和个人经验的判断。阈值设高点还是低点，取决于你的研究目的——如果关注核心菌群，可以严格过滤；如果关注稀有物种，就得适当放宽条件。记住，没有“标准答案”，只有“最合适”的参数。