4. scRNA-seq 数据读取与质量控制

数据读取,听着简单,但真上手的时候,坑真不少。

我第一次拿 10X 数据的时候,习惯性地用了 Read10X 函数。结果呢?跑出来一堆空细胞,指标全飘红。嗯,从那以后,我再也不敢跳过前两步了——先读进来,再仔细看每个 barcode 的质量信号。

4.1 读取 10X 数据

大多数 scRNA-seq 数据都是 10X 格式。我个人习惯用 Seurat::Read10X,或者 DropletUtils::read10xCounts。注意,现在官方给的 H5 文件也越来越常见,用 Read10X_h5 就能直接读。

library(Seurat)
# 读取 raw feature barcode matrix
raw_counts <- Read10X(data.dir = "filtered_feature_bc_matrix/")
# 或者 H5 格式
raw_counts <- Read10X_h5("filtered_feature_bc_matrix.h5")
小提示: 如果用的是 Cell Ranger 输出的 raw_feature_bc_matrix,里面包含所有 barcode(包括空滴和背景)。我建议用 raw 数据做质控,这样你能自己决定过滤阈值,而不是黑盒处理。

读取之后,我一般会检查一下稀疏矩阵的维度、总 UMI 数和基因数。别急着往下跑,先看一眼直方图,心里有个数。

4.2 质量控制的核心指标

质量控制说白了就是三件事:过滤掉质量差的细胞、排除 doublet、保留真实单细胞。我通常按这个顺序来:

  • UMI count 过低 → 可能是空滴或破碎细胞
  • 基因数太少 → 细胞活性低或背景噪音
  • 线粒体基因比例过高 → 细胞膜损伤,cytoplasmic RNA 泄漏
  • doublet 检测 → 两个细胞被标记为一个 barcode

我遇到过一台数据,UMI 中位数有 5000,但线粒体比例高达 40%。一查,原来是细胞冻存时间太长,膜都破了。如果不滤掉这些细胞,下游聚类会多出一堆“应激”簇。

4.2.1 过滤阈值的选择

阈值没有绝对标准,但我通常这样设置:

指标典型阈值原因
min.features200 – 500低于 200 的很可能只是背景 RNA
min.UMI500 – 1000捕集效率低,信噪比差
max.percent.mt10% – 20%人类细胞通常 < 10%,某些组织(肝、肌)可放宽到 20%
max.features5000 – 7000超出范围很可能就是 doublet 或污染

我曾经在肺癌样本上用过 10% 的线粒体阈值,结果滤掉了大量正常的肿瘤细胞——肿瘤细胞代谢率高,线粒体比例天然偏高。后来我改成根据数据分布手动截取:把线粒体比例高的小峰单独拿出来看。

4.3 双细胞检测(doublet detection)

doublet 是 scRNA-seq 里最常见的假象。两个细胞跑进同一个液滴,你看到的是一个“混合”的转录组。如果不处理,聚类里可能会出现中间态的伪簇,特别难解释。

常用的方法有三类:

  • DoubletFinder – 基于人工模拟 doublet 和参数化预测
  • Scrublet – 用 kNN 图计算 doublet 得分
  • cxds / bcds – 基于基因共表达或比例特征

我习惯用 DoubletFinder,因为它不需要设置复杂的参数,而且兼容 Seurat 流程。当然,你要注意它接受的是 UMI counts 矩阵,不是标准化后的。

library(DoubletFinder)
# 先做基本过滤,再归一化、降维
seu <- NormalizeData(seu) %>% FindVariableFeatures() %>% ScaleData() %>% RunPCA()
# 预测 doublet
seu <- doubletFinder_v3(seu, PCs = 1:10, pN = 0.25, pK = 0.09, nExp = nExp_poi, reuse.pANN = FALSE)
# 结果存在 metadata 中,列名类似 "DF.classifications_0.25_0.09_XXX"
doublet_cells <- seu$DF.classifications_0.25_0.09_598
table(doublet_cells)
注意: DoubletFinder 里 nExp 这个参数很关键。我建议不要直接用自动估算,而是结合 cellranger 的 expect doublet rate 来设定。比如你上机了 10000 个细胞,预期 doublet 率 5%,那就设 nExp = 500。太保守会漏,太激进会误杀真正的细胞。

Scrublet 也是一个不错的选择,尤其在处理大规模数据集时,速度更快。用法:

library(scrublet)
scrub <- Scrublet$new(counts_matrix)
results <- scrub$scrub_doublets()
doublet_scores <- results$doublet_scores
predicted_doublets <- results$predicted_doublets

我个人的经验是:不要完全依赖一种方法。我通常会同时跑 DoubletFinder 和 Scrublet,然后把两个方法的预测取交集(都是 doublet 的才过滤掉),或者取并集(只要有一个判为 doublet 就滤掉),取决于你的数据量和容忍度。

4.4 完整的质控流程

下面这张 SVG 图概括了我常用的读取+质控管线。箭头表示处理顺序,每个节点都对应一个关键的决策点。

scRNA-seq 读取与质量控制流程 读取原始数据 UMI / 基因数 / 线粒体过滤 Doublet 检测 标记并移除 Doublet 高质量单细胞矩阵 质控后 cell x gene 矩阵 可进入降维和聚类

4.5 常见坑与避坑指南

  • 空滴过滤不足: 用 EmptyDrops 方法可以判断哪些 barcode 是真的细胞,哪些是背景。我碰到过一次,UMI 中位数只有 300,但用 EmptyDrops 一测,竟然还有 2000 多个真实细胞,说明阈值不能一刀切。
  • 线粒体比例的样本差异: 同一实验的不同样本,线粒体比例可能相差很大。我习惯先分别计算每个样本的 cutoff,而不是统一用 20%。
  • Doublet 检测的假阳性: 如果 doublet 率预期很低(比如捕获了 500 个细胞),Scrublet 的预测可能不准。这时候我会检查 doublet 候选细胞在中高维空间的分布,如果它们聚成一个紧密簇,大概率是真正的 doublet。
  • 过滤顺序: 先过滤低质量细胞,再做 doublet 检测。如果先做 doublet,低质量细胞的异常表达谱可能被误判为 doublet。

我记得有一回,一个客户的数据 doublet 率异常高(15%),我看了半天发现是上机时细胞浓度太高。我后来建议他重新调整稀释比例,果然 doublet 率降到了 3% 以下。所以质控不只是数据处理,也是实验设计的一面镜子。

4.6 质控后的检查清单

在进入下游分析前,我通常会用下面几个指标确认质控效果:

检查项预期结果
细胞数较原始减少通常保留 60%–80%
UMI 分布直方图单峰,没有很长的低右尾
基因数与 UMI 的关系正相关,R² > 0.8
线粒体比例柱状图大部分细胞低于阈值,无异常高值
Doublet 标注散点图doublet 细胞通常分布在 UMI 较高的区域或离群位置

只要这五项都通过,你就可以放心地进入归一化和高变基因选择了。质控虽然繁琐,但省掉这一步,后面的分析全是白费功夫。嗯,数据清洗这件事,花再多时间都不过分。