企业官网建设流程全解析

单细胞转录组分析（scRNA-seq）这几年几乎成了肿瘤免疫、发育生物学、神经科学里最常见的起手式。我最早接触这个技术是做一批肿瘤样本的异质性解析，当时手里只有bulk RNA-seq的结果，肿瘤微环境里的免疫细胞到底由哪些亚群组成、每个亚群的表达谱有什么特征，完全是一团黑箱。换成scRNA-seq之后，每个样本从一张平均表达谱变成了成千上万个细胞各自一张表达谱，信息量大了几个量级，但麻烦也随之而来：稀疏矩阵、批次效应、低质量细胞、双细胞污染……每一步都可能把结果带偏。这篇文章我尽量用一线操作的视角，把从原始数据到生物学结论的完整处理链路、每一步为什么这么做、哪些参数不能照抄，一次讲透。适合刚拿到第一批单细胞数据、正准备跑通全流程的人阅读，也适合已经被质控和下游分析折磨过一轮、想真正搞明白每一步原理的人参考。

1. scRNA-seq到底测到了什么：从组织到稀疏矩阵的完整链路

1.1 bulk表达看不见的东西：细胞异质性才是单细胞的核心价值

很多人第一次接触单细胞的时候，第一反应是“它不就是把RNA-seq做得更精细吗？”这个理解不算错，但很容易让你在后面分析的时候走偏。常规转录组测序的对象是一团组织或者一群细胞，最后得到的是一张平均表达谱，就像一个班里所有人的成绩平均分——你只知道这个班整体水平不错，但完全不知道谁是尖子生、谁在拖后腿。而肿瘤、免疫、发育这些场景，恰恰是异质性在起决定性作用：同一块肿瘤组织里，既有增殖活跃的癌细胞，也有耗竭的T细胞，还有各类髓系细胞和成纤维细胞，它们的转录状态天差地别。平均化之后，很多关键信号直接被抹平了。

scRNA-seq的价值就在于把“平均分”拆回了“每个人的成绩单”。每个细胞单独作为一个样本，我们可以回答三个层面的问题：组织里有哪些细胞类型和状态、每种细胞类型的占比是多少、这些细胞类型之间有什么状态转换或通讯关系。这三个问题对应的就是聚类注释、比例分析、轨迹推断和细胞通讯分析，后续的所有分析流程都是围绕它们展开的。

但必须提前说清楚：单细胞测序不是万能的。它测的是RNA，不是蛋白表达，更不是细胞在组织里的真实空间位置。而且由于每个细胞的RNA总量极低，捕获和扩增过程存在巨大的技术噪音，这直接导致了scRNA-seq数据里一个非常核心的特征——稀疏性。

1.2 一条细胞一个barcode：10x平台的建库与UMI去重逻辑

现在做单细胞转录组，绝大多数人用的都是10x Genomics Chromium平台，也是我后面要讲的操作里默认的平台。10x的核心逻辑可以概括为“给每个细胞贴一个条码”。它先把细胞悬液和带有凝胶微珠的油滴混合，形成成千上万个油包水结构。只要控制好细胞浓度，绝大多数油滴（GEMs）里只包含一个细胞和一个微珠。这个微珠上挂着两样关键东西：一个是Cell Barcode，一段用来标记“这个转录本来自哪个细胞”的序列；另一个是UMI（Unique Molecular Identifier，唯一分子标识符），一段用来标记“这条RNA分子是哪一个原始分子”的随机序列。

建库测序之后，生物信息分析的第一步就是根据Cell Barcode把reads分配到对应的细胞，再根据UMI把同一个转录本的重复扩增reads合并。这个思路很像快递包裹上的双标签：收件人地址是barcode，包裹单号是UMI。只要这两个标签都在，几千条来自同一个转录本的重复reads也能被认为是“同一个原始RNA分子”，从而避免PCR扩增偏差导致定量虚高。这也是10x方案相对于早期Smart-seq2方案的一个巨大优势：Smart-seq2没有UMI，判断基因表达量只能靠比对位置加reads数，PCR偏好性会直接影响基因表达定量的准确性。

有一件事很多人初学时容易忽略：细胞悬液的浓度和质量直接决定了后面数据的质量。如果上样浓度不合适，油滴里要么大部分是空的，要么一个油滴包了两个细胞，后者会在数据里表现为“双细胞（doublet）”。这个问题的处理我后面会专门讲，因为它是质控环节里最隐蔽、也最影响聚类结果的问题之一。

1.3 为什么单细胞矩阵这么稀疏：dropout与真实零表达的边界

