基于图傅里叶变换的空间可变基因识别与分析系统

1.本发明涉及空间转录组学技术领域，特别是涉及基于图傅里叶变换的空间可变基因识别与分析系统。


背景技术：

2.本部分的陈述仅仅是提到了与本发明相关的背景技术，并不必然构成现有技术。
3.相较于单细胞转录组学数据，空间转录组学技术（例如 visium，st，slide-seq等）不需要解离，从而保留了测序位点的空间结构，这些测序位点（spot）包含两种信息：1.基因表达信息，2.空间坐标信息。这项技术为不同分辨率(亚细胞、细胞和微细胞群)的分子组织生物学研究提供了强大工具。在空间转录组学的研究和应用中，一个基础且重要的任务是如何融合基因表达信息和spot的空间位置，进而完成生物发现。
4.识别空间转录组数据中的空间可变基因（spatially variable gene， svg）是空间转录组分析中最基础且重要的一步，它的目的是识别出具有特定空间表达模式的基因，称为空间可变基因，即这些基因的表达在空间各位置上并不是相互独立的，而是与空间位置有密切关联，表现为空间临近的位置有相似的基因表达。
5.现有的空间可变基因预测方法，主要是基于统计框架(如spatialde)或图神经网络(如spagcn)开发的。虽然这些方法也具有一定的svg检测性能，具有严格的统计评价，为不同的生物系统提供了有价值的见解，但仍存在两个局限性，包括：1.这些由前提假设驱动的方法可以有效识别某些明确的模式(如径向热点、曲线带或梯度条纹)，但它们在不规则模式下对异质性肿瘤组织的识别能力出现下降；2.现有的大多数工具都具有较高的svg识别精度，但牺牲了可拓展性，即难以应用到大规模数据集。
6.regev和theis将组织模块(tisue module，tm)定义为包含周期性细胞群落并在不同部位执行特定功能的结构区域，是研究分子组织生物学中关于健康或患病空间转录组数据中分子组成和功能的关键概念。然而，由于组织模块的分子特征和空间结构都是未知的，组织模块的识别还没有一个严格的识别方法。


技术实现要素：

