一种统一细胞类型和状态特征的细胞相似性度量方法

摘要

本申请涉及一种统一细胞类型和状态特征的细胞相似性度量方法，它包括如下步骤：S1：选择特征基因；S2：计算细胞类型差异；S3：计算细胞状态差异；S4：计算细胞间的距离；S5：根据步骤S4中得到的细胞间的距离构建距离矩阵W

说明书

技术领域

本申请涉及系统生物学领域，具体涉及一种统一细胞类型和状态特征的细胞相似性度量方法。

背景技术

细胞是生物学中构成生物体的基本单位，经常被称为生命的积木。多细胞物种包括大量具有不同功能的细胞类型，将细胞分为不同类型为理解组织的结构和功能铺平了道路。随着生物技术的进步，来自scRNA-seq的数据是发现新的假定细胞类型和细胞状态的最新信息来源之一。与传统的bulk测序相比，scRNA-seq分析能够在单细胞分辨率快速确定成千上万个单个细胞的精确基因表达模式，揭示细胞身份和功能，为我们理解细胞的生态提供了新的视角，一个典型应用便是研究肿瘤内异质性。肿瘤组织包括肿瘤细胞与肿瘤基质细胞(包括基质成纤维细胞、内皮细胞和免疫细胞)，它们共同组成了肿瘤微环境。肿瘤微环境的动态相互作用对于刺激肿瘤细胞的异质性、克隆进化至关重要。而肿瘤内异质性是肿瘤治疗和患者生存率的关键决定因素，这种异质性包括时间异质性(Temporalheterogeneity)和空间异质性(Spatial heterogeneity)。肿瘤的时间异质性主要表现为肿瘤在自然或者人为干预下随时间演化的多克隆特性，表现出明显的动态性。肿瘤的空间异质性是指在原发肿瘤内或原发肿瘤与转移瘤之间，遗传信息和细胞形态等特征会存在差异。肿瘤内异质性是肿瘤耐药性形成的主要原因并严重影响了个体化肿瘤治疗的效果和患者的预后。基于单细胞转录组数据的聚类提供了一种从异质系统中区分和识别肿瘤细胞及其状态的直观方法，为揭示和理解肿瘤内异质性铺平了道路。

尽管在scRNA-seq数据中对细胞类型的识别已经付出了巨大的努力，但传统算法的一个关键限制是不能明确定义或区分不同细胞类型和细胞状态，它们通常都将这两个概念视为“同义词”，即它们都是具有某种基因表达模式的、具有生物学意义的细胞群。

事实上，细胞的类型和状态所涉及的基因的表达模式并不相同。通常，控制细胞类型的谱系特征基因具有on/off的表达模式，其表达量大小变化不大。而决定细胞状态的基因，表达量会随着状态的改变呈梯度变化。因此，有必要设计一种能够区分两种表达模式的相似性度量方法。

发明内容

本发明的目的在于，提供一种统一细胞类型和状态特征的细胞相似性度量方法，能够避免表达矩阵中存在大量零表达的影响，有效减少单细胞数据固有的噪声的影响，可应用于单细胞数据集的细胞鉴定和衡量肿瘤内异质性分析。

本发明采取的技术方案是：一种统一细胞类型和状态特征的细胞相似性度量方法，包括如下步骤：

S1：使用平均-离散图选择特征基因，其具体过程为：输入标准化后的基因表达矩阵，计算每个基因在每个细胞中的平均表达和离散度，并将基因分到nbin个箱子中；每个箱子内，计算所有基因的z-scores分数，选择z-scores大于1的基因为特征基因；

S2：利用S1得到的特征基因，计算细胞类型差异，对于细胞i和细胞j，细胞类型差异的计算公式为：

其中，

S3：利用S1得到的特征基因，计算细胞状态差异，对于细胞i和细胞j，细胞状态差异的计算公式为：

其中，‖x

S4：计算细胞间的距离，对于细胞i和细胞j，细胞之间的距离计算公式为：

其中，α为平衡系数；

S5：将平衡系数α以固定步长，从[0,1]的范围内，根据S4得到的细胞间的距离生成数个距离矩阵，所述距离矩阵为对角矩阵；利用最近邻方法构建每个距离矩阵对应的邻接矩阵W