如果你第一次打开单细胞表达矩阵，大概率会被吓一跳：一个几万人乘两万基因的矩阵，非零值占比通常只有5%到15%，剩下全是0。这是scRNA-seq数据最典型的特征——稀疏性。稀疏性的来源有两个，一个是生物学上真实的“这个基因在这个细胞里没表达”，另一个是技术上的“这个基因其实表达了，但因为RNA捕获效率太低，没有被检测到”，后者在专业上叫dropout事件。

需要特别强调的是，低捕获效率是单细胞技术固有的物理限制。一个哺乳动物细胞里大约有10到20万条mRNA分子，但10x平台的实际捕获率通常只有10%到30%，也就是说大部分转录本根本没被反转录和扩增。这导致两个后果：第一，基因表达量在单细胞层面是“数字游戏”，一个基因在某个细胞里检测到5个UMI和检测到10个UMI，很可能只是抽样误差，并不代表表达量差了一倍；第二，很多低表达基因在相当比例的细胞里直接变成了0，这种0和生物学意义上的“不表达”混在一起，很难区分。

这个特性是后面所有分析方法设计的起点。比如归一化、聚类、差异分析，都要考虑到数据的稀疏性，不能用传统bulk RNA-seq的分析逻辑去硬套。你如果看到有人直接把单细胞矩阵当成普通表达矩阵跑edgeR或者DESeq2，基本可以判断这个分析是有问题的。理解了稀疏性，你也就理解了为什么主流的Seurat、Scanpy流程都强调用特定的归一化方法、特定的特征基因筛选策略，因为每一步都在和dropout噪音对抗。

2. 拿到原始测序数据后的第一轮处理：从fastq到基因表达矩阵

2.1 Cell Ranger count之前，先确认这三件事

不管是自己测序还是从公共数据库下载数据，你拿到的原始数据通常是fastq或者BCL文件。如果是公共数据，多数情况下可以直接拿到cellranger count的输入——两个fastq文件或者一个包含所有数据的文件夹。但在跑cellranger count之前，我建议先确认三件事，三件事都搞清楚了再动命令，否则后面会浪费大量时间重跑。

第一，确认样本名称和10x建库类型。10x平台有3'转录组和5'转录组之分，还有针对VDJ、ATAC等不同文库的流程。如果你的数据是3'转录组，就应该用cellranger count；如果是5'转录组加上VDJ，那还要跑cellranger vdj，两者分析的逻辑完全不同，千万别用错参数。第二，确认参考基因组和注释版本。人的话常用GRCh38，小鼠常用GRCm39，但你一定要知道这个数据当初比对的是哪个版本。如果从公共数据库下载已经比好的bam或表达矩阵，那无所谓；但如果要自己跑fastq，参考基因组版本不一致会导致基因注释对不上，尤其是做后续差异分析和富集分析的时候，版本错了基因ID都乱套。第三，确认测序深度和读长。10x平台在大多数情况下推荐的双端reads是R1 28bp、R2 90bp或更长，如果你拿到的fastq读长很奇怪，需要看是不是测序有问题。

2.2 比对与定量：cellranger内部到底在做什么

cellranger count的工作流程可以分成五步：解复用（demultiplexing）、比对、barcode处理、UMI去重和基因定量。第一步，它会根据样本索引把混合的reads拆分成单个样本；第二步，用STAR比对器把reads比对到参考基因组上；第三步，根据Cell Barcode把reads分到不同细胞；第四步，根据UMI去掉重复扩增的reads；第五步，统计每个细胞中每个基因的UMI counts，最终生成表达矩阵。

实际操作中，我最常用来启动分析的命令行是这样的：

cellranger count \ --id=sample1 \ --transcriptome=/path/to/refdata-cellranger-GRCh38-2020-A \ --fastqs=/path/to/fastq_dir \ --sample=sample1 \ --localcores=16 \ --localmem=64

这里有几个参数值得多说几句。--sample参数对应的是fastq文件名里的样本名前缀，10x的fastq文件名中间有一段是样本名，cellranger就是靠这个匹配的，如果你把--sample写错，它会找不到对应的fastq，或者把整个文件夹下的杂信号都读进来。--transcriptome指向参考基因组目录，这个目录名字里通常包含参考基因组版本和cellranger配套注释版本，建议不同物种、不同版本分开存放，不要混用。--localcores和--localmem控制运行资源，16核64G内存跑一个10x样本通常问题不大，但如果细胞数特别多，或者你用的是更加耗资源的参考基因组，内存不够会直接OOM报错。

