CN117457065A - 一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统 - Google Patents

一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统 Download PDF

Info

Publication number
CN117457065A
CN117457065A CN202311080320.8A CN202311080320A CN117457065A CN 117457065 A CN117457065 A CN 117457065A CN 202311080320 A CN202311080320 A CN 202311080320A CN 117457065 A CN117457065 A CN 117457065A
Authority
CN
China
Prior art keywords
cell
data
cell type
phenotype
regulatory
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.)
Pending
Application number
CN202311080320.8A
Other languages
English (en)
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.)
Academy of Mathematics and Systems Science of CAS
Original Assignee
Academy of Mathematics and Systems Science of CAS
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 Academy of Mathematics and Systems Science of CAS filed Critical Academy of Mathematics and Systems Science of CAS
Priority to CN202311080320.8A priority Critical patent/CN117457065A/zh
Publication of CN117457065A publication Critical patent/CN117457065A/zh
Pending legal-status Critical Current

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
    • G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/22—Matching criteria, e.g. proximity measures
    • 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
    • G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • G16B20/40—Population genetics; Linkage disequilibrium
    • 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
    • G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding

Landscapes

  • Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Theoretical Computer Science (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Medical Informatics (AREA)
  • Evolutionary Biology (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Biophysics (AREA)
  • Data Mining & Analysis (AREA)
  • Biotechnology (AREA)
  • Genetics & Genomics (AREA)
  • General Health & Medical Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Chemical & Material Sciences (AREA)
  • Evolutionary Computation (AREA)
  • Analytical Chemistry (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Molecular Biology (AREA)
  • Databases & Information Systems (AREA)
  • Bioethics (AREA)
  • Physiology (AREA)
  • Epidemiology (AREA)
  • Ecology (AREA)
  • Public Health (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本发明涉及一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统。该方法包括：获取人类的单细胞多组学数据；利用获取的单细胞多组学数据建立细胞类型层面的匹配数据；对于每个细胞类型，以建立的匹配数据为输入构建其调控网络；利用每个细胞类型的调控网络构建调控功能区域，在所有细胞类型的调控功能区域上对遗传力进行建模，并定义表型与细胞类型的关联分数；度量表型与细胞类型的关联分数的显著性，利用表型与细胞类型的关联分数的显著性识别表型相关的细胞类型。本发明能够更有效的基于单细胞数据构建细胞类型层面匹配多组学数据，更有效的构建细胞类型的调控网络，并更准确的估计表型和细胞类型之间关联分数及其显著性。

Description

一种基于单细胞多组学数据识别表型相关细胞类型的方法和 系统
技术领域
本发明属于医学技术、信息技术领域，具体涉及一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统。
背景技术
精准医学是一种将个体基因型与环境因素差异考虑在内的疾病预测、预防、干预和治疗的新策略，旨在对疾病进行精确分类及诊断，为患者提供个性化、针对性的预防和治疗措施。个性化医疗的一个重要环节是在正确的细胞环境中，例如器官、组织甚至是更加具体的细胞类型，进行精准干预，进而改善表型或者治疗疾病。但是，对于复杂表型或者疾病，表型相关的细胞环境是不知道的，这使得正确识别表型相关的细胞环境以实现更加精准的预防或者治疗成为精准医学研究和实践的一个重要的研究热点。
表型和疾病是由基因型和环境因素共同的复杂作用形成的。但是，微观分子层面的基因型并不是直接对宏观层面的表型和疾病产生作用的，而是通过多层次分子复杂的调控完成的。具体而言，基因型及其遗传变异会在调控元件上定位，进而对调控网络产生影响，而由于调控网络是在时空特异的细胞环境中才会有活性和功能，因此基因型和遗传变异是在特定的细胞环境中才会对调控产生影响进而影响表型和疾病。这里的细胞环境，既包含相对宏观层次上的器官或者组织，又包含相对微观层次的细胞类型。研究人员已经发现了大量基因型和遗传变异在特定细胞环境中影响表型的实例。特别的，一些表型相关的SNP(Single Nucleotide Polymorphisms，单核苷酸多态性)在特定的发育阶段起作用，例如与智力相关的表型的SNP，在人脑发育的早期的染色质开放区域更加富集。一些遗传变异仅与特定的组织器官有关，例如，因为乳腺组织对性腺激素的敏感性，性腺激素受体(ESR1)基因的突变仅在乳腺作用，突变的ESR1基因会导致对雌激素的反应性发生变化，增加乳腺癌的风险。有很多遗传变异与特定的细胞类型有关，例如EWS-FLI1基因融合是原始神经外胚层瘤的典型遗传变异，该基因的融合只在神经干细胞中发挥作用，并导致该肿瘤的形成。这些实例说明，识别基因型和遗传变异影响表型起作用的细胞环境，尤其是细胞类型，是实现表型干预和疾病治疗的精准医学的重要步骤。
目前，已经有了许多方法去识别表型相关的细胞环境，但大多在组织器官的层次上进行识别。例如，RolyPoly模型是基于各个组织的表达数据构建的，它用特异表达基因在各个组织中的表达值去回归SNP效应分数，进而得到每个组织和该表型的关联性；CoCoNet将特异性扩展到网络层次，它首先利用GETx中各个组织的表达数据识别每个组织特异的基因共表达网络，然后基于基因共表达网络识别表型相关的组织；S-LDSC基于基因组注释区域或者多组织的染色质开放区域的功能区域建立SNP的效应分数和其LD分数之间的回归方程计算表型的遗传力在每个功能区域的富集分数，作为表型和组织的关联性；S-LDSC随后被应用到多组织的表达数据中，LDSC-SEG利用组织的特异表达基因的位置构建功能区域，并计算遗传力富集分数，从而识别表型相关的组织。但是这些方法的缺陷是均基于群体细胞(Bulk)数据，只能识别组织器官层次上的表型相关的细胞环境，例如将糖尿病表型和胰腺组织相关联，但无法关联到胰腺内部多种细胞类型，例如胰岛的alpha细胞类型或者是beta细胞类型。
细胞类型层次上的细胞环境对表型遗传力的的解释能力更强。随着单细胞数据的出现和增长，最近也出现了一些识别表型相关细胞类型的方法。例如LDSC-SEG或者MAGMA可以很直接地扩展到单细胞数据中，通过单细胞数据识别出的细胞类型以及细胞类型特异表达的基因的位置作为调控功能区域，可以识别相关细胞类型。但是此类方法的缺点是只考虑了单细胞基因表达数据。为了考虑染色质开放数据，scLinker被开发出来，它首先利用单细胞基因表达数据识别细胞类型特异的基因程序，随后利用公开数据库中的增强子和基因联系识别其相关的调控元件，随后对特异基因程序及其增强子构成的功能区域计算遗传力富集分数，进而识别表型相关的细胞类型。但是此类方法只使用了公共的染色质开放数据，并没有考虑与单细胞表达数据匹配的单细胞染色质开放数据，利用单细胞多组学数据中蕴含的更丰富的信息，从构建调控网络图谱的角度来集成数据。
目前技术主要缺点：
(1)大部分识别表型相关细胞环境的方法均是基于群体细胞的多组学数据，这种数据大多在组织器官层次获得，因此只能识别组织器官层次上的细胞环境。
(2)单细胞数据的信息挖掘不充分。现阶段虽然有一些识别细胞类型的方法，但是其使用的数据往往只是一个局部范围的单细胞数据，并没有利用当前快速发展的单细胞图谱。
(3)当前对单细胞数据的使用，只使用了单细胞RNA-seq数据，而没有考虑使用同样快速增长的匹配的单细胞ATAC-seq数据，以及整合单细胞多组学数据的调控网络。
(4)缺乏系统性集细胞类型匹配数据构造、细胞类型调控网络构建、识别表型相关相关细胞类型的平台。
发明内容
为了突破上述难点和痛点，本发明基于大量的单细胞RNA-seq和单细胞ATAC-seq数据及其定义出的细胞类型，利用每个细胞类型的单细胞数据构建其调控网络，在细胞类型特异的调控网络上对遗传力富集进行建模并识别表型相关的细胞类型。
本发明采用的技术方案如下：
一种基于单细胞多组学数据识别表型相关细胞类型的方法，包括以下步骤：
获取人类的单细胞多组学数据；
利用获取的单细胞多组学数据建立细胞类型层面的匹配数据；
对于每个细胞类型，以建立的匹配数据为输入构建其调控网络；
利用每个细胞类型的调控网络构建调控功能区域，在所有细胞类型的调控功能区域上对遗传力进行建模，并定义表型与细胞类型的关联分数；
度量表型与细胞类型的关联分数的显著性，利用表型与细胞类型的关联分数的显著性识别表型相关的细胞类型。
进一步地，所述获取人类的单细胞多组学数据，包括：获取人类的多个组织或者器官的匹配的单细胞RNA-seq和单细胞ATAC-seq数据。
进一步地，所述建立细胞类型层面的匹配数据，包括：
对于匹配的单细胞RNA-seq和单细胞ATAC-seq数据，标注其细胞类型；
通过Pseudo-bulk策略构建每个细胞类型的匹配的基因表达和染色质开放数据。
进一步地，所述通过Pseudo-bulk策略构建每个细胞类型的匹配的基因表达和染色质开放数据，包括：
对于一个细胞类型，给定来自其单细胞RNA-seq数据的包含c1个细胞和N个基因的单细胞基因表达矩阵E，以及来自单细胞ATAC-seq数据的包含c2个细胞和P个peak的单细胞染色质开放矩阵O，分别合并c1个和c2个细胞以得到第i个基因的表达PEi和第k个peak的开放性POk：
其中，e是单细胞RNA-seq数据的的RNA片段总数目，o是单细胞ATAC-seq数据的DNA片段总数目，Eic是细胞c中第i个基因的表达，Okc是细胞c中第k个peak的开放性；通过POk≥2的阈值选择开放的peak，从而构建该细胞类型的匹配的基因表达和染色质开放数据。
进一步地，所述调控网络以TF-RE-TG三元组为基本单元，其中TF表示转录因子，TG表示目的基因，RE表示调控元件；将匹配的基因表达和染色质开放数据输入到PECA2模型中构建所述调控网络。
进一步地，所述利用每个细胞类型的调控网络构建调控功能区域，是对于每个细胞类型i，其调控网络中的RE集合组成该细胞类型的调控功能区域Ci；所述表型与细胞类型的关联分数通过整合遗传力富集分数τi及其统计显著性P值pi来定义：
Ri＝τi·(-log pi)
其中，Ri表示表型与第i个细胞类型的关联分数；τi表示表型在Ci中的遗传力富集分数；pi表示遗传力富集分数τi的显著性P值。
进一步地，使用Block jackknife的方法来估计关联分数Ri的标准误差、P值和FDR。
一种基于单细胞多组学数据识别表型相关细胞类型的系统，其包括：
数据获取模块，用于获取人类的单细胞多组学数据；
数学建模模块，用于利用获取的单细胞多组学数据建立细胞类型层面的匹配数据；对于每个细胞类型，以建立的匹配数据为输入构建其调控网络；利用每个细胞类型的调控网络构建调控功能区域，在所有细胞类型的调控功能区域上对遗传力进行建模，并定义表型与细胞类型的关联分数；度量表型与细胞类型的关联分数的显著性，利用表型与细胞类型的关联分数的显著性识别表型相关的细胞类型；
模型测试模块，用于获取给定表型或者疾病的全基因组关联分析GWAS统计汇总数据，对数学建模模块得到的模型进行测试。
本发明公开了一种基于单细胞RNA-seq和单细胞ATAC-seq数据，构建细胞类型层次的调控网络图谱，利用图谱识别表型相关细胞类型的方法，其优点和有益效果包括：
1.提出了更有效的基于单细胞数据构建细胞类型层面匹配多组学数据的方法。
2.提出了更有效的构建细胞类型的调控网络的方法。
3.提出了更准确的估计表型和细胞类型之间关联分数及其显著性的方法。
附图说明
图1是表型相关细胞类型识别的计算系统结构图。
图2是本发明和scLinker在测试集上的遗传力富集比较图。
具体实施方式
为使本发明的上述目的、特征和优点能够更加明显易懂，下面通过具体实施例和附图，对本发明做进一步详细说明。
本发明要解决的技术问题是提供一种识别表型相关细胞类型的方法。该方法基于高通量测序技术获得的大量的单细胞RNA-seq和单细胞ATAC-seq数据，定义细胞类型，构造每个细胞类型的匹配基因表达和染色质开放数据，推断构建其调控网络，形成多细胞类型的调控网络图谱，最终基于调控网络图谱识别给定表型的相关细胞类型。需要突破的难点包括：细胞类型上匹配数据的构造；基于单细胞数据的细胞类型的调控网络的构建；度量遗传力在调控网络中的分布及其显著性等。为突破这些技术瓶颈，本发明提出了基于单细胞数据构建调控网络图谱来识别表型相关细胞类型新方法，方法的核心是定义单个细胞类型的调控网络与表型GWAS的关联性评分，基于此关联性评分，搭建了识别表型相关细胞类型的完整计算流程和系统。
为了克服现有技术的不足，本发明提出了如图1所示的基于调控网络识别表型相关细胞类型的系统。该系统包括三个模块：
(1)数据获取模块
该模块需要收集的数据分为两部分：
(1.1)单细胞数据：作为全面的细胞类型层次的调控网络图谱的数据基础，需要收集人类全面的单细胞RNA-seq和单细胞ATAC-seq数据。在这里，收集的两个层次的数据为匹配的，意思是这两种类型的单细胞数据来自同样的组织、器官或者样本。
(1.2)对于给定表型或者疾病的全基因组关联分析(GWAS)统计汇总数据。
(2)数学建模模块
该模块将以(1)中收集的数据为输入构建模型。该模块分为以下四个步骤：
(2.1)建立细胞类型层面的匹配数据：利用(1.1)中收集的匹配的单细胞RNA-seq和单细胞ATAC-seq数据，利用已有的知识或者计算工具注释细胞类型。然后通过Pseudo-bulk方法构建每个细胞类型的匹配的基因表达和染色质开放数据。
(2.2)细胞类型调控网络的构建：对于每个细胞类型，以(2.1)中构建的匹配数据为输入，构建其调控网络，调控网络以若干TF-RE-TG三元组为基本单元。其中TF表示转录因子，TG表示目的基因，RE表示调控元件。
(2.3)构建细胞类型的调控功能区域并对遗传力进行建模：对于每个细胞类型的调控网络，使用调控网络中的调控元件构建调控功能区域。在所有细胞类型的调控功能区域上，对遗传力进行建模，定义表型与细胞类型的关联分数。
(2.4)度量表型与细胞类型的关联分数的显著性：对于每个细胞类型与表型的关联分数，使用block jackknife方法估计该分数的显著性P值，用来识别细胞类型。
(3)模型测试模块
在(1.2)中收集了六种表型的GWAS进行验证。这六种表型具有其相关细胞类型的先验知识，用以验证并与其他方法进行比较。
本发明所提出的识别表型相关细胞类型的模型分为三个模块，如图1所示，下面针对数据获取、模型构建与求解、模型测试进行详细阐述。
(1)数据获取
从已发表的文献中收集了人类的26个组织或者器官的匹配的单细胞RNA-seq和单细胞ATAC-seq数据，如表1所示，将以此为输入构建细胞类型调控网络图谱。
对于给定的表型，收集该表型的GWAS统计汇总数据。在这里，共收集了六个血液指标的表型(红细胞数目、红细胞分布、单核细胞数目、单核细胞比例、淋巴母细胞数目和淋巴母细胞比例)数据，用于验证和比较。
表1.本发明收集的26个人类组织或者器官的单细胞多组学数据
Organs Stage Organs Stage
Adrenal gland Fetal Kidney Adult
Adrenal gland Adult Liver Fetal
Artery-Aorta Adult Liver Adult
Artery-Tibial Adult Lung Fetal
Cerebellum Fetal Lung Adult
Cerebrum Fetal Muscle Fetal
Colon Adult Pancreas Fetal
Esophagus gastroesophageal Adult Pancreas Adult
Esophagus mucosa Adult Placenta Fetal
Esophagus muscularis Adult Skin Adult
Eye Fetal Small intestine Adult
Frontal cortex Adult Spleen Fetal
Heart Fetal Stomach Fetal
Heart Adult Stomach Adult
Intestine Fetal Thymus Fetal
Kidney Fetal Uterus Adult
(2.1)建立细胞类型层面的匹配数据
从(1)中得到单细胞RNA-seq和单细胞ATAC-seq数据之后，可以借助已有的计算工具对数据进行聚类以标注细胞类型。细胞类型注释可以来自下面三种方式：(a)直接使用生物学家发表数据时定义好的细胞类型；(b)在单细胞RNA-seq和单细胞ATAC-seq分别使用Seurat(软件网址：https://satijalab.org/seurat/)等工具进行聚类得到细胞类型，然后将两层数据的细胞类型对应起来。(c)使用CoupledNMF(Z Duren,X Chen,M Zamanighomi,WZeng,A Satpathy,H Chang,Yong Wang,Wing Hung Wong.Integrative analysis ofsingle cell genomics data by coupled nonnegative matrix factorizations.ProcNatl Acad Sci U S A 115(30),7723-7728,2018.)等方法在两层数据上统一定义细胞类型。这样，每个细胞类型均关联单细胞RNA测序数据的一组细胞和单细胞ATAC测序数据的一组细胞。
细胞类型在两层数据上定义完成之后，接下来使用“Pseudo-bulk”策略建立细胞类型层面的匹配的表达和染色质开放数据。以一个细胞类型为例，给定来自其单细胞RNA-seq数据的包含c1个细胞和N个基因的单细胞基因表达矩阵E，以及来自单细胞ATAC-seq数据的包含c2个细胞和P个peak(峰，代表染色质上开放的区域)的单细胞染色质开放矩阵O，分别合并c1个和c2个细胞以得到“pseudo-bulk”的第i个基因的表达PEi和第k个peak的开放性POk：
这里e是单细胞RNA-seq数据的的RNA片段总数目，o是单细胞ATAC-seq数据的DNA片段总数目。Eic是细胞c中第i个基因的表达，并且Okc是细胞c中第k个peak的开放性。通过POk≥2的阈值选择开放的peak。这样构建了这个细胞类型的匹配的基因表达和染色质开放数据。
(2.2)细胞类型调控网络的构建
然后，同样对于上述的细胞类型，将(2.1)中匹配的“pseudo-bulk”基因表达和染色质开放数据输入到PECA2模型(Duren,Z.et al.(2020)‘Time Course RegulatoryAnalysis based on paired expression and chromatin accessibility data’,GenomeResearch,30(4),pp.622–634.doi:10.1101/gr.257063.119.)中以构建调控网络。PECA2计算反式调控分数来衡量第i个转录因子(TF)对第j个目的基因(TG)的调控强度。它假设TF通过在调控元件(RE)上结合来调节下游TG，通过整合多个由TF结合的RE调节TG来计算反式调控分数。具体的，第i个TF和第j个TG的反式调控分数Rij定义为：
这里PEi和PEj是第i个TF和第j个TG的“pseudo-bulk”基因表达；Bik是第i个TF结合在第k个RE的强度；POk是第k个peak的“pseudo-bulk”开放程度。Ikj表示第k个RE和第j个TG之间的相互作用强度，来自文献“Duren,Z.et al.Integrative analysis of single-cellgenomics data by coupled nonnegative matrix factorizations.Proc Natl Acad SciU S A 115,7723-7728(2018).”；Dij是第i个TF和第j个TG在不同DNA百科全书计划(ENCODE)样本中的表达相关性，也来自该文献。PECA2的输出是TF-TG调控强度矩阵R和所有TF-REs-TG调控三元组。
这样按照上述计算过程，可以得到每个细胞类型的调控网络。即得到每个细胞类型调控网络中的TF，RE和TG，用于后续细胞类型层次的遗传力富集模型构建。汇总所有的细胞类型的调控网络，形成调控网络图谱资源，用于对表型相关组织的系统识别。
(2.3)表型和细胞类型关联分数的计算
假设(2.1)中定义了M个细胞类型，并利用(2.2)的方法得到每个细胞类型的调控网络。对于每个细胞类型i，其调控网络中的RE集合将会组成该细胞类型的调控功能区域Ci。在遗传学的线性加和模型下，分层LD分数拟合方法即S-LDSC“Finucane,H.K.etal.Partitioning heritability by functional annotation using genome-wideassociation summary statistics.Nat Genet 47,1228-1235(2015).”假设表型的因果SNP的效应具有零均值和固定方差的分布，并进行如下建模：
其中，βj表示第j个SNP的效应分数，Var(βj)表示βj的方差。
同时S-LDSC假设对遗传力更富集的功能区域的连锁不平衡分数(LD)的增加比对遗传力富集没有贡献的功能区域的LD的增加对SNP的遗传效应影响更大，S-LDSC对SNP的遗传效应值χ2统计量的期望建模如下：
这里N是GWAS的样本量，Ci表示由第i个细胞类型的调控网络形成的调控功能区域，是来自GWAS的SNP与表型的关联效应统计量，/>是第i个调控功能区域中第j个SNP的LD分数，rjk表示第j个SNP和第k个SNP之间的LD分数，a是偏差，τi表示表型在Ci中的遗传力富集分数，由S-LDSC计算得到。Block jackknife被用于估计τi的标准误差，并使用标准误差来计算遗传力富集分数τi的P值pi。
为了在遗传力富集分数和假设检验产生的P值之间进行平衡，本发明整合遗传力富集分数和统计显著性(P值)来定义该表型与第i个细胞类型的关联分数如下：
Ri＝τi·(-log pi)
(2.4)关联分数显著性的检验
下面使用Block jackknife的方法来估计(2.3)中定义的关联分数Ri的标准误差、P值和假阳性率(FDR)。具体来说，(2.3)计算了表型和第i个细胞类型的关联分数Ri。然后按照以下步骤估计标准误差，P值和FDR：
(a)对于第i个细胞类型，将其调控网络中的RE分成100份。
(b)通过移除100份中的一份来生成一个子样本，可以生成了第i个细胞类型RE的100个子样本。
(c)这个细胞类型的100个RE子样本将形成100个新的调控功能区域，用于拟合S-LDSC。对于每个子样本，通过(2.3)获得了遗传力富集分数，P值和关联分数。
(d)利用100个子样本的100个关联分数，可以估计表型和第i个细胞类型的关联分数Ri的标准误差SDi。
(e)计算第i个细胞类型的关联分数的标准化分数：并估计Ri的P值和FDR值。
表型和第i个细胞类型的关联分数Ri及其FDR值可用于选择表型相关细胞类型。使用阈值(R分数≥100和FDR≤0.01)来选择表型的相关细胞类型。
(3)模型测试
为了测试和验证，利用一个较小的人类骨髓单细胞数据集(BMMC)。从中标注了Erythrocyte、Monocyte、CD4细胞、CD8细胞、NK细胞和B细胞。进一步收集了人类血细胞指标的GWAS统计汇总数据，分别为红细胞数目、红细胞分布、单核细胞数目、单核细胞比例、淋巴母细胞数目和淋巴母细胞比例。根据生物学知识，这两个数据可以提供用于验证的黄金标准：红细胞数目和分布与Erythrocyte有关，单核细胞数目和比例与Monocyte有关，淋巴母细胞数目和比例与CD细胞相关。将基于此黄金标准和scLinker方法进行比较。
本发明的关键点是提出了一种全新的基于单细胞数据构建细胞类型调控网络图谱识别表型相关细胞类型的框架，主要包括：
(1)基于单细胞多组学数据集建立细胞类型层面匹配伪群体(Pseudo-bulk)数据的方法。
(2)构建覆盖多组织器官的人类细胞类型特异的调控网络图谱的方法。
(3)基于调控网络计算遗传力富集识别表型相关细胞类型的方法。
(4)对表型与细胞类型的关联分数进行显著性检验的统计方法。
本发明已经经过实验验证。在模型测试环节，如表2所示，发现本发明可以准确识别表型相关的细胞类型。对于红细胞数目和红细胞分布两个表型，本发明和scLinker均可识别出红细胞与之关联性最强。对于单核母细胞数目和单核母细胞比例两个表型，本发明识别出单核母细胞细胞类型与之关联分数最高，但是scLinker识别出B细胞与之关联分数最高，这说明本发明比scLinker更加准确。对于淋巴母细胞数目和淋巴母细胞比例两个表型，本发明发现CD4+T细胞和CD8+T细胞与之关联分数更强，但是scLinker只能找到CD4+T细胞、NK细胞或者B细胞与之关联分数更强，这些结果也说明本发明的精确性更高。此外，如图2所示，本发明能够比scLinker得到更高的富集分数。例如在本实验中，scLinker得到最高的富集分数仅为19.80，但是本发明在各个细胞类型中的富集分数最低为59.77，这说明通过细胞类型特异的调控网络的限制，本发明可以将调控功能区域定义地更加准确，从而得到更高的遗传力的解释能力。
表2.本发明和scLinker在测试集上识别的表型相关细胞类型的比较
表型 本发明 scLinker
红细胞数目 Erythrocyte Erythrocyte
红细胞分布 Erythrocyte Erythrocyte
单核母细胞数目 Monocyte B
单核母细胞比例 Monocyte B
淋巴母细胞数目 CD4、CD8 CD4、NK
淋巴母细胞比例 CD4、CD8 CD4、B
本发明的另一实施例提供一种计算机设备(计算机、服务器、智能手机等)，其包括存储器和处理器，所述存储器存储计算机程序，所述计算机程序被配置为由所述处理器执行，所述计算机程序包括用于执行本发明方法中各步骤的指令。
本发明的另一实施例提供一种计算机可读存储介质(如ROM/RAM、磁盘、光盘)，所述计算机可读存储介质存储计算机程序，所述计算机程序被计算机执行时，实现本发明方法的各个步骤。
以上公开的本发明的具体实施例，其目的在于帮助理解本发明的内容并据以实施，本领域的普通技术人员可以理解，在不脱离本发明的精神和范围内，各种替换、变化和修改都是可能的。本发明不应局限于本说明书的实施例所公开的内容，本发明的保护范围以权利要求书界定的范围为准。

Claims (10)

1.一种基于单细胞多组学数据识别表型相关细胞类型的方法，其特征在于，包括以下步骤：
获取人类的单细胞多组学数据；
利用获取的单细胞多组学数据建立细胞类型层面的匹配数据；
对于每个细胞类型，以建立的匹配数据为输入构建其调控网络；
利用每个细胞类型的调控网络构建调控功能区域，在所有细胞类型的调控功能区域上对遗传力进行建模，并定义表型与细胞类型的关联分数；
度量表型与细胞类型的关联分数的显著性，利用表型与细胞类型的关联分数的显著性识别表型相关的细胞类型。
2.根据权利要求1所述的方法，其特征在于，所述获取人类的单细胞多组学数据，包括：获取人类的多个组织或者器官的匹配的单细胞RNA-seq和单细胞ATAC-seq数据。
3.根据权利要求1所述的方法，其特征在于，所述建立细胞类型层面的匹配数据，包括：
对于匹配的单细胞RNA-seq和单细胞ATAC-seq数据，标注其细胞类型；
通过Pseudo-bulk策略构建每个细胞类型的匹配的基因表达和染色质开放数据。
4.根据权利要求3所述的方法，其特征在于，所述通过Pseudo-bulk策略构建每个细胞类型的匹配的基因表达和染色质开放数据，包括：
对于一个细胞类型，给定来自其单细胞RNA-seq数据的包含c1个细胞和N个基因的单细胞基因表达矩阵E，以及来自单细胞ATAC-seq数据的包含c2个细胞和P个peak的单细胞染色质开放矩阵O，分别合并c1个和c2个细胞以得到第i个基因的表达PEi和第k个peak的开放性POk：
其中，e是单细胞RNA-seq数据的的RNA片段总数目，o是单细胞ATAC-seq数据的DNA片段总数目，Eic是细胞c中第i个基因的表达，Okc是细胞c中第k个peak的开放性；通过POk≥2的阈值选择开放的peak，从而构建该细胞类型的匹配的基因表达和染色质开放数据。
5.根据权利要求1所述的方法，其特征在于，所述调控网络以TF-RE-TG三元组为基本单元，其中TF表示转录因子，TG表示目的基因，RE表示调控元件；所述调控网络采用以下步骤构建：
将匹配的基因表达和染色质开放数据输入到PECA2模型中；
PECA2模型计算反式调控分数来衡量第i个TF对第j个TG的调控强度，第i个TF和第j个TG的反式调控分数Rij定义为：
其中，PEi和PEj是第i个TF和第j个TG的pseudo-bulk基因表达；Bik是第i个TF结合在第k个RE的强度；POk是第k个peak的pseudo-bulk开放程度；Ikj表示第k个RE和第j个TG之间的相互作用强度；Dij是第i个TF和第j个TG在不同DNA百科全书计划样本中的表达相关性；PECA2模型的输出是TF-TG调控强度矩阵R和所有TF-REs-TG调控三元组。
6.根据权利要求5所述的方法，其特征在于，所述利用每个细胞类型的调控网络构建调控功能区域，是对于每个细胞类型i，其调控网络中的RE集合组成该细胞类型的调控功能区域Ci；所述表型与细胞类型的关联分数通过整合遗传力富集分数τi及其统计显著性P值pi来定义：
Ri＝τi·(-log pi)
其中，Ri表示表型与第i个细胞类型的关联分数；τi表示表型在Ci中的遗传力富集分数；pi表示遗传力富集分数τi的显著性P值。
7.根据权利要求6所述的方法，其特征在于，采用以下步骤估计关联分数Ri的标准误差、P值和FDR：
对于第i个细胞类型，将其调控网络中的RE分成100份；
通过移除100份中的一份来生成一个子样本，生成第i个细胞类型RE的100个子样本；
第i个细胞类型的100个RE子样本形成100个新的调控功能区域，用于拟合S-LDSC；对于每个子样本，获得遗传力富集分数、P值和关联分数；
利用100个子样本的100个关联分数，估计表型和第i个细胞类型的关联分数Ri的标准误差SDi；
计算第i个细胞类型的关联分数的标准化分数用以估计Ri的P值和FDR值；
其中，关联分数Ri及其FDR值用于选择表型相关细胞类型。
8.一种基于单细胞多组学数据识别表型相关细胞类型的系统，其特征在于，包括：
数据获取模块，用于获取人类的单细胞多组学数据；
数学建模模块，用于利用获取的单细胞多组学数据建立细胞类型层面的匹配数据；对于每个细胞类型，以建立的匹配数据为输入构建其调控网络；利用每个细胞类型的调控网络构建调控功能区域，在所有细胞类型的调控功能区域上对遗传力进行建模，并定义表型与细胞类型的关联分数；度量表型与细胞类型的关联分数的显著性，利用表型与细胞类型的关联分数的显著性识别表型相关的细胞类型；
模型测试模块，用于获取给定表型或者疾病的全基因组关联分析GWAS统计汇总数据，对数学建模模块得到的模型进行测试。
9.一种计算机设备，其特征在于，包括存储器和处理器，所述存储器存储计算机程序，所述计算机程序被配置为由所述处理器执行，所述计算机程序包括用于执行权利要求1～7中任一项所述方法的指令。
10.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储计算机程序，所述计算机程序被计算机执行时，实现权利要求1～7中任一项所述的方法。
CN202311080320.8A 2023-08-25 2023-08-25 一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统 Pending CN117457065A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202311080320.8A CN117457065A (zh) 2023-08-25 2023-08-25 一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202311080320.8A CN117457065A (zh) 2023-08-25 2023-08-25 一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统

Publications (1)

Publication Number Publication Date
CN117457065A true CN117457065A (zh) 2024-01-26

Family

ID=89588074

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202311080320.8A Pending CN117457065A (zh) 2023-08-25 2023-08-25 一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统

Country Status (1)

Country Link
CN (1) CN117457065A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117854600A (zh) * 2024-03-07 2024-04-09 北京大学 基于多组学数据的细胞识别方法、装置、设备及存储介质

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117854600A (zh) * 2024-03-07 2024-04-09 北京大学 基于多组学数据的细胞识别方法、装置、设备及存储介质
CN117854600B (zh) * 2024-03-07 2024-05-21 北京大学 基于多组学数据的细胞识别方法、装置、设备及存储介质

Similar Documents

Publication Publication Date Title
CN103797129B (zh) 使用多态计数来解析基因组分数
JP7299169B2 (ja) 体細胞突然変異のクローン性を決定するための方法及びシステム
Li et al. Identification of selection signals by large-scale whole-genome resequencing of cashmere goats
CN108198621B (zh) 一种基于神经网络的数据库数据综合诊疗决策方法
US20140040264A1 (en) Method for estimation of information flow in biological networks
JP6312253B2 (ja) 形質予測モデル作成方法および形質予測方法
US20210151125A1 (en) Methods and systems for decomposition and quantification of dna mixtures from multiple contributors of known or unknown genotypes
CN109192316B (zh) 一种基于基因网络分析的疾病亚型预测系统
CN108335756B (zh) 鼻咽癌数据库及基于所述数据库的综合诊疗决策方法
CN108830045A (zh) 一种基于多组学的生物标记物系统筛选方法
CN108320807A (zh) 一种鼻咽癌人工智能辅助诊疗决策云系统
CN109599157A (zh) 一种精准智能诊疗大数据系统
CN110770839A (zh) 来自未知基因型贡献者的dna混合物的精确计算分解的方法
CN117457065A (zh) 一种基于单细胞多组学数据识别表型相关细胞类型的方法和系统
CN108206056A (zh) 一种鼻咽癌人工智能辅助诊疗决策终端
Dou et al. Single-nucleotide variant calling in single-cell sequencing data with Monopogen
CN108320797B (zh) 一种鼻咽癌数据库及基于所述数据库的综合诊疗决策方法
CN117409860A (zh) 一种基于组织特异调控网络图谱的多基因遗传风险评分计算方法和系统
CN112823391A (zh) 基于检测限的质量控制度量
CN108335748A (zh) 一种鼻咽癌人工智能辅助诊疗决策服务器集群
US20200135300A1 (en) Applying low coverage whole genome sequencing for intelligent genomic routing
Ackerman IV et al. IFPA Meeting 2013 Workshop Report II: Use of ‘omics’ in understanding placental development, bioinformatics tools for gene expression analysis, planning and coordination of a placenta research network, placental imaging, evolutionary approaches to understanding pre-eclampsia
WO2024217010A1 (zh) 基因与疾病关联分析方法、装置、计算机设备及存储介质
US20240153641A1 (en) Methods for genomic identification of phenotype risk
CN118824377A (zh) 基因与疾病关联分析方法、装置、计算机设备及存储介质

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