CITESeQC：首个用于CITE-Seq数据多层级定量质量控制的创新软件工具

本文推荐一款创新软件CITESeQC，作为首个专门针对CITE-Seq数据的系统性质量评估工具。它通过多层级分析（包括RNA、表面蛋白及互作关系），引入香农熵（Shannon entropy）和相关性定量指标，有效解决现有方法在跨模态验证、系统性误差检测和标准化评估方面的不足，显著提升细胞分型准确性与临床应用的可靠性。

　　

背景

单细胞RNA测序（scRNA-seq）技术仅通过捕获和测序RNA分子来分析基因表达，而转录组和表位测序的细胞索引（CITE-Seq）则通过利用独特的DNA条形码抗体（也称为“抗体衍生标签（ADT）”），同时检测RNA分子和感兴趣的表面蛋白。由于细胞表面蛋白作为细胞身份和功能的标记物和通讯分子，CITE-Seq数据不仅能够识别细胞类型特异的基因表达模式，还能识别由特定表面蛋白定义的细胞类型，这些蛋白可能用于进一步的临床应用。例如，尽管某些免疫细胞类型，如γ/δ T细胞、黏膜相关不变T细胞（MAIT细胞）、先天淋巴样细胞（ILC）和中性粒细胞，已显示出显著的临床潜力，但仅凭scRNA-seq数据通常不足以可靠地检测它们。这一限制源于谱系定义转录本可能具有低RNA含量、存在高水平的RNase，或者mRNA表达模式并不总是与蛋白质表达相关。

为确保从CITE-Seq数据中获得高质量发现，第一个关键步骤是控制输入CITE-Seq数据的质量（QC）。对于CITE-Seq数据的QC，以往的研究执行了有限的分析，并且没有独立的方法。为开发一个理想的独立CITE-Seq QC方法，我们认识到当前CITE-Seq研究中的三个局限性。首先，一些研究仅在RNA水平进行QC，例如在转录组文库大小、转录组技术伪影（如RNA污染）或可能的空滴或环境RNA方面。然而，由于CITE-Seq同时检测RNA和细胞表面蛋白数据，CITE-Seq QC必须评估单个RNA质量，还必须评估蛋白质数据及其与RNA数据交互作用的质量。具体来说，i）必须分别控制单个蛋白质和RNA数据的质量，以忠实识别具有某些表面蛋白的细胞类型并捕获细胞的分子谱；ii）必须研究RNA与蛋白质之间的关系，因为如果某些细胞表达一个易于翻译并运输到表面的特定基因，则表面蛋白丰度水平预计与细胞中的基因表达相关。其次，虽然少数其他研究使用表面蛋白信息进行QC，但它们仅检查了部分检测的表面蛋白，因为它们对由表面蛋白标记的特定细胞类型感兴趣。例如，一项研究检查了数据中188个可用标记中的7个蛋白标记（CD3、CD4、CD8、CD14、CD16、CD19和CD56），以区分五种细胞类型（B细胞、CD4 T细胞、CD8 T细胞、经典单核细胞和自然杀伤细胞）；另一项研究检查了17个可用标记中的四个蛋白标记，以区分四种细胞类型（T细胞、单核细胞、B细胞和细胞毒性T淋巴细胞）。然而，为检测影响数据中大多数检测的系统误差，检查大多数RNA和蛋白质而非一小部分至关重要。第三，当上述研究证明基因与相应蛋白质之间的关系时，它们主要依赖于降维空间（如UMAP）的可视化检查，用于基因和蛋白质的丰度水平关系或其细胞类型特异性。然而，需要定量测量来客观评估丰度水平和细胞类型特异性之间的关系。定量测量有助于进一步比较各种CITE-Seq数据集的数据质量，并使QC分析可扩展。

在本研究中，我们介绍了CITESeQC，这是首个专门设计用于提供全面且可解释的定量指标集合的软件包，以评估CITE-Seq数据的质量。CITESeQC不直接执行细胞或特征的过滤或移除，而是作为一个诊断框架，指导用户根据其数据集做出明智的质量控制（QC）决策。CITESeQC通过提供七个用于单独评估RNA或蛋白质数据的模块，以及五个用于评估跨模态关系（如RNA-蛋白质一致性）的额外模块，支持多层级QC。为确保系统覆盖，这12个模块共同评估数据集中的所有基因和表面蛋白，同时使用个体QC指标标记低质量特征。为进行定量评估，CITESeQC计算香农熵以评估细胞类型特异性表达模式，并计算相关系数以捕获基因表达与蛋白质丰度之间的预期关系。CITESeQC设计为具有广泛可用性，通过一个简化的流程引导用户生成完整的markdown报告，包括信息丰富的可视化和解释，且需要最少的用户干预。这种灵活的、用户引导的方法使研究人员能够以细致和生物学知情的方式评估数据质量——支持标准化工作流程和探索性分析——而不依赖于严格的预定义阈值。

