HD数据的细胞注释,我前前后后改了至少四版思路,每次都是被真实数据里那些“不讲道理”的簇逼着改的。HD这里指的是亨廷顿病(Huntington's Disease),这类课题拿到的通常不是常规单细胞悬液,而是人脑组织的单核转录组数据,也就是snRNA-seq。单细胞注释看起来是个标准流程——找一套marker基因库,对着UMAP标签一个个核对就完事了。真到了HD这类神经退行性疾病数据上,你会发现经典marker大量失灵,疾病相关的应激状态、RNA丰度变化、细胞比例重组全都在干扰判断。这篇把我摸出来的坑和更新后的注释思路整理出来,适合正在做脑组织单细胞或单核数据的同学参考,尤其是被神经元亚型注释卡住的朋友。
1. HD数据的特殊性,决定注释流程不能照搬常规
1.1 HD数据到底是什么场景
亨廷顿病是由HTT基因外显子1的CAG三核苷酸重复扩增引起的常染色体显性遗传病,病理核心是纹状体中型多棘神经元(SPN)的选择性丢失,尤其以表达多巴胺D2受体的SPN损伤最早、最重。因为人脑组织难以解离出高质量活细胞,研究通常用死后脑组织做成单核悬液,做snRNA-seq。这就带来第一个现实问题:核转录组里细胞质mRNA占比很低,很多经典marker的检出率远低于普通单细胞数据。
即便有人选择HD小鼠模型,比如Q140或者zQ175品系,数据特征与病人样本类似,仍然绕不开神经组织本身的问题:少突胶质细胞占比高,神经元占比低但转录组复杂度高,疾病中后期还有明显的小胶质细胞激活和星形胶质细胞反应性增生。也就是说,HD数据的细胞组成本身就不是静态的,它天然带着“疾病重构”的痕迹。
另一个容易忽略的点是样本质量和死后间隔(PMI)的影响。神经元对缺氧极其敏感,PMI一长,神经元核RNA降解明显,导致部分神经元簇表现出“低质量”特征。如果你按常规的质控标准把低基因数簇直接去掉,可能就把最脆弱的D2 SPN簇给扔了。这类数据做注释,第一步不是找marker,而是先搞清楚你手里的数据是哪一种质量状态。
1.2 疾病状态让经典marker频繁失灵
我在HD数据上踩得最狠的第一个坑,就是小胶质细胞标记。教科书里小胶质细胞常用P2RY12、CX3CR1、TMEM119,但这些marker在激活状态下会显著下调。HD里大量小胶质细胞处于疾病相关状态,P2RY12表达往下掉,反而TREM2、CD68、APOE、TYROBP一类的基因迅速上升。你要是拿着P2RY12去找小胶质细胞,会漏掉一大片真正被疾病激活的群体,或者更糟——把激活小胶质细胞误判成另一种巨噬细胞亚群。
星形胶质细胞也类似。GFAP是反应性星形胶质细胞的核心标志,但在稳态星形胶质细胞里表达并不高。snRNA-seq的灵敏度本来就低,GFAP检出率更不稳定,单看GFAP极容易把一部分星形胶质细胞漏判。加上反应性星形胶质细胞会上调VIM、SERPINA3、OSMR等基因,这些基因又会干扰聚类结构。如果你把GFAP当成星形胶质细胞的“充分必要条件”,在HD数据里直接看不了完整星形胶质细胞图谱。
还有一个隐性干扰源,就是热休克蛋白和即早基因。HD神经元处于持续应激状态,HSP90AA1、HSPA1A、FOS、JUN这类基因在不少细胞里都高表达。聚类算法是按全转录组相似性分的,应激基因一旦抬高,会把不同细胞类型往一起推,形成“应激团块”。这一团里可能既有神经元也有胶质细胞,marker特征模糊,注释时特别容易糊弄过去。我试过直接对这团做常规注释,结果怎么分都别扭,后来才意识到要先看应激模块,而不是硬套marker。
2. 从“marker驱动”转向“聚类驱动”的注释思路
2.1 旧思路为什么经常翻车
我第一次做HD注释用的还是经典三板斧:聚类完,拿一套已知marker列表,跑DotPlot,挨个cluster看谁像谁,然后贴标签。这套流程用在肿瘤或外周血数据上问题不大,因为那些样本的细胞身份相对清晰,细胞状态差异不会把marker表达完全扭曲。但HD数据不一样,它的问题是很多cluster同时带“身份模糊”和“状态偏移”两个属性。
举个例子,部分星形胶质细胞在反应性状态下会低表达AQP4,但它的转录组总体跟少突胶质细胞前体(OPC)有一部分重叠。你用单个marker看,会觉得这个簇“有一点点星形,也有一点点OPC”,最后注释非常纠结。还有一个更常见的问题:单核数据里DRD2的检出率很低,MSN(中型多棘神经元)里D2亚型的特征被严重弱化。如果只看marker阳性比例来注释SPN亚型,D2簇经常被误并进D1簇或者直接丢掉。
旧思路最大的问题在于“先定身份再找证据”。你心里假设了这个cluster应该是什么,然后去marker列表里找支持证据。HD数据里大多数cluster不是教科书状态,硬套必然出现大量“无法注释”或“注释错误”的中间簇。
2.2 更新后的注释流程骨架
我现在跑HD数据的注释顺序是这样的:质控和标准化之后,先用Harmony或Seurat的标准整合流程处理多批次样本,然后跑一个相对低分辨率的粗聚类(resolution大约0.5到0.8)。在贴任何细胞身份标签之前,先做三件看起来“多余”的事。
第一,把每个簇的高变基因和top差异基因拉出来看。这一步不是为了直接注释,而是判断这个簇到底是“真细胞类型”还是“技术伪影”或“应激反应簇”。如果top基因里出现了大量线粒体基因、核糖体基因,就要警惕低质量簇;如果出现清一色的HSP家族基因,那就是应激模块簇。第二,看每个簇的样本组成和条件组成。某个簇如果全部来自一个样本,那很可能是批次效应而不是生物学差异。第三,看不同分辨率下这个簇的稳定性。如果它在分辨率0.5时存在,到1.0时被拆散混到别的簇里,说明它本身不是一个可靠的身份单元。
做完这三步,才开始大类注释。我习惯先用“谱系marker组合”把神经细胞和非神经细胞分开:神经元看RBFOX3、SNAP25、SYT1,少突胶质细胞看MOBP、PLP1、MBP,星形胶质细胞看AQP4、SLC1A2、GJA1,小胶质细胞看CX3CR1、CSF1R,OPC看PDGFRA、CSPG4,内皮看CLDN5、FLT1,T细胞看CD3D、CD3E。这里特别强调“组合”而不是“单个marker”,因为任何单一基因在HD数据里都可能出错,组合打分能显著提升鲁棒性。
大类注释完,再把感兴趣的谱系subset出来,提高分辨率继续聚类做亚型注释。这个“先划大类、再钻亚型”的顺序,能避免你在全细胞图谱上做精细注释时被跨谱系的相似转录组状态干扰。我试过直接在总聚类上分辨D1和D2 SPN,结果因为簇太大、分辨率不够,边界都是糊的。subset之后重新聚类,清晰度完全不一样。
3. 自动注释工具:起跑线而不是终点线
3.1 几款常用工具的实际体感
讲到细胞注释,绕不开自动工具。我在HD数据上试过的主流工具基本都跑了一遍,简单整理一下感受。
| 工具 | 原理 | 在HD数据上的表现 | 主要坑点 |
|---|---|---|---|
| SingleR | 与参考转录组数据集算相关性打分 | 对大类细胞注释尚可,但疾病细胞容易被映射到相近的正常类型 | 参考库如果来自健康脑组织,激活小胶质细胞会被判为普通巨噬细胞 |
| scType | 基于marker基因集富集打分 | 自定义marker列表时灵活,结果直观 | marker列表质量直接决定上限,需要自己反复调 |
| CellTypist | 预训练逻辑回归模型 | 常见细胞类型识别稳定,但神经元亚型粒度不够 | 对SPN的D1/D2区分基本无能为力 |
| Azimuth | 参考映射到标准图谱 | 在标准脑区数据上效果不错,对疾病状态样本容易“强行归位” | 异常状态簇会被映射到最相近的正常类型,掩盖真实状态 |
| scGate | 逐级marker组合门控 | 透明度高,能自定义多级判断规则,适合过滤杂质簇 | 前期规则构建繁琐,需要熟悉marker逻辑 |
体感很明确:没有任何一个工具能直接输出HD数据里“可信的细胞注释”。它们适合用来做交叉验证,或者在大类层面上快速给一个初筛结果,但疾病相关的异常细胞状态,自动工具普遍理解不了。比如SingleR会把激活小胶质细胞判成巨噬细胞甚至单核细胞样本,Azimuth会把应激状态下转录组偏移的D1 SPN映射成其它皮层神经元亚型。这类错误非常隐蔽,因为你只看返回的标签时会觉得“好像也没错”。
所以我现在的定位是:自动工具先跑一轮,拿结果当“待验证假说”,再用人工marker复核和生物学背景知识做最终裁决。不要一上来就信任自动注释的完整标签,更不要直接拿它做下游差异分析的注释列。
3.2 自动注释结果的四步复核法
每次用自动工具跑出结果后,我固定做四步检查,这一步让我发现了不少自动注释的隐蔽错误。第一步是DotPlot复核重点marker,直接在图上肉眼看每个簇的marker表达是否符合预期。第二步是FeaturePlot看关键marker的空间分布,比如PENK是否特异性出现在你认为的D2簇里,而不是散在到别的神经元簇中。第三步是提取每个簇的top差异基因,手动过一遍,看有没有明显属于另一个谱系的marker混合进来。第四步是用独立打分工具,比如AddModuleScore或AUCell,对一组核心身份marker打分看一致性。
这套四步复核法看起来费时间,实际用熟了以后很快。因为HD数据的核心问题往往集中在少数几个“疑难簇”上,绝大部分簇在大类注释上是一致的,不需要逐簇深究。你要做的就是把那两三个有争议的簇挑出来,集中火力确认它们的真实身份。
比如我遇到过的一个案例:某个cluster在SingleR里被判为“兴奋性神经元”,但它的top差异基因里出现了GAD1和GAD2。我一开始以为只是污染,后来通过四步复核发现,这个簇其实是GABA能中间神经元的一个亚型,因为疾病应激导致兴奋性相关基因异常上调,把SingleR的相关性打分带偏了。如果只看自动注释结果,这个簇就彻头彻尾弄反了。复核到这一步,才能真正体会到“自动工具给出的是起点,不是结论”这句话的分量。
4. HD数据里最难的神经元亚型注释:SPN的D1/D2与中间神经元
4.1 为什么D1/D2分型是HD注释的核心
HD课题绕不开纹状体SPN,这类神经元占纹状体神经元总数的绝大多数,又分成两个功能亚型。D1 SPN主要表达DRD1和TAC1,走直接通路,投射到黑质网状和脚桥核;D2 SPN主要表达DRD2和PENK,走间接通路,投射到外侧苍白球。在HD里D2 SPN比D1 SPN更早出现退行性变化,所以能否准确区分这两个亚型,直接决定你对疾病选择性损伤的解读对不对。
问题在于,snRNA-seq里DRD2的检出率很低。我遇到的情况是,D2 SPN簇里DRD2阳性细胞比例经常只有20%到40%,跟普通scRNA-seq数据完全不在一个量级。如果你严格按“DRD2阳性率”来定义D2 SPN,可能一大半真正的D2细胞会被漏掉。这时候就得换“替补marker”策略,比如用PENK作为D2的稳定标记,用TAC1配合DRD1指定D1,同时用PPP1R1B(DARPP32)先确认SPN的身份,再细分亚型。
4.2 实际操作怎么把SPN亚型分开
我的固定做法是先把所有神经元簇subset出来,排除胶质细胞干扰,然后把这个子集重新跑一遍聚类,分辨率直接放到1.0到1.2。之所以要重新聚类,是因为在全局聚类里SPN往往被压成一个大簇或两三个大簇,亚型信号全被淹没。subset之后再聚类,D1和D2的差异基因才有机会主导分群。
分型时我会建一个组合得分,而不是单看一个基因。D1模块用DRD1、TAC1、C1QL1,D2模块用DRD2、PENK。如果某个簇在两个模块的得分都比较模糊,先不要硬归,看它的FOXP1、FOXP2表达以及它在新聚类图中的位置。D1和D2 SPN在UMAP上通常会形成一个连续的过渡带,中间少数细胞可能是双阳性或处于转换状态,把它们单独分成一个“SPN transitional”注释,也比硬塞进D1或D2更诚实。下游分析时你可以再决定是保留还是排除这个过渡簇。
还有一个容易混淆的点:纹状体中间神经元。GABA能中间神经元表达PVALB、SST或NOS1,胆碱能中间神经元表达CHAT。在HD数据里,SPN数量减少后,相对比例上中间神经元会显得变多,如果不小心把一部分中间神经元并进了SPN簇,后续差异表达分析会引入大量无关信号。我在注释时用CHAT、NOS1、PVALB、SST做排除条件,凡是这些marker明显阳性的簇绝不注释成SPN。
4.3 做亚型注释时的不确定性处理
HD后期的D2 SPN可能因为大量丢失,聚类中只挤成一团很小的簇,甚至完全消失。这时候不要强行找D2,而是要回到样本比例和病理学数据交叉验证。如果你用的是人死后脑组织样本,尾状核背侧D2 SPN丢失严重,那聚类里D2簇变弱本身就是疾病信号,不是注释错误。反过来,如果D2比例在晚期样本里还跟健康对照一样高,那才要警惕注释可能把其它细胞误认成了D2。
对于这种“该缺失却不缺失”的反常情况,我的排查思路是先看免疫组化或文献数据里这个样本的病理阶段,再看其它marker是否支持D2身份,比如PENK、DRD2、ADARB2这类相对特异的基因是否确实表达。如果支持证据不足,我宁可标为“SPN_unresolved”,也不会硬贴D2标签。带着不确定性的注释在后续分析里不会坑你,硬贴的标签大概率会。
5. 注释结果的验证:不能被UMAP骗过去
5.1 三层验证法
注释完成后我习惯做三层验证,第一层是marker特异性验证,看目标marker是否只在对应簇里高表达,平均表达量和阳性比例是不是显著高于其它簇。第二层是跨样本一致性验证,把注释好的结果按样本拆分,检查同一个细胞类型是否在每个样本或大部分样本里都存在。如果某个亚型只在单个样本里出现,优先怀疑它是个体差异或批次效应,而不是真实疾病相关变化。第三层是生物学合理性验证,检查细胞类型比例是否符合已知病理规律。
这层验证听着基础,但我在实际项目中见过很多次“Cluster全对,比例全错”的情况。原因是自动注释时把某个疾病激活状态的亚群单独分了出来,导致原本的细胞类型被拆成两个注释,比例自然就偏了。比如激活小胶质细胞和稳态小胶质细胞如果被标成两个不同身份,那“小胶质细胞”总比例是正确的,但“稳态小胶质细胞”比例异常下降,下游比较就会得出错误的结论。这种问题marker特异性验证看不出毛病,只有回到比例分布才能暴露。
5.2 混合簇的判定与处理
还有一种常见情况是某些cluster同时高表达两个谱系的核心marker,比如GFAP和MBP同时阳性。这在HD数据里需要特别小心。它可能来自真实的双核污染,也可能是两个谱系在UMAP上因为细胞状态相似而黏在一起。我的判断方法是先把这个cluster单独提取出来,重新跑一次PCA和聚类,看能不能稳定再分。能分开就说明是原始分辨率不够,属于注释粒度问题;分不开再看这个cluster的双细胞评分,比如Scrublet或DoubletFinder的得分,拉高的话直接删掉。
如果既不是分辨率问题也不是双细胞问题,那有可能是发育中间态或者疾病异常重编程状态。比如反应性星形胶质细胞和部分OPC在某些条件下会共享一批基因表达,造成谱系归属模糊。这种簇我用“transitional”或“disease-perturbed”这类描述性注释,而不会强行归到某一个大类。先保留不确定性,再去看这个簇的细胞组成与空间位置,往往比硬分更有价值。
5.3 用模块打分把身份与疾病状态拆开
这一步是我在HD数据上最大的思路升级,就是把“细胞身份”和“疾病状态”拆成两个维度分开打分。具体做法是设计两套独立的基因模块,一套是身份marker,比如SPN_D1、SPN_D2、小胶质稳态、星形稳态;另一套是状态marker,比如小胶质激活、星形反应性、神经元应激。用AddModuleScore分别打分,然后再看这些分数在cluster里的分布。
library(Seurat) id_modules <- list( SPN_D1 = c("DRD1", "TAC1", "C1QL1"), SPN_D2 = c("DRD2", "PENK"), Microglia_Homeostatic = c("P2RY12", "CX3CR1", "TMEM119") ) state_modules <- list( Microglia_Activated = c("TREM2", "CD68", "TYROBP", "APOE"), Astro_Reactive = c("GFAP", "VIM", "SERPINA3", "OSMR"), Neuron_Stress = c("FOS", "JUN", "HSP90AA1", "HSPA1A") ) obj <- AddModuleScore(obj, features = id_modules, name = "ID") obj <- AddModuleScore(obj, features = state_modules, name = "STATE")这样打完之后,同一个cluster的身份得分和状态得分是正交的信息。比如某个簇身份模块显示它是星形胶质细胞,同时状态模块显示反应性得分非常高,那它就是一个“反应性星形胶质细胞”,而不是一种新细胞类型。如果不拆分维度,这种簇很容易被误认为“疾病特异的新亚群”,下游分析就麻烦了。
6. HD数据细胞注释的坑位清单与排查实战
6.1 同一个基因在两个细胞类型里都表达怎么办
这不是特例,而是常态。比如GJA1(Connexin 43)在星形胶质细胞里高表达,但在部分内皮亚群和周细胞里也有。S100B在星形胶质细胞里高,但髓系细胞也可能有低表达。如果你因为一个marker同时出现在两个不同类型的簇里,就想当然给其中一个簇贴错标签,那后面所有分析都会带偏。
我的办法是构建“组合marker规则”,而不是依赖单一基因。比如星形胶质细胞的判定条件是AQP4和SLC1A2至少有一个明显阳性,同时GFAP或GJA1作为辅助证据;内皮细胞必须CLDN5和FLT1共同支持。任何单基因阳性都不足以完成注释。把这种组合规则整理成一个小型评分函数,每次拿到新数据集后直接复用,能省大量时间。
6.2 疾病效应和细胞类型身份混淆是HD数据的核心陷阱
HD样本中,疾病状态会引起某些细胞类型的转录组发生系统性变化,这种变化有时候比细胞类型之间的差异还要大。这会导致同一细胞类型在健康组和疾病组分别聚类,看起来就像两个独立的细胞类型。第一次遇到这种情况时,我差点把一个疾病组的D1 SPN亚群注释成“新型神经元亚型”,后来一查发现它的核心身份marker完全没有变,只是多了一堆应激基因。
处理这类问题的关键是在做差异表达之前,先把身份注释锚定在“保守的身份marker”上,把疾病引起的可变基因排除在身份判定之外。然后用Harmony或类似整合算法拉齐不同条件之间的批次/状态差异,再重新聚类和注释。这样注释出来的“类型”是跨条件可比的,下游你才有资格说某个细胞类型的比例变化或表达变化是疾病引起的。
另外还要养成一个好习惯:每次注释完都留一份“marker证据表”,记录每个cluster的判定依据,包括用了哪些marker、阳性比例多少、自动工具的原始注释是什么、你是否做了修改以及修改原因。这个表看起来麻烦,但在审稿或被合作方质疑注释时,它就是最强的底牌。HD数据里太多注释争议来自“这个cluster为什么叫这个名字”,有了证据表,这种争议至少能建立在事实基础上。
6.3 快速排查表
| 现象 | 可能原因 | 处理方法 |
|---|---|---|
| 一个簇同时表达两类marker | 分辨率不足或双核污染 | 提取子集重新聚类,看能否分开;检查双细胞评分 |
| 经典marker阳性率极低 | snRNA-seq灵敏度限制 | 换成核内高表达的替代marker,或使用组合打分 |
| 多个簇都高表达应激基因 | 疾病相关应激反应 | 用模块评分把状态和身份拆开,不强行注释状态簇 |
| 某一细胞类型只在单个样本出现 | 个体差异或批次效应 | 检查样本贡献度,用整合算法验证 |
| 疾病组和健康组同一类型分成两簇 | 疾病效应大于类型效应 | 基于保守marker重新注释,整合后再比较 |
| 自动注释结果与生物学背景不符 | 参考库或模型不适用疾病样本 | 以人工marker复核为准,自动结果降级为参考 |
写在后面
回头看我最初那版“机械套marker”的注释流程,现在这套方案最大的改变就是不再把细胞身份当成本质,而是把它当作一个需要在特定背景下持续检验的判断。HD数据尤其如此,它所有细胞都处于被疾病重塑的动态过程里,身份边界是模糊的。我自己现在的固定动作,是把疾病状态相关基因单独列一个模块打分,跟细胞身份完全分开看,这一招帮我避掉了很多假注释。每次遇到一个不按教科书出牌的cluster,先别急着找它“是什么”,退一步问一句“这个簇的存在是不是被疾病或技术影响扭曲出来的”,往往就有答案了。HD数据不会给你一条清晰的分界线,注释本质上是跟不确定性打交道,把每一步检查做得扎实,比迷信任何工具都更可靠。