CN114121158A

CN114121158A - 一种基于深度网络自适应的scRNA-seq细胞类型识别方法

Info

Publication number: CN114121158A
Application number: CN202111471768.3A
Authority: CN
Inventors: 王树林; 刘孟林
Original assignee: Hunan University
Current assignee: Hunan University
Priority date: 2021-12-01
Filing date: 2021-12-01
Publication date: 2022-03-01

Abstract

本发明涉及生物信息学中的数据挖掘，特别是涉及对scRNA‑seq数据的挖掘。具体涉及一种基于深度网络自适应的scRNA‑seq细胞类型识别方法。本发明的方法包括对scRNA‑seq数据的处理；构建神经网络训练scRNA‑seq数据；加入自适应层优化神经网络架构以克服不同批次数据集之间的差异性；对类型信息未知的scRNA‑seq数据集中的细胞类型进行准确识别。本发明可用于识别未知scRNA‑seq数据集的细胞类型，可有效克服类型信息已知的数据集和类型信息未知的数据集之间的技术差异和批次效应影响。

Description

一种基于深度网络自适应的scRNA-seq细胞类型识别方法

技术领域

本发明涉及生物信息学中的数据挖掘，特别是涉及对scRNA-seq数据的挖掘，具体涉及一种基于深度网络自适应的scRNA-seq细胞类型识别方法。

背景技术

细胞被认为是生物体基本的结构和功能单位。人体细胞中包含大约2万个基因，每个细胞存在自身特异的基因表达模式，仅对部分基因进行表达。这导致了细胞特异性的蛋白质成分和生物功能。scRNA-seq以单个细胞为单位，通过全基因组或转录组扩增进行高通量测序，能够揭示单个细胞的基因结构和基因表达状态并反映细胞间的异质性。近十年scRNA-seq技术发展迅速，测序数据的规模从几十个增加到几千个乃至几百万细胞，出现了大量新的测序平台，如10xGenomics Chromium、inDrop和Drop-seq。细胞类型的识别在scRNA-seq数据的分析中起着重要作用，良好注释的scRNA-seq数据使生物学家能够进行进一步的下游分析并提高我们对疾病细胞机制的理解。

目前通过生物信息学识别scRNA-seq数据类型的方法主要分为三类：第一类方法首先对细胞群进行聚类，然后通过差异表达分析找到每个簇特有的标记基因，最后根据其基因的本体功能对细胞进行注释。然而这类方法的泛化性能通常较差。此外，随着测序数据规模的增加，通过寻找标记基因来注释细胞的任务变得越来越繁重和耗时。第二类方法是利用注释良好的参考数据集的信息来辅助新数据中的细胞类型识别。其中有代表性的方法通过将目标数据集中的细胞投影到一个由注释良好的源数据集中选择的高信息性基因决定的空间，然后根据目标数据中的细胞与源数据中平均细胞类型特异性基因表达的相关性，为它们分配细胞身份。然而这类方法只能利用参考数据中的细胞类型信息，忽略了目标数据中的有用信息。最后一类方法主要是通过深度神经网络克服大规模scRNA-seq数据类型识别的负担，这些方法使用非线性自动编码器将测序数据嵌入到低维空间中进行后续的聚类和分类任务。同样这些方法也没有考虑技术变异和批次效应可能导致的性能崩溃，尤其是当目标数据和参考数据来自不同的测序平台时，细胞分类的准确度会大大降低。

综上所述，现有的方法没有充分考虑到不同测序平台、不同组织、不同物种数据集之间的差异性，极少地充分利用注释良好的参考数据集和未知数据集的基因表达信息和数据分布信息，因此如何设计一种鲁棒性的方法对未知scRNA-seq细胞类型进行准确识别仍然是一个挑战。

发明内容

本发明针对以上方法存在的问题和scRNA-seq细胞类型准确识别的重要性，提出了一种基于深度网络自适应的scRNA-seq细胞类型识别方法。本发明的方法采用深度网络自适应对注释良好的参考数据集和未知的目标数据集进行基因表达信息提取和数据分布对齐，是一种对不同批次scRNA-seq数据集进行细胞类型识别的方法。所叙述方法步骤包括：

1.数据收集阶段

本方法从多个数据平台收集了多种情况的数据集。第一类是通用的基准数据集，通过两种测序方式生成，分别是10x和CelSeq2；第二类是采用不同测序方式生成的人类胰腺组织数据集，通过五种测序方式生成，分别是CelSeq、CelSeq2、SmartSeq2、Fluidigmc1、inDrop；第三类是同一物种内不同组织的数据集，从Figshare上下载的小鼠衰老细胞图谱(Tabula Muris Senis)数据集，包括来自96307个细胞的23341个基因表达信息，包含22种组织。综合上述数据集可评估本方法在多个物种下不同组织的细胞类型的识别准确性。