cellranger count跑完以后，会在输出目录里生成一大堆文件。你最需要关注的是outs/目录下的几个：filtered_feature_bc_matrix是过滤空液滴后的表达矩阵，绝大多数下游分析都用这个；raw_feature_bc_matrix是未过滤的矩阵，主要在评估空液滴、做背景去除的时候用；web_summary.html是一个可视化报告，里面包含了测序饱和度、细胞数、比对比率等质控指标；cloupe.cloupe文件是给10x官方可视化软件用的，我基本不用，但如果你想快速浏览数据，打开它看看也方便。

2.3 产出文件里哪些该看、哪些容易误读

cellranger的web_summary.html是很多人在拿到结果后第一个打开的东西，但这里有两个指标特别容易误导人。一个是“Estimated Number of Cells”，这个数字是根据barcode rank plot的拐点自动估计的，它只是一个估算值，不是金标准。如果样本质量不好或者细胞类型比较特殊（比如细胞RNA含量差异很大），这个估计值可能偏低或偏高，后面还需要结合质控重新判断。另一个是“Fraction Reads in Cells”，它代表有多少reads被分配到了“细胞”的barcode里面，比例越高说明背景噪音越低，但如果你用的是一个本身表达量就很低的样本，这个比例低不一定代表测序失败。

另外还有一个很多人误读的指标叫“Valid Barcodes”和“Q30 Bases”。Valid Barcodes高说明barcode质量好，Q30高说明碱基测序准确度高，但这两个指标再好，也弥补不了上游细胞悬液质量差带来的问题。细胞活率低、细胞碎片多、细胞粘连，这些在测序报告里不会直接显示，却会在表达矩阵里制造出大量低质量细胞和双细胞信号。所以我的习惯是，cellranger跑完以后，web_summary.html先大致看一眼，真正决定过滤策略的一定是在R或Python里对表达矩阵做进一步质控时，结合细胞本身的表达特征来判断。

3. 质控不是筛掉低质量细胞那么简单：阈值选择背后的逻辑

3.1 三层质控指标：UMI数、基因数、线粒体比例的联动关系

拿到filtered_feature_bc_matrix之后，大部分教程会让你画三个指标的分布图：每个细胞的总UMI数、每个细胞的检出基因数、每个细胞的线粒体基因表达比例。这三个指标通常被称为单细胞质控的“三件套”，但很多人只是机械地照着阈值切，没有真正理解它们之间的联动关系。

总UMI数衡量的是这个细胞的转录本捕获量。捕获量过低，可能是细胞状态不好、RNA降解严重，或者压根就是一个空的液滴残留。基因数衡量的是转录组复杂度，它和总UMI数有正相关关系，但又不是完全同步——一个基因表达量很高、其余基因都沉默的细胞，UMI数可能不低，但基因数很低。线粒体基因表达比例衡量的是细胞质完整性：正常情况下，线粒体基因占比在5%到20%之间；如果一个细胞受损、细胞质内的mRNA流失，线粒体转录本因为在线粒体内相对稳定，占比会异常升高。所以这三个指标是联动的，不能只看某一个。一个细胞如果基因数低、线粒体比例高，基本可以定为低质量细胞；但一个细胞基因数低、线粒体比例正常，未必是坏细胞，它可能是某种本来就转录沉默的细胞类型，比如静息期的T细胞或者红细胞。

实际处理时，我不会在一开始就设置一个绝对阈值，而是先看所有样本的分布图，找到明显的拐点。比如UMI数分布图上如果有一个独立的低值峰和主峰之间有明显间隔，那这个低值峰就是明显的低质量细胞群。如果主峰很平滑、没有明显拐点，就需要用相对宽松的阈值，并记录自己切掉的比例。一般我会控制在总细胞数的5%到15%以内，如果一批样本切掉了超过25%的细胞，那上游的细胞制备流程一定有问题，这时候应该回到湿实验去排查，而不是在分析上不断调低阈值自欺欺人。

3.2 空液滴和双细胞：容易被忽略的两类污染

除了低质量细胞，还有两类污染是新手最容易忽略的：空液滴和双细胞。空液滴是指一个油滴里没有细胞，但仍然带有一小段微珠上的barcode，里面可能捕获了环境中游离的RNA，表现为一个“细胞”有极低的UMI数和基因数。cellranger自带的EmptyDrops算法就是专门用来区分真实细胞和空液滴的，它通过比较每个barcode的RNA组成和背景噪音的差异来做判断。用cellranger count生成的filtered_feature_bc_matrix已经做过一轮空液滴过滤，但你仍然需要警惕那些勉强通过过滤、信号极弱的barcode。