结果

CITESeQC量化CITE-Seq质量的各个方面

CITESeQC提供12个R模块来评估RNA、表面蛋白及其交互作用在多个方面的质量，以及一个R模块来定义细胞簇或导入细胞簇定义。这些模块在可能的情况下还提供定量测量，以测试关于质量的特定假设。

1. 1.

   “RNA_read_corr()”生成一个散点图，将分子/基因数量与转录组中识别的基因数量相关联。由于高质量细胞的 cutoff 将作为函数的参数传递，用户可以根据其数据修改它们。默认值来自Seurat引导的聚类教程。计算斯皮尔曼相关系数以允许用户测试假设：转录组中检测到的基因数量随着基因总数的增加而增加。

2. 2.

   “ADT_read_corr()”生成一个散点图，将检测到的ADT数量与细胞表面识别的ADT分子总数相关联。由于识别高质量细胞的 cutoff 在图上注释并作为函数的参数传递，用户可以根据其数据修改它们。默认值来自Seurat引导的聚类教程。计算斯皮尔曼相关系数以允许用户测试假设：细胞表面检测到的ADT数量随着ADT分子总数的增加而增加。

3. 3.

   “RNA_mt_read_corr()”生成一个散点图，将转录组中识别的基因数量与线粒体基因的百分比相关联。计算斯皮尔曼相关系数以允许用户测试假设：线粒体百分比保持不变，无论识别出的分子数量如何。

4. 4.

   “def_clust()”根据输入的基因表达矩阵定义细胞簇或导入定义。为定义细胞簇，它使用Seurat包和输入的聚类分辨率。对于每个细胞簇，无论是内部定义还是导入，此函数都会识别标记基因以供后续使用。

5. 5.

   “RNA_dist()”可视化输入基因表达在定义或导入的细胞簇中的特异性。为量化和比较，它计算跨簇表达分布的香农熵，定义如下：H_normalized = -1/log₂(N) ∑_i=1ⁿ p_i log₂(p_i)，其中N是簇的数量（字母表大小）。香农熵值越低表示基因跨簇的表达越特异。

6. 6.

   “multiRNA_hist()”是def_clust()中识别的标记基因的香农熵值的直方图。直方图显示标记基因跨簇的特异性。用户可以修改标记基因的数量。直方图在高熵值处峰值表明标记基因缺乏特异性。

7. 7.

   “ADT_dist()”可视化输入ADT丰度在细胞簇中的特异性。具体来说，它计算跨簇表达分布的归一化香农熵。注意，簇是基于基因表达定义的，除非由用户提供。

8. 8.

   “multiADT_hist()”是所有ADT针对细胞簇的归一化香农熵值的直方图。直方图显示ADT标记跨簇的特异性。注意，簇是基于基因表达定义的，除非由用户提供。直方图在高熵值处峰值表明标记基因缺乏特异性。

9. 9.

   “RNA_ADT_read_corr()”生成一个散点图，显示跨细胞中检测到的转录组基因数量与检测到的细胞表面蛋白数量之间的相关性。计算斯皮尔曼相关系数以允许用户测试假设：检测到的蛋白数量随着检测到的基因数量的增加而增加。

10. 10.

    “RNA_ADT_UMAP_corr()”生成UMAP图对和散点图。每个UMAP图对分别针对输入ADT的丰度和相应基因表达绘制。散点图绘制ADT的丰度和输入基因的RNA表达。

11. 11.

    “RNA_ADT_cluster_corr()”是一组散点图，每个为每个细胞簇绘制，显示该簇输入ADT丰度与相应基因表达之间的相关性。

12. 12.

    “RNA_ADT_hist()”是所有ADT与相应基因表达对中相关系数的直方图。

13. 13.

    “RNA_ADT_cluster_hist()”是一组直方图，每个显示每个细胞簇所有ADT与相应基因对中相关系数的分布。

CITESeQC诊断质量指标的解释

