SCTransform与Harmony参数设置详解:单细胞转录组分析避坑指南
我做了四五年单细胞转录组分析见过太多人在SCTransform和Harmony这两个环节翻车。明明数据质量挺好跑完聚类却乱七八糟最后发现不是数据问题是参数设置埋了雷。这篇文章专门聊SCTransform和Harmony的参数设置把每个参数背后到底在干什么、为什么这么设、设错了会出什么问题一条条掰开讲清楚。如果你正在处理多样本单细胞数据被批次效应搞得头疼或者发现SCTransform跑完后的结果和预想差很多这篇文章就是给你写的。新手可以把里面的参数组合直接抄走已经跑过流程的老人也能对照检查一下自己有没有踩进某些隐蔽的坑里。1. 先搞清楚SCTransform到底帮你做了什么很多教程会告诉你“用SCTransform替代LogNormalize”但没讲清楚它和传统标准化到底差在哪里。不理解原理就去调参数基本等于盲人摸象。这一节先把SCTransform的核心逻辑讲明白后面聊参数才有意义。1.1 为什么大家开始放弃log-normalize传统Seurat流程里NormalizeData用的是log1p(CPM/10000)这类操作本质是先把UMI数转成相对丰度再取对数压缩动态范围。这个方案的问题在于它假设所有基因的方差只和表达量均值有关但单细胞数据里还存在技术噪声比如测序深度差异和生物学噪声比如细胞周期、线粒体应激这些因素会让高表达基因的方差被系统性放大。你去做ScaleData的时候通常会顺手做一步z-score标准化但z-score的均值和方差估计都会被技术因素污染。结果就是你挑出来的高变基因可能不是因为生物学差异大而是因为某些细胞测序深度特别深。传统流程里有人用“回归掉UMI数和线粒体比例”来补救但本质上还是在用线性回归处理非线性关系效果有限。SCTransform换了个思路它直接用负二项回归建模UMI计数把测序深度、线粒体比例这些技术因素作为协变量放进模型里然后取残差作为校正后的表达量。这一步同时完成了标准化、方差稳定化和技术因素校正比“NormalizeData ScaleData regress out”这条老路干净得多。1.2 SCTransform的“回归”和“标准化”是一次完成的很多人不知道SCTransform跑完以后SCT这个assay里其实包含了三类数据counts和原始RNA assay的counts一样但只保留了参与建模的基因。data对数化后的表达量log1p用于差异表达分析和可视化。scale.data模型残差标准化后的结果主要用于PCA和聚类。关键理解是后续PCA用的是scale.data也就是“去掉技术因素后的残差”。所以如果你在SCTransform之后又画了FeaturePlot去看某个基因的表达量默认用的其实是data里的log1p值而不是残差。这会导致一个经典困惑为什么基因表达图看起来没问题但聚类结果里这个基因却没有贡献因为聚类看的是scale.dataFeaturePlot看的是data两者不是一回事。这个信息在处理“某个marker基因在umap上明明表达但找差异基因时居然不显著”这类问题时特别有用。真实场景中我见过有人因为这个差异去反复调整FindMarkers参数最后才发现是自己对数据的理解偏了。1.3 SCTransform核心参数逐个拆解先给出一份我常用的SCTransform参数清单后面逐个讲为什么这么设。obj - SCTransform( obj, assay RNA, new.assay.name SCT, vars.to.regress c(percent.mt), variable.features.n 3000, return.only.var.genes TRUE, conserve.memory TRUE, seed.use 42, verbose FALSE )assay指定从哪个assay取原始counts默认是RNA。如果你之前跑过cellranger、STARsolo等流程需要确认导入Seurat后counts确实在RNA assay里否则会报错或者拿到空值。new.assay.name默认叫SCT。如果你想保留多套标准化结果做对比可以改成不同的名字。这是我比较推荐的做法尤其是想比较SCT和LogNormalize结果差异的时候。vars.to.regress需要回归掉的协变量。最常见的做法是回归percent.mt线粒体比例因为这个指标通常代表细胞应激或破损程度。但这里有一个很大的坑并不是说“能回归的都回归掉”就是好事。细胞周期分数、样本ID这类变量别随便放进去理由稍后单独细说。variable.features.n要保留的高变基因数量默认3000。如果样本类型复杂度高比如肿瘤组织包含多种微环境细胞我会调到4000-5000。如果数据比较简单纯细胞系、纯T细胞亚群2000-3000就够了。基因太少会丢失稀有亚群的信号太多则把噪声也当成高变特征PCA前几个主成分可能被无关基因主导。return.only.var.genes默认为TRUE意思是SCTransform后SCT assay里的scale.data只保留高变基因。这能显著省内存也不影响后续PCA因为PCA默认就是在高变基因上做的。但如果你之后想用某个不在高变基因列表里的基因去做热图或者打分就需要小心了——这时候你可能需要临时把return.only.var.genes设成FALSE或者用RNA assay的数据去提取表达量。conserve.memory强烈建议设为TRUE。SCTransform对内存的消耗非常大尤其是在10万细胞以上的数据集里。这个参数开启后会让SCTransform分块处理内存占用能下降不少。缺点是速度会慢一点但比起内存溢出来说完全值得。seed.use随机数种子保证结果可复现。写论文或者做项目时这一项必须设置否则每次运行结果可能在局部细节上有差异。verbose设为FALSE可以省去大量中间日志不然控制台刷屏刷得很难受而且SCTransform本身日志量巨大容易刷掉前面的报错信息。1.4 vars.to.regress最容易被误用的坑vars.to.regress这个参数坑特别多我单独拎出来讲。首先说percent.mt。常规认知里线粒体比例高代表细胞状态差或者细胞破碎所以把它回归掉是合理的。但我遇到过一个案例处理某个带有线粒体基因高表达特征的细胞类型比如某些代谢活跃的细胞时强行回归percent.mt后这个细胞类型的特征基因信号被大幅削弱聚类时这个群体直接消失了。所以一个稳健的做法是先用默认设置回归percent.mt跑一遍看看主要细胞类型能不能分出来。如果某个你预期的细胞类型消失了再尝试不回归percent.mt对比两个结果。千万不要机械执行“别人都回归percent.mt所以我也回归”。再说细胞周期。很多人喜欢用CellCycleScoring算出S.score和G2M.score然后放进vars.to.regress。这个操作对增殖能力差异很大的组织比如肿瘤有时是必要的因为细胞周期信号太强会把其他生物学信号盖掉。但如果你研究的就是干/祖细胞的增殖状态或者你的目标细胞类型恰好和细胞周期高度相关那这一回归就把你最关心的信号给干掉了。我的经验是先别急着回归细胞周期。跑完聚类看结果如果发现聚类主要被细胞周期相关基因分开而不是按细胞类型分开再回归不迟。这个“先跑后决策”的思路在单细胞分析里适用性极广。另外千万不要把样本ID或者批次信息放进vars.to.regress。SCTransform是每个细胞独立建模的不是用来做批次整合的。如果把样本ID放进去SCTransform会尝试把样本间的表达差异全部当作技术因素去掉这可能把真实的生物学差异也一并抹掉。批次整合请交给Harmony、CCA或者scVI去处理各司其职别混用。2. Harmony不是无脑跑的参数背后的统计学含义Harmony是目前整合多样本/多批次数据最常用的工具之一。它的效果确实好但也正因为效果好很多人不管三七二十一直接拿默认参数一顿跑结果跑出来的整合结果是“看起来都混在一起了”但生物学差异也被磨平了。要避免这个问题必须理解Harmony在数学上做了什么。2.1 Harmony在单细胞流程里的位置Harmony不是替代SCTransform的它是在PCA之后运行的。标准流程是SCTransform完成标准化和技术因素校正。RunPCA在SCT assay的scale.data上跑PCA得到降维后的细胞嵌入。RunHarmony以PCA的嵌入结果作为输入根据group.by.vars指定的批次变量做整合。后续的FindNeighbors、FindClusters、RunUMAP全部基于Harmony校正后的嵌入。这里有个很多人搞不清楚的点Harmony到底改了什么东西它不是改表达矩阵而是在PCA降维后的低维空间里对细胞的坐标做迭代校正。它先在每个聚类簇里估计“这个簇里不同批次的代表细胞”然后通过混合多个批次的代表细胞来重新放置每个细胞的位置最终让来自不同批次的同类细胞在低维空间里重叠。这意味着Harmony不会改变你的基因表达值只改变细胞在低维空间里的相对位置。如果你后续要做差异表达分析用的还是SCTransform校正后的数据而不是Harmony的输出。搞清楚这一点就不会纠结“Harmony跑完后表达矩阵在哪里”这种问题了。2.2 group.by.vars怎么选别把生物学变量当批次group.by.vars是一个字符串向量表示“你认为哪些列代表了技术批次”。最常见的选择是sample ID或者library batch。这个参数是Harmony整合效果的核心。设对了整合效果好设错了后果往往很严重。举几个反面例子第一把生物学分组比如疾病组vs对照组放进group.by.vars。这样Harmony会强行把疾病组和对照组的同类细胞往一起拉最终导致疾病和对照之间的差异被磨平后续找差异基因基本找不到或者找到的都是假阳性。这个错误出现的频率高得惊人尤其是有代码习惯的人复制粘贴上一轮脚本时最容易犯。第二把太多变量塞进group.by.vars比如同时放sample ID、library batch、date、plate号。Harmony会尝试在所有维度上整合结果就是每个维度都被弱化而且计算量剧增。正确的做法是只放“你确信代表技术批次”的1到2个变量。如果你不确定某个变量是不是技术批次可以在整合前先画个不需要整合的UMAP看细胞是否按这个变量分堆。如果确实分堆了再考虑放进去。第三样本数量太少时用Harmony要谨慎。我见过只有两个样本的数据集其中一个样本细胞数特别少Harmony跑完以后基本等于把所有细胞硬拉到一起聚类结果完全丧失生物学意义。这种情况下先用BBKNN、merge后直接聚类等方式做预实验或考虑直接放弃整合改用严格的样本间差异分析策略。2.3 theta、lambda、sigma这些参数影响什么Harmony的默认参数在多数数据集上表现都还行但你如果想追求更好的效果或者遇到了某些特殊场景就需要调整下面几个核心参数。先看默认参数是什么obj - RunHarmony( obj, group.by.vars sample, theta 2, lambda 1, sigma 0.1, max.iter.cluster 20, epsilon.cluster 1e-4, epsilon.harmony 1e-4, plot_convergence FALSE, verbose FALSE )theta多样性惩罚参数默认2。theta越大Harmony就越强烈地推动不同批次在同一个聚类里混合theta越小混合力度越弱。如果你发现整合不足不同批次还是明显分离可以适当调大theta比如3或4。如果你发现整合过度不同细胞类型被强行拉到一起可以调小theta比如1或0.5。lambda岭回归正则化参数默认1。它控制批次代表细胞估计的平滑程度。lambda越大估计的批次代表细胞越平滑对异常值的容忍度越高lambda越小代表细胞估计越容易受到少数细胞的影响。多数情况下默认值1就可以不用动。只有在样本量极不平衡时才需要考虑调大比如一个样本有5万细胞另一个只有1000适当增大lambda可以让小样本的代表细胞估计更稳定一些。sigma高斯核的带宽参数默认0.1。它影响细胞距离相似度的衰减速度。sigma越小只有非常近的细胞才被认为是相似的sigma越大较远的细胞也会相互影响。大多数情况下0.1是一个稳健的默认值。如果细胞类型很多且离散度高可以试0.2如果细胞类型之间边界很模糊可以试0.05。max.iter.cluster聚类迭代的最大次数默认20。如果你发现Harmony没有收敛日志里一直显示还有迭代在跑可以适当提高这个值。但如果数据量很大这个值也不要设太高否则运行时间会显著拉长。实际处理20万细胞时我会把max.iter.cluster降到10因为20次迭代的运行成本太高而且通常在10次左右就已经收敛了。epsilon.cluster和epsilon.harmony收敛阈值默认都是1e-4。一般不用调。如果你跑大数据集时觉得太慢可以放宽到1e-3能显著减少迭代次数但整合精度会略微下降。我的建议是先默认跑一遍观察日志里的收敛情况再决定要不要放宽。plot_convergence设为TRUE可以画出收敛曲线。第一次跑某个数据集时建议打开看一下确认Harmony确实在第几轮迭代后收敛了免得白白浪费算力。2.4 维数选择的蝴蝶效应RunHarmony的输入来自RunPCA得到的PC嵌入。默认情况下RunPCA会计算50个PC但RunHarmony默认只取前30个PC作为输入。这个“dims”的选择对整合效果影响很大经常被忽略。如果dims设得太小比如只取5-10个PC你等于只保留了表达变化最剧烈的前几个方向这通常是批次效应主导的方向而真正的生物学信号可能排在后面直接就被截掉了。结果是Harmony整合出来的空间里细胞类型分不开。如果dims设得太大比如取到50个那你把大量噪声方向也放进来了Harmony会花很多精力去整合噪声整合效果反而变差还增加运行时间。我的经验是先用ElbowPlot画出每个PC的方差占比曲线找到曲线变平的拐点。一般取拐点前后浮动5个PC。如果拐点不太清楚保守起见取20-30个PC。对于大多数10x平台产生的单细胞数据20-30是比较舒服的范围。你也可以用两种维度设置跑一遍UMAP做对比看哪个结果里细胞类型分得更清楚且批次混合更均匀拿实际结果说话。另外注意RunPCA和RunHarmony时的dims必须保持一致。千万别在RunPCA时用了npcs50然后到RunHarmony时只取dims10这样浪费了PCA后半段的信息还有可能在后续FindNeighbors时维度对不上。建议把这些都写在同一段代码里用同一个变量控制减少出错概率。3. 完整实操流程与参数组合建议理论部分讲完了这里给出一套可以直接复制到RStudio跑的完整代码。我以10x Genomics的PBMC数据为例数据格式是cellranger输出的三个文件barcodes.tsv.gz、features.tsv.gz、matrix.mtx.gz。实际项目中你只需要把数据路径和metadata列名替换成自己的即可。3.1 环境准备和必要包检查需要确保R版本不低于4.1Seurat版本不低于4.0。seurat-disk和harmony这两个包缺一不可。Harmony的安装比较简单但从GitHub安装时偶尔会遇到依赖包编译问题建议优先用CRAN版本。library(Seurat) library(harmony) library(tidyverse) library(RColorBrewer) # 检查版本 packageVersion(Seurat) packageVersion(harmony)如果你的数据是cellranger输出建议用Read10X读取然后创建Seurat对象时加上min.cells和min.features过滤条件。这里min.cells3表示基因至少在3个细胞中表达才保留min.features200表示细胞至少检测到200个基因才保留。这两个过滤条件能去掉绝大部分空液滴和破碎细胞但是别把过滤阈值设得太狠否则稀有细胞类型可能被直接过滤掉。data_dir - path/to/your/cellranger/outs/filtered_feature_bc_matrix counts - Read10X(data.dir data_dir) obj - CreateSeuratObject( counts counts, project my_project, min.cells 3, min.features 500 )3.2 SCTransformHarmony完整流程下面这一段是我日常项目里用得最多的一套参数组合。注意看我在读取数据后先计算了percent.mt然后传给SCTransform的vars.to.regress。# 添加质量指标 obj[[percent.mt]] - PercentageFeatureSet(obj, pattern ^MT-) # 可视化QC确认过滤阈值 VlnPlot(obj, features c(nFeature_RNA, nCount_RNA, percent.mt), ncol 3) plot1 - FeatureScatter(obj, feature1 nCount_RNA, feature2 percent.mt) plot2 - FeatureScatter(obj, feature1 nCount_RNA, feature2 nFeature_RNA) plot1 plot2 # 过滤细胞 obj - subset(obj, subset nFeature_RNA 500 nFeature_RNA 5000 percent.mt 20)QC这块有个细节percent.mt的阈值不要一上来就拍脑袋定成5%或者10%。先画VlnPlot看分布再决定阈值。如果大多数细胞都集中在2%-5%那阈值定在10%就很宽松如果大多数细胞都有10%左右的线粒体比例贸然过滤到5%会把一大半细胞删掉。要基于数据本身做判断不要套模板。# SCTransform obj - SCTransform( obj, vars.to.regress percent.mt, variable.features.n 3000, conserve.memory TRUE, seed.use 42, verbose FALSE ) # PCA obj - RunPCA(obj, npcs 30, verbose FALSE) # 确认主成分 ElbowPlot(obj, ndims 30)如果你之前没有跑过标准化SCTransform会自动完成NormalizeData和ScaleData的操作所以不需要再用NormalizeData。但如果你的数据里还保留了RNA assay后续做FeaturePlot时默认绘制的可能是RNA assay的log1p表达量所以没关系不用额外处理。# Harmony整合 obj - RunHarmony( obj, group.by.vars sample, theta 2, lambda 1, sigma 0.1, max.iter.cluster 20, epsilon.cluster 1e-4, epsilon.harmony 1e-4, plot_convergence TRUE, verbose TRUE )这里group.by.vars“sample”是假设你的metadata里有一列叫sample表示每个细胞的样本来源。如果你的批次信息在别的列比如batch、library就改成对应的列名。# 聚类与UMAP obj - FindNeighbors(obj, reduction harmony, dims 1:30) obj - FindClusters(obj, resolution 0.8) obj - RunUMAP(obj, reduction harmony, dims 1:30, seed.use 42) # 可视化 DimPlot(obj, reduction umap, group.by sample) ggtitle(By Sample) DimPlot(obj, reduction umap, group.by seurat_clusters) ggtitle(By Cluster)我特别提醒一下FindNeighbors里的reduction必须指定成“harmony”很多人漏掉这一步导致后续聚类用的还是未经整合的PCA结果。如果你在日志里看到类似“Using PCA as reduction”的提示那就是没写上reduction参数。这种情况尤其常见在只需要跑单样本数据的脚本上因为单样本根本不需要Harmony所以没问题但一换成多样本整合漏掉reduction参数结果就完全不同了。3.3 参数组合速查表下面是我在不同场景下的参数组合建议可以直接当作参考配置抄走。使用场景SCTransform核心设置Harmony核心设置注意事项常规PBMC/血液样本vars.to.regress percent.mt, variable.features.n 3000theta2, lambda1, dims1:20通用首选稳定但缺乏针对性肿瘤组织高复杂度variable.features.n 5000theta3, dims1:30高变基因多一点否则稀有亚群信号容易丢脑组织/复杂组织vars.to.regress c(percent.mt), variable.features.n 4000theta2, lambda1, dims1:30注意细胞类型注释难度极小样本量2-3个样本return.only.var.genes TRUEtheta1theta调低过度整合风险更大大样本量20万细胞conserve.memory TRUE, variable.features.n 3000max.iter.cluster 10保证服务器内存注意运行时间细胞周期影响明显vars.to.regress c(percent.mt, S.score , G2M.score)theta2先试试不回归的区别确认细胞周期确实是混淆因素再回归3.4 可视化检查批次整合效果跑完Harmony之后第一件事不是急着去看聚类注释而是检查整合效果到底好不好。两个图必看一个按样本上色的UMAP一个按聚类上色的UMAP。按样本上色的UMAP用来检查批次混合程度。理想情况是不同样本的细胞在同一个细胞类型区域里均匀混合而不是各自聚成独立的一堆。但这里有个容易误判的点如果某种细胞类型只出现在某一个样本中那这个区域的细胞就只能来自这个样本看起来“没有混合”其实是正常的生物学现象。所以检查整合效果时一定要结合细胞类型注释来看不要只看“是否完全混合”。按聚类上色的UMAP用来检查聚类分辨率是否合适。如果用了0.8的resolution但聚类数量还是很少比如只有4-5个说明数据复杂度低或者之前的基因筛选/整合过度了。如果聚类数量特别多比如30个可以考虑降低resolution到0.4-0.5或者先用SingleR做初步注释再看是否需要细分。我见过一个更直观的检查方法对某一个已知marker基因做FeaturePlot看这个marker是否只在预期细胞类型里表达。如果marker表达模式良好说明整合没有破坏生物学信号如果marker表达被打散了说明Harmony的强度过大需要调低theta。4. 常见问题与排查技巧实录这一节把我在实际项目中遇到的高频问题整理出来每一个都是真实场景排查思路和解决方案可以直接套用。4.1 问题一SCTransform之后找不到counts数据这种情况最常见于使用旧版Seurat4.x之前或者没有正确设置new.assay.name。如果你在SCTransform后运行FeaturePlot(obj, features CD3D)有可能会报错说找不到这个基因。这是因为SCT assay的scale.data只包含了高变基因return.only.var.genesTRUE。非高变基因的CD3D并没有进入scale.data。但FeaturePlot默认读取的是data槽而不是scale.data所以如果报错多半是你把DefaultAssay切换得不对。排查思路DefaultAssay(obj) - SCT FeaturePlot(obj, features CD3D)如果还报错看下这个基因是否被过滤掉了。在小细胞量的数据集中很多基因会因为min.cells过滤条件被删除。这时你可以临时把DefaultAssay切换成RNA再试因为RNA assay保留了全部检测到的基因。4.2 问题二Harmony跑了但聚类还是按样本堆我在一个六样本的PBMC数据集上遇到过这种情况跑了Harmony但UMAP上样本还是各占一块区域没有混合。排查后发现两个原因一个是在FindNeighbors时忘了指定reductionharmony导致聚类用的还是PCA嵌入另一个是theta设成了默认的2但样本间差异实在太大混合力度不够。解决办法是两步走检查FindNeighbors/FindClusters/RunUMAP是否都用了reductionharmony”。把theta调大从2调到4同时适当增大max.iter.cluster比如从20调到50让Harmony有更多迭代次数去混合。这里要注意调大theta可能带来过度整合的风险所以每次只调一个参数观察一下UMAP变化再决定下一步。有时候问题也不在Harmony本身而在上游的PCA维度选择。如果dims1:10而主要批次效应恰好集中在前10个PC之外那么Harmony根本没有把那个PC纳入进来自然整合不到。这种情况的解决办法是增加dims比如改成1:30看看是否有改善。4.3 问题三回归掉线粒体比例后细胞类型消失了之前我处理一个肿瘤浸润免疫细胞的数据集在SCTransform的vars.to.regress里回归了percent.mt结果发现一群高表达线粒体基因的巨噬细胞亚群完全消失了。后来仔细查了一下这群细胞高表达mitochondrial相关的基因是它们的生物学特征而不纯粹是技术噪声。强行回归掉percent.mt等于把这群细胞的身份特征给抹掉了。解决方案是换一套参数重新跑SCTransform不回归percent.mt对比两套结果。如果主要细胞类型在两种设置下都能分出来说明不回归也没关系。如果少了某个细胞类型再用回归percent.mt的结果来做后续分析但需要在论文里说明这个处理对特定细胞类型可能的影响。还有一种更稳妥的做法把percent.mt作为“已知混杂因素”在Harmony整合时处理而不是在SCTransform里回归。Real world数据里批次效应、样本间差异、线粒体比例经常交织在一起你需要在不同环节用不同工具处理不同因素不能指望一个参数解决所有问题。4.4 问题四数据量大后内存爆掉SCTransform在大数据集上非常吃内存。我之前处理一个20万细胞的数据集16G内存的服务器直接OOM。后来发现是return.only.var.genes和conserve.memory这两个参数的组合出了问题。如果你把return.only.var.genes设为FALSE比如后续想用非高变基因做分析SCT assay里会保存全部基因的scale.data内存占用会爆炸式增长。更优的策略是SCTransform时保持return.only.var.genesTRUE等需要用非高变基因时再随时提取RNA assay的原始counts做分析。这样可以大幅降低内存占用。另外如果你不需要保留原始RNA assay的全部信息可以在一开始读入数据后先用NormalizeData或者直接删除掉RNA assay但一般不建议删因为后续找marker基因时经常要用到counts数据。最稳妥的方案是在创建Seurat对象时就只保留必要的细胞和基因减少数据量。4.5 避坑清单速查把上面提到的问题和关键点整理成一个速查表方便你保存对照。检查项检查内容推荐做法SCTransform回归变量是否把细胞周期、样本ID、批次信息放入了vars.to.regress只用percent.mt必要时再用细胞周期样本ID和批次留给Harmony高变基因数量是否设置了合理的高变基因数量标准数据3000复杂组织4000-5000return.only.var.genes是否与你的后续分析冲突保持默认TRUE需要时再临时修改conserve.memory大数据集是否开启数据集超过5万细胞时强烈建议开启dims选择PCA和Harmony的维度是否一致用ElbowPlot选拐点常用1:20或1:30reduction指定FindNeighbors/RunUMAP是否用了harmony检查日志确认没有用PCAtheta值是否根据整合效果动态调整混合不足调大过度整合调小group.by.vars是否把生物学分组误当批次变量只放技术批次变量不放疾病组/对照组收敛性检查是否确认Harmony收敛用plot_convergenceTRUE查看曲线4.6 一个容易被忽略的细节细胞注释和Harmony的先后顺序最后说一个很多教程都不提前讲清楚的事情如果你打算用SingleR或CellTypist做自动注释注释结果是基于原始表达矩阵的和Harmony整合后的低维嵌入没有直接关系。所以先用注释工具给细胞打上类型标签再基于Harmony整合后的嵌入去看不同类型细胞的分布是一种非常高效的工作流。这样不但可以检查整合效果还能快速发现“某些细胞类型是否在不同样本间有偏移”。反过来如果你先跑聚类再注释聚类结果完全取决于你设置的resolution和Harmony参数如果这些参数设得不好注释结果也会被带偏。先注释再聚类至少能保证你对细胞类型的判定是独立于整合参数的后续调参时也有一个稳定的参照系。我在实际项目中几乎总是先做一轮初步注释再回头调整Harmony参数这样能少走很多弯路。5. 写在最后的个人体会这套SCTransformHarmony的组合我前前后后跑了不下几百个数据集最大的感受就是没有一套参数能吃遍天下但理解了参数背后的原理之后你会知道在什么场景下该动哪个旋钮。调参不是玄学它是在“保留真实生物学信号”和“去除技术批次效应”之间找平衡而找平衡的依据永远来自你对自己数据的理解和一次次比较实验。我自己的习惯是每一个新数据集都会先让SCTransform和Harmony用默认参数跑一遍把QC图、elbow图、整合前后的UMAP图全部存下来留档然后再根据这些图决定要不要调整参数。很多时候问题不是参数不够好而是你根本没来得及看默认结果就急着往下走了。多花十分钟看看那些图后面能省下好几个小时的Debug时间。最后分享一个小技巧给每次分析建立独立的输出目录把不同参数组合跑出来的UMAP图和聚类统计保存成带标签的文件比如umap_theta2_dims30.png、umap_theta4_dims30.png这样你可以随时翻回去对比哪个参数组合在你的数据上表现最好也能避免复现实验时搞混了版本。这一点听起来很平常但我见过太多人栽在这上面了宁可多做一步存档也不要事后抓狂。