双细胞则是两个细胞被包在同一个油滴里，它会表现为一个barcode下有混合的转录组信号。双细胞的危害在于它可能被聚类成一个虚假的中间态细胞群，这在做分化轨迹分析的时候尤其致命——你以为发现了一个过渡态的细胞，其实只是两个不同细胞类型的混合信号。处理双细胞目前的主流工具是DoubletFinder（R语言）和Scrublet（Python），原理都是通过模拟人工双细胞、比较每个细胞与模拟双细胞的特征相似度来给每个细胞打分。我的建议是：DoubletFinder可以在Seurat流程里方便地跑，参数一般用默认的pN=0.25，pK根据细胞数和PC数自动选择，不要随便改；另外DropletUtils包里的doubletCells也可以用来做替代验证。

需要注意的是，双细胞的比例和上样细胞浓度直接相关，10x官方一般推荐捕获5000个细胞时双细胞率约0.4%，捕获20000个细胞时约1.6%，如果你在高通量上样方案下拿到极高的双细胞率，这是正常现象，但下游分析时一定要把双细胞去除掉，否则聚类注释会变得一团糟。

3.3 我实践中的质控顺序与Seurat实操代码

质控的顺序很关键，我通常的做法是：先做基础的细胞过滤（UMI、基因数、线粒体比例），再去除双细胞，最后才进入归一化和降维。如果顺序颠倒，先跑了PCA和聚类再回头看质控指标，那些异常群会让你的聚类受到污染，特征基因选择也会被带偏。

下面是我最常用的Seurat质控代码，你可以直接复制使用：

library(Seurat) library(DoubletFinder) # 读入cellranger结果 data_dir <- "path/to/filtered_feature_bc_matrix" counts <- Read10X(data_dir) obj <- CreateSeuratObject(counts = counts, project = "sample1", min.cells = 3, min.features = 200) # 质控指标 obj[["percent.mt"]] <- PercentageFeatureSet(obj, pattern = "^MT-") # 实际阈值需要结合分布图来确定，这里只是示例 obj <- subset(obj, subset = nFeature_RNA > 200 & nFeature_RNA < 6000 & percent.mt < 10) # 先做标准归一化，用于DoubletFinder obj <- NormalizeData(obj) obj <- FindVariableFeatures(obj, selection.method = "vst", nfeatures = 2000) obj <- ScaleData(obj) obj <- RunPCA(obj) obj <- FindNeighbors(obj, dims = 1:20) obj <- FindClusters(obj, resolution = 0.5) # 双细胞检测 sweep.res <- paramSweep(obj, PCs = 1:20, sct = FALSE) sweep.stats <- summarizeSweep(sweep.res, GT = FALSE) bcmvn <- find.pK(sweep.stats) pK <- as.numeric(as.character(bcmvn$pK[bcmvn$BCmetric == max(bcmvn$BCmetric)])) annotations <- obj@meta.data$seurat_clusters homotypic.prop <- modelHomotypic(annotations) nExp_poi <- round(0.05 * nrow(obj@meta.data)) nExp_poi.adj <- round(nExp_poi * (1 - homotypic.prop)) obj <- doubletFinder(obj, PCs = 1:20, pN = 0.25, pK = pK, nExp = nExp_poi.adj, sct = FALSE) # 列名需要根据实际输出调整 Doublet_col <- grep("Doublet", colnames(obj@meta.data), value = TRUE) obj <- subset(obj, cells = rownames(obj@meta.data)[obj@meta.data[[Doublet_col]] == "Singlet"])

这段代码有几个地方要提醒你：第一个是pattern = "^MT-"，这是人类基因名的写法，如果是小鼠，要改成"^mt-"；第二个是min.cells = 3表示一个基因至少在3个细胞里检出来才保留，这个参数会影响后续的基因总数，如果不是特别缺数据，我建议保留这个默认值；第三个是nFeature_RNA < 6000这个上限不是所有物种都适用，它在免疫细胞、肿瘤细胞和神经元之间的适用性差异很大，如果你做的是干细胞或者大而复杂的细胞，基因数上限可以放宽到8000甚至10000，但一定要先看分布。

4. 归一化、降维与聚类：参数是把双刃剑，选“差不多”而不是“最优”

4.1 归一化为什么选LogNormalize，什么时候考虑SCTransform

