基于自编码器变异聚类的单细胞RNA测序数据癌症亚克隆推断新方法CluVar

《Briefings in Bioinformatics》：CluVar: clustering of variants using autoencoder for inferring cancer subclones from single cell RNA sequencing data

本研究针对单细胞RNA测序(scRNA-seq)数据中高缺失率和噪声干扰的难题，开发了基于自编码器的CluVar框架。该研究通过定制掩码损失函数和贝叶斯高斯混合模型聚类，成功实现了癌症亚克隆的系统发育树重建，在模拟和真实数据中均表现出优越性能，为解析肿瘤进化轨迹提供了新工具。

肿瘤并非均质的细胞群体，而是由遗传背景各异的亚克隆组成的复杂生态系统。这些亚克隆在进化过程中不断积累突变，形成独特的进化轨迹，进而导致肿瘤异质性、药物耐药性和复发风险增加。传统上，研究人员主要通过单细胞DNA测序(scDNA-seq)来解析肿瘤进化树，但该方法通量有限，通常只能分析数百至数千个细胞，难以捕捉肿瘤群体的全貌。相比之下，单细胞RNA测序(scRNA-seq)能够分析更大规模的细胞群体，且公共数据库中资源丰富，为研究肿瘤进化提供了新的机遇。

然而，利用scRNA-seq数据进行突变分析面临严峻挑战。由于mRNA表达量低导致的高缺失率和技术噪声，使得从转录组数据中准确识别遗传变异变得异常困难。尽管已有深度学习技术如自编码器被用于处理scRNA-seq数据的噪声问题，但它们往往未能有效处理缺失值，且通常仅利用编码器的潜在表示，忽略了解码器层蕴含的丰富信息。

为了解决这些瓶颈问题，韩国秀实大学的研究团队在《Briefings in Bioinformatics》上发表了题为"CluVar: clustering of variants using autoencoder for inferring cancer subclones from single cell RNA sequencing data"的研究论文，开发了CluVar这一创新性计算框架。

关键技术方法

研究团队建立了从原始FASTQ文件到系统发育树重建的完整分析流程。关键技术包括：使用CellRanger和minimap2进行序列比对和重映射；采用GATK最佳实践流程进行变异检测；构建细胞-突变矩阵（值定义为：1表示变异存在，0表示参考型，-1表示缺失）。核心创新是专门设计的自编码器架构，其采用掩码损失函数忽略缺失值，并融合潜在向量与解码器第一隐藏层输出进行贝叶斯高斯混合模型(Bayesian Gaussian Mixture Model, BGMM)聚类。最后通过多数投票法确定亚克隆基因型，利用简约法构建系统发育树。分析使用的乳腺癌细胞系数据和胰腺癌患者数据均来自公共数据库。

CluVar框架概述

CluVar的整体设计思路清晰而高效。如图1所示，该框架从10x Genomics的scRNA-seq数据出发，经过变异检测、过滤等预处理步骤，生成基因型矩阵。这个矩阵作为自编码器的输入，通过编码器-解码器结构学习数据的潜在表示。研究团队创新性地将编码器的潜在嵌入与解码器第一隐藏层的输出相结合，形成综合特征向量，再通过BGMM进行聚类分析。每个聚类代表一个独特的癌症亚克隆，最终通过简约算法重建出系统发育树。

CluVar在模拟数据集上的优异表现

为了验证CluVar的可靠性，研究团队进行了全面的模拟实验。结果显示，自编码器训练损失与聚类性能（以调整兰德指数ARI衡量）呈强负相关，表明在缺乏真实标签的情况下，选择低训练损失的模型可以获得可靠的聚类结果。

在聚类算法选择方面，BGMM表现出色，其性能优于标准高斯混合模型(Gaussian Mixture Model, GMM)和k-means聚类（图2b）。更重要的是，研究团队发现将潜在向量z与解码器第一隐藏层输出ψ:1(z)结合，能够在高达5%的错误率下保持较高的聚类精度，显著优于仅使用潜在向量的配置（图2d-f）。这一发现凸显了利用解码器表征的重要性，为后续真实数据分析奠定了基础。

乳腺癌细胞系中揭示主干-分支亚克隆结构

将CluVar应用于ER阳性乳腺癌细胞系数据集（包含药物敏感群体和经ribociclib长期暴露产生的耐药群体），研究团队成功识别出7个亚克隆，呈现出清晰的树干-分支架构（图3a）。亚克隆0和1突变负荷较低，位于系统发育树的根部，构成主干；而亚克隆2-6具有较高的突变数量，形成分支。