2.数据预处理阶段

将不同的scRNA-seq数据集随机分为源域和目标域，源域的类型信息已知，目标域的类型信息未知。处理对象为scRNA-seq数据的基因表达矩阵，其中，行为细胞名称、列为基因名称。附加列为细胞类型信息。数据预处理包括质量控制、数据标准化，细胞类型转换三个步骤。质量控制是检查原始数据集中是否存在异常值并设置阈值去除，数据标准化是过滤少于5000个reads和500个基因的低质量细胞，以及少于10个细胞表达的基因。然后使用SCANPY将每个细胞正则化为10000个读数计数；最后将数据集进行对数处理和归一化处理。细胞类型转换是将数据集的细胞类型信息转换为数字编号便于后续的细胞分类。

3.搭建神经网络架构阶段

本方法使用的神经网络由一个输入层和两个全连接层组成，输入层的神经元数量为数据预处理之后的基因数，全连接层的第一层使用1000个神经元，第二层使用100个神经元。全连接层中神经元的活动通过归一化层Layer Normalization进行标准化。LayerNormalization定义为：

全连接层中的非线性激活函数使用SELU,定义为：

SELU(x)＝scale*(max(0,x)+min(0,α*(exp(x)-1)))

在预训练阶段，使用神经网络的镜像作为解码器，整体构成自编码器对目标域进行预训练，使用均方误差MSE作为自编码器的重建损失函数；在正式训练阶段，源域和目标域均采用上述神经网络作为基本网络结构，源域网络还包含一个分类层，分类层的神经元数量为细胞类型的数量，使用交叉熵(cross-entropy)作为源域网络的分类损失函数，定义为：

其中y表示细胞的真实类型标签，如果细胞属于第j个细胞类型则y[j]定义为1，而y中其他位置定义为0。y′表示输出的类型标签，y′[j]表示该细胞是第j个细胞类型的后验概率。

4.优化神经网络架构阶段

在源域和目标域的网络结构的第二个全连接层后均加入域自适应层，自适应层能够使得源域和目标域的数据分布更加接近，减少批次效应对最终分类结果的影响。自适应度量方法采用多核MMD(MK-MMD)，MK-MMD的平方公式定义为：

其中p,q分别表示源域和目标域的概率分布，H_k表示具有特征核k的再生希尔伯特空间RKHS，d_k(p,q)表示p和q的平均嵌入之间的RKHS距离。重要的性质是如果

那么p＝q。与特征映射φ相关联的特征核定义为：

k(X^s,X^t)＝<φ(X^s),φ(X^t)>

其多核表示形式为多个PSD核{k_u}的凸组合：

其中对系数{β_u}施加约束，以保证导出的多核k具有特征性。MK-MMD通过多个不同的核进行加权最终得到的表征能力强于只有一个核的MMD。

最终网络模型的优化目标由两部分组成：分类损失函数和自适应损失。优化目标通过最小化分类损失和MK-MMD实现，总的损失函数定义为：

其中Θ表示网络的所有权重和bias参数，是需要学习的目标；λ>0是惩罚参数；l₁到l₂表示需要进行自适应的层数；

和

分别表示源域和目标域的第l层隐藏表示；x_a和n_a分别表示源域和目标域中所有包含类型信息的数据的集合；J(·)表示分类损失函数。学习阶段关键任务是学习网络参数Θ和MK-MMD的β。

5.准确识别未知scRNA-seq数据集的细胞类型

使用小批量随机梯度下降法(mini-batch SGD)对源域网络和目标网络进行参数更新和迭代优化。源数据集和目标数据集通过PyTorch自带的DataLoder分为多个mini-batch作为输入训练和优化网络，最终训练完成的目标域网络作为分类器对类型信息未知的目标数据集进行准确类型识别。

附图说明

图1：深度网络自适应模型框架图

图2：目标数据预训练流程图

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实验，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。

硬件环境主要是一台PC主机。其中，PC主机的CPU为Intel(R)Core(TM)i5-6400，2.70GHz，内存为16GB RAM，64位操作系统。软件以Windows 10为平台，在PyCharm环境下用Python语言实现，PyCharm版本为2021.1.3，Python语言版本为3.7.0。

1.数据收集与整理阶段

本方法所用数据包括三大类，第一类是基准数据集，通过10x和CelSeq2两种测序方式生成；第二类是采用不同测序方式生成的人类胰腺组织数据集，通过五种测序方式生成，分别是CelSeq、CelSeq2、SmartSeq2、Fluidigmc1、inDrop；第三类是从Figshare上下载的小鼠衰老细胞图谱(Tabula Muris Senis)数据集。包括来自96307个细胞的23341个基因表达信息，包含22种组织，测序数据较为完备。通过这三类数据可以对本方法的实用性进行较为完善的评估。所有数据对象均通过AnnData进行存储。scRNA-seq初始数据由表1所示的几个部分组成。