经过质控以后，你会得到一个过滤后的表达矩阵。接下来面临的问题是归一化。为什么要归一化？因为不同细胞捕获到的总RNA量不同，导致同样的基因在不同细胞里的counts没有直接可比性。比如一个细胞里某个基因是100 UMI、总UMI是10000；另一个细胞里同一个基因是200 UMI、总UMI是50000，从绝对数值看后者比前者高，但从比例看前者其实更高。归一化就是为了消除这个“测序深度”差异。

Seurat最常用的LogNormalize方法的计算公式是：每个基因的counts除以该细胞总counts，乘以一个缩放因子（默认10000），再取log1p。这个方法的优点是简单、稳定、可解释性好，适合绝大多数场景。但它的缺点也很明显：对高表达基因和低表达基因的方差处理做得不够好，而且它假设所有细胞的RNA组成是一致的，这个假设在包含分化程度差异极大的细胞群里并不成立。

如果数据里有明显的细胞状态差异、扩增批次复杂，或者你准备做跨样本整合，我更推荐用SCTransform。SCTransform用一个正则化负二项回归模型来同时建模测序深度和基因表达的关系，然后再做标准化。它的好处是可以减少测序深度带来的虚假技术差异，同时能更好地平衡高表达和低表达基因的权重。实际操作中，SCTransform跑起来比LogNormalize慢不少，尤其是大样本（几万细胞）的时候，内存和时间的开销要心里有数。我的经验是：先跑LogNormalize做探索性分析，如果聚类结果里出现了明显由测序深度驱动的细胞群（比如所有UMI特别高的细胞聚在一起），再换SCTransform验证一把，效果确实有提升就说明换得对。

4.2 主成分数、UMAP困惑度、聚类分辨率怎么定更稳妥

降维是单细胞分析中最容易被“玄学化”的一步。很多人直接照抄教程里的dims = 1:20，然后就看图。但PCA到底取多少个主成分，直接影响后面聚类和UMAP的结果，这里其实有一个相对客观的判断方法：看PC的P值分布。Seurat的JackStraw工具就是做这个用的，它通过置换检验来判断每个PC是否显著。另一种更简便的方法是用ElbowPlot，画出每个PC解释的方差比例，找到“肘部”位置。如果肘部不明显，我一般会多取一些PC，比如20到30之间，宁可多也不能少，因为信息冗余总比信息丢失好。但要注意，取太多PC会把部分技术噪音也纳入聚类，所以也不是越多越好。

UMAP的两个核心参数是n.neighbors和min.dist。这两个参数决定了UMAP嵌入的局部和全局结构展示。n.neighbors越大，看到的整体结构越平滑，小群落越容易被糊掉；min.dist越小，细胞在图上挤得越紧，群之间的边界越清晰。我通常用n.neighbors = 30, min.dist = 0.3作为初始参数，但一定要明白，UMAP只是可视化工具，聚类的结果不是由UMAP决定的，而是由cluster算法决定的。实际聚类用的是基于SNN图的社区发现算法（Louvain或Leiden），你可能在UMAP图上看到两群细胞分得很开，但如果聚类分辨率不够，它们可能还是被划成一个cluster。

聚类分辨率（resolution）是个磨人的参数。分辨率越高，聚类得到的簇越多，反之越少。我没有一个万能的分辨率值，但我的策略是：从resolution = 0.1开始，逐渐调到0.5甚至1.0，观察每个簇是否表达已知的细胞类型marker基因。如果你在某个分辨率下能得到一群清晰的、有明确marker的细胞类型，而不是把两个已知的细胞类型混在一起，那这个分辨率就是合适的。切记，聚类不是越细越好，超过真实生物学状态的划分只会让后续注释变成一场灾难。

4.3 注释细胞类型：marker基因是参考，不是判决书

聚类完成之后，就到了我认为整个单细胞分析中最需要经验、也最容易翻车的环节：细胞类型注释。几乎所有教程都会给你一张marker基因表，比如T细胞用CD3D、B细胞用MS4A1/CD79A、巨噬细胞用LYZ/CD68、NK细胞用NKG7……但实际数据往往不会这么干净。这是因为marker基因本身有组织特异性，同一类细胞在不同组织里的marker谱会有差异，而且肿瘤微环境中的细胞经常处于异常状态，比如肿瘤相关巨噬细胞（TAM）可能同时表达M1和M2的marker。

