📚 课程目录第2章差异表达分析概述

2. 单细胞数据读取与质控：从10X格式到可用对象

拿到10X的原始数据，第一件事是什么？不是跑差异分析，也不是画图。先读进来，再做一遍质控。这一步做不好，后面的分析全是白搭。我见过太多人直接跳过质控，结果下游差异基因找出来全是线粒体基因，嗯，那画面太美我不敢看。

今天咱们就把这个流程走一遍。从原始的10X格式数据，到创建Seurat对象，再到过滤掉那些“坏细胞”。每一步我都会把实际项目里的坑指出来。

2.1 读取10X格式数据

10X Genomics 官方输出通常是三个文件：barcodes.tsv.gz、features.tsv.gz（以前叫genes.tsv）、matrix.mtx.gz。Seurat 提供了一个函数 Read10X，直接搞定。

library(Seurat)

# 假设数据放在 data/pbmc3k/ 目录下
data_dir <- "data/pbmc3k/"
counts <- Read10X(data_dir)

这个函数返回一个稀疏矩阵。为什么要用稀疏矩阵？因为10X数据里大部分是零，稀疏矩阵能省内存。你想想看，单细胞数据动辄几万个细胞，两万个基因，如果全用稠密矩阵存，内存早炸了。

个人习惯：我通常会把读取后的 counts 对象快速看一眼维度：dim(counts)。确认基因数和细胞数跟预期一致。有一次我读错目录，拿到了空矩阵，差点多花两小时debug。

2.2 创建Seurat对象

拿到counts矩阵后，下一步就是创建Seurat对象。这是Seurat分析的核心数据结构。里面不仅存表达矩阵，还存元数据、降维结果、差异标记等。创建时我们只需要传入原始计数。

pbmc <- CreateSeuratObject(counts = counts,
                           project = "pbmc3k",
                           min.cells = 3,
                           min.features = 200)

这里有两个参数很关键：min.cells 和 min.features。

min.cells = 3：基因至少在3个细胞中有表达，否则过滤掉。这能去掉那些几乎不表达的“噪音基因”。
min.features = 200：细胞至少检测到200个基因，否则过滤掉。这能去掉空液滴或死细胞。

这两个参数是经验值。我记得最早做10X PBMC数据时，官方教程就是这两个数。后来我试过调大min.features到500，对于某些质量差的数据反而更合理。但刚上手时，先按标准值来。

2.3 质控指标计算

Seurat对象创建好之后，我们需要算几个质控指标：每个细胞的UMI总数、检测到的基因数、以及线粒体基因占比。这些指标能帮我们识别出异常细胞。

# 线粒体基因 —— 以 MT- 开头的人类基因
pbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")

为什么要算线粒体比例？说白了，如果一个细胞大量表达线粒体基因，说明它的细胞膜可能破了，细胞质里的mRNA跑掉了，线粒体RNA反而被富集。这种细胞属于“濒死细胞”，必须剔除。

其他常用质控指标还有：核糖体基因比例、血红蛋白基因比例（红细胞污染）。不过初期我们只需要关注UMI总数、基因数和线粒体比例这三个就够了。

注意：人类线粒体基因前缀是 ^MT-，小鼠是 ^mt-。我有个同事没注意物种，在小鼠数据上用 ^MT-，结果线粒体比例全是0。嗯，排查了好久。

2.4 质控过滤阈值

拿到指标后，怎么设阈值？没有绝对标准，但有几个常用范围。我一般先画小提琴图或散点图，看看分布，再定阈值。

下面是一个典型过滤策略：

指标	常用阈值	过滤原因
nCount_RNA（UMI总数）	500 ~ 10000	过少是空液滴，过多可能是双细胞
nFeature_RNA（基因数）	200 ~ 5000	过少是死细胞/空液滴，过多可能是双细胞
percent.mt（线粒体%）	< 20%（免疫细胞可更严）	线粒体高 -> 细胞受损

比如PBMC数据，我常用：

pbmc <- subset(pbmc, subset = nFeature_RNA > 200 &
                              nFeature_RNA < 2500 &
                              percent.mt < 5)

注意，这里的 percent.mt < 5 比较严格，因为我做的PBMC数据质量很好。如果是肿瘤组织或冷冻组织，线粒体比例可能到20%甚至更高。千万别死板套用数字，先看数据。

核心原则：质控阈值要根据数据的实际分布来定。先可视化，再决策。我习惯用 VlnPlot(pbmc, features = c("nFeature_RNA", "nCount_RNA", "percent.mt")) 快速扫一眼。

2.5 完整流程梳理

把上面几步串起来，就是标准的数据读取+质控流程。我画了一张流程图，方便你对照操作。

从图里可以看到，整个流程清晰：原始文件 → 读取 → 构建对象 → 计算指标 → 过滤。每一步都有对应的Seurat函数，千万别跳步。

2.6 常见问题与避坑

灵长类线粒体前缀：人类 ^MT-，小鼠 ^mt-，斑马鱼是 ^mt- 但名称略有不同。不确定时先看基因名列表。
过滤阈值不固定：我曾在某肿瘤数据中因为设了 percent.mt < 5% 导致过滤掉60%细胞，后来放宽到20%反而结果更好。先看数据的指标分布图。
空液滴 vs 真实细胞：有些协议会包含空液滴检测（如DropletUtils的emptyDrops），如果数据来自10X v3以上且没有细胞过滤，最好额外做一步。不过对于标准PBMC数据，min.features=200 已经足够。

一个实用小技巧：过滤完成后，及时查看Seurat对象的元数据表 head(pbmc@meta.data)。我会确保三个关键列（nCount_RNA, nFeature_RNA, percent.mt）都在，并且数值范围合理。这一步只需要10秒，但能省下后续排查的大把时间。

2.7 本节小结

这一章我们完成了单细胞分析的“地基工作”：读取10X格式、创建Seurat对象、计算质控指标、按阈值过滤。别小看这几行代码，实际项目中80%的数据问题都能在质控阶段暴露出来。

质控完成后，你的Seurat对象里只有那些“干净”的细胞。后面做归一化、找高变基因、降维聚类，结果才可靠。记住：垃圾进，垃圾出。先把数据质量守住，后面的分析才能站得住脚。