数据分析方法和装置

1.本发明涉及计算机技术应用领域，尤其涉及一种数据分析方法和装置。


背景技术：

2.单细胞测序技术能够在单细胞分辨率上定量测量基因表达、染色质可及性、甲基化等模态的数据。单细胞测序技术的快速发展促使全世界收集了大量的多模态、跨尺度、高维度的单细胞数据集。上述全面的数据集与单细胞数据分析的计算方法极大地促进了人们对与性状和疾病相关的细胞异质性、细胞发育、编码区和非编码区遗传变异的理解。然而，每种模态的数据都具有自身的优劣势。每种模态的数据都提供了不同的细胞生物学视角，单模态数据的分析可能会因为其特定于某一模态的技术噪声而丢失生物信号。例如，仅通过分析单细胞rna测序(single cell rna sequencing，简称scrna-seq)数据很难分离t细胞的亚细胞类型。因此，前沿的一些单细胞技术(例如，细胞索引转录组和抗原决定簇测序(cellular indexing of transcriptomes and epitopes by sequencing，简称cite-seq)、单细胞测序技术(single-nucleus chromatin accessibility and mrna expression sequencing，简称snare-seq)、10xmultiome、甲基化组和转录组测序(single-cell nucleosome,methylation,and transcription sequencing,简称scnmtseq)、sci-car)已被开发用于同时测量来自同一细胞的两种或多种模态的数据。由于不同模态的数据质量和信息等上都有较大差异，整合这些新数据的任务仍然具有挑战性。解决这一问题，还需要开发高效的数据整合方法，通过联合分析的方法准确定义不同组织的多模态数据细胞图谱，以揭示正常情况下和病理情况下的细胞异质性、功能调节和细胞状态转变。
3.随着单细胞多模态测序技术的出现，一些用于联合分析同一细胞中同时测量的多模态组学数据的方法被提出：seurat-v4方法使用加权最近邻策略(weighted-nearest neighbor，简称wnn)学习每个模态特定于该细胞的模态权重，该权重反映每种模态的信息相对重要性；mofa 是基于变分推断和分组自动相关性确定先验(automatic relevance determination，简称ard)技术，开发的一种推断潜在因子和相关的特征权重矩阵，从而实现多模态和多样本组的同时集成的统计方法；最近提出的方法包括变分自动编码器(variational auto encoder，简称vae)的变体(例如cobolt、scmm和totalvi)，它们通过多模态编码器学习非线性联合嵌入。然而，大多数现有方法都有其自己的局限性。具体来说，seurat-v4在没有统计模型的情况下推断出一个集成的网络来对细胞状态进行联合定义，无法解释细胞在每种模态下的样本方差的来源；具有复杂参数模型的方法(例如mofa )很难扩展到大型数据集上；一些深度学习方法(如cobolt和totalvi)尽管可以扩展到大型数据集，但由于对低维细胞表示的预测不准确，它们可能无法揭示稀有或小的细胞亚群。
4.针对目前相关技术中由于现有技术的局限性导致的无法对细胞的多模态进行分析的问题，目前尚未得到有效的解决。


技术实现要素：