我的注释流程通常是“两步法”。第一步，用一套基础的“谱系marker”把所有细胞大致分到几个大类：免疫（PTPRC/CD45）、上皮（EPCAM/KRT）、间质（COL1A1/PDGFRB）、内皮（PECAM1/VWF）等。第二步，在每个大类内部，再根据细分的marker做亚群注释。比如免疫细胞内部先区分髓系（LYZ/CD68/FCGR3A）和淋巴系（CD3D/CD79A），再在T细胞内部看CD4、CD8A/CD8B、FOXP3（调节性T）、CCL5（效应记忆T）等。这种分层注释的思路比从头到尾用一个超长marker列表做一次聚类省力得多，也更容易发现异常群体。

在使用marker基因注释时，有几个坑需要特别提醒。第一，不要只看一个marker，至少用3到5个marker联合判断，因为单个基因容易出现dropout，某个细胞型里CD3D丢了一个也不奇怪；第二，要结合RNA的细胞来源背景，比如肺组织里的免疫细胞和肝组织里的免疫细胞marker谱会有差异；第三，遇到一些跨类别的marker，比如LILRB4在单核细胞/树突状细胞/巨噬细胞里都可能有表达，需要结合其他marker和聚类位置综合判断。如果实在无法确定某个簇的身份，我会先标记为“Unknown”，不要硬套名字，等后续做差异分析或者更多功能验证后再决定。

5. 差异分析、富集与批次效应：哪些分析必须做，哪些容易翻车

5.1 差异基因：从细胞水平到样本水平，pseudobulk为什么更可靠

聚类注释完成以后，最常见的下一步是想知道组间（比如疾病组和对照组）或者细胞亚群之间有哪些差异基因。这里有一个在单细胞领域反复被强调但又反复被违反的原则：如果你的实验设计里有多个生物学重复样本，做组间差异分析时不要用单个细胞作为独立样本。原因非常简单，同一只小鼠或同一个病人的细胞之间存在很强的相关性，如果把每个细胞都当成独立的样本去统计，你的样本量会被严重虚增，最终导致p值极度膨胀，产生大量假阳性差异基因。这在统计上叫“伪重复（pseudoreplication）”。

正确的做法是先在样本层面汇总，再进行差异比较，这就是所谓的“pseudobulk”策略。具体来说，就是把每个样本内所有某类细胞的表达量加起来（或者取中位数），这样每个样本每个细胞亚群就得到一个表达谱，然后用DESeq2或edgeR这类传统的bulk差异分析方法去做组间比较。Seurat里有一个非常方便的AggregateExpression函数，可以直接做这个汇总：

bulk <- AggregateExpression(obj, group.by = c("sample_id", "cell_type"), return.seurat = TRUE) # 然后提取表达矩阵，用DESeq2或edgeR做差异分析

如果你的实验设计本身没有生物学重复（比如一个样本里的两个亚群比较，或者公共数据里的一个case和一个control各只有一个样本），那pseudobulk也做不了，这时候用Seurat自带的FindMarkers（默认Wilcoxon秩和检验）做探索性分析可以接受，但在写文章和解读结果时要非常谨慎，因为这种探索性结果除非后续做实验验证，否则说服力有限。我见过不少初学者用FindMarkers跑出了几百个差异基因，画了一大堆火山图，最后发现只要换一种差异检验方法，这些基因就剩不下几个了。所以我给的建议是：先想清楚你的“生物学重复单位”是什么，再决定分析方法，这才是设计层面的问题，不是代码层面的问题。

5.2 富集分析的背景基因问题：为什么显著不等于有意义

差异基因列表出来后，绝大多数人都会去做GO/KEGG富集分析。这个步骤看着简单，但里面的一个核心陷阱是背景基因的选择。富集分析的逻辑是：你拿到一个差异基因列表，看这个列表里的基因在哪些通路或GO条目中显著富集。这里所谓的“显著”，是在和背景基因比较后得到的——差异基因列表里落在某个通路的基因数，相对于随机抽取同样数量的背景基因时落在这个通路的期望值，是否显著偏高。问题就出在“背景基因”应该是你实际检测到的所有基因，而不是整个基因组的全部基因。

如果我的差异分析是基于pseudobulk得到的，那么背景基因应该取DESeq2结果里所有做过检验的基因；如果差异分析是用Seurat的FindMarkers，背景基因默认就是矩阵里的全部基因。但很多人不检查这些细节，直接把差异基因列表丢进DAVID或Metascape，用默认的“全基因组”作为背景，这样富集到的结果往往会偏向那些基因数特别多的通路，而且会漏掉真正的关键通路。我自己踩过这个坑：有一次用全基因组背景跑出来一堆“核糖体”“氧化磷酸化”这种通用通路，把真正关注的炎症和T细胞激活通路完全淹没了，换成表达矩阵的背景基因后，整个富集结果焕然一新。