表1：AnnData的主要组成部分

在本方法中，矩阵数据为细胞*基因的矩阵；观察值数据包含细胞类型信息以及测序数据的批次信息等。

2.数据预处理阶段

将不同的scRNA-seq数据集随机分为源数据集和目标数据集，源数据集的类型信息已知，目标数据集的类型信息未知。首先确定源数据集和目标数据集中均检测到的基因，然后基于共同基因将源数据集和目标数据集合并为一个矩阵。接着开始对scRNA-seq初始数据进行预处理，预处理包括质量控制、数据标准化，细胞类型转换三个步骤。质量控制主要是删除细胞中无细胞类型信息的数据，这部分数据的类型信息为‘nan’或‘NA’；数据标准化通过SCANPY包实现，过滤少于5000个reads和500个基因的低质量细胞，以及少于10个细胞表达的基因，然后将每个细胞正则化为10000个读数计数，最后可以根据实际情况选择是否对数据集进行对数处理和归一化。本方法在实验过程中均对所有数据集进行了对数处理和归一化。细胞类型转换是将字符形式的细胞类型信息映射为数字编号，方便后续对细胞进行分类以及评价方法的类型识别准确度。经过预处理之后，将数据重新拆分为源数据集和目标数据集。源域和目标域的示例数据信息如表2所示。

表2：源域和目标域的数据信息

3.搭建神经网络架构阶段

本方法使用的神经网络由一个输入层和两个全连接层组成，输入层的神经元数量为数据预处理之后的基因数，全连接层的第一层使用1000个神经元，第二层使用100个神经元。其中每个连接层中包含四个步骤，分别是：(1)对输入数据应用线性变换；(2)对神经元的活动通过归一化层Layer Normalization进行标准化；(3)使用激活函数SELU对神经元的活动进行非线性变换；(4)使用dropout实现正则化。

在预训练阶段，使用神经网络的镜像作为解码器，整体构成一个自编码器对目标域进行预训练，使用均方误差MSE作为自编码器的重建损失函数；本方法默认使用的参数pretrain_epochs为10。可根据数据集的实际情况及训练结果选择是否开启预训练步骤。

在正式训练阶段，源域和目标域均采用上述神经网络作为基本网络结构，源域网络还包含一个分类层，分类层的神经元数量为细胞类型的数量。神经网络中的主要参数设置如下：起始学习率为0.001，学习率呈阶梯指数衰减，衰减步长设置为20。这意味着每20个epochs后，学习率将是原始学习率乘以0.95，使用50个epochs来训练神经网络，mini-batch大小为32，这是在每个epochs中使用的细胞数量。

4.优化及训练阶段

在源域和目标域的网络结构的第二个全连接层后均加入域自适应层，源域和目标域之间的自适应损失使用多核MMD(MK-MMD)度量。在训练过程中计算源域和目标域在自适应层的MMD距离，具体是将源域和目标域映射到具有特征核k的再生希尔伯特空间RKHS中，然后计算两者在高维空间中的数据分布距离。本方法使用40个n_iter_per_epoch训练自适应层，这意味着在每一次全局训练步骤中还包含40个迭代训练mmd过程。

最终网络模型的优化目标由两部分组成：分类损失函数和自适应损失。优化目标通过最小化分类损失和MK-MMD实现。总损失函数中自适应损失部分的惩罚参数默认设置为10。使用小批量随机梯度下降法(mini-batch SGD)对源域网络和目标网络进行参数更新和迭代优化。源数据集和目标数据集通过PyTorch自带的DataLoder分为多个mini-batch作为输入训练和优化网络，最终训练完成的源域网络作为分类器对类型信息未知的目标数据集进行准确类型识别。

5.结果分析验证

将基准数据集中由10x和CelSeq2测序生成的scRNA-seq数据分别作为源数据集和目标数据集，统计相应的准确率，然后将两者交换，统计相应的准确率，结果如表3所示。

表3：基准数据集中细胞类型识别准确率

对于不同测序平台生成的人类胰腺组织数据集，本方法使用由CelSeq和CelSeq2分别作为源数据集和目标数据集，统计相应的准确率，然后将两者交换，统计相应的准确率，实验结果如表4所示。

表4：人类胰腺组织数据集中细胞类型识别准确率

从表3和表4中可以看出，虽然源数据集和目标数据集由不同测序平台生成，但是本方法对基准数据集和人类胰腺组织数据集中无类型信息的目标数据集的类型识别准确率是比较高的，前者接近于100％，后者也达到了92％。一定程度上表明本方法可以很好地克服基准数据集和人类胰腺组织中不同批次数据之间的差异，利用已有数据对未知数据的类型进行准确识别。