5.本发明实施例提供了一种数据分析方法和装置，以至少解决相关技术中由于现有技术的局限性导致的无法对细胞的多模态进行分析的问题。
6.根据本发明实施例的一个方面，提供了一种数据分析方法，包括：获取待分析数据；从待分析数据中提取第一个模态的低维特征；分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据。
7.可选的，从待分析数据中提取第一个模态的低维特征包括：通过第一个特征提取器从待分析数据中提取第一个模态的低维特征。
8.进一步地，可选的，第一个特征提取器为非线性模型，第一个特征提取器的结构为编码器。
9.可选的，分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据包括：通过第二个非线性模型对低维特征进行重构，得到第一个模态的数据；通过第三个非线性模型对低维特征进行重构，得到第二个模态的数据。
10.进一步地，可选的，第二个非线性模型的结构为解码器。
11.可选的，第三个非线性模型的结构为解码器。
12.可选的，在待分析数据为单细胞数据的情况下，第一模态的数据和第二个模态的数据为单细胞数据的模态分析数据，其中，第一模态的数据和第二个模态的数据的模态不同。
13.根据本发明实施例的另一个方面，提供了一种数据分析装置，包括：获取模块，用于获取待分析数据；提取模块，用于从待分析数据中提取第一个模态的低维特征；分析模块，用于分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据。
14.可选的，提取模块包括：提取单元，用于通过第一个特征提取器从待分析数据中提取第一个模态的低维特征。
15.进一步地，可选的，第一个特征提取器为非线性模型，第一个特征提取器的结构为编码器。
16.本发明实施例中，基于获取待分析数据；从待分析数据中提取第一个模态的低维特征；分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据。也就是说，本发明实施例能够解决了由于现有技术的局限性导致的无法对细胞的多模态进行分析的问题，从而达到了从成对的多模态组学数据中识别细胞类型，解析其细胞功能的技术效果。
附图说明
17.此处所说明的附图用来提供对本发明的进一步理解，构成本技术的一部分，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
18.图1为本发明实施例提供的一种数据分析方法的流程示意图；
19.图2为本发明实施例提供的一种数据分析方法中基于多视图潜在空间学习和典型关联分析相结合的单细胞多模态数据联合分析方法vimcca的流程示意图；
20.图3为本发明实施例提供的一种数据分析方法中分布使用两个算法的细胞类型向量与具有金标准的bulk rna-seq数据集中的对应细胞类型向量计算了皮尔逊相关系数的示意图；
21.图4为本发明实施例提供的一种数据分析方法中三种算法的ari系数的对比示意图；
22.图5为本发明实施例提供的一种数据分析方法中数据集的对比示意图；
23.图6为本发明实施例提供的一种数据分析装置的示意图。
具体实施方式
24.为了使本技术领域的人员更好地理解本发明方案，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分的实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都应当属于本发明保护的范围。
25.需要说明的是，本发明的说明书和权利要求书及附图中的术语“第一”、“第二”等是用于区别不同对象，而不是用于限定特定顺序。
26.根据本发明实施例的一个方面，提供了一种数据分析方法，图1为本发明实施例提供的一种数据分析方法的流程示意图。如图1所示，本技术实施例提供数据分析方法包括：
27.步骤s102，获取待分析数据；
28.步骤s104，从待分析数据中提取第一个模态的低维特征；
29.可选的，步骤s104中从待分析数据中提取第一个模态的低维特征包括：通过第一个特征提取器从待分析数据中提取第一个模态的低维特征。
30.进一步地，可选的，第一个特征提取器为非线性模型，第一个特征提取器的结构为编码器。
31.具体的，本技术实施例中基于多视图潜在空间学习和典型关联分析相结合的单细胞多模态数据联合分析的方法，简称vimcca(variational inference multi-view canonical correlation analysis，变分推理多视角典型相关分析)，用于从成对的多模态组学数据中识别细胞类型，解析其细胞功能。vimcca计算框架包括：第一个特征提取器用于提取来自第一个模态的低维特征，它的结构为编码器。此外，vimcca计算框架还包括：第二个非线性模型和第三个非线性模型。
32.步骤s106，分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据。
33.可选的，步骤s106中分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据包括：通过第二个非线性模型对低维特征进行重构，得到第一个模态的数据；通过第三个非线性模型对低维特征进行重构，得到第二个模态的数据。
34.进一步地，可选的，第二个非线性模型的结构为解码器。
35.可选的，第三个非线性模型的结构为解码器。
36.可选的，在待分析数据为单细胞数据的情况下，第一模态的数据和第二个模态的数据为单细胞数据的模态分析数据，其中，第一模态的数据和第二个模态的数据的模态不同。
37.具体的，第二个非线性模型用于通过低维特征重构出第一个模态的数据，它的结构是解码器；第三个非线性模型用于通过低维特征重构第二个模态的数据，它的结构是解
码器。vimcca通过一个低维嵌入特征来重构两个模态的数据，促使神经网络的低维瓶颈层学习到两个模态的互补低维表达。通过vimcca学习到的多模态数据的互补低维特征可以用于多项下游分析任务，例如聚类、细胞类型、细胞亚型识别等。μσε
38.其中，图2为本发明实施例提供的一种数据分析方法中基于多视图潜在空间学习和典型关联分析相结合的单细胞多模态数据联合分析方法vimcca的流程示意图。如图2所示，基于深度神经网络，包含一个用于提取来自第一个模态的低维特征的非线性模型(黑色部分)、一个用于重构第一个模态数据的非线性模型(白色部分)和一个用于重构第二个模态数据的非线性模型(条纹部分)。
39.综上，结合步骤s102至步骤s105中模型的具体内容包括以下部分：
40.(a)三个非线性模型：
41.给定一个单细胞多模态数据集从n个细胞上采集两个模态的数据。其中xi∈x和yi∈y表示第i个细胞中的两个不同的模态数据。vimcca的目标是整合x和y以找到多于一个模态的多模态信息来更详细分析细胞的真实生物学状态。以单细胞基因表达加蛋白质表达数据为例，x是一个维度为n
×
p的矩阵，表示n个细胞在一组p个基因维上的表达；y是一个维度为n
×
q的矩阵，表示n个细胞在一组q个蛋白质表达矩阵。还可以将模型扩展到单细胞多模态数据的其他组合，例如转录组数据与染色质可及性的组合，则y表示染色质可及区域的读取计数。
42.为了对多模态数据进行联合分析，涉及到为n个细胞找到一个维度为n
×
d的表达矩阵，其中d＜＜p和d＜＜q。低维矩阵z应该反映细胞的真实生物学状态，可用于许多下游分析，例如细胞亚群的识别、轨迹推断、聚类和可视化。为了估计z，vimcca将基因表达和蛋白质计数建模为非线性模型，该模型将潜在变量z分别转换为基因表达空间和蛋白质表达空间，具体公式如下：
43.xi＝f
x
(zi；θ
x
) ei；
ꢀꢀꢀ
(1)
44.yi＝fy(zi；θy) ei；
ꢀꢀꢀ
(2)
45.zi～n(0,i)，ei～n(0,σ2i)；
ꢀꢀꢀ
(3)
46.其中，xi为单细胞基因表达向量，yi为单细胞蛋白质表达向量，f
x
(
·
)，fy(
·
)为非线性回归函数，分别使用一组参数θ
x
和θy将zi从一个潜在的d维空间转换为两个可观察空间，ei为残差向量。
47.变分推断假设隐变量zi服从标准正态分布。在这个假设下，n个细胞的潜在变量都位于相同的d维上。非线性函数f
x
(
·
)和fy(
·
)是利用两个深度神经网络(称为解码器)构建的，它们都以zi为输入，输出重构为xi和yi。对于解码器f
x
(
·
)，zi被以d
→
32
→
64
→
128
→
p的形式连接到输出层。fy(
·
)的结构与f
x
(
·
)的区别仅在于输出层的维数为q。vimcca的中间层网络使用relu激活函数激活，然后连接一个bn层和一个dropout层。
48.(b)kl散度：
49.给定两个随机变量(x,y)表示单个细胞中的两个不同模态数据，vimcca基于典型相关分析(cca)的概率潜在模型，该模型定义了以下联合分布：
50.p(x,y,z)＝p(z)p(x∣z)p(y∣z)
ꢀꢀꢀ
(4)
51.p(x,y)＝∫p(x,y,z)dz
ꢀꢀꢀ
(5)
52.其中x和y是可观察变量，是无法观察到的隐变量。cca将x和y视为z的两种
不同的外部表现形式。并假设这两个视图在以隐变量z为条件的情况下是独立的。对于单细胞多模态数据联合分析而言，x和y分别代表单细胞中的两个模态，例如基因表达值和蛋白质计数。vimcca的目标是整合x和y以便找到可以揭示细胞深层信息并改善许多下游分析的z。
53.在vimcca中，非线性观测模型p
θ
(x∣z；θ
x
)和p
θ
(x∣z；θy)分别由θ
x
和θy参数化，它们是深度神经网络权重的集合，用条件密度q
φ
(z∣x；φz)逼近p
θ
(z∣x)，其中φz是另一个深度神经网络的参数集合。可以使用q
φ
(z∣x)推导出边缘数据对数似然的下界：
[0054][0055]
其中d
kl
(q
φ
(z∣x)||p(z))表示潜在变量的近似后验q
φ
(z∣x)和先验q(z)之间的kl散度。
[0056]
vimcca在训练集上最大化数据对数似然的变分下界：
[0057][0058]
当正确选择参数化q(z|x)时，该项可以精确地以封闭形式计算。让变分近似后验为具有对角协方差的多元高斯。也就是说，对于一个样本对(xi,yi)，有：
[0059][0060][0061]
其中均值μi，和协方差σi是编码器f的输出，则kl散度公式如下：
[0062][0063]
同样，为了梯度反向传播的问题，使用了重参数化技巧来替代蒙特卡洛采样，的生成公式为：
[0064][0065]
这样在梯度反向传播的时候就可以正常传播梯度。
[0066]
(c)目标函数
[0067]
vimcca使用如下高斯观测模型：
[0068][0069][0070]
p
θ
(x∣z)和p
θ
(y∣z)分别使用两个深度神经网络测量来自样本z的每个视图输入的重建误差。综上所述，模型的目标函数如下：
[0071][0072]
式中λ为超参数，用于调节两个模态误差的权重。
[0073]
下面结合具体实施方式，举例详细说明本发明的应用效果。
[0074]
(a)分析人类pbmc(外周血单核细胞(peripheral blood monoculear cell,简称pbmc))基因 atac多模态数据
[0075]
这一应用中采用的数据集是来自健康人类供体的10,412个pbmc的atac 基因表达数据集，该数据集通过细胞分选去除粒细胞。它可以从10x genomics单细胞门户网站免费下载。该数据集同时在36,601个基因和106,056个开放染色质上对每个细胞进行联合测序。
[0076]
首先使用本发明提出的方法vimcca训练得到两种模态的互补降维特征。然后在细胞低维互补特征上使用leiden聚类算法对细胞进行聚类。通过对聚类结果的基因差异表达分析，vimcca找出了每个簇的标记基因。然后参考人类蛋白质图谱(human protein atlas)(参考文献thul pj et al.,a subcellular map of the humanproteome.science.(2017).)，根据聚类结果注释了20种pmbc细胞类型。为了验证vimcca标记的细胞类型的合理性，可以通过绘制标记基因在注释的细胞类型上的热图得到。
[0077]
接下来，使用umap算法对原始rna表达数据进行可视化，并用vimcca标记的细胞类型标签进行染色。在rna表达数据中，cd4阳性中央型效应型记忆t细胞亚群(cd4 effector memory t，简称cd4 tem)和cd4阳性效应型记忆t细胞亚群(cd4 central memory t，简称cd4 tcm)不能很好地被区分开。同样使用umap可视化算法对原始atac表达数据进行可视化，并用标记的细胞类型标签进行染色。在atac表达数据的可视化结果中，嗜碱性球细胞(basophil)与cd14阳性单核细胞(cd14 monocyte)不能被很好的区分开。然后，使用umap可视化算法对vimcca捕获的的互补低维特征进行可视化。上述细胞类型在vimcca的可视化结果中得到了很好的区分。同时为了进一步验证细胞类型注释的合理性以及算法的有效性，将vimcca细胞类型注释结果与seurat 4.0的注释结果进行了全面的对比。首先，计算了vimcca细胞类型注释与seurat 4.0细胞类型注释之间的相关性。
[0078]
然后，通过绘制vimcca细胞类型注释与seurat 4.0细胞类型注释之间的桑基能量分流图，可以得到每个细胞类型对应于seurat 4.0结果的流向。同时可以得到对于造血干和祖细胞(haematopoietic stem and progenitor cell，简称hspc)，vimcca与seurat 4.0识别出的完全一致。
[0079]
此外，从geo数据库(gene expression omnibus)中下载编号为gse07011的bulk rna-seq数据，对vimcca和seurat 4.0注释的细胞类型对应的细胞的向量求均值，作为对应的细胞类型向量。然后分布使用两个算法的细胞类型向量与具有金标准的bulk rna-seq数据集中的对应细胞类型向量计算了皮尔逊相关系数。从图3可以看出，对于12种重叠的细胞类型，vimcca注释的细胞类型向量与bulk rna-seq细胞类型向量之间的皮尔逊相关系数更高。
[0080]
除了与seurat 4.0重叠的细胞类型以外，还注释了三个细胞亚型，包括嗜碱性粒细胞(basophil)、cd4阳性干细胞记忆t细胞(cd4 memory stem t cells,简称cd4 tscm)、分化前b细胞(pre-b)。对于这三个细胞，找到了对应的标记基因，并通过绘制这三种亚型相对于seurat 4.0中确定的细胞类型的差异基因表达的小提琴图，根据人类蛋白质图谱，rab31基因和mtrnr2l2基因是basophil的标记基因。由此可知，它们在basophil上明显更加富集。
[0081]
ding,c等人的研究表明，stat3基因和cd44基因在pre-b细胞中具有更高的表达
值。从stat3和cd44标记基因提琴图可以看出stat3基因和cd44基因在vimcca定义的pre-b细胞中具有更高的表达值，从而验证了该细胞亚型注释的合理性。
[0082]
除了中央记忆型t细胞(central memory t cell，简称tcm)和tscm的共同标记外，lef1在tscm细胞中富集。从lef1标记基因提琴图可以看出lef1标记基因在注释的cd4tscm上比在cd4 tcm上更加富集。
[0083]
通过以上综合对比实验，验证了vimcca可以学习多模态数据的互补特征，从而更准确地注释细胞和挖掘子细胞类型。
[0084]
(b)分析人类pbmc基因 蛋白质多模态数据
[0085]
这一应用中采用的数据集是来自健康人类供体的8,258个pbmc的蛋白质 基因表达数据集。它可以从10x genomics单细胞门户网站免费下载。该数据集同时在33,538个基因和17个细胞表面蛋白上对每个细胞进行联合测序。
[0086]
为了验证本发明提出的方法vimcca在基因 蛋白质的多模态数据集上的可扩展性，将vimcca应用于来自健康人体供体的8,258个pbmc。对于每个细胞，在基因表达数据集中测量了33,538个基因，在细胞表面蛋白数据集中测量了17个细胞表面蛋白。
[0087]
首先使用vimcca训练得到基因表达模态和蛋白质表达模态的互补低维特征。然后在细胞低维特征上使用leiden算法进行聚类。通过对聚类结果进行基因差异表达分析，注释了21种细胞类型。为了验证vimcca注释的细胞类型的合理性，可以通过绘制细胞簇和标记基因的表达热图。然后使用umap可视化算法得到原始rna表达的可视化结果，然后根据vimcca注释的细胞类型染色，可以得到基因表达的umap降维结果中不能很好地区分igg阳性cd14阳性单核细胞(igg cd14 monocyte)和初始cd8阳性初始t细胞(cd8 naive)。同样，对原始蛋白质表达数据进行了umap可视化，可以得到蛋白质表达数据的可视化结果中的嗜碱性粒细胞(basophil)和浆细胞(plasma)和浆细胞样树突状细胞(plasmacytoid dendritic cells，pdc)没有很好的被区分出来。然而，上述细胞类型在vimcca的可视化结果中得到了很好的区分。
[0088]
此外，将vimcca与两个多模态联合分析算法seurat 4.0和totalvi进行了比较。通过绘制seurat 4.0的umap可视化结果和totalvi的umap可视化结果。seurat 4.0的umap可视化结果无法将中间型单核细胞(intemidiate monoocyte)与cd14 monocyte分离，也无法将cd45ra阳性cd4阳性效应记忆t细胞(cd4 tem(ra ))与cd45ro阳性cd4阳性效应记忆t细胞(cd4 tem(ro ))分离。totalvi的umap可视化结果显示该算法无法区分cd8阳性初始t细胞(cd8 naive)和gamma delta t细胞(gdt)。
[0089]
然后进一步计算了三种算法的ari系数，绘制了对比图，如图4。可以看出vimcca在所有三种算法中取得最高的ari值。这表明vimcca在保留细胞异质性和识别细胞亚型上有更好的效果。
[0090]
为了比较三种算法的运行效率，在配备16gb ram的macbook pro上测试了三种算法的运行效率。首先对数据进行下采样和上采样，生成细胞数量分别为3,000、6,000、9,000、12,000和15,000的数据集。然后在这五个数据集上运行三个算法并统计每个算法的运行时间，绘制了对比图，如图5。可以从图中看到vimcca和seurat 4.0的运行速度比totalvi快很多，因此vimcca和seurat 4.0都可以适用于大型数据集。由于vimcca基于神经网络并在训练期间使用小批量梯度下降方法，因此它对内存的要求非常低。同时，vimcca还可以利
用gpu加速，从而大幅提高运行速度。
[0091]
为了进一步说明vimcca注释的细胞亚型的合理性，可以通过绘制这四种亚型相对于上层细胞类型的差异基因表达基因与蛋白质的小提琴图，蛋白质igg2a，igg2b，igg1在vimcca识别的igg cd14 monocyte细胞类型中明显比在传统cd14 monocyte中更加富集。
[0092]
蛋白质cd45ra在vimcca定义的cd45ra cd14 monocyte中比在传统cd14 monocyte中更加富集。
[0093]
vimcca还识别出来cd4阳性效应记忆t细胞的两个亚型，cd45ro cd4 tem和cd45ra cd4 tem细胞。蛋白质cd45ra和cd45ro分别在cd45ro cd4 tem和cd45ra cd4 tem细胞中富集。
[0094]
vimcca还识别出了b细胞的亚型，cd3 b细胞。cd3e基因和cd3g基因在cd3 b细胞亚型中明显富集。
[0095]
以上综合实验结果表明，验证了vimcca可以学习基因表达和蛋白质表达的多模态数据的互斥特征，从而更精确地标注细胞和挖掘细胞亚型。并且在所有三种算法中，vimcca在ari系数上的表现也是最好的。同时，vimcca具有很高的运行效率。
[0096]
本发明实施例中，基于获取待分析数据；从待分析数据中提取第一个模态的低维特征；分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据。也就是说，本发明实施例能够解决了由于现有技术的局限性导致的无法对细胞的多模态进行分析的问题，从而达到了从成对的多模态组学数据中识别细胞类型，解析其细胞功能的技术效果。
[0097]
根据本发明实施例的另一个方面，提供了一种数据分析装置，图6为本发明实施例提供的一种数据分析装置的示意图。如图6所示，本技术实施例提供数据分析装置包括：获取模块62，用于获取待分析数据；提取模块64，用于从待分析数据中提取第一个模态的低维特征；分析模块66，用于分别对低维特征进行重构，得到第一个模态的数据和第二个模态的数据。
[0098]
可选的，提取模块64包括：提取单元，用于通过第一个特征提取器从待分析数据中提取第一个模态的低维特征。
[0099]
进一步地，可选的，第一个特征提取器为非线性模型，第一个特征提取器的结构为编码器。
[0100]
以上所述，仅为本发明的较佳实施例而已，并非用于限定本发明的保护范围。