我们使用来自健康供体的两个示例CITE-Seq数据集展示了CITESeQC的适用性。第一个包括外周血单核细胞（PBMC），第二个包括脐带血单核细胞（CBMC）。在数据集上，三个以“RNA”或“ADT”开头并以“read_corr”结尾的函数检查跨细胞的总读数与与RNA或蛋白质对齐的读数之间的相关性，使用户能够测试对齐过程是否对质量有贡献。CITESeQC计算斯皮尔曼相关系数和基于排列的p值作为定量测量。我们对PBMC和CBMC数据集的分析证实，有效的对齐应产生正相关。函数RNA_dist()和ADT_dist()使用香农熵计算单个标记基因或表面蛋白跨细胞簇的分布以量化靶标特异性。为说明其效用，我们检查了PBMC中的CCR7和CST7——分别是幼稚T细胞和细胞毒性淋巴细胞的经典标记。尽管两者都是公认的标记，但Seurat的内置模块缺乏区分它们跨簇相对特异性的分辨率。相比之下，我们基于熵的量化提供了清晰、可解释的特异性测量。例如，CCR7的特异性低于CST7（熵值分别为2.53和2.34），使研究人员能够在下游分析（如细胞类型注释、差异表达和实验验证）中优先选择CST7而非CCR7。这种额外的可解释性层代表了相对于现有方法的关键优势。我们还展示了CCR7在CBMC中的特异性以及CD14 ADT在PBMC和CBMC中的特异性。CD14在经典和中间单核细胞中稳健表达，在PBMC和CBMC细胞簇中也显示出强特异性，香农熵值分别为2.39和3.83。“multiRNA_hist()”和“multiADT_hist()”分别可视化标记基因和表面蛋白的香农熵值分布。在我们的分析中，我们使用了每个簇的前10个标记基因以及PBMC和CBMC中识别的所有表面蛋白。此外，三个以“RNA_ADT”开头并以“corr”结尾的函数允许从业者量化RNA与表面蛋白之间的相关性。我们对PBMC上的CD14和CBMC上的CCR7的分析在UMAP上和使用相关性视觉展示了它们跨细胞簇的特异性。最后，两个以“RNA_ADT”开头并以“hist”结尾的函数可视化跨所有簇或每个簇的相关性分布。在CCR7和ADT14上运行这些函数显示了标记的簇特异性行为。在运行需要细胞簇定义的函数（如RNA_dist()）之前，应调用def_clust()来定义或导入它们。

CITESeQC对CITE-Seq数据中技术噪声敏感性的系统评估

为展示CITESeQC如何检测系统误差，我们使用PBMC数据集中随机选择的10%细胞进行了两个受控噪声注入实验。首先，为模拟由特征计数关系中的系统干扰引入的噪声，我们对5%、10%和20%的RNA特征以及10%、20%和30%的ADT特征 shuffling 表达值。我们为ADT数据选择更高百分比以确保尽管其特征集较小（33,538个RNA对17个ADT）也能产生明显效果。对于RNA，每个条件重复10次；对于ADT，重复50次以获得统计显著性和计算效率。为量化噪声效应，CITESeQC计算nFeature（细胞中检测到的唯一基因或蛋白质数量）与nCount（每个细胞的总计数）之间的斯皮尔曼相关性。在高质量数据中，这些指标预计显示强正相关——具有更多检测特征的细胞往往具有更高的总计数。我们的 shuffling 策略是保留细胞水平关系同时破坏基因或蛋白质水平关系。在结果中，我们观察到随着噪声水平的增加，RNA和ADT的相关性值一致下降。RNA模态由于其更多的特征显示出更宽的降解动态范围。这些结果证实CITESeQC基于相关性的指标对全局干扰敏感，并能有效捕获系统质量问题。其次，我们评估了增加的随机性如何影响跨簇的基因/蛋白特异性，这是下游分析的关键步骤。我们分别随机 shuffling 10%、20%和30%的RNA和ADT特征，并使用函数def_clust()定义簇。为效率，我们通过根据基因跨细胞表达的标准差排名并保留变异最高的基因，选择了10,000个RNA特征。使用定义的簇，我们运行multiRNA_hist()和multiADT_hist()函数来计算所有 shuffling 特征的香农熵。在高质量数据中，具有特异性的标记应显示低熵。随着噪声水平的增加，熵值呈现系统增加，整体分布向更高值移动（即右移）。对于RNA特征，我们观察到香农熵从10%到20%以及从20%到30%的显著移动（p值：分别为0.04和0.05），表明簇特异性表达模式的丢失。对于ADT特征发现了类似的移动，尽管不显著（p值：在10%–20%和20%–30%均为0.2），可能由于测量的ADT数量有限（n = 17）。这些发现表明CITESeQC中基于熵的指标有效捕获了由于随机噪声导致的生物信号侵蚀。两个实验共同验证了CITESeQC在多个水平——全局结构和簇特异性——检测质量问题的敏感性，使其成为跨应用和平台的CITE-Seq数据QC的有价值工具。

CITESeQC促进标记特异性分析

