Nature Methods |单细胞转录组的深度生成建模

作者 | 戴迟迟 编辑 | 戴迟迟 校对 | 李仲深

今天给大家介绍加利福尼亚大学的Nir Yosef教授等人发表在Nature Methods上的一篇文章 “Deep generative modeling for single-cell transcriptomics” 。单细胞转录组测量可以揭示未开发的生物多样性，但它们受到技术噪音和偏差的影响，必须建模以解释下游分析中产生的不确定性。本文介绍了single-cell variational inference (scVI)，一个现成的可扩展框架，用于概率表示和分析单细胞中的基因表达。scVI使用随机优化和深度神经网络来聚合相似细胞和基因的信息，并近似观察到的表达值的分布，同时考虑批次效应和有限的灵敏度。本文将scVI用于一系列基本的分析任务，包括批处理校正、可视化、聚类和差异性表达，并为每个任务实现了较高的精度。

一、研究背景

单细胞RNA测序 (scRNA-seq) 是一个强大的工具，它开始对不同的研究领域做出重要贡献，如自身免疫和癌症。然而，解读scRNA-seq数据仍然具有挑战性，因为它受到一些干扰因素的影响，如灵敏度有限且可变、批次效应和转录噪声。最近的几项研究通过为每个细胞中的每个基因测量拟合一个概率模型，模拟了scRNA-seq偏差和不确定性，这代表了一个较低的和潜在的较少噪声维度的数据。一旦这些模型被拟合，它们就可以用于各种任务，如聚类、插补和差异性表达分析。

本文开发了scVI，一种用于归一化和分析scRNA-seq数据的全概率方法。scVI基于层次贝叶斯模型，该模型具有由深度神经网络指定的条件分布，即使对于非常大的数据集也可以非常有效地训练。每个细胞的转录组是通过非线性转化为正常随机变量的低维潜在向量进行编码的。这个潜在的表示随后被另一个非线性转换解码，以产生每个细胞中每个基因分布参数的后验估计。这种变换被假定为零膨胀的负二项分布。为了演示它的灵活性，我们进行批次删除、正常化、降维、聚类、和差异性表达。我们在这里展示，对于这些任务中的每一个，scVI都优于当前最先进的方法。

二、模型与方法

本文将每个细胞n中每个基因g的观察表达xng建模为一个样本，该样本取自一个零膨胀的负二项 (ZINB) 分布

，该分布以每个细胞的批次注释sn以及另外两个未观察到的随机变量为条件。第一个变量ln是一维高斯，它表示由于捕获效率和测序深度的差异而厌变，并且它被看做特异性比例因子。第二个变量zn是一个高斯的低维向量 (这里设为10维) 表示更能反映细胞间的生物学差异的保留的变异。本文使用它将每个细胞表示为在用于可视化和聚类的低维潜在空间中的一个点。在scVI模型中，神经网络将潜在变量映射到ZINB分布的参数 (图1，神经网络5和6) 。这种映射通过中间值­

，它为在每个细胞n中源自每个基因g的转录的百分比提供批次修正的且规范化的估计。这些估计用于差异性表达分析及其扩展版本 (­

乘以估计的库大小ln­) 来插补。最后通过使用变分推理和一个可扩展的随机优化过程训练一个神经网络并得到了一个潜在变量的后验分布

的近似 (图1a，神经网络1 - 4)。

三、实验结果

3.1 模型评估

本实验使用一套基准方法对scVI进行了评估，这些基准方法用于概率建模和使用了一系列已发布的数据集对scRNA-seq数据进行插补。

首先，为了评估训练的可扩展性，随机抽取了包含130万个小鼠脑细胞的数据集 (BRAIN-LARGE)。为了便于与目前最先进的算法进行比较，对scRNA-seq数据进行概率建模和降维，本实验将分析限制在所有细胞中s.d. 最大的720个基因 (图1b) 。结果表明，大多数方法在耗尽内存 (使用32 GB RAM) 之前能够处理多达50,000个细胞。相比之下，scVI通常更快，并扩大到100万个细胞，这要归功于它在迭代随机优化的每次迭代中依赖固定数量的细胞。除此之外，在DCA上观察到类似的可伸缩性，DCA是一种去噪的自编码器，也使用了随机优化。值得注意的是，当数据集大小接近100万个单元时，训练迭代的次数更少，因此“停止启发式”的策略更节省时间。的确，标准scVI (使用固定数量的迭代次数) 比DCA (默认使用“停止启发式”策略) 要慢，但scVI的早期停止选项大大提高了速度 (在1小时内训练)，而不影响数据拟合 (图2)。

图1. scVI概述以及不同方法运行时间随数据集大小变化

图2. 在BRAIN-LARGE数据集的100万个细胞样本上“停止启发式”策略的稳定性

其次，本实验通过评估它们精确地插补缺失值的能力来评估这些方法对数据的拟合程度。在五个不同大小的数据集 (BRAIN-LARGE，CORTEX，PBMC，RETINA，HEMATO；3 - 27000细胞)，本实验将9%的非零项设置为零，然后测试每个方法恢复的这些值的能力。在大多数情况下，基于ZINB分布的方法——即scVI、DCA和ZINB-wave——比使用替代分布的方法表现更好 (例如，ZIFA中的log normal)。此外，在HEMATO造血分化数据集上，scVI的表现被MAGIC (通过在细胞-细胞相似性图中传播) ，该数据集包含的细胞 (4,016) 少于基因 (7,397) 。在这种情况下，scVI预计会对数据拟合不足，可能导致更差的插补精度。图3和4显示随机选择9%的非零项的结果，图5和6显示按照偏好低值的原则选择9%的非零项的结果。

作为模型拟合的额外评估，本实验还测试了在训练期间提供的数据的似然值，并获得了与综合“dropout”测试一致的结果，详细结果被作者放置在文章的补充材料中。