另外还有一个问题：富集分析的结果是“相关性”而不是“因果性”。富集到的通路只能说明你的差异基因里有很多成员参与了某个已知通路，但这个通路是否真的在你的细胞状态转变中起驱动作用，需要后续的实验验证。我建议除了看p值排序，一定要看富集分数（比如enrichment score）和基因在这个通路中的具体成员。如果一个通路富集非常显著，但实际落在里面的差异基因就两三个，这个富集结果的意义也要打个问号。

5.3 多个样本整合：批次效应的识别与Harmony/CCA的选择时机

只要你的数据是从多个样本、多个批次、甚至多个平台来的，批次效应就是你绕不开的问题。批次效应指的是由于实验处理时间、建库批次、测序深度等因素造成的技术性差异，它的特点是：这种差异和生物学差异混在一起，如果不处理，聚类时细胞可能先按照“样本来源”而不是“细胞类型”聚在一起。判断有没有批次效应的直观方法是用UMAP图按样本着色：如果同一细胞类型的细胞聚在一起但每个样本有自己的小分支，问题不大；如果同一细胞类型的细胞在不同样本里被分到了完全不同的聚类簇，就需要做整合了。

目前处理批次效应最主流的方法是Harmony和Seurat的CCA整合。我个人的选择逻辑是这样的：如果你的核心任务是比较几个样本之间的细胞类型组成，而且样本之间细胞类型本身就比较接近，那么用Harmony更快、更稳，因为它直接在低维空间（PCA）上做迭代校正，不改变原始表达矩阵，跑起来也很快。如果你的重点是寻找跨样本保守的细胞状态、希望在整合后做亚群比较，而且样本间的细胞类型构成差异非常大（比如肿瘤样本里某些样本特有某类细胞），那么Seurat的CCA整合通常会更合适，因为它通过识别跨样本共享的“锚点”来对齐数据，能更好地保留样本特异性的真实生物学。

需要特别警惕的是过度整合。有些人一上来不分青红皂白就把所有样本全部整合，导致本来真实存在的样本特异细胞群或状态被磨平了。我的习惯是先不整合，看原始聚类里批次效应到底有多严重；再根据前面说的问题类型选择整合方法；整合完后一定要验证关键marker是否仍然清晰分离（比如CD8+T细胞的marker不能因为整合糊掉）。一句话总结：整合是为了去技术垃圾，但绝不能把生物学信号一起扔掉。

6. 进阶分析方向和一份可以直接抄的排错清单

6.1 拟时序、细胞通讯、CNV推断：什么数据条件才配做

基础聚类注释做完之后，很多人的目标是上“进阶分析三项”：拟时序分析（轨迹推断）、细胞通讯分析和拷贝数变异（CNV）推断。这三项分析每个都很强大，但每个都有严格的使用前提，不适合的数据硬跑只会得到一堆没办法解释的图。

拟时序分析解决的核心问题是“细胞状态如何从A过渡到B”。它最适用于胚胎发育、干细胞分化、免疫细胞活化这些存在连续状态转换的场景。如果我们做的是肿瘤组织，细胞类型之间的分化关系可能并不明确，某些“过渡态”在数据里根本没有采样到，轨迹算法就可能强行构建出一条虚假的轨迹。我在实际分析中遇到过这种情况：由于某个中间态细胞比例太低、捕获数量不足，Monocle3自动推断出的轨迹直接跨过了一个真实的细胞群，把两个不相关的状态连成一条线性通路。所以在跑拟时序之前，一定要先问自己：我的数据里有没有足够的过渡态细胞？我关心的起点和终点是否明确？如果没有，宁可不做。

细胞通讯分析（CellChat、CellPhoneDB等）基于配体-受体相互作用来推断细胞之间的通讯网络。这个分析的问题在于，它本质上用的是“这个细胞表达了配体，另一个细胞表达了受体”的共表达信息，不是真的看到了蛋白水平的互作。因此它的结果更偏向于“预测”而不是“验证”。即使如此，它在我们筛选关键互作、解释细胞相互作用时依然很有价值，尤其是配合空间转录组数据的时候。跑细胞通讯分析时，一个容易被忽略的细节是：低表达基因的配体-受体对会因为dropout而获得较低的检测灵敏度，所以建议至少在细胞亚群层面做pseudobulk聚合后再计算，而不是直接用单细胞矩阵。

CNV推断（inferCNV等）则是利用基因在染色体上的位置信息，根据表达量推断拷贝数变异。它在肿瘤研究里被广泛用来区分恶性肿瘤细胞和正常细胞。但要注意，inferCNV的结果是相对推断，不能当成真实的CNV证据，而且它对基因表达水平、参考细胞的选择非常敏感。我见过有人直接用inferCNV的结果写进文章当作拷贝数变异的结论，这是很不严谨的，如果没有WES或者WGS等直接基因组证据，建议在文内明确说明这只是“基于转录组的推断”。

6.2 三个高频故障的完整排错思路

实际分析中90%的时间都花在排错上。这里分享三个我反复遇到、八成的人也会遇到的问题，每个问题给一个完整的排查链路，而不是直接贴修复代码。

第一个问题是“细胞数远少于预期”。如果cellranger报告里估计的细胞数只有你预期的一半，先不要调--expect-cells参数硬撑。先看web_summary.html里的barcode rank plot，判断是“没有明显拐点”（说明背景RNA污染严重）还是“拐点位置偏低”（说明细胞悬液浓度太低或细胞活力不足）。前者需要在过滤环节更严格地设置阈值，后者则要考虑上游实验流程。如果证实是细胞悬液质量问题，最好的做法是重跑实验，不要在生信上硬救，因为后续的聚类注释都会被低质量信号干扰。

第二个问题是“UMAP图上所有细胞叠成一团，没有任何分离”。这个现象通常有三个原因：一是主成分数取太少，信息没进到聚类里；二是特征基因选择不当，比如大部分特征基因都是线粒体或核糖体基因；三是归一化方法选择不匹配数据。我的排查顺序是：先看VariableFeatures里都是什么基因，如果核糖体基因占了半壁江山，用SCTransform或手动排除^RP[SL]基因；再看PC的显著性；最后检查是不是所有细胞都来自同一个样本且样本本身确实没有明显异质性。如果是最后一种情况，那“分不开”反而是符合预期的。

第三个问题是“细胞分群图看起来很好，但注释困难”。这个时候先别急着怀疑自己的背景知识，很可能是聚类分辨率选得过高，同一类型的细胞被硬生生切成了好几个簇。我的处理办法是把分辨率调低一档，看这几个簇是否合并，合并后marker表达是否更清晰。如果分辨率已经很低了，注释仍然困难，你就需要扩大marker列表，用SingleR或CellTypist这类自动化注释工具做交叉验证，但记住自动注释结果只是参考，不能直接照单全收。

6.3 一份实用的分析决策清单

基于上面的经验，我整理了一份单细胞分析从上到下的决策清单，每次跑新数据的时候都会过一遍，省掉了很多回头返工的时间。

• 拿到数据先确认物种、参考基因组版本、10x建库类型，这三个信息错了后边全白做。
• cellranger count跑完后，先看web_summary.html的barcode rank plot和测序饱和度，再决定是否继续后期分析。
• 质控阈值不要照抄，先画图观察分布，线粒体比例上限建议在5%到20%之间根据样本实际情况调整。
• 双细胞检测建议放在质控后、归一化前，用DoubletFinder或Scrublet至少做一个。
• 归一化方法优先用SCTransform做复杂样本，LogNormalize适合快速探索和小样本。
• 降维时PC数根据JackStraw或ElbowPlot判断，不要无限多加。
• 聚类分辨率从低到高多调几档，选择能分离出已知生物学细胞类型的最低分辨率。
• 细胞注释采用“谱系分大类、再分小类”的两步法，marker至少联合3个以上判断。
• 组间差异分析优先pseudobulk策略，避免伪重复。
• 做富集分析前确认背景基因与检测基因一致，不要直接用全基因组做背景。
• 多样本先画UMAP看批次效应，根据“共享细胞状态”还是“特异细胞状态”选择整合方法。
• 拟时序、细胞通讯、CNV推断都有严格的适用前提，不满足条件宁可只展示基础聚类结果。

这份清单看起来条目很多，但实际操作起来都很短，大部分就是看一眼图表就确定的。真正花时间的还是在“理解数据”上。

我自己的习惯是每个项目固定创建一个分析SOP文件，把每一步的阈值、参数和判断依据记录下来。单细胞分析里变量太多，你今天改一个阈值，明天换一个归一化方法，结果可能就完全不一样。没有记录的话，三个月后回看你只会记得“大概处理过”，但完全不记得当时为什么这么选。把决策过程记录成文件，既方便自己复盘，也给合作的人和审稿人留一份可追溯的分析脉络。毕竟单细胞转录组分析最讲究的，不是你会跑多少个流程，而是每一步做决策时，你心里知道自己在做什么、为什么这么做。