基于单细胞RNA测序数据的系统发育树推断：SCITE-RNA方法

《Genome Biology》：Phylogenetic tree inference from single-cell RNA sequencing data with SCITE-RNA

　　

研究人员提出SCITE-RNA，一种专为单细胞RNA测序（scRNA-seq）数据设计的系统发育树推断新方法，其输入为单核苷酸变异（Single-Nucleotide Variants, SNV）的参考等位基因与突变等位基因读数。该方法采用最大似然随机扫描贪婪搜索策略，在细胞谱系树与突变树两种表征之间交替迭代，以逃离局部最优解，直至在两种空间中均达到收敛。研究人员在模拟数据上验证了该方法的优越性能，与现有方法相比表现更优。此外，研究展示了该方法在癌症单细胞RNA测序数据中的适用性，能够将细胞的进化轨迹与其基因表达谱相联系。

基于单细胞RNA测序数据的系统发育树推断：SCITE-RNA方法

近年来，测序技术的进步推动了单细胞测序数据的可获得性大幅增加，同时也促进了新型计算工具的开发。与传统方法如 bulk 测序（批量测序）相比，单细胞测序技术为分析细胞群体的异质性提供了更高的分辨率。在单细胞RNA测序（scRNA-seq）中，通过逆转录获得的扩增互补DNA（cDNA）被用于测序。由于转录组由基因组转录而来，scRNA-seq 能够间接提供关于表达基因组区域的遗传信息，同时揭示细胞的基因表达谱特征，因此在原则上也可用于推断细胞间的进化关系。系统发育推断旨在从可遗传特征（如遗传变异）中重建这些进化历史。

目前，针对单细胞DNA测序数据已开发了多种系统发育树推断方法。然而，由于成本高昂，单细胞DNA测序的应用仍受限。更为常用的是批量DNA测序和 scRNA-seq。scRNA-seq 技术包括全长度转录本测序方案（如 Smart-seq2）和高通量3'端测序方法（如10x Genomics平台）。全长度方法通常在基因检测灵敏度、较低的技术噪音、较高覆盖度以及较少缺失事件方面优于3'端方案，后者 additionally 存在3'端覆盖偏倚问题。尽管如此，3'端方法因具有更高的细胞通量，在检测稀有细胞群体方面具有优势。

现有的基于 scRNA-seq 的系统发育树重建方法大多依赖单核苷酸变异（SNV），因此全长度技术（如 Smart-seq2）相比3'端方法具有上述优势。现有方法包括 DENDRO、SClineager、PhylinSic、Canopy2 和 PhylEx 等。DENDRO 和 SClineager 均使用转录SNV的突变与参考等位基因读数作为输入：DENDRO 计算细胞间的成对遗传分歧矩阵并基于此聚类克隆；SClineager 则通过贝叶斯层次模型与马尔可夫链蒙特卡洛（Markov Chain Monte Carlo, MCMC）采样推断真实变异等位基因频率，再采用层次聚类构建克隆树。PhylinSic 在树推断前先进行SNV calling 并从突变位点读数推断核苷酸序列，随后使用 BEAST2 模型推断细胞谱系树。PhylEx 和 Canopy2 则整合批量DNA数据与 scRNA-seq 数据重建克隆树，分别采用树结构 stick breaking 过程和贝叶斯层次模型配合MCMC采样。

现有方法存在若干局限性：SClineager 和 DENDRO 依赖将细胞聚类为有限数量的克隆，可能缺乏分辨率且无法准确反映真实系统发育结构；PhylinSic 在系统发育重建前从参考和突变读数中 calling 基因型，未能充分利用读数信息进行树推断；PhylEx 和 Canopy2 对批量DNA数据的依赖在仅有 scRNA-seq 数据时构成限制。此外，树优化方法常见的局部最优陷阱问题在 PhylinSic、PhylEx 和 Canopy2 中同样存在。

为此，研究人员开发了 Single-Cell Inference of Tumor Evolution from RNA sequencing data（SCITE-RNA），一种针对 scRNA-seq 数据的树推断方法。该方法不依赖额外数据源（如批量DNA测序），不预先假设克隆数量，专门解决局部最优陷阱问题，并能直接从数据中估计关键参数（如缺失概率）。研究人员开展了模拟数据验证和真实癌症数据集应用，证明其在克隆数量众多的困难重建问题上的优越性，并发表于《Genome Biology》杂志。

研究采用的主要关键技术方法包括：基于SNV参考与突变等位基因读数的输入处理；最大似然随机扫描贪婪搜索优化策略；细胞谱系树与突变树空间交替切换以避免局部最优；基于Beta-二项分布的基因型似然建模，包含等位基因缺失（Allelic Dropout）的概率模型；SNV特异性参数优化（采用有限内存BFGS算法）；以及自助法（Bootstrapping）生成一致性树以评估不确定性。样本队列来源包括：多发性骨髓瘤患者MM16（46细胞，含化疗前后样本）、MM34（127细胞，含原发肿瘤和转移灶样本）以及胶质母细胞瘤数据集BT_S2（1,169细胞，来源于肿瘤核心和瘤周组织）。

