ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

单细胞测序t-SNE聚类与marker基因筛选实战指南

单细胞测序t-SNE聚类与marker基因筛选实战指南 做单细胞测序分析的朋友走到t-SNE聚类分析这一步多半已经熬过了前面最枯燥的流程从FASTQ比对、定量到创建Seurat对象再到质控、标准化、高变基因筛选。你以为接下来就是跑两行代码出张图的事我第一次做的时候也是这么想的结果图是出来了但分群乱得像一锅粥marker基因怎么找都对不上已知的细胞类型。后来一步步排查才发现问题出在降维参数的选取和聚类分辨率上。这篇就把我这几年跑单细胞测序流程五的经验整理出来从t-SNE聚类的底层逻辑讲到marker基因筛选的完整链路顺便把踩过的坑也一并交代清楚。1. 走到t-SNE这一步之前你的数据应该长什么样1.1 前置流程做没做扎实直接决定聚类质量很多人上来就跳过中间步骤把原始count矩阵直接扔进RunTSNE()这基本等于拿毛坯房当精装房住。单细胞测序流程走到t-SNE这一步至少需要完成下面这条链路CreateSeuratObject→NormalizeData→FindVariableFeatures→ScaleData→RunPCA→FindNeighbors→FindClusters→RunTSNE。每条命令都在做一件不可替代的事缺一步后面的结果都会变形。我把这些步骤的作用和常见陷阱整理成了下面这张表方便你在跑之前对照检查自己的数据步骤核心作用常见错误NormalizeData消除测序深度差异把每个细胞的count归一化到可比尺度不归一化直接聚类高表达基因会主导分群FindVariableFeatures挑选在细胞间变异最大的2000个基因用于后续降维选了太多或太少基因导致降维捕捉的是技术噪声ScaleData基因表达量标准化到均值为0、方差为1避免高表达基因权重过大忘记设置features只缩放高变基因内存暴涨RunPCA将高维表达矩阵压缩到几十个主成分去掉冗余和噪声PCA维度数没有根据拐点图或方差贡献率合理选择FindNeighbors基于PCA结果构建KNN图和SNN图定义细胞之间的邻居关系使用的PCA维度数和后续不一致导致图形结构错位FindClusters用Louvain/Leiden算法在图结构上找社区得到初步细胞群resolution设置不会调分群过粗或过碎RunTSNE把高维结构投影到二维供人眼视觉检查整体分群格局没设seed导致每次跑图结果不一样或perplexity不合适记住一点t-SNE聚类分析中的聚类结果实际上在FindClusters()这一步就已经确定了t-SNE只是把已经分好的群画出来给你看。很多人在这一步犯迷糊后面我专门用一节讲清楚这个关系。1.2 质控不合格的数据t-SNE会给你加倍奉还如果你的数据里有大量低质量细胞——比如线粒体基因占比超过20%的濒死细胞、或者双细胞doublet混在其中——t-SNE并不会好心地把它们单独挑出来变成一群而是会硬塞进某些群或者形成一团无法解释的过渡态细胞。最典型的表现是聚类热图上几个群之间的marker表达差异模模糊糊或者某个cluster里同时出现两种完全不同的细胞类型的marker。所以我强烈建议在跑进t-SNE之前先用subset()把不符合条件的细胞剔除干净。我的常用阈值是nFeature_RNA 200排除空液滴和碎片、nFeature_RNA 6000排除双细胞或高复杂度异常值具体上限根据建库方式调整、percent.mt 20小鼠和人的经验值不太一样脑组织可以适当放宽到25%。做这一步时可以用VlnPlot()先看一下分布再定阈值不要无脑套参数。质控这关省下的时间会在后面找marker时十倍地还给你。2. 先分清一件事t-SNE是可视化方法不是聚类算法本身2.1 为什么说t-SNE聚类这个叫法不太严谨你会在各种教程里看到t-SNE聚类分析这种说法包括很多高分文章的methods部分也这么写但严格来讲这个叫法是有问题的。聚类clustering指的是把细胞划分到不同群里的计算过程这个任务在Seurat标准流程里由FindClusters()完成它用的是基于图的社区发现算法Louvain或Leiden跟t-SNE没有关系。t-SNEt-distributed Stochastic Neighbor Embedding做的事情是把高维空间中细胞与细胞之间的相似性映射到二维或三维坐标上让相似的细胞在图上靠得近、不相似的细胞离得远本质是一个非线性降维可视化工具。它在单细胞流程中的作用是展示不是划分。那为什么大家习惯叫t-SNE聚类因为实际项目里这套动作是连在一起做的聚类算法先算好cluster标签再用t-SNE画图观察整体格局整个过程对使用者来说是一气呵成的所以口语上经常合并叫。这个可以理解但你在理解结果、排查问题的时候头脑里一定要有这根弦——如果分群不清楚问题可能出在FindClusters()的输入图结构上而如果只是图不好看、但群标签是合理的那问题可能出在t-SNE的投影参数上。这两类问题的修法完全不同。2.2 t-SNE和UMAP到底选哪个现在单细胞领域早就不是t-SNE一枝独秀了UMAP在绝大多数新文章里反而更常见。两者的核心区别在于对全局结构的保留程度t-SNE倾向于只保留局部结构距离远近在图上会失真群与群之间的间隙大小没有实际意义而UMAP在保留局部结构的同时还会尽量维持群与群之间的相对距离全局关系更可靠。我用一个生活化的类比来帮助你理解t-SNE更像一张地铁线路图站与站之间的相对顺序是准的但实际距离全被扭曲了UMAP更像一张普通地图整体比例尺更接近真实距离。具体到实际项目的选择我给的建议是如果为了审稿人、合作方或老板能一眼看清分离效果用UMAP做主图t-SNE作为补充图展示。现在大多数期刊对UMAP的接受度更高而t-SNE在部分审稿人眼里已经有点老派了。如果要做精细亚群的视觉判断比如同一个T细胞亚群内不同状态的过渡t-SNE在部分数据集上会把过渡态压缩得更明显看起来更干净UMAP有时会把这类过渡细胞拉成一条连续的桥反而更真实但不够美观。两个都跑一遍不丢人。我在实际项目中通常把两者都跑出来发现结论不一致时优先相信UMAP的全局结构再去检查t-SNE是否因为perplexity设置不当导致局部失真。另外要强调一点无论选哪个真正决定分群的是聚类步骤不是降维结果。所以每次调整完FindClusters()的分辨率后记得重新画一次降维图不要只在旧图上加标签。3. 实操Seurat从PCA到t-SNE聚类的完整参数调优3.1 PCA维度的选择决定下游聚类的输入质量RunPCA()跑完之后紧接的FindNeighbors()需要指定dims参数也就是使用多少个主成分PC。这个参数对聚类结果的影响非常大却常常被人忽视。选择PCA维度的核心逻辑是保留真实生物学变异去掉噪声维度。如果PC选太少会丢失微弱但真实的亚群信号选太多等于把噪声也喂给了聚类算法分群边界会被糊掉。实际操作中我会用两种方式交叉验证第一种看ElbowPlot()的拐点。这个图会展示每个PC解释的方差百分比通常前面几个PC方差占比很高之后迅速下降并趋于平缓你的目标是在下降变缓的位置附近取一个值。但我必须提醒你这个拐点在真实数据里往往没有那么锐利肉眼判断常会有5~10个维度的浮动空间。第二种直接看多个dims下的聚类结果对比。我习惯把dims从10、15、20、30依次跑一遍用DimPlot()看一下分群是否稳定。如果某个dims值下出现一群细胞从聚在一起变成散成一团说明这个维度下噪声已经开始干扰聚类需要退回更保守的值。这个方法的缺点是耗时因为每次都要重新算FindNeighbors和FindClusters但换来的是分群可靠性非常值。有一个经验值供参考对于10X平台约5000~10000个细胞、2000个高变基因的数据集dims取20~30是比较常见的区间。如果你的数据集更复杂比如有多个样本合并、或来自肿瘤微环境这样高度异质的样本dims可以适当取高一点。但记住堆太多PC不会让结果更细只会让聚类失去稳健性。3.2 FindClusters的分辨率参数掌握好粗与细的分寸FindClusters()里的resolution参数决定了聚类算法会把细胞切得多细。我把这个参数理解为一把分辨率旋钮调大每个群的内部差异会被进一步挖掘产生更多亚群调小亚群会被合并成大群更容易看出大的细胞类别。这个参数没有绝对正确值它取决于你的生物学问题。如果你想回答的问题是这个组织里都有哪些免疫细胞类型那resolution设0.1~0.5就足够了如果你想找的是CD8阳性T细胞里是否存在一个衰竭亚群那至少需要0.8~1.5甚至更高。我常用的探查法是这样从resolution 0.1开始逐步往上加0.2、0.5、0.8、1.0、1.5每跑一次都配合FindAllMarkers()看新增的cluster有没有可靠的marker支持。如果增加resolution之后出现的新群无法用已知marker解释那多半是过拟合了噪声应该回到更低的resolution。如果继续增加分辨率后marker仍然特异说明这个亚群可能是真实存在的。这里有一个我在实际项目中反复体会到的点cell type细胞类型和cell state细胞状态是两个不同层次的东西。低分辨率下聚出来的群往往对应真正的细胞类型——比如T细胞、B细胞、巨噬细胞高分辨率下细分出来的子群有时候不是新的细胞类型而是同一类细胞的不同状态比如静息态、活化态、耗竭态。在注释结果时一定要区分这两种情况否则会把同一个细胞类型拆成多个新细胞类型给下游分析埋雷。3.3 RunTSNE的关键参数与稳定性问题在聚类完成后就该把结果投影到t-SNE图上了。Seurat里的核心命令是sce - RunTSNE(sce, dims 1:20, perplexity 30, seed.use 42) DimPlot(sce, reduction tsne, label TRUE, pt.size 0.5)dims参数必须和你前面FindNeighbors()里用的保持一致这是很容易犯的低级错误但后果很严重——如果前后不一致图上的分群格局会和你实际定义的cluster标签完全不同步。perplexity参数困惑度控制的是t-SNE在计算邻居关系时看多远可以把它理解为每个点在计算相似度时要考虑附近多少个邻居。这个值的推荐范围是5~50Seurat默认30。数据量越大perplexity通常可以适当调大一点但如果你发现t-SNE图上出现大量细胞被挤成一根根细丝或者群内出现明显的小空洞很可能是perplexity过低。反过来如果所有细胞糊成一团分不开而UMAP上明明分得很清楚可以试试把perplexity调低到10左右。还有一个特别实用的细节t-SNE算法带有随机性每次运行结果会有差异。所以你一定要设置seed.use参数并且在论文或报告中注明保证结果的可重复性。我见过不止一次同一位研究者隔天跑同一个脚本发现图变了最后才发现是忘了固定种子。3.4 分群结果不理想先别急着调参数这个我放到这一节稍微提一下详细排查链路在后面的实战踩坑笔记里再展开。当你看到分群图一团乱、群与群之间没有清晰边界时第一反应不应该是无限调perplexity或resolution。你首先要检查的是源头——数据本身干不干净。我排查的标准顺序是打开FeaturePlot看几个已知的重要marker基因看它们的表达模式是否具有空间连续性。如果某个已知marker在所有细胞里都是均匀弱表达先警惕数据质量。查看percent.mt、nCount_RNA、nFeature_RNA在t-SNE图上的分布。如果发现某一群细胞恰好对应于高线粒体比例或极低的基因数那这群细胞很可能不是真正的细胞亚群而是濒死细胞聚成一团。检查是否存在明显的样本来源标签结构——如果两个样本混合后t-SNE图上细胞严格按样本来源分成左右两团而不是按生物学类型混合分布那大概率是批次效应需要用数据整合方法处理这属于另一套流程。4. 找marker基因命令只是一步看懂结果才是关键4.1 FindAllMarkers背后的统计逻辑分群完成后下一个重头戏就是寻找每个cluster的marker基因。Seurat里的核心命令是markers - FindAllMarkers(sce, only.pos TRUE, min.pct 0.25, logfc.threshold 0.25)拆开来看FindAllMarkers实际上是对每个cluster都把它和所有其他cluster的细胞做一次差异表达检验默认是Wilcoxon秩和检验找出在该cluster中显著高表达的基因。only.pos TRUE意味着只返回上调基因对于找marker来说通常够用了因为我们要找的是标记这个群的特征基因而不是这个群里下调的基因。这里三个参数的含义很容易被忽略我逐个说一下min.pct基因必须在多大比例的细胞中被检测到。默认0.1我常用0.25。如果设得太低你会得到一堆只在极少数细胞里表达、却在统计上显著的基因它们往往是噪声。logfc.threshold只保留log2倍数变化大于该值的基因。如果不设这个阈值Wilcoxon检验可能会把表达量差0.1倍的基因也判定为显著但生物学家看了只想打人。only.pos只返回该群中高表达的基因。如果设为FALSE你会同时拿到一堆下调基因对找marker来说往往会造成干扰。跑完之后我通常会给结果加上一行筛选逻辑top_markers - markers %% filter(p_val_adj 0.05, avg_log2FC 1) %% group_by(cluster) %% slice_max(n 10, order_by avg_log2FC)p_val_adj是校正后的p值因为你要对成千上万个基因做几千次检验不做多重假设校正假阳性会非常严重。avg_log2FC 1意味着该基因在目标群中的平均表达量是其他群的2倍以上这是一个比较保守的筛选线如果你发现筛选完没有足够的marker可以放宽到0.58相当于1.5倍。4.2 怎么判断一个marker基因是否真的能用很多初学者拿到FindAllMarkers输出的表格之后第一反应是把每个cluster的top基因拿去做GO富集结果发现一堆看不太懂的通路名。我个人建议在富集分析之前先做一步更朴素但更关键的事——逐个用可视化确认marker的特异性。我判断一个marker基因是否可用的标准有三个第一表达特异性。这个基因在这个cluster里高表达在其他cluster里低表达或几乎不表达。用FeaturePlot()看它在t-SNE/UMAP图上的分布应该呈现出只在某个区域内亮起来的模式而不是全局扩散。第二阳性细胞比例。看结果中的pct.1该群中表达此基因的细胞占比和pct.2其他群中表达此基因的细胞占比。我的习惯是要求pct.1明显高于pct.2比如70%对10%。如果pct.1只有30%而pct.2有25%即使p值显著这个基因也只能算弱marker不适合作为群身份的决定性证据。第三生物学常识。marker基因最终要服务于细胞身份注释因此它一定要对应已知的细胞生物学功能或文献报道。比如你用Cd3d标记T细胞、用Cd79a标记B细胞、用Lyz2标记巨噬细胞这些是教科书级别的marker。如果一个全新基因在统计上高度显著但你在文献里查不到它的谱系关联那它更适合描述为该亚群的潜在标志物而不是直接作为注释依据。4.3 热图和小提琴图marker验证的标准动作选出一批候选marker后我会按照以下流程做最终确认。首先生成每个cluster的top marker热图top10 - markers %% group_by(cluster) %% top_n(n 10, wt avg_log2FC) DoHeatmap(sce, features top10$gene) NoLegend()热图应该呈现明显的方块化结构每个cluster对应一块独立的高表达区域。如果热图上不同cluster之间你中有我、我中有你界限不清说明你选的marker不够特异或者前面的聚类本身就有问题。然后再对每个关键的marker做VlnPlot()VlnPlot(sce, features c(Cd3d, Cd79a, Lyz2), pt.size 0)这一步看的是marker在每个群中的表达量分布。理想情况下某个群的表达量分布应该整体上移而不是只有少数离群点拉高均值。用pt.size 0可以让小提琴图更干净不会被几万个点糊满。这一步做完你手里的marker列表才真正可以用来做细胞类型注释。我见过太多人只跑了一句FindAllMarkers就拿着结果开始注释结果把一个cluster因为Gapdh高表达注释成代谢活跃细胞这种错误完全可以通过热图和小提琴图避免。5. 单细胞聚类实战中的高频翻车场景与排查链路5.1 翻车现场一t-SNE图上所有细胞糊成一团先交代一下背景我处理过一批某疾病模型的组织样本建库质量不算差QC的时候各项指标也都正常但跑完t-SNE之后整张图上所有细胞挤成一大坨几乎看不到任何分群结构。一开始我还以为是perplexity设太低了换来换去都不见好转后来才意识到问题出在PCA维度选择上——我在上一步保守地只取了前10个PC结果把真正用于区分细胞类型的信号给丢掉了许多。排查链路是先看多个dims20、30下的聚类结果变化再看高变基因分析是否合理最后回到标准化那一步检查是不是用了错误的SCT和LogNormalize混用。从这个坑里得到的教训是如果你在UMAP和t-SNE上都看不到任何分群优先怀疑输入信号不足而不是降维参数的问题。我现在的做法是先用一个宽松的dims比如30快速跑一遍看能否推出大致格局如果不行再回过头检查上游的高变基因数和样本异质性。5.2 翻车现场二批次效应让分群按样本而不是按细胞类型走单细胞项目很少有只跑一个样本的一旦合并多个样本批次效应就是绕不开的坎。最典型的症状是t-SNE图上细胞不按T细胞聚在一起、B细胞聚在一起分布而是样本A的所有细胞都在这边、样本B的所有细胞都在那边。如果你在FeaturePlot里用样本ID给细胞上色把样本信息存在metadata里DimPlot(sce, group.by sample_id)发现颜色严格分区几乎可以断定是批次效应。这时候不要指望调聚类参数能解决应该回到数据整合环节。目前主流的解法是Harmony也可以在Seurat里用SCTransform配合IntegrateData做整合。Harmony的使用逻辑很简单RunHarmony替代RunPCA之后的FindNeighbors步骤让聚类基于去批次后的嵌入坐标进行。整合后t-SNE图上细胞应该按生物学类型混合分布这才是可以继续往下走的信号。有一个细节值得提醒你做完整合后后面所有降维和找marker都要基于整合后的对象进行不要再回退到整合前的PCA结果上。5.3 翻车现场三某个cluster同时表达两个谱系的marker这个场景在我第一次处理肿瘤样本时遇到过某一个cluster的top marker里既有T细胞标志物Cd3d又有髓系标志物Lyz2。当时我第一反应是聚类分辨率不够把两种不同细胞硬分到一个群里了于是把resolution一路上调结果分出来的子群依然同时表达两个谱系的marker。之后我意识到问题不在聚类而是数据分析前就混入了双细胞doublet——一个液滴里同时包裹了一个T细胞和一个巨噬细胞测到的基因表达谱天然就是两者叠加。双细胞在单细胞数据里造成的假象非常隐蔽不加处理的情况下它们会形成一个不伦不类的中间群或者被硬挤进已有的群导致marker信号被稀释。解决办法是提前使用DoubletFinder或scDblFinder预测并过滤双细胞。跑完降维聚类后再看一眼每个cluster里的双细胞预测比例——如果某个群的双细胞比例远高于平均水平那这个群的marker需要格外谨慎地解读最好把它排除后再重新聚类。5.4 翻车现场四细胞周期效应把细胞按正在分裂分组这类问题在增殖活跃的组织比如发育中的大脑、肿瘤中特别常见。它的特征是你本来想分细胞类型结果t-SNE图上分成两三个巨大的阴阳格局——一大群细胞表达高水平的G2/M期marker如Mki67、Top2a另一大群表达G1/S期marker如Ube2c第三种是G1期细胞。这时候细胞类型的信息反而被淹没了。一个快速验证方法用CellCycleScoring()给每个细胞打周期分数再用DimPlot按周期分数上色。如果分群方向刚好沿着周期的梯度分布那就是周期效应。处理方式有两种思路。一种是如果细胞周期变异与你的研究问题无关直接通过ScaleData回归掉周期分数的影响另一种是如果周期相关基因本身是研究对象比如肿瘤增殖就不要简单回归而是在注释时区分增殖型亚群和静息型亚群。5.5 排查链路总结按顺序来不要跳步把上面几个场景的排查顺序整理成一张决策链方便你在实际项目中照着走先确认QC阈值是否合理低质量细胞有没有干扰聚类检查是否存在明显批次效应必要时做数据整合验证双细胞比例排除中间状态造成的虚假cluster评估细胞周期效应是否主导了分群方向再回到dims、resolution、perplexity这些参数逐个检查是否设置得合理最后才用marker的生物学意义来验证分群的解释力。我一直认为单细胞分析中90%的聚类效果差都不是调参能解决的真正要做的是回到数据源头找原因。这条链路我在每个新项目里都会完整走一遍虽然费时间但能帮你省下后面注释时的无数烦恼。6. marker基因的生物学验证从统计显著到真凭实据6.1 注释细胞类型时用哪些数据库和工具当你拿到每个cluster的top marker之后下一步就是把它们对应到已知的细胞类型上。这一步我用到的工具包括SingleR基于参考转录组数据集自动注释。优点是快速、省力适合对细胞大类做初步判断缺点是参考数据的选择会影响结果而不同来源的参考集对同一群细胞的注释可能不同。我通常把它当助手而不是裁判。CellMarker、PanglaoDB、CellTypist这些数据库中存有文献中验证过的细胞类型marker基因列表。你可以在里面搜索自己的top marker是否有已知的细胞类型关联。人工查阅文献这是绕不开的最终王道。尤其在做新亚群注释时自动化工具的准确度往往不够最后还是要回到原始文献中确认这群细胞的表型特征。我在实际中倾向于先跑一遍SingleR做初步判断再用已知的经典marker做人工验证。比如我预期某个cluster是CD8阳性T细胞就会看Cd8a、Cd3d、Gzmb等经典marker是否在这个cluster中高表达并且要求多个marker同时支持同一结论。如果只有单个marker支持我会把它标记为待验证不会直接拍板。6.2 区分marker基因与差异表达基因很多初学者把FindAllMarkers输出的所有基因都当作marker基因这是一个很需要纠正的误解。差异表达基因DEG和marker基因之间是包含关系所有的marker基因都是差异表达基因但不是所有差异表达基因都能作为marker。marker基因强调的是类别的指示标志它在目标群中高表达在其他所有群中不表达或低表达且这个模式具有稳定性。而DEG只强调不同群之间有统计差异——一个基因如果在A群中表达量是10在B群中是2它已经算差异了但作为marker远不够。所以我在筛选marker时会刻意提高阈值。除了前面提到的p_val_adj 0.05和avg_log2FC 1还会特别关注pct.1和pct.2之间的差距。一个真正的markerpct.1通常要接近0.8或更高而pct.2应该低于0.2。如果两者都在0.5附近我会怀疑这个基因只是普遍高表达不能用来定义细胞身份。6.3 热图验证时的几个细节最后再补几个DoHeatmap使用中的实操细节。第一DoHeatmap默认显示每个身份类最多50个细胞如果你的cluster细胞数量特别多热图看起来会很拥挤。可以通过DoHeatmap(sce, features top10$gene, size 4)调整字体大小或者先对数据subset成每类细胞抽样再画。第二热图上的表达值默认经过ScaleData标准化所以Z值的正负只能说明相对表达高低不能直接当作绝对表达量理解。第三如果在同一张热图里画了几十上百个基因图像会变得过密难读。更推荐的做法是先选20~30个核心marker基因画一张精简热图再对重点关注的具体基因单独用VlnPlot和FeaturePlot做深入可视化和验证。精简热图给审稿人看全局格局单独feature plot给自己确认细节两不误。写在最后一些操作性很强的建议如果你现在正准备跑单细胞测序流程五我的建议是把注意力放在数据质量和参数的因果关系上而不是机械地执行脚本。梳理一下本篇文章中个人最想强调的几点第一t-SNE只是可视化工具真正决定聚类的是FindNeighbors和FindClusters它们各自依赖的PCA维度、resolution参数永远值得你先想清楚其生物学意义再跑第二找marker不只是跑一句FindAllMarkers要结合热图、小提琴图和生物学常识做交叉验证宁可多花半天验证也不要急着出图发表第三数据质量问题双细胞、批次效应、周期效应如果要靠降维聚类参数修复那一定修不好。最后分享一个我在实际项目中体会很深的小技巧每次调参后我都会把当时的t-SNE/UMAP图、聚类resolution、PCA维度Marker基因列表存档在同一份运行笔记里。这样即使一个月后回头复查也能回忆起当时的画面为何长这样以及改参数后结果是如何演变的。单细胞分析的迭代次数非常多没有这份记录你很容易在调参迷宫里彻底迷失方向。祝你聚类顺利marker一找一个准。
返回列表