CN115240772B - 一种基于图神经网络的解析单细胞通路活性的方法 - Google Patents

一种基于图神经网络的解析单细胞通路活性的方法 Download PDF

Info

Publication number
CN115240772B
CN115240772B CN202211005236.5A CN202211005236A CN115240772B CN 115240772 B CN115240772 B CN 115240772B CN 202211005236 A CN202211005236 A CN 202211005236A CN 115240772 B CN115240772 B CN 115240772B
Authority
CN
China
Prior art keywords
cell
gene
network
cells
genes
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202211005236.5A
Other languages
English (en)
Other versions
CN115240772A (zh
Inventor
郭雪江
李妍
祝辉
韩旭东
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing Medical University
Original Assignee
Nanjing Medical University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing Medical University filed Critical Nanjing Medical University
Priority to CN202211005236.5A priority Critical patent/CN115240772B/zh
Publication of CN115240772A publication Critical patent/CN115240772A/zh
Application granted granted Critical
Publication of CN115240772B publication Critical patent/CN115240772B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/10—Sequence alignment; Homology search
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/08—Learning methods
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • General Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • General Health & Medical Sciences (AREA)
  • Mathematical Analysis (AREA)
  • Biophysics (AREA)
  • Health & Medical Sciences (AREA)
  • Computing Systems (AREA)
  • Databases & Information Systems (AREA)
  • Analytical Chemistry (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Biotechnology (AREA)
  • Algebra (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Medical Informatics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Chemical & Material Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Biomedical Technology (AREA)
  • Computational Linguistics (AREA)
  • Evolutionary Computation (AREA)
  • Molecular Biology (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本发明公开了一种基于图神经网络的解析单细胞通路活性的方法，该方法包括：基于单细胞多组学测序数据，利用深度神经网络自编码器学习基因‑细胞的潜在关联；利用图神经网络自编码器构建基因‑基因关联网络和细胞‑细胞关联网络；合并基因‑细胞的关联关系构建基因‑细胞关联网络；利用网络权重融合的方法来整合单细胞多组学数据产生多组学支持的基因‑细胞关联网络；利用重启动随机游走算法在单细胞单一组学或多组学的基因‑细胞关联网络推断通路在每个细胞中的活性分数和识别细胞表型关联的基因模块。本发明能高效、准确地解析单细胞的活性通路和基因模块，深度挖掘单细胞多组学数据的生物学机制并促进单细胞多组学的未来发展和应用。

Description

一种基于图神经网络的解析单细胞通路活性的方法
技术领域
本发明涉及单细胞多组学数据和图神经网络技术领域，尤其涉及一种基于图神经网络的解析单细胞多组学中活性通路的方法。
背景技术
细胞是生物体的基本结构与功能单位，人体内有200多种细胞，特定种类的细胞还可以呈现不同的亚型或功能状态，具有一定的异质性。单细胞测序技术能够检测单个细胞中的分子活性和调控状态，如单细胞RNA测序技术(scRNA-seq)能检测单个细胞的基因表达，单细胞染色质可及性测序技术(scATAC-seq)是利用Tn5转座酶研究单个细胞的染色质可进入性(可接近性)的高通量测序技术，能反映基因表达相关的染色质开放和调控状态。由多种单细胞测序技术产生的单细胞多组学数据的整合和分析，能够帮助系统地解析细胞的异质性和功能状态，在前所未有的精度下研究复杂的细胞调控事件以及细胞的交互作用，加深对生物系统的理解。然而，由于细胞的异质特性和单细胞测序技术的缺陷，单细胞测序谱图中往往包含大量的零值，其中部分零值是由于测序深度的限制导致的假阴性。这为后续的数据分析和生物学解释造成巨大的挑战。
生物通路是相互作用或功能相似的分子集合，可以对高通量测序数据进行功能注释，以解决固有的噪声并深入解析生物学机制。通路富集分析已广泛应用于不同的生物学领域，如分子机制研究、药物再利用、疾病分类和诊断。单细胞测序技术的发展为在单细胞水平整合不同组学数据来解析通路活性提供前所未有的机会。
几种基于Kolmogorov-Smirnov统计量的通路富集分析方法包括GSEA、ssGSE和GSVA，主要应用于大量细胞的RNA-seq数据。一些基准研究的结果表明，它们并不适用于单细胞数据分析。最近，一些针对scRNA-seq数据的通路富集方法被提出来研究细胞异质性，例如AUCell、Pagoda2和UniPath。AUCell计算通路在每个细胞的基因排序列表中的受试者工作特征曲线(ROC)下面积(AUC)分数作为通路活性分数。Pagoda2拟合一个模型来重新标准化基因表达谱，并使用第一加权主成分来量化通路活性分数。UniPath将基因表达的分布建模为双峰，并将非零的基因表达值转换为p值。通过布朗的方法来合并通路中基因的p值，并使用空背景模型将p值调整为通路富集分数。
虽然这三种单细胞通路富集方法可以估计单个细胞中的通路活性，但它们仍然存在局限性，难以从单细胞数据中深入挖掘潜在的信息。例如，AUCell依赖于基因的排序列表，这使得它一次只能识别与顶端基因相关的少数通路。Pagoda2只关注第一个主成分，会导致数据信息的丢失。UniPath需要针对不同物种构建背景模型，影响方法的可扩展性，并且背景模型的完整性直接影响其性能。对于具有许多零值的scRNA-seq数据，这些方法没有对假阴性的零表达基因进行推断。此外，UniPath提出超几何或二项式检验对scATAC-seq数据进行分析，但仍然需要依赖背景分布。而AUCell和Pagoda2只能对单细胞转录组数据进行通路富集分析。另一方面，这些方法只能分析预定义的通路或基因组，不能自主识别细胞表型相关的基因模块。细胞表型相关基因模块的分析研究有利于了解细胞状态转变过程中的机制和不同细胞表型之间的调控关系。
最近，一些新的单细胞测序方法能在同一细胞中获得不同类型的组学数据，例如CITE-seq、scNMT-seq、SNARE-seq、scNOMeRe-seq和NEAT-seq。在同一细胞内整合多种组学数据能够将不同层次的关系关联起来，从而更加准确地重建了分子过程。目前，已经出现了几种用于联合单细胞多组学数据的计算方法，例如，scJoint、Cobolt和scMVP。它们使用神经网络或自动编码器来压缩和联合来自不同单细胞组学的特征。然而，这些联合分析方法只是将多个单细胞组学的特征整合到一个共同的低维空间中以将细胞映射在一起，这导致它们具有很差的生物学解释性。而在通路和基因模块水平整合多组学数据具有更好的生物学解释，可以挖掘单个组学数据无法发现的潜在机制。但是，目前还没有方法能够整合单细胞多组学数据推断活性通路和细胞表型相关基因模块。
发明内容
本发明的目的在于提供一种基于图神经网络的解析单细胞多组学中活性通路的方法，能够对单细胞多组学数据进行整合分析，高度精确地解析单细胞中的生物通路活性和细胞表型关联的基因模块。
为实现上述目的，本发明采用的技术方案是：一种基于图神经网络的解析单细胞多组学中活性通路的方法，其包括以下步骤：
步骤(1)、对单细胞高通量测序数据预处理去除低质量的细胞，识别高变基因且高变基因被选择来构建基因-细胞矩阵；
步骤(2)、对于包含高变基因的基因-细胞矩阵，利用左截断高斯混合模型提取的转录调控信号正则化的深度神经网络学习细胞和基因在低维空间中的嵌入特征，推断基因-细胞之间的关联关系；
步骤(3)、对于包含高变基因的基因-细胞矩阵，利用皮尔森相关性构建细胞-细胞关联网络和基因-基因关联网络；
步骤(4)、将细胞和基因的低维嵌入特征分别作为细胞-细胞关联网络和基因-基因关联网络的节点特征，利用图卷积神经网络学习和推断细胞-细胞和基因-基因之间的关联关系；整合基因-细胞、细胞-细胞和基因-基因之间的关联关系构建基因-细胞关联网络；
步骤(5)、通过网络权重融合，整合单细胞不同组学的基因-细胞关联网络产生多组学支持的基因-细胞关联网络；
步骤(6)、利用重启动随机游走算法在单细胞单一组学或多组学的基因-细胞关联网络推断通路在每个细胞中的活性分数和识别细胞表型关联的基因模块。
作为本方案的进一步优选，所述步骤(1)中，对于单细胞高通量测序数据进行预处理，包含超过1％非零表达基因的细胞和在超过1％细胞中表达的基因被保留用于进一步分析。
作为本方案的进一步优选，采用基于全局尺度的标准化方法来标准化单细胞测序数据，其公式为：其中，Cij是第i个基因在第j个细胞中的基因读段计数，Cj是第j个细胞的所有基因的读段计数，f是尺度因子、一般设为1×104。
作为本方案的进一步优选，所述步骤(2)中，基于转录调控的动力学关系，左截断高斯混合模型从基因-细胞矩阵中提取转录调控信号；左截断高斯混合模型设置一个阈值Zcut将基因表达划分为两个部分：X＝{x1,…,xM}(X<Zcut)和X＝{x1,…,xN-M}(X≥Zcut)；K个高斯分布对应每个转录调控状态，公式为： 其中Θ是K个高斯分布，aj,μj,和σ分别是混合概率权重、均值和标准差；其中，利用期望最大化算法来评估Θ和Zcut，利用贝叶斯信息准则来确定高斯分布数量K；最终，基因在每个细胞中的表达值根据所属的高斯分布将其离散为转录调控信号。
作为本方案的更进一步优选，转录调控信号正则化的深度神经网络的编码器从基因-细胞矩阵Xmn(m个基因，n个细胞)中学习基因和细胞在低维空间中的特征，操作如下：
其中，Wl是第l层神经网络的可学习权重，和/>分别是基因和细胞的d维压缩特征矩阵,σ是非线性激活函数；
基于因式分解的深度神经网络的解码器来重构和推断基因-细胞关联矩阵，公式为：
转录调控信号正则化的均方误差作为损失函数，公式为： 其a是正则化权重，STRS∈Rm×n是转录调控信息号矩阵。
作为本方案的进一步优选，所述步骤(3)中，以经验的P值用来筛选强相关的基因-基因对，公式为：其中，Pi g是第i各基因与其他基因相关性值的向量，/>是第i个基因与第j个基因的相关性值，设定p-values<0.05来筛选强相关的基因对，并构建基因相关性网络。
作为本方案的进一步优选，所述步骤(4)中，图卷积神经网络定义为： 其中A为基因相关性网络的邻接矩阵，D是邻接矩阵的度矩阵，W0和W1是可学习的权重矩阵；
基于图卷积过程的变分自编码器来聚合节点的邻居特征，并重构基因相关性网络为基因关联网路，其编码器定义为：其中μi是来自μ＝GCNμ(E,A)的均值向量，σi是方差并且logσ＝GCNσ(E,A)，Z是基因关联网络在低维空间的嵌入特征矩阵，GCNμ(E,A)和GCNσ(E,A)共用图卷积的第一层权重矩阵W0；解码器来产生基因关联网络，公式为：/> 损失函数定义为：L＝Eq(Z|E,A)[logp(A|Z)]-KL[q(Z|E,A)||p(Z)]，KL是Kullback-Leibler散度。
作为本方案的进一步优选，所述步骤(5)中，对于来自同一细胞的单细胞多组学数据，通过图神经网络过程分别产生对应的基因-细胞关联网络；且基于Empirical Brown’smethod的网络权重融合的方法来整合来自scRNA-seq和scATAC-seq基因-细胞关联网络，产生一个多组学支持的基因-细胞关联网络。
作为本方案的进一步优选，所述步骤(6)中，对于单一组学或多组学的基因-细胞关联网络，将通路包含的基因对应的网络节点设为重启动节点及种子，利用重启动随机游走算法来计算通路在每个细胞中的活性分数，操作如下：pt+1＝(1-r)Wpt+rp0，其中p0是初始化概率向量，只有种子对应的概率为非零值，t是迭代次数，r是重启概率(r＝0.7)，当pt+1与pt的差异小于10-6，迭代停止，得到种子在每个细胞中的静态概率值。
作为本方案的更进一步优选，随机抽取N次相同数量的基因作为新的种子，重新迭代得到种子在每个细胞中扰动的静态概率值p′，并用来校正原始的静态概率值作为通路活性分数；操作如下：其中PASij是第i个通路在第j个细胞中的活性分数；相同的，属于相同表型的细胞设为种子，通过RWR和扰动分析校正来识别细胞表型关联的基因模块。
本发明的有益效果如下：
本发明是基于单细胞测序数据，构建了一个图神经网络框架，其包含深度神经网络自编码器和图卷积自编码器，来学习基因与细胞之间的潜在关联，将基因-细胞矩阵转换为稳定的基因-细胞关联网络。
本发明通过聚合相邻节点信息，将稀疏的不稳定的单细胞谱数据转换为稳定的基因-细胞关联关系。通过整合节点之间的关联效应来对假阴性的零表达或低表达的基因进行推断。对于单细胞多组学数据，网络融合方法备用来整合来自不同组学的基因-细胞关联网络。RWR算法进一步测量通路活性分数并识别细胞表型相关的基因模块。真实和模拟的单细胞数据集被用来对该方法的性能进行基准测试，并证明它在细胞聚类、活性通路推断和识别细胞表型相关基因模块方面具有优秀的性能。
应用本发明提供的方法或者系统，在模拟和真实的单细胞多组学数据进行性能评估。细胞聚类、通路活性打分和识别细胞表型关联的基因模块的准确性和稳健性方面取得了明显的提升。本发明作为一个便捷、低成本的计算框架，能高效、准确地大规模实解析单细胞多组学数据的活性通路和基因模块，并促进单细胞测序数据的深度解析和应用。
本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。
附图说明
图1为本发明实施提供的一种基于图神经网络的解析单细胞多组学中活性通路的方法的流程图。其中：(A)深度神经网络自编码器(Deep neural network autoencoder)和图神经网络自编码器(Graph autoencoder)学习单细胞数据中的潜在特征并构建基因-细胞关联网络，RWR算法来计算通路活性分数和识别细胞表型关联的基因模块；(B)本发明所具备的功能，包括将单细胞基因表达谱数据转换为通路活性分数矩阵、识别每个通路中关键的基因、细胞表型关联的基因模块和构建细胞群体网络。
图2为本发明实施提供的基于网络权重融合的整合单细胞多组学数据的流程图。
图3为本发明在真实scRNA-seq数据集评估细胞聚类精度。其中：(A)本方法与已报导的单细胞通路活性打分方法(AUCell、Pagoda2和UniPath)在人类多能干细胞的不同发育时间的scRNA-seq数据集中细胞的降维聚类图；(B)应用本方法和其他三种单细胞通路活性打分方法在20个scRNA-seq数据集中，用三种细胞聚类精度指标(平均的兰德指数，AverageARI；平均的标准互信息，Average NMI；平均的轮廓系数，Average SW)来比较本方法与其他方法的通路活性分数对细胞聚类的精度。
图4为本方法与已报导的单细胞通路活性打分方法在识别通路和细胞表型关联的基因模块的性能比较。其中：(A)本方法与其他三种通路打分方法在两个同质(Homogenous)和两个异质(Heterogenous)的scRNA-seq数据集中识别四种细胞(K562、A549、GM12878和ESC)对应的生物标志基因集在通路活性分数降序排列前五的细胞比例；(B)本方法与其他三种通路打分方法在包含B细胞或T细胞的同质和异质的scRNA-seq数据集中，检测B细胞受体信号通路和T细胞受体信号通路的活性分数在通路活性分数降序排列前五的B细胞或T细胞比例；(C)以K562、A549、GM12878和ESC细胞对应的生物标志基因集作为金标准，使用ROC分析本方法在同质和异质的scRNA-seq数据集中识别细胞表型关联的基因模块的准确性。
图5为本方法应用在scATAC-seq数据以及整合单细胞的scRNA-seq和scATAC-seq数据的性能评估。其中：(A)本方法应用在scATAC-seq数据识别各细胞类型对应的生物标志基因集在通路活性分数降序排列前五的细胞比例；(B)本方法应用在scATAC-seq数据识别T细胞受体信号通路在在通路活性分数降序排列前五的T细胞比例；(C)本方法整合细胞的scRNA-seq和scATAC-seq数据的三个细胞聚类指标条形图；(D)本方法在整合细胞的scRNA-seq和scATAC-seq数据后星形胶质细胞(Astrocyte)对应的生物标志基因集在通路活性分数降序排列前五的星形胶质细胞比例。
具体实施方式
下面结合附图和具体实施例对本发明进行详细说明。
如图1所示，为本发明提供的一种实施例：一种基于深度学习的空间蛋白质组学亚细胞定位的预测方法，包括步骤如下：
步骤(1)、单细胞高通量测序数据预处理去除低质量的细胞，识别高变基因；
步骤(2)、对于包含高变基因的基因-细胞矩阵，利用左截断高斯混合模型提取的转录调控信号正则化的深度神经网络学习细胞和基因在低维空间中的嵌入特征，推断基因-细胞之间的关联关系；
步骤(3)、对于包含高变基因的基因-细胞矩阵，利用皮尔森相关性构建细胞-细胞关联网络和基因-基因关联网络；
步骤(4)、将细胞和基因的低维嵌入特征分别作为细胞-细胞关联网络和基因-基因关联网络的节点特征，利用图卷积神经网络学习和推断细胞-细胞和基因-基因之间的关联关系；整合基因-细胞、细胞-细胞和基因-基因之间的关联关系构建基因-细胞关联网络；
步骤(5)、通过网络权重融合，整合单细胞不同组学的基因-细胞关联网络产生多组学支持的基因-细胞关联网络。
步骤(6)、利用重启动随机游走(RWR)算法在单细胞单一组学或多组学的基因-细胞关联网络推断通路在每个细胞中的活性分数和识别细胞表型关联的基因模块。
本实施例的步骤(1)中，对于单细胞测序数据(例如单细胞RNA数据或ATAC数据)进行预处理。包含超过1％非零表达基因的细胞和在超过1％细胞中表达的基因被保留用于进一步分析。一个基于全局尺度的标准化方法被用来标准化单细胞测序数据，其公式为：其中，Cij是第i个基因在第j个细胞中的基因读段计数，Cj是第j个细胞的所有基因的读段计数，f是尺度因子，一般设为1×104。高变基因被选择来构建基因-细胞矩阵。
本实施例的步骤(2)中，基于转录调控的动力学关系，左截短高斯混合模型从基因-细胞矩阵中提取转录调控信号。左截短高斯混合模型设置一个阈值Zcut将基因表达划分为两个部分，X＝{x1,…,xM}(X<Zcut)和X＝{x1,…,xN-M}(X≥Zcut)。K个高斯分布对应每个转录调控状态，公式为其中Θ是K个高斯分布，aj,μj,和σ分别是混合概率权重、均值和标准差。期望最大化(EM)算法用来评估Θ和Zcut。贝叶斯信息准则被用来确定高斯分布数量K。最终，基因在每个细胞中的表达值根据所属的高斯分布将其离散为转录调控信号。转录调控信号正则化的深度神经网络的编码器从基因-细胞矩阵Xmn(m个基因，n个细胞)中学习基因和细胞在低维空间中的特征，操作如下：
其中，Wl是第l层神经网络的可学习权重，和/>分别是基因和细胞的d维压缩特征矩阵,σ是非线性激活函数。基于因式分解的深度神经网络的解码器来重构和推断基因-细胞关联矩阵，公式为/>转录调控信号正则化的均方误差作为损失函数，公式为：/>其a是正则化权重，STRS∈Rm×n是转录调控信息号矩阵。
本实施例的步骤(3)中，基于基因-细胞矩阵，皮尔森相关性用来构建基因-基因相关性网络和细胞-细胞相关性网络。以基因-基因相关性网络的邻接矩阵为例，经验的P值用来筛选强相关的基因-基因对，公式为：/>其中Pi g是第i各基因与其他基因相关性值的向量，/>是第i个基因与第j个基因的相关性值。我们设定p-values<0.05来筛选强相关的基因对，并构建基因相关性网络。相同的过程应用于细胞中，可构建细胞相关性网络。
本实施例的步骤(4)中，将步骤(2)得到的转录调控信号正则化的深度神经网络编码的基因和细胞的低维特征分别作为步骤(3)得到的基因相关性网络和细胞相关性网络的节点特征，通过图卷积神经网络来构建基因关联网络和细胞关联网络。以构建基因相关性网络为例，一个两层的图卷积神经网络被定义： 其中A为基因相关性网络的邻接矩阵，D是邻接矩阵的度矩阵，W0和W1是可学习的权重矩阵。基于图卷积过程的变分自编码器来聚合节点的邻居特征，并重构基因相关性网络为基因关联网路，其编码器被定义为：/>其中μi是来自μ＝GCNμ(E,A)的均值向量，σi是方差并且logσ＝GCNσ(E,A)，Z是基因关联网络在低维空间的嵌入特征矩阵，GCNμ(E,A)和GCNσ(E,A)共用图卷积的第一层权重矩阵W0。解码器来产生基因关联网络，公式为：
损失函数被定义为：L＝Eq(Z|E,A)[logp(A|Z)]-KL[q(Z|E,A)||p(Z)]，KL是Kullback-Leibler散度。基因-细胞、细胞-细胞和基因-基因之间的关联关系被整合为一个基因-细胞关联网络。
本实施例的步骤(5)中，对于来自同一细胞的单细胞多组学数据，例如scRNA-seq和scATAC-seq，通过上述步骤的图神经网络过程分别产生对应的基因-细胞关联网络。如图2所示，基于Empirical Brown’s method的网络权重融合的方法来整合来自scRNA-seq和scATAC-seq基因-细胞关联网络，产生一个多组学支持的基因-细胞关联网络。
本实施例的步骤(6)中，对于单一组学或多组学的基因-细胞关联网络，将通路包含的基因对应的网络节点设为重启动节点及种子，重启动的随机游走用来计算通路在每个细胞中的活性分数，操作如下：pt+1＝(1-r)Wpt+rp0，其中p0是初始化概率向量，只有种子对应的概率为非零值，t是迭代次数，r是重启概率(r＝0.7)。当pt+1于pt的差异小于10-6，迭代停止，得到种子在每个细胞中的静态概率值。为了去除随机噪声的影响，随机抽取N次相同数量的基因作为新的种子，重新迭代得到种子在每个细胞中扰动的静态概率值p′并用来校正原始的静态概率值作为通路活性分数，操作如下：其中PASij是第i个通路在第j个细胞中的活性分数。相同的，属于相同表型的细胞设为种子，通过重启动随机游走和扰动分析校正来识别细胞表型关联的基因模块。
为了评估本发明对细胞聚类的性能，一个包含人类多能干细胞的不同发育时间的scRNA-seq数据集(GSE75748)被用来计算通路活性分数。如图3A所示，相对于三种已报导的单细胞通路活性打分方法(AUCell、Pagoda2和UniPath)，本方法计算的通路活性分数能够将处于不同时间的细胞区分开，相同时间的细胞聚集的更加紧密。接下来，20个scRNA-seq数据集(表1)、10种单细胞聚类方法(表2)和三个细胞聚类精度指标(Average ARI、AverageNMI和Average SW)被用来系统评估本方法的准确性和稳定性。如图3B所示，在多重数据和细胞聚类方法中，相较于其他三种通路活性打分方法，本方法在三个细胞聚类精度指标中表现出更加精确和稳定地对细胞进行聚类。
为了评估本发明识别通路和细胞关联的基因模块的精度，两个同质(Homogenous)和两个异质(Heterogenous)的scRNA-seq数据集被用来计算通路活性分数。将四种细胞(K562、A549、GM12878和ESC)的通路活性分数进行降序排列并统计这四种细胞相应的已知的生物标志基因集排列在前五的细胞比例。如图4A所示，相较于其他三种通路活性打分方法，本发明无论在同质还是异质性数据集中都能够精确地识别生物标志基因集。B细胞受体信号通路和T细胞受体信号通路是已知的分别在B细胞和T细胞中激活的生物通路。因此，在包含B细胞或T细胞的同质和异质性数据集中，检测本方法和其他三种方法对B细胞受体信号通路和T细胞受体信号通路的打分在通路活性分数降序排列中前五的在B细胞和T细胞的比例。如图4B所示，本发明相比其他三种通路活性打分方法能够更精准地识别B细胞和T细胞的激活的生物标志通路。本方法识别K562、A549、GM12878和ESC细胞关联的基因模块，将这四种细胞对应的生物标志基因集作为金标准，ROC分析表明本方法识别基因模块与这四种细胞表型具有紧密的关联关系(如图4C所示)。进一步地，将由Seurat工具在细胞类型(Cell type,GSE84133)、细胞亚型(Cell subtype,GSE76381)和时序(Time series,GSE75748)数据集中识别的生物标志基因作为金标准进行ROC分析，本方法同样地表现出精准地细胞表型关联的基因模块的识别能力。
如图5A和B所示，对于来自10x Genomics测序平台的人类外周血的scATAC-seq数据集，相对于UniPath的两种针对scATAC-seq的通路打分方法(Binomial和Hypergeometeric)，本方法能够更加精确地识别各细胞类型的生物标志基因集和T细胞的T细胞受体信号通路。本方法通过网络权重融合的方法整合来自成年小鼠大脑的scRNA-seq和scATAC-seq数据集(GSE126074)。如图5C所示，三种细胞聚类精度指标表明本方法在整合单细胞多组学数据后能够更好的对细胞进行聚类。同时，如图5D所示，在整合单细胞多组学后依然能够精确地识别星形胶质细胞对应的生物标志基因集。如表3所示，整合单细胞多组学数据后，本方法能够识别到单一组学无法识别到的生物通路，并且已被研究所证实。
表1
表2
表3
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分步骤是可以通过程序来指令相关的硬件来完成，所述的程序可以存储于一计算机可读取存储介质中，所述的存储介质，如ROM/RAM、磁盘、光盘等。
以上显示和描述了本发明的基本原理、主要特征和优点。本领域的普通技术人员应该了解，上述实施例不以任何形式限制本发明的保护范围，凡采用等同替换等方式所获得的技术方案，均落于本发明的保护范围内。
本发明未涉及部分均与现有技术相同或可采用现有技术加以实现。

Claims (5)

1.一种基于图神经网络的解析单细胞通路活性的方法，其特征在于，包括以下步骤：
步骤(1)、对单细胞高通量测序数据预处理去除低质量的细胞，识别高变基因且高变基因被选择来构建基因-细胞矩阵；
步骤(2)、对于包含高变基因的基因-细胞矩阵，利用左截断高斯混合模型提取的转录调控信号正则化的深度神经网络学习细胞和基因在低维空间中的嵌入特征，推断基因-细胞之间的关联关系；具体的：
基于转录调控的动力学关系，左截断高斯混合模型从基因-细胞矩阵中提取转录调控信号；左截断高斯混合模型设置一个阈值Zcut将基因表达划分为两个部分：X＝{x1,…,xM},X<Zcut和X＝{xM,…,xN-M},X≥Zcut；K个高斯分布对应每个转录调控状态，公式为：其中Θ是K个高斯分布，ai,μi,和σi分别是混合概率权重、均值和标准差；其中，利用期望最大化算法来评估Θ和Zcut，利用贝叶斯信息准则来确定高斯分布数量K；最终，基因在每个细胞中的表达值根据所属的高斯分布将其离散为转录调控信号；
转录调控信号正则化的深度神经网络的编码器从基因-细胞矩阵X中学习基因和细胞在低维空间中的特征，操作如下：
其中，m个基因，n个细胞，Wl是第l层神经网络的可学习权重，和/>分别是基因和细胞的d维压缩特征矩阵,σ是非线性激活函数；
基于因式分解的深度神经网络的解码器来重构和推断基因-细胞关联矩阵，公式为：
转录调控信号正则化的均方误差作为损失函数，公式为： 其中a是正则化权重，STRS∈Rm×n是转录调控信息号矩阵；
步骤(3)、对于包含高变基因的基因-细胞矩阵，利用皮尔森相关性构建细胞-细胞关联网络和基因-基因关联网络；
步骤(4)、将步骤(2)得到的转录调控信号正则化的深度神经网络编码的基因和细胞的低维特征分别作为步骤(3)得到的细胞-细胞关联网络和基因-基因关联网络的节点特征，通过图卷积神经网络来构建基因关联网络和细胞关联网络；
步骤(5)、通过网络权重融合，整合单细胞不同组学的基因-细胞关联网络产生多组学支持的基因-细胞关联网络；
步骤(6)、利用重启动随机游走算法在单细胞单一组学或多组学的基因-细胞关联网络推断通路在每个细胞中的活性分数和识别细胞表型关联的基因模块；具体的：
对于单一组学或多组学的基因-细胞关联网络，将通路包含的基因对应的网络节点设为重启动节点及种子，利用重启动随机游走算法来计算通路在每个细胞中的活性分数，操作如下：pt+1＝(1-r)Wpt+rp0，其中p0是初始化概率向量，只有种子对应的概率为非零值，t是迭代次数，r是重启概率r＝0.7，当pt+1与pt的差异小于10-6，迭代停止，得到种子在每个细胞中的静态概率值；
随机抽取N次相同数量的基因作为新的种子，重新迭代得到种子在每个细胞中扰动的静态概率值p′，并用来校正原始的静态概率值作为通路活性分数；操作如下：其中PASij是第i个通路在第j个细胞中的活性分数；相同的，属于相同表型的细胞设为种子，通过RWR和扰动分析校正来识别细胞表型关联的基因模块。
2.根据权利要求1所述的一种基于图神经网络的解析单细胞通路活性的方法，其特征在于，所述步骤(1)中，对于单细胞高通量测序数据进行预处理，包含超过1％非零表达基因的细胞和在超过1％细胞中表达的基因被保留用于进一步分析。
3.根据权利要求1或2所述的一种基于图神经网络的解析单细胞通路活性的方法，其特征在于，采用基于全局尺度的标准化方法来标准化单细胞测序数据，其公式为：其中，Ckj是第k个基因在第j个细胞中的基因读段计数，Cj是第j个细胞的所有基因的读段计数，f是尺度因子、设为1×104。
4.根据权利要求1所述的一种基于图神经网络的解析单细胞通路活性的方法，其特征在于，所述步骤(3)中，以经验的P值用来筛选强相关的基因-基因对，公式为：其中，/>是第k个基因与其他基因相关性值的向量，/>是第k个基因与第z个基因的相关性值，设定p-values<0.05来筛选强相关的基因对，并构建基因相关性网络。
5.根据权利要求1所述的一种基于图神经网络的解析单细胞通路活性的方法，其特征在于，所述步骤(5)中，对于来自同一细胞的单细胞多组学数据，通过图神经网络过程分别产生对应的基因-细胞关联网络；且基于网络权重融合的方法来整合来自scRNA-seq和scATAC-seq基因-细胞关联网络，产生一个多组学支持的基因-细胞关联网络。
CN202211005236.5A 2022-08-22 2022-08-22 一种基于图神经网络的解析单细胞通路活性的方法 Active CN115240772B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211005236.5A CN115240772B (zh) 2022-08-22 2022-08-22 一种基于图神经网络的解析单细胞通路活性的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211005236.5A CN115240772B (zh) 2022-08-22 2022-08-22 一种基于图神经网络的解析单细胞通路活性的方法

Publications (2)

Publication Number Publication Date
CN115240772A CN115240772A (zh) 2022-10-25
CN115240772B true CN115240772B (zh) 2023-08-22

Family

ID=83680859

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211005236.5A Active CN115240772B (zh) 2022-08-22 2022-08-22 一种基于图神经网络的解析单细胞通路活性的方法

Country Status (1)

Country Link
CN (1) CN115240772B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115631799B (zh) * 2022-12-20 2023-03-28 深圳先进技术研究院 样本表型的预测方法、装置、电子设备及存储介质
CN116452910B (zh) * 2023-03-28 2023-11-28 河南科技大学 基于图神经网络的scRNA-seq数据特征表示和细胞类型识别方法
CN116504314B (zh) * 2023-06-27 2023-08-29 华东交通大学 基于细胞动态分化的基因调控网络构建方法
CN117854592A (zh) * 2024-03-04 2024-04-09 中国人民解放军国防科技大学 一种基因调控网络构建方法、装置、设备、存储介质
CN117854600A (zh) * 2024-03-07 2024-04-09 北京大学 基于多组学数据的细胞识别方法、装置、设备及存储介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114022693A (zh) * 2021-09-29 2022-02-08 西安热工研究院有限公司 一种基于双重自监督的单细胞RNA-seq数据聚类方法
CN114783526A (zh) * 2022-05-11 2022-07-22 南开大学 基于高斯混合图变分自编码器的深度无监督单细胞聚类方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10810213B2 (en) * 2016-10-03 2020-10-20 Illumina, Inc. Phenotype/disease specific gene ranking using curated, gene library and network based data structures

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114022693A (zh) * 2021-09-29 2022-02-08 西安热工研究院有限公司 一种基于双重自监督的单细胞RNA-seq数据聚类方法
CN114783526A (zh) * 2022-05-11 2022-07-22 南开大学 基于高斯混合图变分自编码器的深度无监督单细胞聚类方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
scGNN is a novel graph neural network framework for single- cell RNA- Seq analyses;Wang, J. et al;《Nature Communications》;第12卷(第1期);1-11 *

Also Published As

Publication number Publication date
CN115240772A (zh) 2022-10-25

Similar Documents

Publication Publication Date Title
CN115240772B (zh) 一种基于图神经网络的解析单细胞通路活性的方法
Fan et al. An assembly and alignment-free method of phylogeny reconstruction from next-generation sequencing data
Toh et al. Inference of a genetic network by a combined approach of cluster analysis and graphical Gaussian modeling
Celton et al. Comparative analysis of missing value imputation methods to improve clustering and interpretation of microarray experiments
Patruno et al. A review of computational strategies for denoising and imputation of single-cell transcriptomic data
O’Connor The distribution of common-variant effect sizes
Horimoto et al. Statistical estimation of cluster boundaries in gene expression profile data
Maulik et al. Simulated annealing based automatic fuzzy clustering combined with ANN classification for analyzing microarray data
US20050288871A1 (en) Estimating the accuracy of molecular property models and predictions
Rasheed et al. Metagenomic taxonomic classification using extreme learning machines
Lopez et al. Bayesian inference for a generative model of transcriptome profiles from single-cell RNA sequencing
Biswas et al. Mapping gene expression quantitative trait loci by singular value decomposition and independent component analysis
Celik et al. Biological cartography: Building and benchmarking representations of life
Kaski et al. Associative clustering for exploring dependencies between functional genomics data sets
Toh et al. System for automatically inferring a genetic netwerk from expression profiles
Gorin et al. Distinguishing biophysical stochasticity from technical noise in single-cell RNA sequencing using Monod
Shaik et al. A unified framework for finding differentially expressed genes from microarray experiments
Tarazona et al. Variable selection for multifactorial genomic data
Liu et al. Assessing agreement of clustering methods with gene expression microarray data
Zhou et al. Category encoding method to select feature genes for the classification of bulk and single‐cell RNA‐seq data
Bhyratae et al. Design And Development Of Advanced Similarity Measure for Reconstructing GRN using mRNA Expression Profiles
Wang et al. Constructing biological pathways by a two-step counting approach
Xin-guo et al. Gene cluster algorithm based on most similarity tree
Jayanetti Statistical Methods for Meta-Analysis in Large-Scale Genomic Experiments
Deng Algorithms for reconstruction of gene regulatory networks from high-throughput gene expression data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant