尧图网页设计 
1. 内容整体设计与思路拆解1.1 为什么数据拿回来不能直接跑分析很多刚接触单细胞转录组测序scRNA-seq的朋友从测序公司拿到表达矩阵后第一反应就是直接丢进Seurat或者Scanpy里跑聚类。这个习惯我当年也有直到有一次做出了一群“不明所以”的细胞亚群折腾了整整两周才发现是低质量细胞和双细胞搅在一起。从那以后我彻底明白了单细胞分析的第一步从来不是聚类而是质控而且是极其认真、极其较真的质控。单细胞数据分析的特殊之处在于我们处理的数据本质上是一次又一次的“抽样测量”。每个细胞就是一个独立的样本测序仪把每个细胞里的mRNA分子打成片段、读取、计数最终得到一个整数矩阵行是基因列是细胞矩阵里的数字代表某个基因在该细胞中被捕获到的转录本数量。听起来很干净对不对但实际上这个过程里混进了大量杂质。什么样的情况算杂质一个细胞被酶消化时破裂了胞浆里的mRNA释放到溶液中这些游离的mRNA会被周围的活细胞或者空的微滴捕获于是你可能得到一个“细胞”它里面有几百个基因的微弱信号但实际可能是一个死细胞的残骸。反过来两个细胞靠得太近被包进了同一个微滴或孔里测序仪把它们当成一个“细胞”来读取你会看到一个基因数量异常高、总UMI数异常高的假阳性细胞。还有更隐蔽的问题红细胞污染、环境游离RNAambient RNA造成的背景噪声、线粒体基因占比过高反映的细胞应激状态这些都在质控需要处理的范畴内。我之前做过一个样本组织解离完成后处理得很慢结果测序数据里线粒体基因比例普遍偏高细胞活性数据惨不忍睹。当时如果直接往下游分析大概率会把一群“快死但还没死透”的细胞当成一个真实的细胞类型来注释。这种坑只有做了充分质控的人才能避开。1.2 质控在整个单细胞分析流程中的位置要理解质控的重要性得先看清整个单细胞转录组分析的流程全貌。一个完整的标准分析流程大致是这样的测序数据产出FASTQ→ 比对与定量生成表达矩阵→ 数据质控QC→ 标准化与归一化Normalization→ 高变基因识别HVGs→ 降维PCA→ 聚类Clustering→ 差异表达分析 → 细胞类型注释 → 功能富集等下游分析。可以看到数据质控是位于比对定量之后、所有实质性分析之前的核心关卡。它有双重使命第一把不合格的细胞剔除掉避免它们干扰下游聚类时的细胞群体结构第二把数据中系统性的技术偏差暴露出来为后续归一化方法的选择提供依据。换句话说质控做得彻不彻底直接决定了后面所有分析结果的可靠性。我早期遇到过一位合作者她拿到数据后直接跳过质控就跑聚类结果分出了12个亚群。我帮她补做质控后发现有两个“亚群”其实是低质量细胞和双细胞构成的“垃圾堆”真正的生物信号只剩10个亚群。这个例子想说明的是单细胞数据里的技术噪声不是小问题而是会直接篡改生物学结论的大问题。1.3 核心质控指标的群体画像在深入细节之前我先把质控中最常看的几个指标放在一起做一个整体对比方便大家建立全局观。后续再逐个细讲它们背后的原理和坑。质控指标反映的问题指标偏高的可能原因指标偏低的可能原因细胞计数数据规模与细胞捕获效率测序深度高或双细胞多细胞捕获失败或裂解过度总UMI数每个细胞的测序深度与mRNA捕获效率双细胞或多细胞混合细胞裂解或RNA降解检测到的基因数细胞转录组的复杂程度与细胞活性双细胞、多细胞混合死细胞、空液滴或低质量细胞线粒体基因比例细胞活性与应激状态细胞质RNA流失、细胞破碎增殖旺盛或代谢活跃的细胞核糖体基因比例翻译活性状态部分流程关注代谢活跃细胞细胞应激或处于静息期双细胞比例建库过程中的双重捕获上样浓度过高或微滴加载过量上样浓度适中双细胞率低比对率测序文库质量与比对准确性基因组注释不完整样本降解或测序污染这张表看起来简单但每一条实际上都有大量值得展开的内容。接下来我逐个拆解同时说明怎么定阈值。1.4 质控阈值设置的总体思路阈值怎么设置是单细胞质控里争论最多的问题。很多人希望找到一个“标准答案”比如“线粒体比例超过10%就剔除”“基因数低于200就剔除”。但实际项目里这些通用参数有时并不适用生搬硬套反而会误杀真实细胞群体。我个人的经验是阈值不能拍脑袋要从数据本身的分布出发来判断。线粒体比例高的细胞不一定都是死细胞比如心肌细胞、肝细胞、肾小管上皮细胞这类代谢旺盛的细胞线粒体基因比例本来就高低基因数的细胞也不一定都是空液滴有可能是静止期的T细胞或者成熟的粒细胞。所以我的建议是先用数据分布图和生物学背景相互印证设定一个初步阈值跑完质控后再用下游聚类的结果反向验证如果某些已知的稀有细胞类型消失了就该反思阈值是否太严格了。质控阈值是一个可以迭代、可以调整的决策而不是一条写死的规则。2. 核心细节解析与实操要点2.1 总UMI数和基因数解读测序深度与转录组复杂度UMIUnique Molecular Identifier唯一分子标识符是单细胞测序技术的核心设计之一。建库时每个mRNA分子会被随机连上一段8-10 bp的条形码序列测序后通过比对这段条形码来校正PCR扩增引入的偏差从而更准确地计数原始转录本数量。简单理解UMI就是给每个RNA分子发了一张“身份证”让测序仪知道这个分子的计数是真实存在的而不是PCR扩增出来的假象。因此一个细胞的总UMI数可以近似代表该细胞mRNA的捕获量也就是测序深度。而检测到的基因数则是另一个维度它反映的是这个细胞转录组的复杂程度也就是有多少种不同的基因被表达出来了。这两个指标是强相关的测序深度越高能检测到基因的数量也会趋于饱和但两者并不完全相同。两个细胞总UMI数相同一个基因很分散很多基因都有一点表达另一个基因很集中少数基因表达量极高它们的生物学状态显然不同。在实际质控操作中我需要区分几种典型情况。如果一个细胞的UMI数和基因数同时偏低说明这个细胞的质量很可能有问题可能是细胞即将死亡、mRNA降解、或者捕获效率低如果UMI数和基因数同时异常高则有可能是双细胞污染。但需要注意有些正常细胞天然转录组就比较“丰富”比如巨噬细胞它们的基因检出数天然就高也有一些细胞的基因检出数天然偏低所以在做决定前最好先看一下这个细胞类型在真实生物学背景下的预期特征。2.2 线粒体基因比例细胞活性最灵敏的“晴雨表”线粒体基因比例是质控中最常用的指标之一我个人的经验是这个指标包含了比大家想象中更多的生物信息。线粒体基因主要是MT-ND1、MT-ND2、MT-CO1等与核基因组不同它们有自己的转录调控机制在细胞质中的mRNA稳定性也更高。关键机制在于当细胞膜破碎或细胞即将死亡时细胞质中的mRNA会迅速降解但线粒体因为结构相对独立线粒体中的mRNA降解速度明显慢得多。于是在测序数据上就会表现为线粒体基因占比“虚高”。这就像一间房间着火时里面的纸质文件烧得很快但保险箱里放着的文件保存得更完整于是保险箱文件的相对占比就上升了。这个指标特别适合用来剔除濒死细胞dying cells和破碎细胞broken cells。但问题出在哪里呢不同组织来源的细胞线粒体基因比例的基准线差异很大。胚胎干细胞和诱导多能干细胞通常线粒体基因比例很低可能只有5%左右而肝细胞、心肌细胞这种线粒体丰富的细胞正常状态下线粒体基因比例就可能超过20%。如果拿一个统一的阈值对所有样本进行质控注定会造成生物信息的丢失。我自己在分析公共数据时遇到过这样一个案例一批人胰岛样本α细胞和β细胞的线粒体基因比例普遍偏高接近30%。如果机械地把阈值设置在20%这些真实的内分泌细胞就会被全部剔除剩下的大多是增殖能力较弱、转录活动较低的细胞。所以处理线粒体阈值时一定要结合细胞类型和组织来源来判断。2.3 红细胞与核糖体基因污染两个容易忽略的指标很多人在做质控时只盯着线粒体比例却忽略了另外两个同样重要的指标红细胞基因比例和核糖体基因比例。红细胞基因如HBA1、HBA2、HBB是红细胞特有的标记基因。在组织解离过程中如果血液冲洗不充分红细胞包括有核和无核的会混入细胞悬液。这些红细胞的转录组信号会被一起捕获造成污染。一方面红细胞本身会占据一定的测序深度稀释掉其他真实细胞的信息另一方面红细胞也能形成假细胞这些细胞的特征是基因数低但有大量的血红蛋白基因表达聚类时会聚成一个独立的小群体。核糖体基因如RPS系列和RPL系列的问题则更微妙。核糖体基因在所有细胞中都有较高表达在单细胞数据中它们往往占据大量的UMI计数。更麻烦的是它们对计算细胞周期评分、代谢评分等功能性分析时会造成干扰因为那些分析算法可能把高核糖体表达误判为活跃增殖。不同的分析流程对核糖体基因的处理方式不一致有的流程建议直接去除有的流程建议保留但做回归校正具体采取什么策略取决于你要解决的科学问题。在实操层面如果做的是血液样本我几乎一定会在质控阶段计算红细胞基因比例设定阈值如果做的是实体组织则更关注线粒体比例和基因数这两个核心指标红细胞污染按需排查。2.4 双细胞Doublet检测质控中最容易被低估的环节双细胞污染是指两个或多个细胞被液滴系统或孔板捕获到一个反应单元中在测序后表现为单个“细胞”。它的处理难度比低质量细胞更大因为双细胞往往会通过基因数和UMI数较高的特征被初步识别但仍然有相当一部分双细胞在指标上与正常细胞重叠——比如双细胞是两个低转录活性的细胞合并后的总UMI数可能落在了正常单细胞范围内。双细胞不仅会形成不存在的“伪细胞类型”比如一个巨噬细胞加一个T细胞会同时表达两群基因还会定义出假的细胞亚群或干扰差异表达分析的结果。特别是稀有细胞类型存在时双细胞的比例可能高达总细胞数的5%-10%这对下游分析的干扰非常明显。检测双细胞的算法目前主流的有DoubletFinder、Scrublet、DoubletDecon等。它们的基本思路都是基于一种假设双细胞的基因表达谱是两种细胞类型的“混合体”。算法会模拟生成人工双细胞将真实细胞与仿真双细胞进行对比识别特征与模拟双细胞更相似的真实细胞。我个人的实操经验是不要指望任何一种方法能100%识别出所有双细胞双细胞检测本质上是一个概率判断。推荐的做法是对不同方法识别的双细胞取交集保留高置信度的双细胞集合进行剔除在质控流程中保留“可疑双细胞”的细胞ID在下游聚类后回头查看它们是否聚在一起如果聚成独立群体这个群体很可能就是双细胞群。3. 实操过程与核心环节实现3.1 环境准备与数据读取我日常做单细胞分析用的是R语言加Seurat包这套组合是单细胞转录组分析领域最主流的工具之一文档完善、社区活跃遇到问题通常都能在网上找到解决方案。下面我用Seurat的流程来演示从读取数据到完成质控的全过程。在开始之前需要准备好表达矩阵。如果测序公司返回的是10x Genomics的标准输出通常会有三个文件barcodes.tsv.gz、features.tsv.gz和matrix.mtx.gz。这三个文件分别记录了细胞条形码、基因信息和表达矩阵数据。首先加载必要的包library(Seurat) library(ggplot2) library(dplyr) library(DoubletFinder)如果数据是Cell Ranger输出的标准格式直接用Read10X函数读取counts - Read10X(data.dir path/to/your/filtered_feature_bc_matrix)这里需要注意Cell Ranger的输出分为raw_feature_bc_matrix和filtered_feature_bc_matrix两个文件夹。前者包含所有检测到信号的液滴包括大量的空液滴后者则已经按照Cell Ranger的默认规则做过一次初步过滤。质控阶段我们通常从filtered_feature_bc_matrix开始因为空液滴已经被初步剔除但即便经过这一步数据中仍然可能有死细胞残骸、双细胞以及环境RNA污染。如果条件允许保留原始数据并用raw_feature_bc_matrix重新跑一遍Cell Ranger的细胞鉴定流程如cellranger count中的默认细胞检测可以获得更灵活的把控空间。读取完成后创建Seurat对象并且导入一个重要的元数据字段——线粒体基因比例。线粒体基因在人类数据中的标准写法是以MT-为前缀在小鼠数据中则是mt-前缀。在创建对象时有一个细节我特别要提醒如果你用的是Ensembl基因ID而不是symbol线粒体基因的识别方式会不同需要先通过注释文件转换基因ID。我建议直接用symbol操作这样后续写代码时更直观。seurat_obj - CreateSeuratObject(counts counts, project QC_demo, min.cells 3, min.features 200) seurat_obj[[percent.mt]] - PercentageFeatureSet(seurat_obj, pattern ^MT-)min.cells 3的意思是一个基因至少在3个细胞中有表达才被保留。这个参数过滤掉在极少数细胞中出现的基因能减少极端稀疏表达基因造成的噪声。min.features 200的意思是一个细胞至少检测到200个基因才被保留低于这个值的细胞大概率是空液滴或死亡细胞。3.2 质控指标的分布可视化与阈值设定创建对象后的第一步不是直接跑过滤函数而是先看数据分布。这一步很多新手会跳过但我强烈建议不要省。只有看清你的数据长什么样才能决定阈值设在哪里。用VlnPlot查看关键指标的整体分布VlnPlot(seurat_obj, features c(nFeature_RNA, nCount_RNA, percent.mt), ncol 3, pt.size 0.1)这个图画出来之后你会看到每个样本的基因数分布、UMI数分布和线粒体比例的箱线图。注意观察几个维度基因数是否有大量细胞集中在200左右如果是可能说明细胞活性较差线粒体比例是否有明显的“尾巴”也就是有很大比例的细胞线粒体比例非常高不同样本之间的分布差异是否很大接下来用散点图看双变量关系FeatureScatter(seurat_obj, feature1 nCount_RNA, feature2 nFeature_RNA)这张图非常直观通常nCount_RNA和nFeature_RNA呈正相关。如果出现一部分细胞nCount_RNA极高但nFeature_RNA没有相应提高这些细胞可能是测序深度极高的单个细胞如果两者同时飙升则有双细胞的嫌疑。再画一个线粒体比例与基因数的散点图FeatureScatter(seurat_obj, feature1 nFeature_RNA, feature2 percent.mt)如果看到线粒体比例高的细胞群体基因数低基本可以确定那部分细胞是濒死或已死的细胞如果线粒体比例高的细胞群体基因数也高则需要小心可能是某个代谢旺盛的细胞亚群需要结合生物学背景判断。根据这些图形我就能对阈值有一个初步判断。假设数据分布显示大多数细胞的基因数在800-3000之间线粒体比例在5%以下那我可能会把基因数下限设置为500上限设置为5000线粒体比例设置为10%。3.3 核心过滤操作从阈值设定到细胞筛选设定好阈值后用Seurat的subset函数进行过滤seurat_obj - subset(seurat_obj, subset nFeature_RNA 500 nFeature_RNA 5000 percent.mt 10)这几个参数含义分别是nFeature_RNA 500剔除基因数过少的低质量细胞nFeature_RNA 5000剔除基因数异常高的潜在双细胞percent.mt 10剔除线粒体基因比例过高的濒死细胞。我在这里特别想强调一个容易被忽视的操作点不要只看过滤前后的细胞变化数量一定要逐个检查被过滤的细胞的特征。操作上可以新建一个QC失败的metadata列记录每个细胞被过滤的原因seurat_obj$qc_pass - ifelse( seurat_obj$nFeature_RNA 500 seurat_obj$nFeature_RNA 5000 seurat_obj$percent.mt 10, pass, fail )这样在后续的可视化中你可以随时把QC失败细胞标记出来观察它们在组织中的空间分布或与其他条件的关系。如果在某个特定样本中QC失败的细胞比例特别高这往往意味着该样本的解离过程或保存条件有问题而不是单纯的测序噪声。过滤后还需要做一轮指标复查把过滤前后的指标分布放在一起对比确保没有把某个重要细胞群体“误杀”。如果发现基因数高、线粒体比例正常的细胞被大量剔除就可能意味着阈值太保守了。3.4 双细胞检测实操双细胞检测我通常在完成基本质控后、进入标准化之前进行。操作上推荐使用DoubletFinder这个包在Seurat对象上的集成比较方便。首先需要做一次初步的PCA和聚类DoubletFinder需要利用这些结果来估计双细胞比例seurat_obj - NormalizeData(seurat_obj) seurat_obj - FindVariableFeatures(seurat_obj, selection.method vst, nfeatures 2000) seurat_obj - ScaleData(seurat_obj) seurat_obj - RunPCA(seurat_obj, npcs 30) seurat_obj - RunUMAP(seurat_obj, dims 1:30) seurat_obj - FindNeighbors(seurat_obj, dims 1:30) seurat_obj - FindClusters(seurat_obj, resolution 0.8) # 使用DoubletFinder sweep.res - paramSweep(seurat_obj, PCs 1:30, sct FALSE) sweep.stats - summarizeSweep(sweep.res, GT FALSE) bcmvn - find.pK(sweep.stats) # 根据bcmvn结果选择最优pK值 optimal_pk - as.numeric(as.character(bcmvn$pK[which.max(bcmvn$BCmetric)])) # 估计双细胞比例一般按每1000个细胞约0.8%计算 nExp_poi - round(0.008 * nrow(seurat_objmeta.data)) seurat_obj - doubletFinder(seurat_obj, PCs 1:30, pN 0.25, pK optimal_pk, nExp nExp_poi, sct FALSE)这里有一个常被问到的点双细胞比例估计的总体验一般是每增加1000个细胞双细胞率大约增加0.8%但这是基于10x Genomics的标准上样浓度计算的如果实际实验中的上样浓度比标准浓度高双细胞率会相应升高。对于高浓度上样的文库建议把这个比例上调到1%至1.5%之间。DoubletFinder运行完成后会在metadata里生成一列标注细胞是Singlet还是Doublet。我个人的习惯是不直接删掉所有Doublet而是先看看这些Doublet在聚类图谱上的分布情况——如果它们没有聚成独立群体而是散在各个细胞类型中说明双细胞混杂比较均匀可以放心删除如果它们集中聚成一个或几个群体那么这些群体大概率是从真实细胞类型里“混”出来的需要重新审视之前设定的基因数和UMI数上限是否合理。完成这一步后删掉双细胞seurat_obj - subset(seurat_obj, cells rownames(seurat_objmeta.data[seurat_obj$DoubletFinder_classification Singlet, ]))3.5 质控后的初检标准化与高变基因探索质控完成后在正式进入下游分析之前我通常会再做一遍快速的“血压检测”确保数据质量可靠。流程是标准化、识别高变基因、PCA、UMAP和聚类然后检查几个关键维度。seurat_obj - NormalizeData(seurat_obj, normalization.method LogNormalize, scale.factor 10000) seurat_obj - FindVariableFeatures(seurat_obj, selection.method vst, nfeatures 2000) seurat_obj - ScaleData(seurat_obj) seurat_obj - RunPCA(seurat_obj, npcs 30) seurat_obj - RunUMAP(seurat_obj, dims 1:30) seurat_obj - FindNeighbors(seurat_obj, dims 1:30) seurat_obj - FindClusters(seurat_obj, resolution 0.8)跑完这些之后画出UMAP图然后用已知的细胞类型标记基因进行检测marker_genes - c(PTPRC, CD3D, CD68, MZB1, KLRB1) FeaturePlot(seurat_obj, features marker_genes, ncol 3)如果各个已知标记基因能够富集到对应的独立细胞群体中说明质控后的数据在进行细胞分群时能反映出正确的生物学信号如果标记基因在UMAP上呈弥散分布没有清晰的群体结构那就要警惕了可能是数据本身的质量仍然不好也可能是PCA主成分数选择不当或基因过滤过严。另外质控后的指标复查也包括观察细胞比例结构。比如如果某个样本的细胞存活率太低QC后细胞数不足原有的一半这个样本在下游分析中可能会因为细胞数过少而缺乏统计效力。我在项目中通常会记录一个QC前后的对照表把每个样本的细胞数变化记录下来作为项目交付的一部分。这个记录看起来简单但对于后续分析中出现的任何异常都可以快速回溯到源头上排查是数据本身的问题还是流程操作的问题。4. 常见问题与排查技巧实录4.1 阈值怎么定才对分布驱动与生物背景双轨并进质控阈值这个问题可能是大家问得最多的。动不动就有人问“线粒体基因比例到底应该卡在10%还是20%”“基因数下限是200还是500”。每次遇到这种问题我都很理解问的人的焦虑但答案真的不是固定的。我总结了一套“两看一定”的方法分享给大家。第一看看指标分布图。把nFeature_RNA、nCount_RNA、percent.mt的分布画出来找出明显的“拐点”或“断裂”位置。比如线粒体比例的分布如果呈现双峰那波峰和波谷的边界通常就是设定阈值的好位置如果分布是单峰的那就得从生物背景出发找依据。第二看看与已知生物学知识的匹配度。举个例子如果是人的外周血单个核细胞PBMC数据成年哺乳动物外周血淋巴细胞、单核细胞这些群体的线粒体比例通常不高但如果是心肌组织数据那线粒体比例的基准会高很多。不了解组织的线粒体密度就去设定阈值一定会出问题。一定定阈值之后务必回到UMAP图谱检查已知的细胞类型标记基因是否清晰分离。如果T细胞、B细胞、单核细胞这些本该分得很开的群体出现了混叠但在质控前它们分明是分开的那就说明阈值太严格把决定群体差异的高表达基因数据也一并过滤了。值得警惕的是一些特殊细胞类型非常容易被“误伤”。浆细胞plasma cell的核糖体基因比例极高同时总UMI数和基因数也不低肝细胞和心肌细胞的线粒体比例远超常规阈值中性粒细胞因为RNA含量低基因数经常落在“低质量细胞”的范围内。如果做的是这些组织或细胞类型的研究阈值更需要单独设定。4.2 过滤后细胞数太少怎么办有的样本在质控后细胞数会大幅缩水从1万个细胞变成3000个。这种情况在低质量样本中很常见。处理这个问题的关键是判断“为什么少”。如果是因为解离过程中细胞大量死亡导致线粒体比例普遍升高、基因数普遍降低那么这种样本就算勉强保留后续分析也容易出问题建议直接考虑退回到实验端重新解离建库。如果是因为阈值设置过于严格则可以通过调整阈值来保留更多细胞。这时候需要特别注意区分“真实细胞但表达量低”和“垃圾细胞”的边界。一个有效的策略是用比较宽松的阈值先做一轮聚类观察各类群的标记基因表达和线粒体比例然后针对性地剔除那些确实属于低质量的特征群体。这种迭代式质控在实践中比一次性设阈值更好用。另外还有一个小技巧做同类型多组样本时阈值尽量保持一致。不同样本使用不同阈值会给下游的组间比较引入系统性偏差。比如样本A用10%的线粒体阈值样本B用20%下游差异分析的结果会更难解释。如果各样本的分布差异过大可以先按统一阈值做一版再看各样本的QC后细胞数差异必要时协商调整。4.3 线粒体比例高但基因数正常的细胞是什么这是我被问过最多的一种情况。按常规理解线粒体比例高说明细胞状态不好但基因数正常又说明细胞内mRNA还没有大规模降解。这种组合其实经常出现在一些高代谢活性的细胞中比如激活状态的T细胞、肝细胞、心肌细胞等。在这些细胞中线粒体生物合成增加线粒体转录本数量上升导致比例升高。遇到这种情况我一般先做两件事一是看这些细胞能否聚成独立的群体二是查看群体内是否有特异的生物学标记基因。如果一群细胞特异表达某个已知的与代谢或应激相关的标记基因那么它们很可能是真实的细胞类型不应被当作低质量细胞剔除。反之如果这些细胞分散在各个群体中没有统一的标记基因很可能只是随机性的转录组退化可以按低质量细胞处理。4.4 不同方法或不同批次间的质控标准如何统一多批次数据的质控是所有单细胞分析中最让人头疼的问题之一。每个样本的操作人员可能不同解离时间可能不同建库批次可能不同这些都会导致质控指标分布存在差异。如果不做标准化调整批次效应会混入聚类结果形成以批次为边界的伪群体。在实际处理中我会在单样本质控完成后把各样本的质控指标放在一起做对比。如果发现某些样本的线粒体比例分布的“中位数”明显高于其他样本就需要考虑这个样本的质量对下游分析的潜在影响。处理方式有两种一是为这个样本单独设定更合适的阈值二是在后续分析中加入样本ID作为协变量进行回归校正。需要注意的是虽然Harmony、Seurat的CCA整合等方法可以校正部分批次效应但它们都无法“挽回”严重低质量样本带来的信息损失。质控阶段淘汰低质量样本有时比任何下游校正手段都更有效。我的经验法则是如果一个样本在QC后细胞数不足其它样本的三分之一无论后续用什么方法整合都建议先考虑剔除该样本或重新建库。4.5 质控后仍然有异常聚类群体怎么办有时即便质控做得再仔细UMAP图上还是会出现一个奇怪的群体——它们可能表达多种在生物学上不太可能共现的标记基因或者整体基因表达水平异常低下。这种群体出现时我的第一反应不是直接标注为“未知类型”而是回头去看这个群体的原始QC指标。打开metadata查看这个异常群体的nFeature_RNA、nCount_RNA、percent.mt分布。大概率会出现两种情况一种情况是这些细胞仍然具有低nFeature_RNA和高percent.mt的特征说明这一群体并没有被完全剔除干净需要调整阈值后重新过滤。另一种情况是这些细胞的QC指标看起来正常但无法用已知的生物学知识解释。这时我会先怀疑它们的来源比如是不是红细胞污染、是不是环境RNA的噪声形成了“虚胖”的细胞群或者是两种细胞融合形成的异源双细胞。处理这类群体时可以用FindMarkers找出它们的高表达基因如果高表达基因中同时出现两个不同谱系的特征基因比如同时有T细胞的CD3D和B细胞的MS4A1双细胞的可能性就比较大了。如果所有办法都用完了还是无法解释那就先给这个群体打上一个暂时性的标签如Unknown_Cluster把它们保留在数据中继续往下走而不是强行删除。因为任何异常群体都可能携带尚未被认识的生物学信号过早丢弃反而可能丢掉新的发现。等后续加入更多数据或做空间转录组验证后再作判断是更稳妥的路径。5. 实操心得与习惯养成回到最开始的话题——单细胞分析的第一面为什么一定是数据质控因为单细胞技术的本质决定了它的测量精度受到大量技术噪声的干扰。如果不把噪声和真实信号分离后续所有的分析都会建立在一个不牢固的基础上。我在多次项目实践中逐渐养成了一个习惯每次拿到新数据不管合作者怎么催我一定会先花一整天的时间把质控指标从头到尾看一遍。不只是看统计数据而是逐个样本、逐个批次地检查分布、检查散点图、检查过滤后的细胞数量变化。这个过程看起来很繁琐但恰恰是发现数据问题的关键窗口。最后分享一个小技巧把质控参数和过滤日志记录下来形成一份标准的QC报告。每次分析都按照同样的模板记录输入细胞数、过滤后的细胞数、设定的阈值和过滤理由。这份报告不仅方便自己回溯也是与合作者沟通、写方法学段落时的珍贵素材。我见过太多人在写论文的Methods部分时回忆起自己当时设的质控参数只能给出一个模糊的“with standard QC criteria”这种表达在审稿人眼里是很减分的。数据质控看似是单细胞分析流程中最“没有技术含量”的一步但它恰恰是最能区分一位分析者是“流水线操作工”还是“真正的数据科学家”的分水岭。能踏踏实实把质控做深、做细的人在下游分析中往往也会少踩很多坑。希望这篇文章能帮你跨过单细胞分析的第一道门槛把基础打牢后面的路才能真正走稳。 单细胞转录组数据质控全解析：从UMI到双细胞检测的实操指南
1. 内容整体设计与思路拆解1.1 为什么数据拿回来不能直接跑分析很多刚接触单细胞转录组测序scRNA-seq的朋友从测序公司拿到表达矩阵后第一反应就是直接丢进Seurat或者Scanpy里跑聚类。这个习惯我当年也有直到有一次做出了一群“不明所以”的细胞亚群折腾了整整两周才发现是低质量细胞和双细胞搅在一起。从那以后我彻底明白了单细胞分析的第一步从来不是聚类而是质控而且是极其认真、极其较真的质控。单细胞数据分析的特殊之处在于我们处理的数据本质上是一次又一次的“抽样测量”。每个细胞就是一个独立的样本测序仪把每个细胞里的mRNA分子打成片段、读取、计数最终得到一个整数矩阵行是基因列是细胞矩阵里的数字代表某个基因在该细胞中被捕获到的转录本数量。听起来很干净对不对但实际上这个过程里混进了大量杂质。什么样的情况算杂质一个细胞被酶消化时破裂了胞浆里的mRNA释放到溶液中这些游离的mRNA会被周围的活细胞或者空的微滴捕获于是你可能得到一个“细胞”它里面有几百个基因的微弱信号但实际可能是一个死细胞的残骸。反过来两个细胞靠得太近被包进了同一个微滴或孔里测序仪把它们当成一个“细胞”来读取你会看到一个基因数量异常高、总UMI数异常高的假阳性细胞。还有更隐蔽的问题红细胞污染、环境游离RNAambient RNA造成的背景噪声、线粒体基因占比过高反映的细胞应激状态这些都在质控需要处理的范畴内。我之前做过一个样本组织解离完成后处理得很慢结果测序数据里线粒体基因比例普遍偏高细胞活性数据惨不忍睹。当时如果直接往下游分析大概率会把一群“快死但还没死透”的细胞当成一个真实的细胞类型来注释。这种坑只有做了充分质控的人才能避开。1.2 质控在整个单细胞分析流程中的位置要理解质控的重要性得先看清整个单细胞转录组分析的流程全貌。一个完整的标准分析流程大致是这样的测序数据产出FASTQ→ 比对与定量生成表达矩阵→ 数据质控QC→ 标准化与归一化Normalization→ 高变基因识别HVGs→ 降维PCA→ 聚类Clustering→ 差异表达分析 → 细胞类型注释 → 功能富集等下游分析。可以看到数据质控是位于比对定量之后、所有实质性分析之前的核心关卡。它有双重使命第一把不合格的细胞剔除掉避免它们干扰下游聚类时的细胞群体结构第二把数据中系统性的技术偏差暴露出来为后续归一化方法的选择提供依据。换句话说质控做得彻不彻底直接决定了后面所有分析结果的可靠性。我早期遇到过一位合作者她拿到数据后直接跳过质控就跑聚类结果分出了12个亚群。我帮她补做质控后发现有两个“亚群”其实是低质量细胞和双细胞构成的“垃圾堆”真正的生物信号只剩10个亚群。这个例子想说明的是单细胞数据里的技术噪声不是小问题而是会直接篡改生物学结论的大问题。1.3 核心质控指标的群体画像在深入细节之前我先把质控中最常看的几个指标放在一起做一个整体对比方便大家建立全局观。后续再逐个细讲它们背后的原理和坑。质控指标反映的问题指标偏高的可能原因指标偏低的可能原因细胞计数数据规模与细胞捕获效率测序深度高或双细胞多细胞捕获失败或裂解过度总UMI数每个细胞的测序深度与mRNA捕获效率双细胞或多细胞混合细胞裂解或RNA降解检测到的基因数细胞转录组的复杂程度与细胞活性双细胞、多细胞混合死细胞、空液滴或低质量细胞线粒体基因比例细胞活性与应激状态细胞质RNA流失、细胞破碎增殖旺盛或代谢活跃的细胞核糖体基因比例翻译活性状态部分流程关注代谢活跃细胞细胞应激或处于静息期双细胞比例建库过程中的双重捕获上样浓度过高或微滴加载过量上样浓度适中双细胞率低比对率测序文库质量与比对准确性基因组注释不完整样本降解或测序污染这张表看起来简单但每一条实际上都有大量值得展开的内容。接下来我逐个拆解同时说明怎么定阈值。1.4 质控阈值设置的总体思路阈值怎么设置是单细胞质控里争论最多的问题。很多人希望找到一个“标准答案”比如“线粒体比例超过10%就剔除”“基因数低于200就剔除”。但实际项目里这些通用参数有时并不适用生搬硬套反而会误杀真实细胞群体。我个人的经验是阈值不能拍脑袋要从数据本身的分布出发来判断。线粒体比例高的细胞不一定都是死细胞比如心肌细胞、肝细胞、肾小管上皮细胞这类代谢旺盛的细胞线粒体基因比例本来就高低基因数的细胞也不一定都是空液滴有可能是静止期的T细胞或者成熟的粒细胞。所以我的建议是先用数据分布图和生物学背景相互印证设定一个初步阈值跑完质控后再用下游聚类的结果反向验证如果某些已知的稀有细胞类型消失了就该反思阈值是否太严格了。质控阈值是一个可以迭代、可以调整的决策而不是一条写死的规则。2. 核心细节解析与实操要点2.1 总UMI数和基因数解读测序深度与转录组复杂度UMIUnique Molecular Identifier唯一分子标识符是单细胞测序技术的核心设计之一。建库时每个mRNA分子会被随机连上一段8-10 bp的条形码序列测序后通过比对这段条形码来校正PCR扩增引入的偏差从而更准确地计数原始转录本数量。简单理解UMI就是给每个RNA分子发了一张“身份证”让测序仪知道这个分子的计数是真实存在的而不是PCR扩增出来的假象。因此一个细胞的总UMI数可以近似代表该细胞mRNA的捕获量也就是测序深度。而检测到的基因数则是另一个维度它反映的是这个细胞转录组的复杂程度也就是有多少种不同的基因被表达出来了。这两个指标是强相关的测序深度越高能检测到基因的数量也会趋于饱和但两者并不完全相同。两个细胞总UMI数相同一个基因很分散很多基因都有一点表达另一个基因很集中少数基因表达量极高它们的生物学状态显然不同。在实际质控操作中我需要区分几种典型情况。如果一个细胞的UMI数和基因数同时偏低说明这个细胞的质量很可能有问题可能是细胞即将死亡、mRNA降解、或者捕获效率低如果UMI数和基因数同时异常高则有可能是双细胞污染。但需要注意有些正常细胞天然转录组就比较“丰富”比如巨噬细胞它们的基因检出数天然就高也有一些细胞的基因检出数天然偏低所以在做决定前最好先看一下这个细胞类型在真实生物学背景下的预期特征。2.2 线粒体基因比例细胞活性最灵敏的“晴雨表”线粒体基因比例是质控中最常用的指标之一我个人的经验是这个指标包含了比大家想象中更多的生物信息。线粒体基因主要是MT-ND1、MT-ND2、MT-CO1等与核基因组不同它们有自己的转录调控机制在细胞质中的mRNA稳定性也更高。关键机制在于当细胞膜破碎或细胞即将死亡时细胞质中的mRNA会迅速降解但线粒体因为结构相对独立线粒体中的mRNA降解速度明显慢得多。于是在测序数据上就会表现为线粒体基因占比“虚高”。这就像一间房间着火时里面的纸质文件烧得很快但保险箱里放着的文件保存得更完整于是保险箱文件的相对占比就上升了。这个指标特别适合用来剔除濒死细胞dying cells和破碎细胞broken cells。但问题出在哪里呢不同组织来源的细胞线粒体基因比例的基准线差异很大。胚胎干细胞和诱导多能干细胞通常线粒体基因比例很低可能只有5%左右而肝细胞、心肌细胞这种线粒体丰富的细胞正常状态下线粒体基因比例就可能超过20%。如果拿一个统一的阈值对所有样本进行质控注定会造成生物信息的丢失。我自己在分析公共数据时遇到过这样一个案例一批人胰岛样本α细胞和β细胞的线粒体基因比例普遍偏高接近30%。如果机械地把阈值设置在20%这些真实的内分泌细胞就会被全部剔除剩下的大多是增殖能力较弱、转录活动较低的细胞。所以处理线粒体阈值时一定要结合细胞类型和组织来源来判断。2.3 红细胞与核糖体基因污染两个容易忽略的指标很多人在做质控时只盯着线粒体比例却忽略了另外两个同样重要的指标红细胞基因比例和核糖体基因比例。红细胞基因如HBA1、HBA2、HBB是红细胞特有的标记基因。在组织解离过程中如果血液冲洗不充分红细胞包括有核和无核的会混入细胞悬液。这些红细胞的转录组信号会被一起捕获造成污染。一方面红细胞本身会占据一定的测序深度稀释掉其他真实细胞的信息另一方面红细胞也能形成假细胞这些细胞的特征是基因数低但有大量的血红蛋白基因表达聚类时会聚成一个独立的小群体。核糖体基因如RPS系列和RPL系列的问题则更微妙。核糖体基因在所有细胞中都有较高表达在单细胞数据中它们往往占据大量的UMI计数。更麻烦的是它们对计算细胞周期评分、代谢评分等功能性分析时会造成干扰因为那些分析算法可能把高核糖体表达误判为活跃增殖。不同的分析流程对核糖体基因的处理方式不一致有的流程建议直接去除有的流程建议保留但做回归校正具体采取什么策略取决于你要解决的科学问题。在实操层面如果做的是血液样本我几乎一定会在质控阶段计算红细胞基因比例设定阈值如果做的是实体组织则更关注线粒体比例和基因数这两个核心指标红细胞污染按需排查。2.4 双细胞Doublet检测质控中最容易被低估的环节双细胞污染是指两个或多个细胞被液滴系统或孔板捕获到一个反应单元中在测序后表现为单个“细胞”。它的处理难度比低质量细胞更大因为双细胞往往会通过基因数和UMI数较高的特征被初步识别但仍然有相当一部分双细胞在指标上与正常细胞重叠——比如双细胞是两个低转录活性的细胞合并后的总UMI数可能落在了正常单细胞范围内。双细胞不仅会形成不存在的“伪细胞类型”比如一个巨噬细胞加一个T细胞会同时表达两群基因还会定义出假的细胞亚群或干扰差异表达分析的结果。特别是稀有细胞类型存在时双细胞的比例可能高达总细胞数的5%-10%这对下游分析的干扰非常明显。检测双细胞的算法目前主流的有DoubletFinder、Scrublet、DoubletDecon等。它们的基本思路都是基于一种假设双细胞的基因表达谱是两种细胞类型的“混合体”。算法会模拟生成人工双细胞将真实细胞与仿真双细胞进行对比识别特征与模拟双细胞更相似的真实细胞。我个人的实操经验是不要指望任何一种方法能100%识别出所有双细胞双细胞检测本质上是一个概率判断。推荐的做法是对不同方法识别的双细胞取交集保留高置信度的双细胞集合进行剔除在质控流程中保留“可疑双细胞”的细胞ID在下游聚类后回头查看它们是否聚在一起如果聚成独立群体这个群体很可能就是双细胞群。3. 实操过程与核心环节实现3.1 环境准备与数据读取我日常做单细胞分析用的是R语言加Seurat包这套组合是单细胞转录组分析领域最主流的工具之一文档完善、社区活跃遇到问题通常都能在网上找到解决方案。下面我用Seurat的流程来演示从读取数据到完成质控的全过程。在开始之前需要准备好表达矩阵。如果测序公司返回的是10x Genomics的标准输出通常会有三个文件barcodes.tsv.gz、features.tsv.gz和matrix.mtx.gz。这三个文件分别记录了细胞条形码、基因信息和表达矩阵数据。首先加载必要的包library(Seurat) library(ggplot2) library(dplyr) library(DoubletFinder)如果数据是Cell Ranger输出的标准格式直接用Read10X函数读取counts - Read10X(data.dir path/to/your/filtered_feature_bc_matrix)这里需要注意Cell Ranger的输出分为raw_feature_bc_matrix和filtered_feature_bc_matrix两个文件夹。前者包含所有检测到信号的液滴包括大量的空液滴后者则已经按照Cell Ranger的默认规则做过一次初步过滤。质控阶段我们通常从filtered_feature_bc_matrix开始因为空液滴已经被初步剔除但即便经过这一步数据中仍然可能有死细胞残骸、双细胞以及环境RNA污染。如果条件允许保留原始数据并用raw_feature_bc_matrix重新跑一遍Cell Ranger的细胞鉴定流程如cellranger count中的默认细胞检测可以获得更灵活的把控空间。读取完成后创建Seurat对象并且导入一个重要的元数据字段——线粒体基因比例。线粒体基因在人类数据中的标准写法是以MT-为前缀在小鼠数据中则是mt-前缀。在创建对象时有一个细节我特别要提醒如果你用的是Ensembl基因ID而不是symbol线粒体基因的识别方式会不同需要先通过注释文件转换基因ID。我建议直接用symbol操作这样后续写代码时更直观。seurat_obj - CreateSeuratObject(counts counts, project QC_demo, min.cells 3, min.features 200) seurat_obj[[percent.mt]] - PercentageFeatureSet(seurat_obj, pattern ^MT-)min.cells 3的意思是一个基因至少在3个细胞中有表达才被保留。这个参数过滤掉在极少数细胞中出现的基因能减少极端稀疏表达基因造成的噪声。min.features 200的意思是一个细胞至少检测到200个基因才被保留低于这个值的细胞大概率是空液滴或死亡细胞。3.2 质控指标的分布可视化与阈值设定创建对象后的第一步不是直接跑过滤函数而是先看数据分布。这一步很多新手会跳过但我强烈建议不要省。只有看清你的数据长什么样才能决定阈值设在哪里。用VlnPlot查看关键指标的整体分布VlnPlot(seurat_obj, features c(nFeature_RNA, nCount_RNA, percent.mt), ncol 3, pt.size 0.1)这个图画出来之后你会看到每个样本的基因数分布、UMI数分布和线粒体比例的箱线图。注意观察几个维度基因数是否有大量细胞集中在200左右如果是可能说明细胞活性较差线粒体比例是否有明显的“尾巴”也就是有很大比例的细胞线粒体比例非常高不同样本之间的分布差异是否很大接下来用散点图看双变量关系FeatureScatter(seurat_obj, feature1 nCount_RNA, feature2 nFeature_RNA)这张图非常直观通常nCount_RNA和nFeature_RNA呈正相关。如果出现一部分细胞nCount_RNA极高但nFeature_RNA没有相应提高这些细胞可能是测序深度极高的单个细胞如果两者同时飙升则有双细胞的嫌疑。再画一个线粒体比例与基因数的散点图FeatureScatter(seurat_obj, feature1 nFeature_RNA, feature2 percent.mt)如果看到线粒体比例高的细胞群体基因数低基本可以确定那部分细胞是濒死或已死的细胞如果线粒体比例高的细胞群体基因数也高则需要小心可能是某个代谢旺盛的细胞亚群需要结合生物学背景判断。根据这些图形我就能对阈值有一个初步判断。假设数据分布显示大多数细胞的基因数在800-3000之间线粒体比例在5%以下那我可能会把基因数下限设置为500上限设置为5000线粒体比例设置为10%。3.3 核心过滤操作从阈值设定到细胞筛选设定好阈值后用Seurat的subset函数进行过滤seurat_obj - subset(seurat_obj, subset nFeature_RNA 500 nFeature_RNA 5000 percent.mt 10)这几个参数含义分别是nFeature_RNA 500剔除基因数过少的低质量细胞nFeature_RNA 5000剔除基因数异常高的潜在双细胞percent.mt 10剔除线粒体基因比例过高的濒死细胞。我在这里特别想强调一个容易被忽视的操作点不要只看过滤前后的细胞变化数量一定要逐个检查被过滤的细胞的特征。操作上可以新建一个QC失败的metadata列记录每个细胞被过滤的原因seurat_obj$qc_pass - ifelse( seurat_obj$nFeature_RNA 500 seurat_obj$nFeature_RNA 5000 seurat_obj$percent.mt 10, pass, fail )这样在后续的可视化中你可以随时把QC失败细胞标记出来观察它们在组织中的空间分布或与其他条件的关系。如果在某个特定样本中QC失败的细胞比例特别高这往往意味着该样本的解离过程或保存条件有问题而不是单纯的测序噪声。过滤后还需要做一轮指标复查把过滤前后的指标分布放在一起对比确保没有把某个重要细胞群体“误杀”。如果发现基因数高、线粒体比例正常的细胞被大量剔除就可能意味着阈值太保守了。3.4 双细胞检测实操双细胞检测我通常在完成基本质控后、进入标准化之前进行。操作上推荐使用DoubletFinder这个包在Seurat对象上的集成比较方便。首先需要做一次初步的PCA和聚类DoubletFinder需要利用这些结果来估计双细胞比例seurat_obj - NormalizeData(seurat_obj) seurat_obj - FindVariableFeatures(seurat_obj, selection.method vst, nfeatures 2000) seurat_obj - ScaleData(seurat_obj) seurat_obj - RunPCA(seurat_obj, npcs 30) seurat_obj - RunUMAP(seurat_obj, dims 1:30) seurat_obj - FindNeighbors(seurat_obj, dims 1:30) seurat_obj - FindClusters(seurat_obj, resolution 0.8) # 使用DoubletFinder sweep.res - paramSweep(seurat_obj, PCs 1:30, sct FALSE) sweep.stats - summarizeSweep(sweep.res, GT FALSE) bcmvn - find.pK(sweep.stats) # 根据bcmvn结果选择最优pK值 optimal_pk - as.numeric(as.character(bcmvn$pK[which.max(bcmvn$BCmetric)])) # 估计双细胞比例一般按每1000个细胞约0.8%计算 nExp_poi - round(0.008 * nrow(seurat_objmeta.data)) seurat_obj - doubletFinder(seurat_obj, PCs 1:30, pN 0.25, pK optimal_pk, nExp nExp_poi, sct FALSE)这里有一个常被问到的点双细胞比例估计的总体验一般是每增加1000个细胞双细胞率大约增加0.8%但这是基于10x Genomics的标准上样浓度计算的如果实际实验中的上样浓度比标准浓度高双细胞率会相应升高。对于高浓度上样的文库建议把这个比例上调到1%至1.5%之间。DoubletFinder运行完成后会在metadata里生成一列标注细胞是Singlet还是Doublet。我个人的习惯是不直接删掉所有Doublet而是先看看这些Doublet在聚类图谱上的分布情况——如果它们没有聚成独立群体而是散在各个细胞类型中说明双细胞混杂比较均匀可以放心删除如果它们集中聚成一个或几个群体那么这些群体大概率是从真实细胞类型里“混”出来的需要重新审视之前设定的基因数和UMI数上限是否合理。完成这一步后删掉双细胞seurat_obj - subset(seurat_obj, cells rownames(seurat_objmeta.data[seurat_obj$DoubletFinder_classification Singlet, ]))3.5 质控后的初检标准化与高变基因探索质控完成后在正式进入下游分析之前我通常会再做一遍快速的“血压检测”确保数据质量可靠。流程是标准化、识别高变基因、PCA、UMAP和聚类然后检查几个关键维度。seurat_obj - NormalizeData(seurat_obj, normalization.method LogNormalize, scale.factor 10000) seurat_obj - FindVariableFeatures(seurat_obj, selection.method vst, nfeatures 2000) seurat_obj - ScaleData(seurat_obj) seurat_obj - RunPCA(seurat_obj, npcs 30) seurat_obj - RunUMAP(seurat_obj, dims 1:30) seurat_obj - FindNeighbors(seurat_obj, dims 1:30) seurat_obj - FindClusters(seurat_obj, resolution 0.8)跑完这些之后画出UMAP图然后用已知的细胞类型标记基因进行检测marker_genes - c(PTPRC, CD3D, CD68, MZB1, KLRB1) FeaturePlot(seurat_obj, features marker_genes, ncol 3)如果各个已知标记基因能够富集到对应的独立细胞群体中说明质控后的数据在进行细胞分群时能反映出正确的生物学信号如果标记基因在UMAP上呈弥散分布没有清晰的群体结构那就要警惕了可能是数据本身的质量仍然不好也可能是PCA主成分数选择不当或基因过滤过严。另外质控后的指标复查也包括观察细胞比例结构。比如如果某个样本的细胞存活率太低QC后细胞数不足原有的一半这个样本在下游分析中可能会因为细胞数过少而缺乏统计效力。我在项目中通常会记录一个QC前后的对照表把每个样本的细胞数变化记录下来作为项目交付的一部分。这个记录看起来简单但对于后续分析中出现的任何异常都可以快速回溯到源头上排查是数据本身的问题还是流程操作的问题。4. 常见问题与排查技巧实录4.1 阈值怎么定才对分布驱动与生物背景双轨并进质控阈值这个问题可能是大家问得最多的。动不动就有人问“线粒体基因比例到底应该卡在10%还是20%”“基因数下限是200还是500”。每次遇到这种问题我都很理解问的人的焦虑但答案真的不是固定的。我总结了一套“两看一定”的方法分享给大家。第一看看指标分布图。把nFeature_RNA、nCount_RNA、percent.mt的分布画出来找出明显的“拐点”或“断裂”位置。比如线粒体比例的分布如果呈现双峰那波峰和波谷的边界通常就是设定阈值的好位置如果分布是单峰的那就得从生物背景出发找依据。第二看看与已知生物学知识的匹配度。举个例子如果是人的外周血单个核细胞PBMC数据成年哺乳动物外周血淋巴细胞、单核细胞这些群体的线粒体比例通常不高但如果是心肌组织数据那线粒体比例的基准会高很多。不了解组织的线粒体密度就去设定阈值一定会出问题。一定定阈值之后务必回到UMAP图谱检查已知的细胞类型标记基因是否清晰分离。如果T细胞、B细胞、单核细胞这些本该分得很开的群体出现了混叠但在质控前它们分明是分开的那就说明阈值太严格把决定群体差异的高表达基因数据也一并过滤了。值得警惕的是一些特殊细胞类型非常容易被“误伤”。浆细胞plasma cell的核糖体基因比例极高同时总UMI数和基因数也不低肝细胞和心肌细胞的线粒体比例远超常规阈值中性粒细胞因为RNA含量低基因数经常落在“低质量细胞”的范围内。如果做的是这些组织或细胞类型的研究阈值更需要单独设定。4.2 过滤后细胞数太少怎么办有的样本在质控后细胞数会大幅缩水从1万个细胞变成3000个。这种情况在低质量样本中很常见。处理这个问题的关键是判断“为什么少”。如果是因为解离过程中细胞大量死亡导致线粒体比例普遍升高、基因数普遍降低那么这种样本就算勉强保留后续分析也容易出问题建议直接考虑退回到实验端重新解离建库。如果是因为阈值设置过于严格则可以通过调整阈值来保留更多细胞。这时候需要特别注意区分“真实细胞但表达量低”和“垃圾细胞”的边界。一个有效的策略是用比较宽松的阈值先做一轮聚类观察各类群的标记基因表达和线粒体比例然后针对性地剔除那些确实属于低质量的特征群体。这种迭代式质控在实践中比一次性设阈值更好用。另外还有一个小技巧做同类型多组样本时阈值尽量保持一致。不同样本使用不同阈值会给下游的组间比较引入系统性偏差。比如样本A用10%的线粒体阈值样本B用20%下游差异分析的结果会更难解释。如果各样本的分布差异过大可以先按统一阈值做一版再看各样本的QC后细胞数差异必要时协商调整。4.3 线粒体比例高但基因数正常的细胞是什么这是我被问过最多的一种情况。按常规理解线粒体比例高说明细胞状态不好但基因数正常又说明细胞内mRNA还没有大规模降解。这种组合其实经常出现在一些高代谢活性的细胞中比如激活状态的T细胞、肝细胞、心肌细胞等。在这些细胞中线粒体生物合成增加线粒体转录本数量上升导致比例升高。遇到这种情况我一般先做两件事一是看这些细胞能否聚成独立的群体二是查看群体内是否有特异的生物学标记基因。如果一群细胞特异表达某个已知的与代谢或应激相关的标记基因那么它们很可能是真实的细胞类型不应被当作低质量细胞剔除。反之如果这些细胞分散在各个群体中没有统一的标记基因很可能只是随机性的转录组退化可以按低质量细胞处理。4.4 不同方法或不同批次间的质控标准如何统一多批次数据的质控是所有单细胞分析中最让人头疼的问题之一。每个样本的操作人员可能不同解离时间可能不同建库批次可能不同这些都会导致质控指标分布存在差异。如果不做标准化调整批次效应会混入聚类结果形成以批次为边界的伪群体。在实际处理中我会在单样本质控完成后把各样本的质控指标放在一起做对比。如果发现某些样本的线粒体比例分布的“中位数”明显高于其他样本就需要考虑这个样本的质量对下游分析的潜在影响。处理方式有两种一是为这个样本单独设定更合适的阈值二是在后续分析中加入样本ID作为协变量进行回归校正。需要注意的是虽然Harmony、Seurat的CCA整合等方法可以校正部分批次效应但它们都无法“挽回”严重低质量样本带来的信息损失。质控阶段淘汰低质量样本有时比任何下游校正手段都更有效。我的经验法则是如果一个样本在QC后细胞数不足其它样本的三分之一无论后续用什么方法整合都建议先考虑剔除该样本或重新建库。4.5 质控后仍然有异常聚类群体怎么办有时即便质控做得再仔细UMAP图上还是会出现一个奇怪的群体——它们可能表达多种在生物学上不太可能共现的标记基因或者整体基因表达水平异常低下。这种群体出现时我的第一反应不是直接标注为“未知类型”而是回头去看这个群体的原始QC指标。打开metadata查看这个异常群体的nFeature_RNA、nCount_RNA、percent.mt分布。大概率会出现两种情况一种情况是这些细胞仍然具有低nFeature_RNA和高percent.mt的特征说明这一群体并没有被完全剔除干净需要调整阈值后重新过滤。另一种情况是这些细胞的QC指标看起来正常但无法用已知的生物学知识解释。这时我会先怀疑它们的来源比如是不是红细胞污染、是不是环境RNA的噪声形成了“虚胖”的细胞群或者是两种细胞融合形成的异源双细胞。处理这类群体时可以用FindMarkers找出它们的高表达基因如果高表达基因中同时出现两个不同谱系的特征基因比如同时有T细胞的CD3D和B细胞的MS4A1双细胞的可能性就比较大了。如果所有办法都用完了还是无法解释那就先给这个群体打上一个暂时性的标签如Unknown_Cluster把它们保留在数据中继续往下走而不是强行删除。因为任何异常群体都可能携带尚未被认识的生物学信号过早丢弃反而可能丢掉新的发现。等后续加入更多数据或做空间转录组验证后再作判断是更稳妥的路径。5. 实操心得与习惯养成回到最开始的话题——单细胞分析的第一面为什么一定是数据质控因为单细胞技术的本质决定了它的测量精度受到大量技术噪声的干扰。如果不把噪声和真实信号分离后续所有的分析都会建立在一个不牢固的基础上。我在多次项目实践中逐渐养成了一个习惯每次拿到新数据不管合作者怎么催我一定会先花一整天的时间把质控指标从头到尾看一遍。不只是看统计数据而是逐个样本、逐个批次地检查分布、检查散点图、检查过滤后的细胞数量变化。这个过程看起来很繁琐但恰恰是发现数据问题的关键窗口。最后分享一个小技巧把质控参数和过滤日志记录下来形成一份标准的QC报告。每次分析都按照同样的模板记录输入细胞数、过滤后的细胞数、设定的阈值和过滤理由。这份报告不仅方便自己回溯也是与合作者沟通、写方法学段落时的珍贵素材。我见过太多人在写论文的Methods部分时回忆起自己当时设的质控参数只能给出一个模糊的“with standard QC criteria”这种表达在审稿人眼里是很减分的。数据质控看似是单细胞分析流程中最“没有技术含量”的一步但它恰恰是最能区分一位分析者是“流水线操作工”还是“真正的数据科学家”的分水岭。能踏踏实实把质控做深、做细的人在下游分析中往往也会少踩很多坑。希望这篇文章能帮你跨过单细胞分析的第一道门槛把基础打牢后面的路才能真正走稳。