转录组分析显示，CluVar定义的亚克隆在UMAP图上呈现空间分离的非随机分布（图3c）。敏感细胞主要由主干亚克隆主导，而耐药细胞则富含具有更高突变负荷的分支亚克隆（图3d）。Monocle3伪时间分析再现了这种连续排列，强调了突变层次与转录动态之间的联系（图3e）。

通路富集分析揭示了生物学程序的时间性转变（图3f，g）。早期状态富含雌激素反应基因和涉及一氧化氮/cGMP信号传导的通路，而晚期状态则显示出侵袭性肿瘤特征的显著上调，如缺氧、上皮-间质转化(Epithelial-Mesenchymal Transition, EMT)以及细胞周期程序。值得注意的是，亚克隆1、4和6在耐药细胞中特异性富集，并逐步表现出翻译、核糖体和MYC靶标通路活性的增强，这表明肿瘤细胞通过激活CDK2介导的MYC通路来规避CDK4/6抑制剂的治疗作用。

胰腺癌组织样本中的亚克隆进化

在更具异质性的胰腺癌样本中，CluVar再次展现出强大能力。对三名患者的原发肿瘤和肝转移标本分析后，研究团队识别出7个亚克隆，形成明显的树干-分支架构（图4a）。亚克隆0和1构成主干，而其余亚克隆分裂为两个主要谱系：一个包含亚克隆2、4和5，另一个包含亚克隆3和6。

特别值得注意的是，亚克隆6主要局限于患者2，并在转移性病变中富集（图4d）。通路富集分析显示，早期状态（簇0）富含DNA损伤反应通路，表现出p53通路和UV反应特征的显著激活；中间状态（簇1-4）呈现应激适应/炎症特征，表现为HSF1激活、KRAS、IL-17信号传导和缺氧通路的显著上调；晚期状态则展现出两种不同的表型：簇5表现出缺氧、金属解毒、EMT和金黄色葡萄球菌相关炎症通路的强烈富集，表明其处于高度侵袭性、转移性和应激适应状态；而簇6保留上皮连接基因表达，同时启动强大的I/II型干扰素和抗病毒反应，反映其免疫参与但非侵袭性的表型特征。

与现有方法的比较

研究团队将CluVar与SCClone和bmVAE等现有方法进行了全面比较。在三种模拟场景下（10000个细胞5个簇、10000个细胞7个簇、5000个细胞7个簇），CluVar在聚类准确性（ARI）和系统发育树重建（Robinson-Foulds距离）方面均表现出优越性能（图5a-f）。特别是在错误率增加的情况下，CluVar的稳健性更加明显。此外，CluVar的运行时间也显著短于比较方法，显示出良好的可扩展性（图5g）。

研究结论与意义

CluVar通过整合变异信息和转录组谱，为癌症细胞的功能亚分型提供了创新性解决方案。该研究的核心创新在于：针对scRNA-seq数据中的缺失值问题设计了定制掩码损失函数；发现并利用了解码器隐藏层输出的信息价值；将突变背景与转录状态相结合，全面解析肿瘤异质性。

在乳腺癌细胞系中，CluVar揭示了EIF2S2和NDUFS8等基因的非同义外显子突变如何区分不同的代谢和增殖表型。在胰腺癌中，该技术识别出两个分叉的突变谱系，分别与不同的表型特征相关，其中GSTP1、LGALS1、GSDMD、NEAT1和NCL等基因的突变在免疫调节和应激反应中发挥关键作用。

尽管CluVar存在一些局限性，如变异数量受限于模型稳定性考虑、对高缺失率数据的敏感性以及scRNA-seq数据低读深带来的突变检测准确性挑战，但其通过将稀疏突变数据嵌入可扩展的深度学习框架，为解析遗传改变和转录程序如何共同塑造癌症进化提供了强大工具。

该研究的实际价值还体现在其可重复性和易用性上。研究团队将整个工作流程实现为基于snakemake的管道，利用singularity容器技术，使从预处理到聚类的端到端执行变得简单高效，无需手动集成多个工具，极大促进了方法的推广和应用。

总体而言，CluVar代表了单细胞癌症进化研究领域的重要进展，不仅提供了强大的分析工具，更深化了我们对肿瘤异质性和进化动力学的理解，为未来开发针对特定亚克隆的精准治疗策略奠定了坚实基础。