为展示CITESeQC的定量测量如何改进下游生物学分析，我们系统性地确定了一个香农熵 cutoff 来评估标记基因的特异性。具体来说，我们专注于定义一个经验阈值，以区分真正的簇特异性标记与背景非特异性基因。为建立此阈值，我们首先随机选择1,000个在PBMC数据集中任何簇中未差异表达的表达RNA（平均表达>5）作为阴性对照。然后我们计算这些非标记基因跨预定义簇的香农熵。因为这些基因预计广泛且非特异性地表达，它们的熵分布反映了非特异性表达的零分布。我们将标记特异性 cutoff 定义为此分布的5th百分位数（即左尾），将低于此阈值的熵值识别为统计上特异的。然后我们应用此经验 cutoff 来评估我们分析中识别的前10、20和30个RNA标记（按差异表达p值排名）。尽管具有更多RNA标记的集合表现出熵值的异质分布，但 cutoff 清晰地将显著特异性标记与非特异性标记区分开来。例如，在PBMC中，低于1.45的熵值被视为特异，分别有26（20%）、39（16%）和41（12%）的前10、20和30个标记符合此标准。在CBMC中， cutoff 为0.75，观察到类似趋势。此分析定量验证了哪些标记真正对每个簇特异。通过基于CITESeQC熵特异性选择簇特异性标记，用户可以增强单细胞数据分析的生物学可解释性和临床实用性。这尤其重要，因为高特异性标记对于稳健的细胞类型分类、生物标志物发现、治疗靶向以及确保跨数据集的可重复性至关重要。

讨论

CITESeQC包是首个在单个RNA、表面蛋白及其交互作用方面评估CITE-Seq数据质量的软件包。为进行定量评估，CITESeQC计算香农熵和RNA-ADT相关系数——两个生物学知情的指标。尽管熵本身旨在量化表达分布而不是技术质量的直接指标，但当应用于标记基因或蛋白质时，它变得对数据质量有信息性。在高质量CITE-Seq数据中，公认的细胞类型标记——如T细胞的CD3或B细胞的CD19——应表现出低熵，表达局限于预期簇。如果这些经典标记反而表现出意外的高熵——即广泛或随机分布的表达式——可能表明技术问题，如导致标记表达渗入无关簇的环境RNA污染、反映转录组信号不足的差聚类分辨率，或ADT层中的抗体非特异性或背景染色。类似地，对于一部分特征明确、高表达的表面标记，在生物学一致且技术健全的CITE-Seq数据中，预计mRNA和蛋白质水平之间存在中度至强正相关。当已知一致的标记表现出意外低或 erratic 相关性时，可能表明技术伪影，如抗体脱落或 mislabeling、液滴条形码或环境标签错误分配，或批效应或样品降解。CITESeQC不使用这些指标施加严格阈值或自动丢弃特征；相反，它将它们作为诊断工具提供，以允许用户区分有意义的生物学异质性和技术噪声。总之，我们为CITE-Seq数据提供了一套全面的计算QC措施，评估和量化数据质量在单个RNA和蛋白质水平及其交互作用中的各个方面。

为使用CITESeQC提供的定量测量确定CITE-Seq数据集的质量，下一步是为每个测量确定适当的 cutoff 值。然而，建立一些 cutoff 值并不简单。例如，将RNA与其相应表面蛋白相关联的测量不仅取决于数据质量，还取决于RNA的翻译效率。即使对于相同质量的数据集，由于转录后调控过程（如选择性多聚腺苷酸化和竞争性内源RNA），翻译效率在生物背景中可能不同。因此，为使用相关性测量评估质量，我们建议将值与用户具有先验数据质量知识的其他CITE-Seq数据集的值进行比较。未来，为执行无参考数据集的QC分析，我们计划收集多个高质量和低质量的CITE-Seq数据集，并直接从数据中确定 cutoff 值。

方法

用户友好R markdown中的CITESeQC

CITESeQC（版本0.9.1）是一个具有最小先决条件的R包，公开可用。它使用基础R包——graphics、stats和utils——使用户易于安装。可下载源代码和带有示例数据集的教程。该工具可在R脚本或R Markdown文件中使用。此设计的优点是它可以在单个文档中集成代码、可视化和解释，这有助于数据分析工作流程的可重复性和文档化。此外，R markdown文件不需要熟悉基于Linux环境的许多软件程序所需的命令行语法。

实验数据

PBMC具有单轮核，包括从外周血分离的淋巴细胞（T细胞、B细胞和NK细胞）和单核细胞。我们从网站下载了数据集，CBMC源自脐带血。它们包括造血干细胞/祖细胞和比成人PBMC更幼稚的免疫细胞，使它们对研究免疫发育有价值。我们从NCBI GEO下载了数据集。