方法概述方面，SCITE-RNA 计算每个细胞-突变对基于观测读数、推断基因型和模型参数（如缺失率）的似然值，进而获得任意给定突变树或细胞谱系树的似然值。细胞谱系树中节点代表单个细胞，二叉树结构捕捉细胞间的进化关系；突变树中节点对应突变，编码突变发生顺序。该方法采用贪婪优化策略，通过随机剪枝和重附着操作迭代探索树空间，每次优化后将树转换为另一种表征（细胞谱系树到突变树或反之），直至在两种空间中似然值均达到最优。随后基于推断基因型进一步优化缺失率和过度离散参数，提升模型对RNA特异性现象（如单等位基因表达和不平衡双等位基因表达）的拟合能力，并进行第二轮树推断。

模拟数据性能评估方面，研究首先比较了不同树空间优化策略。通过设置三种细胞和SNV配置（更多SNV：100细胞/500 SNV；平衡：500细胞/500 SNV；更多细胞：500细胞/100 SNV），比较仅优化细胞谱系树空间（c）、仅优化突变树空间（m）、以及从细胞树空间（cm）或突变树空间（mc）开始交替切换的策略。结果显示，除在SNV数量较大时仅优化突变树空间表现不佳外，SCITE-RNA在所有策略下均能找到高似然树。当细胞数与SNV数相近时，空间切换的优势最为明显；空间切换策略在平均多次运行后表现至少不劣于、且通常优于任一单独空间优化的最佳结果。

与现有方法的比较基于50细胞和500 SNV、克隆数可变的模拟数据集，采用预测变异等位基因频率的平均绝对误差（MAE VAF）、路径差异（Path Difference）和标准化 Robinson-Foulds 距离作为评估指标。结果表明，SCITE-RNA在所有设置中明显且持续优于 DENDRO、SClineager 和 PhylinSic，尤其在克隆数较高的数据集上树重建优势更加显著。对于较小克隆数，由于SCITE-RNA假设SNV/细胞独立放置于树中，该假设在克隆数较少时不成立，结合数据噪音会导致过拟合。通过层次聚类或自助法一致性树可改善低克隆复杂度下的性能。在包含更多突变事件（如趋同进化和平行进 copy数变异）的模拟中，SCITE-RNA 在更高噪音、更低覆盖度和更多额外突变事件下表现略有下降，但始终优于替代方法。在高达百万细胞-突变条目的更大规模数据集中，SCITE-RNA 保持优越性能，且运行时间在大规模数据上仅次于 DENDRO，远快于基于MCMC的方法。

真实癌症数据应用方面，研究应用于两个多发性骨髓瘤数据集（MM16和MM34）和一个胶质母细胞瘤数据集（BT_S2）。在MM16数据集中，推断树显示化疗前后细胞部分分离，部分治疗后细胞与原发肿瘤细胞关系更近，与原始研究基于拷贝数变异鉴定的五个存活肿瘤细胞一致；而 DENDRO 将样本分为几乎完全不同的两支，SClineager 和 PhylinSic 则未能反映这五个细胞的预期位置。在MM34数据集中，方法清晰识别原发与转移灶肿瘤细胞的分裂，支持线性进化模型，与原始研究一致，调整兰德指数（Adjusted Rand Index, ARI）达0.938（一致性树为0.969）。在BT_S2胶质母细胞瘤数据集中，肿瘤细胞主要形成独立分支（含85%肿瘤细胞，纯度93%），而健康细胞占剩余部分的95%（纯度89%）；尽管未能区分具体健康细胞类型，但在区分肿瘤与健康细胞的二元分类中 ARI 达0.615，优于其他方法。此外，研究将基因表达映射到推断的系统发育结构上：MM16中 RGS1 基因在治疗前相关分支中高表达；MM34中 CTSB 基因在转移分支中显著高表达；BT_S2中 EGFR 和 SOX9 在肿瘤分支上调、CD86 下调，均与已知生物学特征一致。

讨论部分，研究人员指出 SCITE-RNA 通过空间切换有效解决了局部最优问题，相比聚焦有限克隆数的 DENDRO 和 SClineager 能捕捉更细粒度的生物学异质性。虽然 PhylinSic 表现略优于后两者，但计算更慢且总体仍不及 SCITE-RNA。在真实数据中，推断的系统发育更能反映不同样本来源细胞间的合理进化关系。不过，过拟合仍是关注点，可通过惩罚边数或自助法缓解，但后者在大数据集上因共享克隆稀少而受限。研究还指出，当前模型未整合拷贝数变异或基因表达进行树优化，是未来发展方向。模型假设每位点仅两个拷贝、两种等位基因、无限位点等，虽提升计算效率但也限制了对更复杂遗传变异的处理能力。

结论部分，研究人员总结道：通过在树空间之间交替优化是高效的策略，可扩展至其他组学模态（如 scDNA-seq）。快速优化与数据驱动的参数估计相结合，使得 scRNA-seq 数据不仅可用于基因表达分析，还能用于重建细胞的进化历史，从而将基于表达的分析与系统发育信息相结合，为肿瘤进化研究提供新的视角。