7.为了解决现有技术的不足，本发明提供了基于图傅里叶变换的空间可变基因识别与分析系统；基于图傅里叶变换的空间可变基因识别与分析系统，包括：获取模块，其被配置为：获取待识别的基因表达矩阵；预处理模块，其被配置为：基于待识别的基因表达矩阵，生成无向图；获取无向图对应的邻接矩阵；获取无向图对应的度矩阵；基于邻接矩阵和度矩阵，确定无向图的拉普拉斯矩阵；对拉普拉斯矩阵进行谱分解，得到拉普拉斯矩阵的特征值和特征值对应的特征向量；
转换模块，其被配置为：将特征值作为图傅里叶变换中的频率，将特征值对应的特征向量作为图傅里叶变换中的傅里叶模；采用图傅里叶变换，将当前基因在每个位点的表达值从空间域转换到谱域，得到谱域信号；空间可变基因识别模块，其被配置为：基于拉普拉斯矩阵的特征值和谱域信号，计算出待识别基因表达矩阵在空间域分布的随机性得分；根据随机性得分，确定出空间可变基因。
8.与现有技术相比，本发明的有益效果是：本发明（spagft）是为空间转录组学（spatial transcriptomics）设计的基于图傅里叶变换（graph fourier transform）的空间可变基因识别、组织模块识别、基因表达信号增强算法。
9.本发明利用空间转录组数据的图结构特性，首次将基因表达数据和图拓扑结构有机地结合在一起，通过图傅里叶变换将基因表达数据从空间域(spatial domain)映射到谱域(frequency domain)，得到基因表达数据的新可解释性表示。在谱域中，本发明通过设计空间可变基因分数gftscore，实现了对空间可变基因的高效、快速识别。
10.本发明通过对所有空间可变基因的谱域信号进行聚类，首次用计算方法识别出执行特定功能的组织模块。
11.本发明设计了基因表达数据低通滤波器，实现了基因表达数据有效降噪和信号增强。
12.相较于当前那些仅能识别空间可变基因的方法，本发明无任何对数据分布假设的约束，这意味着本发明可以应用在各类空间转录组技术中，包括基于成像的空间转录组技术（osmfish，seqfish， merfish等）和基于测序的空间转录组技术（st，visium，slide-seq等）。本发明开创性地构建起经典的信号处理领域和新兴的空间转录组学的桥梁，可以应用在生命科学研究中包括肿瘤微环境和胚胎发育在内的各个领域。
附图说明
13.构成本发明的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。
14.图1为实施例一的系统功能模块图。
具体实施方式
15.应该指出，以下详细说明都是示例性的，旨在对本发明提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本发明所属技术领域的普通技术人员通常理解的相同含义。
16.需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本发明的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，术语“包括”和“具有”以及他们的任何变形，意图在于覆盖不排他的包含，例如，包含了一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元，而是可包括没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。
17.在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。
18.本实施例所有数据的获取都在符合法律法规和用户同意的基础上，对数据的合法应用。
19.实施例一本实施例提供了基于图傅里叶变换的空间可变基因识别与分析系统；如图1所示，基于图傅里叶变换的空间可变基因识别与分析系统，包括：获取模块，其被配置为：获取待识别的基因表达矩阵；预处理模块，其被配置为：基于待识别的基因表达矩阵，生成无向图；获取无向图对应的邻接矩阵；获取无向图对应的度矩阵；基于邻接矩阵和度矩阵，确定无向图的拉普拉斯矩阵；对拉普拉斯矩阵进行谱分解，得到拉普拉斯矩阵的特征值和特征值对应的特征向量；转换模块，其被配置为：将特征值作为图傅里叶变换中的频率，将特征值对应的特征向量作为图傅里叶变换中的傅里叶模；采用图傅里叶变换，将当前基因在每个位点的表达值从空间域转换到谱域，得到谱域信号；空间可变基因识别模块，其被配置为：基于拉普拉斯矩阵的特征值和谱域信号，计算出待识别基因表达矩阵在空间域分布的随机性得分；根据随机性得分，确定出空间可变基因。
20.进一步地，获取待识别的基因表达矩阵；其中矩阵的纬度为；其中，表示位点的个数；表示基因的个数。
21.进一步地，所述基于待识别的基因表达矩阵，生成无向图，具体包括：先计算出所有两两位点之间的空间距离；将每一个位点视为无向图的一个节点；判断一个节点是否是另外一个节点的个最近的邻居，如果是，则表示这两个节点之间存在连接边，如果否，则表示这两个位点之间不存在连接边；进而得到无向图。
22.示例性地，生成无向图，其中，对应于个位点；是边集，满足：当且仅当是的近邻或者是的近邻。
23.进一步地，所述获取无向图对应的邻接矩阵，具体包括：对于无向图，它对应的邻接矩阵(adjacent matrix)定义为：进一步地，所述获取无向图对应的度矩阵，具体包括：对于无向图，它对应的度矩阵(degree matrix)为对角矩阵，定义为：进一步地，所述基于邻接矩阵和度矩阵，确定无向图的拉普拉斯矩阵，具体包括：
图对应的拉普拉斯矩阵(laplacian matrix)为：。
24.进一步地，所述对拉普拉斯矩阵进行谱分解，得到拉普拉斯矩阵的特征值和特征值对应的特征向量，具体包括：对拉普拉斯矩阵进行谱分解得到：；；；其中，是对角矩阵，其对角元素为拉普拉斯矩阵的特征值（eigenvalue），是特征值对应的特征向量（eigenvector），也称为傅里叶模（fourier mode，fm）。
25.进一步地，所述将特征值作为图傅里叶变换中的频率，将特征值对应的特征向量作为图傅里叶变换中的傅里叶模；采用图傅里叶变换，将当前基因在每个位点的表达值从空间域转换到谱域，得到谱域信号，具体包括：对于傅里叶模，表示在位点上的表达值。
26.的平滑度（smoothness）反映第个fm在所有相互相邻的点上的总差异性变化，的平滑度公式表示为：的平滑度公式表示为：的平滑度公式，改写为：其中，是的转置；根据平滑度的定义，如果一个特征向量对应一个小的特征值，说明fm值在相邻节点上的变化较小。特征值的增大趋势对应于特征向量振荡强度的增加趋势；因此，用拉普拉斯矩阵的特征值和特征向量分别作为我们的图傅里叶变换中的频率和傅里叶模（fm）。直观上看，一个小的特征值对应一个低频fm，而一个大的特征值对应一个高频fm。
27.对于基因，令表示该基因在个位点上的表达值，是一个维向量；称为基因的空间域信号(signal in spatial domain)。接下来，对该信号进行图傅里叶变换，将其从空间域（spatial domain）的转换到谱域(frequency domain)
的，得：由此，可以看作在傅里叶模上的投影，反映了对的贡献，。这种傅里叶变换统筹了基因表达及其空间分布，得到对空间转录组学数据中基因的新表示方法。
28.进一步地，所述基于拉普拉斯矩阵的特征值和谱域信号，计算出待识别基因表达矩阵在空间域分布的随机性得分，具体包括：采用打分算法gftscore，用来定量测量基因表达在空间域分布的随机性，打分算法gftscore定义为：法gftscore定义为：其中，是拉普拉斯的第个特征值，和为分别为经过标准化处理后的和原始的谱域信号中第个分量，表示该基因在个位点上的表达值。
29.进一步地，所述根据随机性得分，确定出空间可变基因，具体包括：根据随机性得分，采用kneedle算法识别出所有基因得分中的拐点并将拐点对应的得分设为设定阈值；将得分高于设定阈值的基因作为空间可变基因。
30.应理解地，高gftscore的基因在空间域上倾向于非随机分布基因，反之亦然。因此这些高gftscore的基因就是空间可变基因（svg）。至此，已经得到所有基因的gftscore。为了筛选出这些svg，本发明采用kneedle算法识别出具有显著性高的gftscore的基因。
31.进一步地，所述根据随机性得分，确定出空间可变基因；还包括：通过单端秩和检验（wicoxonone-sided rank sum test ）计算各基因的值，然后对值进行fdr（false discovery rate）校正来降低假阳性，从而实现更高的识别精度降低错误率，最终将值小于0.05的基因作为空间可变基因。
32.应理解地，除了采用kneedle算法，本发明也加入非参数检验方法来检验低频信号与高频信号中值的差值。其零假设是，svg的低频信号的中值等于或低于高频元素的中值。通过wilcoxon单端秩和检验计算各基因的值，然后使用错误发现率(fdr)方法进行调整。最终，根据以下两个准则确定空间可变基因，1.采用kneedle算法识别出的高得分的基因，2. wilcoxon检验中adjustedvalue 小于0.05的基因。
33.进一步地，所述系统还包括：组织模块识别模块；其中，组织模块识别模块，其被配置为：采用聚类算法对空间可变基因进行聚类，每一个聚类簇内包含若干个空间可变基因，将每一个聚类簇作为识别的组织模块输出。
34.示例性地，具有相似模式的空间可变基因在频域也具有相似的低频信号，这为聚
类提供了理论可行性。本发明采用了louvain算法对得到的谱域信号进行聚类。根据上述的图傅里叶变换方法得到这些空间可变基因前低频信号作为聚类特征，输入louvain算法得到空间可变基因的聚类结果。
35.为了可视化组织模块，对于在同一个簇（cluster）的空间可变基因，它们支撑了一个组织模块tm。此外，设计了来表示该组织模块在位点上的虚拟表达值：；进一步地，所述系统还包括：空间可变基因信号增强模块；其中，空间可变基因信号增强模块，其被配置为：通过低通滤波器对空间可变基因的谱域信号进行增强；所述低通滤波器设有目标函数；采用凸优化的方式对目标函数进行求解，得到增强后的谱域信号；对增强后的谱域信号，进行傅里叶逆变换，得到增强的基因表达信号。
36.进一步地，所述低通滤波器的目标函数，包括：选取最优化问题的最优解作为的良估计：其中，是2-范数，为解空间（solution space）中的变量，，是一个常数，用来平衡相似性和平滑性。
37.其中，上述目标函数具有两个特点，包括：第一，增强后的信号与原始谱域信号具有相似性；第二，相邻位点的基因表达值具有相似性。
38.某些tm的svg可能会出现低表达或脱落（dropout）问题，表现为检测到的基因表达信号具有较高的噪声。为了解决这个问题，本发明spagft采用了低通过滤器（low-pass filter）来增强空间可变基因的谱域信号，降低基因表达信号中的噪声。对于带有噪声的基因表达信号，设为真实的基因表达信号，并且有，其中，表示噪声。
39.本发明将通过低通滤波器得到，是对一个良估计（well estimation）。一个良估计应该满足两个条件:(1)增强信号(估计的基因表达)应与原始信号相似；(2)在估计的基因表达在保持相邻的位点之间具有低差异。
40.根据凸优化（convex optimization），最优化问题的最优解满足：
其中，，为单位矩阵。为低通滤波器(low-pass filter)，表示增强的谱域信号。
41.最终代表经过傅里叶逆变换（inverse graph fourier transform，igft）的增强基因表达信号。
42.本发明是为空间转录组学（spatial transcriptomics）设计的基于图傅里叶变换（graph fourier transform）的空间可变基因识别、组织模块识别、基因表达数据降噪算法。
43.本发明利用空间转录组数据的图结构特性，首次将基因表达数据和图拓扑结构有机地结合在一起，通过图傅里叶变换将基因表达数据从空间域映射到谱域, 得到基因表达数据的可解释性新表示。
44.在谱域中，本发明通过设计空间可变基因平滑分数，实现了对空间可变基因的高效、快速识别。
45.本发明通过对所有空间可变基因的谱域信号进行聚类，首次用计算方法识别出具有生物学功能的组织模块。
46.本发明设计了基因表达数据低通滤波器，实现了基因表达数据有效降噪，从而大幅提升了功能区域检测算法的性能。
47.相较于当前仅能识别空间可变基因的方法，本发明无任何数据分布假设，这意味着本发明可以应用在各类空间转录组技术中，包括基于成像的空间转录组技术（osmfish，seqfish，merfish等）和基于测序的空间转录组技术（st，visium，slide-seq等）。
48.本发明开创性地构建起经典信号处理领域和新兴空间转录组学的桥梁，可以应用在生命科学研究中包括肿瘤微环境和胚胎发育在内的各个领域。
49.1，识别空间可变基因。空间可变基因可以反映出复杂组织中结构特征和分子功能差异，通过识别和分析这些空间可变基因，进而实现对宏观结构的微端解释。
50.2，识别空间特异性组织模块。组织模块是研究分子组织生物学中关于健康或患病空间转录组数据中分子组成和功能的关键概念，本发明通过对空间可变基因的谱域信号聚类，找出单样本或多样本中的空间特性组织模块，有助于生物功能的分子解释。
51.3，基因表达信号增强。基于测序的空间转录组数据具有往往具有极高的噪声，本发明设计出低通滤波器，有效地降低生物数据中的噪声，有利于生物结构性区域的识别。
52.我们假设一个可靠的组织模块应该由一组空间可变基因(spatially variable gene)决定，这些基因在空间位置上共享一个重复出现的空间表达模式。本发明首次为组织模块识别提供了一种新的计算框架。
53.本发明首次提出了一个无数据分布假设的图傅里叶变换(gtf)框架，命名为spagft，用于从空间转录组准确地检测空间可变基因和组织模块。由于它在分子组织生物学、前沿空间转录组技术和高效的信号处理方法之间建立了桥梁，因此它具有高效的计算效率和良好的可解释性。
54.本发明创造性地利用图傅里叶变换在频率/谱域上得到基因的新表示方法。这样的新表示同时反映了基因在所有位点上的表达信息和由这些位点所构建的空间结构信息。
55.因此，在多组学中，通过分析这些新表示，可以实现三个重要功能，包括：1.空间可变基因的识别，2.组织模块的识别，3.基因表达信号的增强。
56.总的来说，本发明spagft构造了一个最近邻图(knn)来利用spot的空间信息。然后，根据图的拓扑结构得到不同的傅里叶模 (fourier modes，fm)。
57.该方法将一个基因表达作为图上的图信号，采用图傅里叶变换(gft)获得该基因的频域表示。在频域，通过强调频域信号在低频区域上的强度，可以实现空间可变基因的识别；通过在频域上对相似频域信号的基因进行聚类来寻找组织模块；通过低通滤波器（low-pass filter）实现基因的降噪。
58.本发明可以分为两个流程，包括将基因表达信号从空间域转换到频域的图傅里叶变换，和对这些基因在频域的表示进行分析实现若干功能。在第一个流程中，需要实现三个步骤，包括构造knn图（graph），计算傅里叶模(fourier mode)和图的傅里叶变换(gft)。第二个流程，通过对谱域上的信号进行处理和分析，从而实现识别空间可变基因、识别组织模块、增强基因表达信号增强的功能。spagft优越的准确性、可扩展性和可解释性表明，它是一种新颖而强大的工具，可以获得对各种生物学问题的新见解。
59.以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。