小鼠衰老细胞图谱(Tabula Muris Senis)数据集包含22种组织，本方法在这里选取了细胞类型丰富的4种组织(Heart、Limb_Muscle、Brain_Non-Myeloid、Liver)进行了实验。将10x Genomics测序的数据作为源数据集，SmartSeq2测序的数据作为目标数据集，统计相应的准确率，实验结果如表5所示。

表5：小鼠衰老细胞图谱数据集中细胞类型识别准确率

从表5中可以看出，在数据量较为庞大、细胞类型丰富的小鼠衰老细胞图谱(Tabula Muris Senis)中，本方法对多个不同组织中无类型信息的目标数据集的类型识别准确率依然是比较高的，进一步证实了本方法对于具有批次效应的scRNA-seq细胞类型识别的可靠性。

Claims

1.一种基于深度网络自适应的scRNA-seq细胞类型识别方法，其特征在于实施步骤：

(1)收集数据，包括通用的基准数据集，不同测序方式生成的人类胰腺组织数据集，同一物种内不同组织的数据集；

(2)对scRNA-seq数据进行预处理，不同的scRNA-seq数据集随机分为源数据集和目标数据集，源数据集的类型信息已知，目标数据集的类型信息未知，预处理包括质量控制、数据标准化、细胞类型转换三个步骤；

(3)搭建神经网络架构，首先目标域可通过自编码器初始化神经网络参数，然后源域和目标域均采用相同的神经网络作为其基本网络结构；

(4)优化框架，在源域和目标域的网络结构中均加入域自适应层，自适应层能够使得源域和目标域的数据分布更加接近，减少批次效应对最终分类结果的影响；

(5)准确识别未知scRNA-seq数据集的细胞类型，使用小批量随机梯度下降法(mini-batch SGD)对源域网络和目标网络进行参数更新和迭代优化，最终模型拥有对类型信息未知的目标数据集进行准确类型识别的能力。

2.根据权利要求1所属的基于深度网络自适应的scRNA-seq细胞类型识别方法，其特征在于数据收集阶段：

(1)基准数据集通过两种测序方式生成，分别是10x和CelSeq2；

(2)人类胰腺组织数据集通过五种测序方式生成，分别是CelSeq、CelSeq2、SmartSeq2、Fluidigmc1、inDrop；

(3)从Figshare上下载的小鼠衰老细胞图谱(Tabula Muris Senis)数据集，包括来自96307个细胞的23341个基因表达信息，包含22种组织。

3.根据权利要求1所属的基于深度网络自适应的scRNA-seq细胞类型识别方法，其特征在于数据的预处理阶段：

(1)检查原始数据集中是否存在异常值并设置阈值去除；

(2)过滤少于5000个reads和500个基因的低质量细胞，以及少于10个细胞表达的基因，使用SCANPY将每个细胞正则化为10000个读数计数，最后将数据集进行对数处理和归一化；

(3)将数据集的细胞类型信息转换为数字编号便于后续的细胞分类。

4.根据权利要求1所属的基于深度网络自适应的scRNA-seq细胞类型识别方法，其特征在于搭建神经网络架构阶段：

(1)神经网络由一个输入层和两个全连接层组成，输入层的神经元数量为数据预处理之后的基因数，全连接层的第一层使用1000个神经元，第二层使用100个神经元；

(2)在预训练阶段，使用神经网络的镜像作为解码器，整体构成自编码器对目标域进行预训练，使用均方误差MSE作为自编码器的重建损失函数；

(3)在正式训练阶段，源域和目标域均采用上述神经网络作为基本网络结构，源域网络还包含一个分类层，分类层的神经元数量为细胞类型的数量，使用交叉熵(cross-entropy)作为源域网络的分类损失函数。

5.根据权利要求1所属的基于深度网络自适应的scRNA-seq细胞类型识别方法，其特征在于优化网络框架阶段：

(1)在源域网络和目标网络的第二个全连接层后加入自适应层；

(2)自适应度量方法采用多核MMD(MK-MMD)，它度量源域和目标域在再生希尔伯特空间RKHS中数据分布的距离，MK-MMD的平方公式定义为：

与特征映射φ相关联的特征核定义为：

k(X^s,X^t)＝<φ(X^s),φ(X^t)>

其多核表示形式为多个PSD核{k_u}的凸组合：

(3)通过最小化MK-MMD实现源域和目标域的数据分布对齐。

6.根据权利要求1所属的基于深度网络自适应的scRNA-seq细胞类型识别方法，其特征在于能够准确识别未知scRNA-seq数据集的细胞类型，将源数据集和目标数据集通过PyTorch自带的DataLoder分为多个mini-batch作为输入训练和优化网络，优化目标由两部分组成：分类损失和自适应损失，最终训练好的目标域网络作为分类器对目标数据集的类型进行准确识别。