其中，w

S6：使用谱聚类，输入邻接矩阵W

W

进一步地，还包括步骤S7：对两种肿瘤亚型中包含的数十位患者分别计算恶性细胞细胞间距离的四分位距，并作为异质性评分ITH

进一步地，步骤S1中获取标准化后的基因表达矩阵的具体方法为：

S101：从GEO数据库，即Gene Expression Omnibus数据库中下载基因表达数据；

S102：使用ComBat-seq算法移除不同患者间的批次效应，其使用负二项回归模型来估计计数矩阵中批次效应，可保留计数数据的整数性质；

S103：使用Seurat提供的LogNormalize方法标准化表达矩阵，消除细胞文库大小的差异性。

进一步地，步骤S1中的基因分类方法为：计算每个基因在每个细胞中的平均表达和离散度，根据平均表达将基因平均分为nbin个箱子，nbin的取值为两倍目标类簇数量。

本发明的有益效果在于：本发明有效挖掘了基因表达模式的内在特性，将细胞类型的差异和细胞状态的差异有效的结合在一起；与现有的技术相比，能够避免表达矩阵中存在大量零表达的影响，有效减少单细胞数据固有的噪声的影响；在肿瘤数据上，可以根据计算出的距离矩阵，衡量不同的肿瘤组织或样本之间的差异是更多的体现在细胞类型还是细胞状态上，为生物学家进行生物学研究提供可靠依据。

附图说明

为了更清楚地说明本发明实施例中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其它的附图。

图1为本发明实施例采用的UCRSI方法与六种现有的聚类方法在10个不同的数据集上进行基准测试的比较图；

图2为Tirosh数据集与原始UMAP的可视化结果对比；

图3为Grubman数据集与原始UMAP的可视化结果对比；

图4为LUAD和LUSC病人

图5为LUAD和LUSC病人

具体实施方式

为了能够更清楚地理解本发明的上述目的、特征和优点，下面结合附图和具体实施方式对本发明进行进一步的详细描述。在下面的描述中阐述了很多具体细节以便于充分理解本发明，但是，本发明还可以采用其他不同于在此描述的其他方式来实施，因此，本发明并不限于下面公开的具体实施例的限制。

除非另作定义，此处使用的技术术语或者科学术语应当为本申请所述领域内具有一般技能的人士所理解的通常意义。本专利申请说明书以及权利要求书中使用的“第一”、“第二”以及类似的词语并不表示任何顺序、数量或者重要性，而只是用来区分不同的组成部分。同样，“一个”或者“一”等类似词语也不表示数量限制，而是表示存在至少一个。“连接”或者“相连”等类似的词语并非限定于物理的或者机械的连接，而是可以包括电性的连接，不管是直接的还是间接的。“上”、“下”、“左”、“右”等仅用于表示相对位置关系，当被描述对象的绝对位置改变后，则该相对位置关系也相应地改变。

本发明实施例从GEO数据库，即Gene Expression Omnibus数据库中下载基因表达数据包GSE148071，GSE148071基因表达数据包含42位患者的88197个细胞，29527个基因，肿瘤细胞类型包括LUAD、LUSC和NSCLC。

本发明实施例采用了一种统一细胞类型和状态特征的细胞相似性度量方法对GSE148071基因表达数据包进行分析，具体包括如下步骤：

S1：使用平均-离散图选择特征基因，具体方法为：

S101：从GEO数据库，即Gene Expression Omnibus数据库中下载基因表达数据；

S102：使用ComBat-seq算法移除不同患者间的批次效应，其使用负二项回归模型来估计计数矩阵中批次效应，可保留计数数据的整数性质；

S103：使用Seurat提供的LogNormalize方法标准化表达矩阵，消除细胞文库大小的差异性。

输入标准化后的基因表达矩阵，计算每个基因在每个细胞中的平均表达和离散度，基于平均表达将基因分到nbin个箱子中；每个箱子内，计算所有基因的z-scores分数，选择z-scores大于1的基因为特征基因。从而实现在选取特征基因的同时，控制基因表达量的可变性和平均表达之间的关系，在本发明实施例中，nbin的取值为两倍目标类簇数量。

S2：利用S1得到的特征基因，计算细胞类型差异，对于细胞i和细胞j，细胞类型差异的计算公式为：

其中，

S3：利用S1得到的特征基因，计算细胞状态差异，对于细胞i和细胞j，细胞状态差异的计算公式为：

其中，‖x

S4：计算细胞间的距离，对于细胞i和细胞j，细胞之间的距离计算公式为：

其中，α为平衡系数。细胞之间的距离distance

S5：在单细胞数据分析流程中，识别细胞类型是一个非常关键的步骤，而聚类是识别细胞类型的关键；对此，利用S4得到的细胞间距离进行聚类分析，将相似的细胞分为一组，为随后的类型识别打下基础，具体方法为：将平衡系数α以固定步长，从[0,1]的范围内，根据S4得到的细胞间的距离生成数个距离矩阵，所述距离矩阵为对角矩阵，每个距离矩阵代表了在某个α设定值下成对细胞间距离；利用最近邻方法构建每个距离矩阵对应的邻接矩阵W

其中，w

S6：使用谱聚类，输入邻接矩阵W

W

在单细胞数据分析流程中，识别细胞类型是一个非常关键的步骤，而聚类是识别细胞类型的关键。本发明实施例将

为了衡量肿瘤内异质性ITH，即Intra-tumor heterogeneity，本发明实施例还包括步骤S7：对存在两种亚型的肿瘤细胞，每种亚型中包含的数十位患者，分别计算每位患者恶性细胞细胞间距离的四分位距，并作为异质性评分ITH

图1显示了本发明实施例采用的UCRSI方法与六种现有的聚类方法在10个数据集上的基准测试结果，评价指标为NMI、AMI和ARI。图1中每个数据集结果比较图中的柱体从左向右依次为GiniClust3方法、K-means方法、NMF方法、SC3方法、SCCAF方法、Seurat方法和本发明实施例采用的UCRSI方法。图2和图3显示了其中两个数据集Tirosh和Grubman的可视化结果。本发明实施例采用的UCRSI方法在3个数据集上与其它方法有相似的性能，在4个数据集上则明显优于其他方法。具体来说，对于简单的数据集，如GSE138852数据集，含有13214个细胞，6种细胞类型，包括Microglia、astrocytes、neurons、oligodendrocyteprogenitor cells、oligodendrocytes和endothelial cells。从图1很明显可以看出，Seurat方法、SCCAF方法和本发明实施例采用的UCRSI方法的结果明显优于其它方法，三者在三个评价指标上都大于0.9并且差异是细微的。对于异质性较强的数据集，如小鼠(Mouse)数据，包含13种细胞类型，本发明实施例采用的UCRSI方法比其他的要精确得多。对于肿瘤数据Tirosh和Puram，本发明实施例采用的UCRSI方法对肿瘤细胞的分型最为准确，其次是SC3方法。然而K-means方法和SC3方法对于其它类型细胞的分型与数据集的原始研究结果相比有较大差异。对于大型数据集GSE102827，包含48266个细胞，加上亚型共33种细胞类型。本发明实施例采用的UCRSI方法与GiniClust3方法有相似的性能，但与SCCAF方法相比有一定的差距，说明SCCAF方法在识别稀有细胞类型上具有一定的优势。

本发明实施例使用双侧非配对Wilcoxon检验来检测本发明实施例采用的UCRSI方法是否优于现有方法。统计检验结果表明，本发明实施例采用的UCRSI方法在三个评价指标NMI、AMI和ARI上显著优于除了SCCAF方法之外的其他五种方法，p-value最高为0.0195。基准测试的结果显示出了本发明实施例采用的UCRSI方法和SCCAF方法之间的互补性，本发明实施例采用的UCRSI方法在4个数据集上好于SCCAF方法。检验结果表明，对于高度异质性的组织中混合细胞的分型，本发明实施例采用的UCRSI方法对细胞类型和状态的区分可较现有聚类算法效果更好。

图4和图5显示了本发明在GSE148071数据集上采用的UCRSI方法在衡量肿瘤内异质性上的效益。对于专注于细胞类型的细胞距离UCRSI

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。