图3. scVI插补结果1

图4. scVI插补结果2

图5. scVI插补结果3

图6. scVI插补结果4

3.2 捕捉潜在空间中的生物结构

本实验评估了scVI推断的潜在空间在多大程度上反映了细胞之间的生物可变性。评估这一点的一种方法是依赖于预先将细胞分成具有生物学意义的亚群，这通常是通过无监督的聚类和人工检查以及标注来完成的。通过在潜在空间应用k - means聚类并测试带标签的亚种群与其的一致性 (使用与带注释的数据相同的k，即亚群数)，或者通过比较同一族群细胞的邻近性以及不同族群细胞的邻近性来评估精度 (CORTEX和PBMC数据集)。使用一个除了mRNA (CBMC) 之外，还包括单细胞蛋白质的测量数据集作为基准，来评估mRNA潜伏空间中细胞之间的相似性与它们在蛋白质水平上的相似性的一致程度。

本实验使用SIMLR对scVI进行了基准测试，SIMLR是一种将聚类与细胞-细胞相似性矩阵和各自的低维潜在表示结合起来的方法，且性能优于scVI。在蛋白质与mRNA的对比测试中，scVI和DCA表现最好，尽管差距不大 (图7c)。scVI还能更准确地捕捉细胞亚群之间的层次结构，如在CORTEX数据集中。来自相关亚群的29个细胞在scVI的潜在空间中彼此更接近 (图7e-g)。另一个重要的是变异是连续的，而不是离散的，在HEMATO数据集中，SIMLR识别了几个离散的簇，并没有像scVI或PCA那样反映这个系统的连续性 (图8)。其他更加详细的结果作者放在了文章的补充材料中。综上所述， scVI的潜在空间是灵活的，即使数据不适合离散细胞状态的简单结构，也能很好地描述数据。

图7. 利用scVI潜在空间聚类，并针对最新的方法对数据集进行基准测试

图8. 生物信号被保留在scVI潜伏空间中

3.3 解释技术可变性

为了评估scVI纠正批次效应的能力，本实验使用了两批小鼠视网膜双极神经元数据集 (RETINA)。定义了一个熵来评估来自不同批次的细胞在潜在空间的任何局部邻域的混合 (使用k-means)。在本数据集中，scVI比ComBat和最近一种基于相互最近邻匹配的方法更好地解决了批次效应，同时仍然保持预先标注的子群体的聚集性 (图8)。尽管SIMLR和DCA能够很好地聚类每个批次中的细胞，但每个批次中的各自的聚类仍然很大程度上显示分离状态。当本实验应用一个没有批次变量的简化版本的scVI时，得到了类似的结果，从而支持了建模选择。

本文发现，在相对同质的群体中，有scVI推断的库大小与每个细胞的观测深度密切相关 (例如在PBMC中)。但是由于有限的mRNA捕获效率和测序深度，造成了灵敏度较低，这增多了零条目的数量，并可能扭曲同质细胞之间的相似性。此外，本实验发现，大多数零项可以用负二项分量 (图9a, b) 来解释。与其一致的是，零的产生与在每个细胞中随机采样基因的过程一致，该过程与期望频率成正比，并且没有额外的偏差 (图10b)。本文发现，负二项分布的概率与库大小相关的细胞特异性质量因子更加相关，而来自Bernoulli的零概率与显示校准误差的质量因子更加相关 (图9c, d; 图10c, d)。综上所述，这些结果证实了，至少在本文探索的数据集中，大多数零值可以解释为各自转录本的低 (或零)“生物学”丰度，并因有限的采样，这种情况被加剧。

图9. scVI的生成分布

图10. 使用scVI捕获技术可变性

3.4 差异性表达

为了评估scVI与其他差异性表达分析方法的性能，本实验使用了一个来自健康人体供体的拥有12039个PBMCs的数据集 (PBMC)，并对B细胞和树突状细胞簇、CD4+和CD8+ T细胞簇进行了两组的比较。将批量测序级别的结果作为评估标准。为了进行评估，本实验首先将基因定义为真阳性 (在批量测序数据中 BH-adjusted P value < 0.05)，然后根据scVI的贝叶斯因子或其他方法的BH-adjusted P value计算ROC曲线下面积 (AUROC) 。由于真阳性的定义需要一个阈值，本实验还使用了第二个评分来评估基因排序的再现性，使用不可再现发现率 (IDR)。scVI在T细胞比较中具有最高的AUROC，而edgeR在B细胞与树突状细胞比较中优于scVI。在两种比较中，scVI在IDR方面表现最好 (图11; 图12a-e)。与直接应用DESeq2相比，使用DCA后再使用DESeq2性能提升较大，后者是为批量测序数据设计的，因此支持了对单细胞适应模型的需求。

图11. 差异性表达分析的基准评估

图12. 在PBMC数据集上使用scVI进行差异性表达分析

四、总结

因为提供了基因表达的一般概率表示，scVI可以实现本研究未探讨的其他scRNA-seq分析形式，如谱系推断和细胞状态标注等。此外，因为它只需要潜在空间和模型规范来生成任何感兴趣的数据点 (细胞×基因)，scVI可以用作可扩展和交互式可视化工具的有效基线。最后，scVI可以扩展到来自给定组织的多个数据集的并集，同时整合先前的细胞类型的生物学标签。最后，除了上述描述的实验之外，本文的作者还将许多详细的实验结果放置到了补充材料中，如果读者感兴趣，可以自行搜索scVI的原文进行阅读。

代码

https://github.com/YosefLab/scVI

参考文献

Romain, Lopez, Jeffrey, et al. Deep generative modeling for single-cell transcriptomics[J]. Nature Methods, 2018.