CN114864003A - 基于混合实验组和对照组单细胞样本的差异分析方法及系统 - Google Patents
基于混合实验组和对照组单细胞样本的差异分析方法及系统 Download PDFInfo
- Publication number
- CN114864003A CN114864003A CN202210262567.0A CN202210262567A CN114864003A CN 114864003 A CN114864003 A CN 114864003A CN 202210262567 A CN202210262567 A CN 202210262567A CN 114864003 A CN114864003 A CN 114864003A
- Authority
- CN
- China
- Prior art keywords
- data
- cell
- cells
- analysis
- 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.)
- Pending
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 62
- 108090000623 proteins and genes Proteins 0.000 claims abstract description 68
- 238000012216 screening Methods 0.000 claims abstract description 15
- 238000012545 processing Methods 0.000 claims abstract description 12
- 238000003908 quality control method Methods 0.000 claims abstract description 12
- 230000035772 mutation Effects 0.000 claims abstract description 3
- 230000014509 gene expression Effects 0.000 claims description 54
- 238000000034 method Methods 0.000 claims description 38
- 238000007621 cluster analysis Methods 0.000 claims description 23
- 239000011159 matrix material Substances 0.000 claims description 9
- 238000004422 calculation algorithm Methods 0.000 claims description 8
- 238000000585 Mann–Whitney U test Methods 0.000 claims description 6
- 238000004364 calculation method Methods 0.000 claims description 6
- 238000005457 optimization Methods 0.000 claims description 6
- 238000000513 principal component analysis Methods 0.000 claims description 6
- 108020005196 Mitochondrial DNA Proteins 0.000 claims description 4
- 238000002474 experimental method Methods 0.000 claims description 4
- 238000010195 expression analysis Methods 0.000 claims description 2
- 210000004027 cell Anatomy 0.000 description 232
- 201000010099 disease Diseases 0.000 description 13
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 13
- 239000000523 sample Substances 0.000 description 10
- 238000005516 engineering process Methods 0.000 description 8
- 210000002540 macrophage Anatomy 0.000 description 7
- 201000009794 Idiopathic Pulmonary Fibrosis Diseases 0.000 description 6
- 208000036971 interstitial lung disease 2 Diseases 0.000 description 6
- 238000012163 sequencing technique Methods 0.000 description 6
- 238000007405 data analysis Methods 0.000 description 5
- 210000004072 lung Anatomy 0.000 description 5
- 230000001105 regulatory effect Effects 0.000 description 5
- 238000011160 research Methods 0.000 description 5
- 210000001744 T-lymphocyte Anatomy 0.000 description 4
- 210000003719 b-lymphocyte Anatomy 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 238000012854 evaluation process Methods 0.000 description 3
- 238000002360 preparation method Methods 0.000 description 3
- 208000005069 pulmonary fibrosis Diseases 0.000 description 3
- 239000013068 control sample Substances 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000002888 effect on disease Effects 0.000 description 2
- 230000008030 elimination Effects 0.000 description 2
- 238000003379 elimination reaction Methods 0.000 description 2
- 210000003958 hematopoietic stem cell Anatomy 0.000 description 2
- 238000002372 labelling Methods 0.000 description 2
- 238000007781 pre-processing Methods 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 238000011144 upstream manufacturing Methods 0.000 description 2
- 206010061818 Disease progression Diseases 0.000 description 1
- 235000012571 Ficus glomerata Nutrition 0.000 description 1
- 244000153665 Ficus glomerata Species 0.000 description 1
- 206010058467 Lung neoplasm malignant Diseases 0.000 description 1
- 206010028980 Neoplasm Diseases 0.000 description 1
- 238000003559 RNA-seq method Methods 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 201000011510 cancer Diseases 0.000 description 1
- 230000015556 catabolic process Effects 0.000 description 1
- 230000005750 disease progression Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 230000003902 lesion Effects 0.000 description 1
- 201000005202 lung cancer Diseases 0.000 description 1
- 208000020816 lung neoplasm Diseases 0.000 description 1
- 238000002156 mixing Methods 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 108020004707 nucleic acids Proteins 0.000 description 1
- 150000007523 nucleic acids Chemical class 0.000 description 1
- 102000039446 nucleic acids Human genes 0.000 description 1
- 230000008506 pathogenesis Effects 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
Images
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
- G16B20/40—Population genetics; Linkage disequilibrium
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
- G06F18/2135—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on approximation criteria, e.g. principal component analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H50/00—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics
- G16H50/20—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics for computer-aided diagnosis, e.g. based on medical expert systems
Landscapes
- Engineering & Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Health & Medical Sciences (AREA)
- Theoretical Computer Science (AREA)
- Life Sciences & Earth Sciences (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Physics & Mathematics (AREA)
- Medical Informatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- General Health & Medical Sciences (AREA)
- General Engineering & Computer Science (AREA)
- Biomedical Technology (AREA)
- Genetics & Genomics (AREA)
- Artificial Intelligence (AREA)
- General Physics & Mathematics (AREA)
- Public Health (AREA)
- Ecology (AREA)
- Analytical Chemistry (AREA)
- Primary Health Care (AREA)
- Pathology (AREA)
- Databases & Information Systems (AREA)
- Physiology (AREA)
- Chemical & Material Sciences (AREA)
- Epidemiology (AREA)
- Biophysics (AREA)
- Molecular Biology (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Biotechnology (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明涉及一种基于混合实验组和对照组单细胞样本的差异分析方法，包括：a.获取数据，并对所述数据进行质量控制及标准化处理，所述数据包括实验数据和开源数据；b.从标准化处理后的数据中筛选出高变基因的数据，并对筛选出的数据进行降维；c.根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果；d.对细胞亚群进行差异分析。本发明还涉及一种基于混合实验组和对照组单细胞样本的差异分析系统。本发明能够更好地发现由特定组构成的细胞亚群，从而实现在单细胞层面的实验组和对照组的差异分析。
Description
技术领域
本发明涉及一种基于混合实验组和对照组单细胞样本的差异分析方法及系统。
背景技术
高通量的单细胞技术目前是生物信息领域的热门研究方向之一，它带领了基因组学领域的转变，即从整体组织分析转向对单个细胞的详细和全面的研究，提供了对不同生物体组织中细胞多样性的洞察力。当前许多单细胞研究的目的主要集中在鉴定细胞类型(如B细胞，T细胞，以及表达独特基因的细胞亚群)，很少聚焦于在单细胞层面分析实验组(即疾病组)和对照组的单细胞表达差异，且直接将全部单细胞按照实验组和对照组进行整体差异分析捕捉到的差异信号较弱，效果不佳。
在癌症等疾病中，由于组织中某一群细胞的病变就可能导致整个有机体的崩溃，而基于组织样本的分析仅仅获取了组织中所有细胞的平均表达信息，很难分辨细胞间的差异，也很难识别在疾病进展中起到重要作用的稀有细胞群。单细胞技术的出现有效地解决了这一问题，单细胞技术是在组织中分离单细胞并对其进行测序以获得单个细胞的基因表达的技术，它的发展使基因表达的研究达到了更高的分辨率。单细胞技术分为上游和下游两个层面，上游层面指的是细胞样本的制备和测序技术，包括单细胞分离、核酸提取和扩增、测序库制备、RNA测序。而下游层面则是针对单细胞表达矩阵的数据分析技术，包括数据质量控制、数据降维、聚类、差异分析等一系列数据分析。
目前单细胞研究的主要方向是鉴定新的细胞亚群类型和对疾病发展有重要作用的特殊细胞亚群。单细胞研究的取样策略通常有两种：一是通过实验获取疾病组织样本(如肺癌，肺纤维化的肺部组织样本)的单细胞，只针对疾病单细胞样本进行分析；二是选取疾病组织和健康组织的单细胞样本共同分析。两种方式都是通过从组织中分离单细胞进行测序，并鉴定细胞亚群。第一种取样方式更关注于鉴定疾病进程中的关键细胞亚群，即关键基因表达上调的细胞亚群；而第二种方法更关注细胞类型的鉴定，即分析实验组和对照组单细胞样本鉴定组织中不同的细胞类型(如B细胞、T细胞、巨噬细胞及其亚型等)。
上述提到的单细胞分析的主要任务是鉴定关键细胞亚群，所述关键细胞亚群可能表达某些独特基因并对疾病的发展产生重要的影响。但有一个问题一直是单细胞研究中很少关注到的，即是否能够从疾病样本和对照样本单细胞中找到有着共同表达模式(即共同表达某些基因)，且单独由疾病样本或对照样本细胞构成的细胞亚群。
寻找上述细胞亚群的意义在于可以更精细化地进行差异分析。在之前的分析中，发现许多实验组和对照组的细胞亚群共享某种基因表达模式，即在多个细胞亚群中实验组和对照组的单细胞数目相近，这也意味着这些细胞不能体现实验组和对照组的表达差异。而细胞占比差距悬殊，尤其是一个基本由疾病样本细胞构成的细胞亚群较为罕见，而往往这类罕见的细胞群体高度表达了疾病相关的重要基因。如果能够在混合两组样本的细胞中找到这样的细胞亚群，就可以为解释疾病的发病机制提供独特的见解。
在组织层面的研究中，研究人员会直接使用实验组和对照组的组织样本的基因表达数据进行差异分析，但这种分析方法在单细胞层面是行不通的。单细胞测序数据的特点是高维稀疏，仅一个组织样本就能分离出几千个单细胞，并且一个单细胞只表达一部分基因。在这样的数据背景下直接将实验组和对照组的细胞样本整体进行对比会使得到的差异信号很弱。
发明内容
有鉴于此，有必要提供一种基于混合实验组和对照组单细胞样本的差异分析方法及系统。
本发明提供一种基于混合实验组和对照组单细胞样本的差异分析方法，该方法包括如下步骤：a.获取数据，并对所述数据进行质量控制及标准化处理，所述数据包括实验数据和开源数据；b.从标准化处理后的数据中筛选出高变基因的数据，并对筛选出的数据进行降维；c.根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果；d.对细胞亚群进行差异分析。
优选地，所述的步骤a包括：
对原始表达矩阵进行质量控制，包括剔除低质量的细胞和基因；其中，低质量的细胞是指：表达基因数量较少，线粒体基因表达比例较高的细胞；低质量的基因是指：在细胞中表达比例较低的基因。
优选地，所述的步骤a还包括：
通过对唯一分子识别标签计数归一化缩小数值差异，采用LogNormalization方法，公式如下：
其中，i表示细胞，j表示基因；
在对数据归一化之后，继续对数据进行缩放，公式如下：
优选地，所述的步骤b包括：
通过对基因在细胞中的表达均值和方差进行建模，筛选高变基因，对单细胞在高变基因中的表达数据进行主成分分析，并可视化每个主成分的方差百分比，选择排名靠前的主成分用于后续聚类分析。
优选地，所述的步骤c具体包括：
采用基于图的聚类方法：首先计算每个细胞之间的相互距离以确定每个细胞的k最近邻，然后计算每个细胞与其k最近邻之间的邻域重叠构造共享最近邻图；公式如下，其中A，B代表两个细胞的“邻居”细胞数：
优选地，所述的步骤c还包括：
使用Louvain算法作为模块化优化函数确定细胞亚群的数量：通过扫描节点及其邻居节点，对每个节点计算模块度并衡量其加入模块后的模块收益，选择最大收益的邻居节点加入模块，重复迭代，最终形成节点簇；对于加权图，模块性计算公式如下，其中，m是图中边总数，i，j代表两个节点(即细胞)，Aij代表两个节点的权重，ki，kj代表节点i，j所有边缘权重之和，ci，cj代表节点i，j的模块，δ是Kronecker delta函数：当x＝y时，δ(x，y)＝1；否则为0：
优选地，所述的步骤d具体包括：
对每群细胞亚群计算差异基因：比较本群细胞和除本群细胞之外所有细胞的基因表达差异，使用Wilcoxon秩和检验方法进行差异表达分析。
本发明提供一种基于混合实验组和对照组单细胞样本的差异分析系统，该系统包括获取模块、筛选模块、聚类分析模块、差异分析模块，其中：所述获取模块用于获取数据，并对所述数据进行质量控制及标准化处理，所述数据包括实验数据和开源数据；所述筛选模块用于从标准化处理后的数据中筛选出高变基因的数据，并对筛选出的数据进行降维；所述聚类分析模块用于根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果；所述差异分析模块用于对细胞亚群进行差异分析。
本发明改变了现有的单细胞数据分析任务多数聚焦于鉴定表达关键基因的细胞亚群或新的细胞类型的传统做法，在鉴定细胞亚群的基础上，更加关注实验组和对照组单细胞样本的差异研究，即通过本申请能够更好地发现由特定组构成的细胞亚群，从而实现在单细胞层面的实验组和对照组的差异分析。
附图说明
图1为本发明基于混合实验组和对照组单细胞样本的差异分析方法的流程图；
图2为本发明实施例基于混合实验组和对照组单细胞样本的差异的流程示意图；
图3为本发明实施例提供的代表细胞亚群走向的聚类树(以GSE36831数据集中的巨噬细胞为例)；
图4为本发明实施例提供的细胞亚群构成图(以GSE36831数据集中的巨噬细胞为例)；
图5为本发明实施例提供的最佳聚类结果评估流程图；
图6为本发明基于混合实验组和对照组单细胞样本的差异分析系统的硬件架构图。
具体实施方式
下面结合附图及具体实施例对本发明作进一步详细的说明。
参阅图1、图2所示，是本发明基于混合实验组和对照组单细胞样本的差异分析方法较佳实施例的作业流程图。
步骤S1，获取数据，并对所述数据进行质量控制及标准化处理。其中，所述数据包括实验数据和开源数据。具体而言：
本实施例分析了两个数据集：第一个数据集是小鼠肺纤维化造血细胞(Cd45+)数据集，该数据集由实验获得，包含了两个不同时间点的实验组和对照组，共四组单细胞表达数据，每组检测了约5000+单细胞在30000+基因中的表达量；第二个数据集是从GEO数据库下载的人特发性肺纤维化肺组织单细胞数据集(GSE36831)，包括28个健康个体(对照组)和31个特发性肺纤维化(IPF)患者(实验组)的肺组织单细胞样本，共检测了380000+单细胞在40000+基因中的表达量。
首先，对原始表达矩阵进行质量控制，包括剔除低质量的细胞(表达基因数量较少，线粒体基因表达比例较高)和基因(基因在细胞中表达比例较低)。由于每个细胞在不同基因的表达读数差异非常大，影响后续数据分析，因此通过对UMIs计数归一化缩小数值差异。所述表达读数即是唯一分子识别标签(unique molecular identifier，UMI)。这里本实施例采用LogNormalization方法，公式如下，其中i表示细胞，j表示基因：
在对数据归一化之后，继续对数据进行缩放(即线性变换)，这是在执行降维和聚类之前的标准预处理步骤，公式如下：
步骤S2，从标准化处理后的数据中，筛选高度变化特征(Highly VariableFeatures)，即在不同细胞中表达差异较大的基因，并对数据降维。具体而言：
由于单细胞表达矩阵是稀疏矩阵，多数基因在细胞中的表达差异较小，且将所有基因用于分析会占用大量计算资源。因此，选择高度变化的特征(即高变基因)进行后续的降维和聚类分析会节省很多计算资源，并且有助于突出单细胞数据集中的生物信号。通过对基因在细胞中的表达均值和方差进行建模，筛选高变基因，对单细胞在高变基因中的表达数据进行主成分分析(PCA)，并可视化每个主成分(PC)的方差百分比(SD)，选择排名靠前的主成分用于后续聚类分析。
步骤S3，根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果：
聚类分析用于两个任务：第一个是细胞分型，即从组织分离的细胞样本中标记每个细胞的细胞类型(如巨噬细胞，B细胞，T细胞等)；第二个是对实验组和对照组的差异分析。这一步在已经确定每个细胞的细胞类型之后进行，通过提取同一细胞类型的实验组和对照组的单细胞样本，调整参数进行多次聚类，评估每次聚类结果，选择最佳聚类结果。
具体而言：
在本步骤中，本实施例采用基于图的聚类方法。该方法首先计算每个细胞之间的相互距离以确定每个细胞的k最近邻(k-nearest neighbors，k默认为20)，然后计算每个细胞与其k最近邻之间的邻域重叠(Jaccard指数)构造共享最近邻图(Shared NearestNeighbor，SNN)。Jaccard指数计算的是两个细胞共有的“邻居”和所有“邻居”的比值，比值越大，两个细胞越相似，公式如下，其中A，B代表两个细胞的“邻居”细胞数：
然后，使用Louvain算法作为模块化优化函数确定簇(细胞亚群)的数量。Louvain算法是一种从网络中提取模块的方法，是一种贪婪优化方法。Louvain算法通过扫描节点及其邻居节点，对每个节点计算模块度并衡量其加入模块后的模块收益，选择最大收益的邻居节点加入模块，重复迭代，最终形成节点簇。对于加权图，模块性计算公式如下，其中，m是图中边总数，i，j代表两个节点(即细胞)，Aij代表两个节点的权重，ki，kj代表节点i，j所有边缘权重之和，ci，cj代表节点i，j的模块，δ是Kronecker delta函数(当x＝y时，δ(x，y)＝1；否则为0)：
本实施例利用聚类分析完成两个任务：
第一个任务是细胞分型，即对组织样本分离的单细胞鉴定细胞类型。这一步通过对所有细胞进行聚类分析，划分细胞亚群，计算每个亚群表达上调的基因。通过在线资源(文献、数据库)人工分析所述上调的基因是否包含某一细胞类型的markers，从而判断细胞亚群所属细胞类型。
第二个任务是本实施例的关键部分，即对混合实验组和对照组特定细胞类型的单细胞样本聚类。首先，选取上一步分析的特定细胞类型的两组细胞(如巨噬细胞)，调整参数进行多次聚类，构建聚类树，可视化每次调整参数后细胞分群走向(如图3)，并对聚类后每个细胞亚群计算实验组和对照组构成比例(如图4)，选择最佳聚类结果(最佳聚类结果评估流程如图5所示)。
所述最佳聚类结果需综合考虑两组样本构成比例差距悬殊的细胞亚群数目以及所述细胞亚群的细胞数量。如果两者都是最多，则该参数下的聚类结果为最佳聚类结果；否则，如果比例差距悬殊的细胞亚群数目不是最多，但所述细胞亚群的细胞数足够多，则优先考虑这样的聚类结果。细胞数量太少的亚群不能有力的证明其对疾病的影响。
步骤S4，在所述聚类分析完成后，对细胞亚群进行差异分析。具体而言：
对每群细胞亚群计算差异基因。本实施例的计算策略是比较本群细胞和除本群细胞之外所有细胞的基因表达差异，使用了Wilcoxon秩和检验(Wilcoxon Rank Sum test)方法进行差异表达分析。本实施例限制两个细胞群中基因的平均|log fold change|>0.25，并且限制被比较的两个细胞群中任意一群细胞的各基因表达占比>0.1，且差异基因的FDR<0.05。
参阅图6所示，是本发明基于混合实验组和对照组单细胞样本的差异分析系统10的硬件架构图。该系统包括：获取模块101、筛选模块102、聚类分析模块103、差异分析模块104。其中：
所述获取模块101用于获取数据，并对所述数据进行质量控制及标准化处理。其中，所述数据包括实验数据和开源数据。具体而言：
本实施例分析了两个数据集：第一个数据集是小鼠肺纤维化造血细胞(Cd45+)数据集，该数据集由实验获得，包含了两个不同时间点的实验组和对照组，共四组单细胞表达数据，每组检测了约5000+单细胞在30000+基因中的表达量；第二个数据集是从GEO数据库下载的人特发性肺纤维化肺组织单细胞数据集(GSE36831)，包括28个健康个体(对照组)和31个特发性肺纤维化(IPF)患者(实验组)的肺组织单细胞样本，共检测了380000+单细胞在40000+基因中的表达量。
首先，所述获取模块101对原始表达矩阵进行质量控制，包括剔除低质量的细胞(表达基因数量较少，线粒体基因表达比例较高)和基因(基因在细胞中表达比例较低)。由于每个细胞在不同基因的表达读数差异非常大，影响后续数据分析，因此通过对UMIs计数归一化缩小数值差异。所述表达读数即是唯一分子识别标签(unique molecularidentifier，UMI)。这里本实施例采用LogNormalization方法，公式如下，其中i表示细胞，j表示基因：
在对数据归一化之后，所述获取模块101继续对数据进行缩放(即线性变换)，这是在执行降维和聚类之前的标准预处理步骤，公式如下：
所述筛选模块102用于从标准化处理后的数据中，筛选高度变化特征(HighlyVariable Features)，即在不同细胞中表达差异较大的基因，并对数据降维。具体而言：
由于单细胞表达矩阵是稀疏矩阵，多数基因在细胞中的表达差异较小，且将所有基因用于分析会占用大量计算资源。因此，选择高度变化的特征(即高变基因)进行后续的降维和聚类分析会节省很多计算资源，并且有助于突出单细胞数据集中的生物信号。所述筛选模块102通过对基因在细胞中的表达均值和方差进行建模，筛选高变基因，对单细胞在高变基因中的表达数据进行主成分分析(PCA)，并可视化每个主成分(PC)的方差百分比(SD)，选择排名靠前的主成分用于后续聚类分析。
所述聚类分析模块103用于根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果：
聚类分析用于两个任务：第一个是细胞分型，即从组织分离的细胞样本中标记每个细胞的细胞类型(如巨噬细胞，B细胞，T细胞等)；第二个是对实验组和对照组的差异分析。这一步在已经确定每个细胞的细胞类型之后进行，通过提取同一细胞类型的实验组和对照组的单细胞样本，调整参数进行多次聚类，评估每次聚类结果，选择最佳聚类结果。
具体而言：
所述聚类分析模块103采用基于图的聚类方法。首先计算每个细胞之间的相互距离以确定每个细胞的k最近邻(k-nearest neighbors，k默认为20)，然后计算每个细胞与其k最近邻之间的邻域重叠(Jaccard指数)构造共享最近邻图(Shared Nearest Neighbor，SNN)。Jaccard指数计算的是两个细胞共有的“邻居”和所有“邻居”的比值，比值越大，两个细胞越相似，公式如下，其中A，B代表两个细胞的“邻居”细胞数：
然后，使用Louvain算法作为模块化优化函数确定簇(细胞亚群)的数量。Louvain算法是一种从网络中提取模块的方法，是一种贪婪优化方法。Louvain算法通过扫描节点及其邻居节点，对每个节点计算模块度并衡量其加入模块后的模块收益，选择最大收益的邻居节点加入模块，重复迭代，最终形成节点簇。对于加权图，模块性计算公式如下，其中，m是图中边总数，i，j代表两个节点(即细胞)，Aij代表两个节点的权重，ki，kj代表节点i，j所有边缘权重之和，ci，cj代表节点i，j的模块，δ是Kronecker delta函数(当x＝y时，δ(x，y)＝1；否则为0)：
本实施例利用聚类分析完成两个任务：
第一个任务是细胞分型，即对组织样本分离的单细胞鉴定细胞类型。这一步通过对所有细胞进行聚类分析，划分细胞亚群，计算每个亚群表达上调的基因。通过在线资源(文献、数据库)人工分析所述上调的基因是否包含某一细胞类型的markers，从而判断细胞亚群所属细胞类型。
第二个任务是本实施例的关键部分，即对混合实验组和对照组特定细胞类型的单细胞样本聚类。首先，选取上一步分析的特定细胞类型的两组细胞(如巨噬细胞)，调整参数进行多次聚类，构建聚类树，可视化每次调整参数后细胞分群走向(如图3)，并对聚类后每个细胞亚群计算实验组和对照组构成比例(如图4)，选择最佳聚类结果(最佳聚类结果评估流程如图5所示)。
所述最佳聚类结果需综合考虑两组样本构成比例差距悬殊的细胞亚群数目以及所述细胞亚群的细胞数量。如果两者都是最多，则该参数下的聚类结果为最佳聚类结果；否则，如果比例差距悬殊的细胞亚群数目不是最多，但所述细胞亚群的细胞数足够多，则优先考虑这样的聚类结果。细胞数量太少的亚群不能有力的证明其对疾病的影响。
所述差异分析模块104用于在所述聚类分析完成后，对细胞亚群进行差异分析。具体而言：
所述差异分析模块104对每群细胞亚群计算差异基因。本实施例的计算策略是比较本群细胞和除本群细胞之外所有细胞的基因表达差异，使用了Wilcoxon秩和检验(Wilcoxon Rank Sum test)方法进行差异表达分析。本实施例限制两个细胞群中基因的平均|log fold change|>0.25，并且限制被比较的两个细胞群中任意一群细胞的各基因表达占比>0.1，且差异基因的FDR<0.05。
本申请通过混合特定细胞类型的实验组和对照组单细胞样本进行多轮聚类分析，评估每轮聚类结果，选择最佳聚类结果并进行差异研究。与基于组织测序数据的差异分析方法不同，本申请首先通过聚类分析找到实验组和对照组细胞占比悬殊的细胞亚群，然后分析上述细胞亚群的表达模式，比直接分析实验组和对照组的单细胞具备更强的差异基因信号。
虽然本发明参照当前的较佳实施方式进行了描述，但本领域的技术人员应能理解，上述较佳实施方式仅用来说明本发明，并非用来限定本发明的保护范围，任何在本发明的精神和原则范围之内，所做的任何修饰、等效替换、改进等，均应包含在本发明的权利保护范围之内。
Claims (8)
1.一种基于混合实验组和对照组单细胞样本的差异分析方法，其特征在于，该方法包括如下步骤：
a.获取数据，并对所述数据进行质量控制及标准化处理，所述数据包括实验数据和开源数据；
b.从标准化处理后的数据中筛选出高变基因的数据，并对筛选出的数据进行降维；
c.根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果；
d.对细胞亚群进行差异分析。
2.如权利要求1所述的方法，其特征在于，所述的步骤a包括：
对原始表达矩阵进行质量控制，包括剔除低质量的细胞和基因；其中，低质量的细胞是指：表达基因数量较少，线粒体基因表达比例较高的细胞；低质量的基因是指：在细胞中表达比例较低的基因。
3.如权利要求2所述的方法，其特征在于，所述的步骤a还包括：
通过对唯一分子识别标签计数UMI归一化缩小数值差异，采用LogNormalization方法，公式如下：
其中，i表示细胞，j表示基因；
在对数据归一化之后，继续对数据进行缩放，公式如下：
4.如权利要求3所述的方法，其特征在于，所述的步骤b包括：
通过对基因在细胞中的表达均值和方差进行建模，筛选高变基因，对单细胞在高变基因中的表达数据进行主成分分析，并可视化每个主成分的方差百分比，选择排名靠前的主成分用于后续聚类分析。
5.如权利要求4所述的方法，其特征在于，所述的步骤c具体包括：
采用基于图的聚类方法：首先计算每个细胞之间的相互距离以确定每个细胞的k最近邻，然后计算每个细胞与其k最近邻之间的邻域重叠构造共享最近邻图；公式如下，其中A，B代表两个细胞的“邻居”细胞数：
6.如权利要求5所述的方法，其特征在于，所述的步骤c还包括：
使用Louvain算法作为模块化优化函数确定细胞亚群的数量：通过扫描节点及其邻居节点，对每个节点计算模块度并衡量其加入模块后的模块收益，选择最大收益的邻居节点加入模块，重复迭代，最终形成节点簇；对于加权图，模块性计算公式如下，其中，m是图中边总数，i,j代表两个节点(即细胞)，Aij代表两个节点的权重，ki,kj代表节点i,j所有边缘权重之和，ci,cj代表节点i,j的模块，δ是Kronecker delta函数：当x＝y时，δ(x,y)＝1；否则为0：
7.如权利要求6所述的方法，其特征在于，所述的步骤d具体包括：
对每群细胞亚群计算差异基因：比较本群细胞和除本群细胞之外所有细胞的基因表达差异，使用Wilcoxon秩和检验方法进行差异表达分析。
8.一种基于混合实验组和对照组单细胞样本的差异分析系统，其特征在于，该系统包括获取模块、筛选模块、聚类分析模块、差异分析模块，其中：
所述获取模块用于获取数据，并对所述数据进行质量控制及标准化处理，所述数据包括实验数据和开源数据；
所述筛选模块用于从标准化处理后的数据中筛选出高变基因的数据，并对筛选出的数据进行降维；
所述聚类分析模块用于根据降维后的数据，对单细胞进行多轮不同参数的聚类分析，计算每轮聚类分析的细胞亚群实验组和对照组的构成比例，以得到最优结果；
所述差异分析模块用于对细胞亚群进行差异分析。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210262567.0A CN114864003A (zh) | 2022-03-17 | 2022-03-17 | 基于混合实验组和对照组单细胞样本的差异分析方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210262567.0A CN114864003A (zh) | 2022-03-17 | 2022-03-17 | 基于混合实验组和对照组单细胞样本的差异分析方法及系统 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN114864003A true CN114864003A (zh) | 2022-08-05 |
Family
ID=82627886
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210262567.0A Pending CN114864003A (zh) | 2022-03-17 | 2022-03-17 | 基于混合实验组和对照组单细胞样本的差异分析方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114864003A (zh) |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115527610A (zh) * | 2022-11-09 | 2022-12-27 | 上海交通大学 | 一种单细胞组学数据的聚类分析方法 |
| CN115881218A (zh) * | 2022-12-15 | 2023-03-31 | 哈尔滨星云医学检验所有限公司 | 用于全基因组关联分析的基因自动选择方法 |
| CN116825206A (zh) * | 2023-08-30 | 2023-09-29 | 四川大学华西医院 | 探究fh-缺陷型肾癌关键细胞亚群的方法、装置和设备 |
| CN117409965A (zh) * | 2023-09-28 | 2024-01-16 | 江苏先声医学诊断有限公司 | 适用于亚洲her2阳性乳腺癌患者风险预测系统 |
| CN117747093A (zh) * | 2024-02-20 | 2024-03-22 | 神州医疗科技股份有限公司 | 一种特发性肺纤维化诊断模型的构建方法及诊断系统 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109979538A (zh) * | 2019-03-28 | 2019-07-05 | 广州基迪奥生物科技有限公司 | 一种基于10x单细胞转录组测序数据的分析方法 |
| CN111951892A (zh) * | 2020-08-04 | 2020-11-17 | 荣联科技集团股份有限公司 | 一种基于单细胞测序数据分析细胞轨迹的方法及电子设备 |
| CN112700820A (zh) * | 2021-01-07 | 2021-04-23 | 广州华银健康医疗集团股份有限公司 | 一种基于单细胞转录组测序的细胞亚群注释方法 |
-
2022
- 2022-03-17 CN CN202210262567.0A patent/CN114864003A/zh active Pending
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109979538A (zh) * | 2019-03-28 | 2019-07-05 | 广州基迪奥生物科技有限公司 | 一种基于10x单细胞转录组测序数据的分析方法 |
| CN111951892A (zh) * | 2020-08-04 | 2020-11-17 | 荣联科技集团股份有限公司 | 一种基于单细胞测序数据分析细胞轨迹的方法及电子设备 |
| CN112700820A (zh) * | 2021-01-07 | 2021-04-23 | 广州华银健康医疗集团股份有限公司 | 一种基于单细胞转录组测序的细胞亚群注释方法 |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115527610A (zh) * | 2022-11-09 | 2022-12-27 | 上海交通大学 | 一种单细胞组学数据的聚类分析方法 |
| CN115527610B (zh) * | 2022-11-09 | 2023-11-24 | 上海交通大学 | 一种单细胞组学数据的聚类分析方法 |
| CN115881218A (zh) * | 2022-12-15 | 2023-03-31 | 哈尔滨星云医学检验所有限公司 | 用于全基因组关联分析的基因自动选择方法 |
| CN116825206A (zh) * | 2023-08-30 | 2023-09-29 | 四川大学华西医院 | 探究fh-缺陷型肾癌关键细胞亚群的方法、装置和设备 |
| CN116825206B (zh) * | 2023-08-30 | 2023-11-03 | 四川大学华西医院 | 探究fh-缺陷型肾癌关键细胞亚群的方法、装置和设备 |
| CN117409965A (zh) * | 2023-09-28 | 2024-01-16 | 江苏先声医学诊断有限公司 | 适用于亚洲her2阳性乳腺癌患者风险预测系统 |
| CN117747093A (zh) * | 2024-02-20 | 2024-03-22 | 神州医疗科技股份有限公司 | 一种特发性肺纤维化诊断模型的构建方法及诊断系统 |
| CN117747093B (zh) * | 2024-02-20 | 2024-06-07 | 神州医疗科技股份有限公司 | 一种特发性肺纤维化诊断模型的构建方法及诊断系统 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN114864003A (zh) | 基于混合实验组和对照组单细胞样本的差异分析方法及系统 | |
| Yao et al. | An integrated transcriptomic and epigenomic atlas of mouse primary motor cortex cell types | |
| Diggins et al. | Methods for discovery and characterization of cell subsets in high dimensional mass cytometry data | |
| Wirth et al. | Expression cartography of human tissues using self organizing maps | |
| Jiang et al. | Cluster analysis for gene expression data: a survey | |
| Melsen et al. | A comprehensive workflow for applying single-cell clustering and pseudotime analysis to flow cytometry data | |
| Fuda et al. | Artificial intelligence in clinical multiparameter flow cytometry and mass cytometry–key tools and progress | |
| Lan et al. | scIAC: clustering scATAC-seq data based on Student’s t-distribution similarity imputation and denoising autoencoder | |
| JP2023546645A (ja) | シングルセルゲノミクスデータセットからの細胞のサブサンプリング方法及びシステム | |
| Rossin et al. | A framework for analytical characterization of monoclonal antibodies based on reactivity profiles in different tissues | |
| Zhen et al. | A review and performance evaluation of clustering frameworks for single-cell Hi-C data | |
| CN117423391A (zh) | 一种基因调控网络数据库的建立方法、系统及设备 | |
| CN116842996A (zh) | 一种基于深度压缩感知的空间转录组方法及装置 | |
| CN115527610A (zh) | 一种单细胞组学数据的聚类分析方法 | |
| Tasoulis et al. | Unsupervised clustering of bioinformatics data | |
| CN111755074A (zh) | 一种酿酒酵母菌中dna复制起点的预测方法 | |
| Shaik et al. | Fuzzy-adaptive-subspace-iteration-based two-way clustering of microarray data | |
| Padma et al. | A modified algorithm for clustering based on particle swarm optimization and K-means | |
| Cortal | Development of bioinformatics methods for high-dimensional single-cell data analysis and their application to the study of cell heterogeneity | |
| Lopez | Charting Cellular States, One Cell at a Time: Computational, Inferential and Modeling Perspectives | |
| McLoughlin | Scaling the Shared Identified Differential Expression (SIDEseq) Measure for Massive Scale Single Cell RNA Sequencing Data and Exploring Extension of the Measure to Chromatin Accessibility Co-Assays | |
| Pinello et al. | SIMBA: SIngle-cell eMBedding Along with features | |
| Lee et al. | Finite mixture models in biostatistics | |
| Aouf et al. | Gene Expression Data For Gene Selection Using Ensemble Based Feature Selection | |
| Sun | Novel statistical methods in analyzing single cell sequencing 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 |