CN110379459A - 一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统 - Google Patents
一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统 Download PDFInfo
- Publication number
- CN110379459A CN110379459A CN201910742777.8A CN201910742777A CN110379459A CN 110379459 A CN110379459 A CN 110379459A CN 201910742777 A CN201910742777 A CN 201910742777A CN 110379459 A CN110379459 A CN 110379459A
- Authority
- CN
- China
- Prior art keywords
- gene
- expression
- difference
- function
- transcript profile
- 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.)
- Granted
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
- G16B20/20—Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
-
- 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/30—Detection of binding sites or motifs
-
- 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
- G16B25/00—ICT specially adapted for hybridisation; ICT specially adapted for gene or protein expression
- G16B25/10—Gene or protein expression profiling; Expression-ratio estimation or normalisation
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Engineering & Computer Science (AREA)
- Genetics & Genomics (AREA)
- Biotechnology (AREA)
- Biophysics (AREA)
- Molecular Biology (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Theoretical Computer Science (AREA)
- Chemical & Material Sciences (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Analytical Chemistry (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明公开了一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统,通过对样本的多时序采样和高通量转录组测序技术,获得目标病理生理过程、对照病理生理过程的动态时序转录组变化;通过基因功能关联分析,评估差异共表达集与目标病理生理过程基因集之间的关联强度,从而发现组学数据中稳健可靠、同时具有病理生理意义的标志物。本发明充分利用时序关联和组学指标的内在关联,标志物发现的过程经过不同类型信息的相互支撑和验证,鲁棒性高,因此对数据的测量精度和方法的参数优化要求较低,没有丰富实验经验和分析经验的人员也能有效使用。
Description
技术领域
本发明涉及利用组学技术寻找生物标志物的领域,包括但不限于诊断试剂研发,具体地说,是一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统。
背景技术
基础医学研究的一项重要内容是用分子手段表征病理生理的机制并找到关键的基因。关键基因可以是生物标志物,用于疾病的预警、诊断、预后;也可以是靶标,用于疾病的治疗。
标志物与疾病状态仅有相关性,而靶标要求在相关性的基础上,对其干预能够改变疾病的病理生理进程。因此从一定意义上说,靶标是标志物的一个子集,可能从标志物中进一步筛选与病理生理进程有因果关系的基因获得。基础医学研究的一项重要任务,就是要发现与特定病理生理相关的生物标志物。
传统的生物标志物发现手段是假说验证。首先依据领域知识,猜测哪些基因或其他分子可能与病理生理变化相关,然后组织样本进行验证。这一策略存在一定局限:第一,依赖假说。之所以能想到这些假说,主要是因为相关文献的提示。而如果文献有提示,研究人员再去做,那么所发现的标志物的创新性就可能不足,缺乏科学价值和经济价值。另一方面,如果没有足够的前期文献支撑,验证风险大,最终获得标志的成本高,难以获得立项或投资。第二,在大量样本中表征一组候选关键基因的工作,实验量大,周期长,不同的实验人员,不同批次,可能引入较高的测量误差,导致评价不准确。三,在现有知识的基础上提示的候选关键基因与病理进程间的关联可能不够紧密,导致发现的标志物在灵敏度和特异性上都存在缺陷,往往需要多个标志物联合使用,且解读复杂,临床转化价值有限。
为了弥补由假说驱动验证的标志物发现方法的不足,现代生物学发展了一系列组学研究技术,通过全面表征细胞生理,产生数据进行系统筛选,以获得与表型相关的标志物。迄今为止,医学界已经建立了多种组学技术,包括基因组,转录组,蛋白组等,其中以转录组与病理生理特征最为相关。随着组学技术的发展,特别是测序技术的发展,通过转录组测序系统表征细胞状态,无偏倚地在全基因组范围内发现基因标志物的方法,逐渐成为主流。
这一策略的优势是数据驱动而不是假说驱动,对领域知识依赖小,容易发现全新、高敏感度、高特异性的标志物。由于组学分析测试已经形成了完善的科研外包体系,所以组学标志物发现的策略,对于研究人员而言,反而实验工作量更小,数据质量更为可控。
当前,转录组标志物研究的一般流程是首先设计合理的对照,采集足够的样本,通过统计检验发现差异,锁定候选标志物,再组织验证,报告结果。这一转录组标志物发现的策略也有其固有的局限。一是由于转录组方法通量大,一次检测成千上万个基因表达,而使用的样本量却相对很小,所以其单个基因表达量的检测误差和稳定性不高。二,由于转录组测试通量很大,在进行差异的统计显著性评估时,需要用到多重统计矫正。高通量实验的多重统计矫正往往非常严格。只有那些最大的差异才能通过。因此转录组数据中通过多重统计矫正的差异往往是那些信号通路下游的差异,差异会比较大。而关键的信号通路上游的基因,往往不需要过大的变化就能够有效地启动下游通路。这些自身变化较小,但实际上更为重要的上游基因,其变化很难通过多重统计矫正的显著性检验。三是转录组表征了全面的细胞生理,而不仅仅是与病理进程有关的细胞生理。两个对照组在样本量不够大的情况下,有大概率存在与病理进展无关的细胞生理差异,比如说一组恰巧都是胖子的样本而一组恰巧都是瘦子的样本,转录组研究会发现与胖瘦相关的标志物,而不仅仅是与疾病相关的标志物。这些标志物在扩大样本验证时就会失败。
[1]Watson M.CoXpress:differential co-expression in gene expressiondata[J].BMC Bioinformatics.2006,7:509.
[2]Zhou X,Chen P,Wei Q,et al.Human interactome resource and gene setlinkage analysis for the functional interpretation of biologically meaningfulgene sets[J].Bioinformatics.2013,29(16):2024-2031.
[3]Pertea M,Kim D,Pertea G M,et al.Transcript-level expressionanalysis of RNA-seq experiments with HISAT,StringTie and Ballgown[J].NatProtoc.2016,11(9):1650-1667.
[4]Kim D,Langmead B,Salzberg S L.HISAT:a fast spliced aligner withlow memory requirements[J].Nat Methods.2015,12(4):357-360.
[5]Pertea M,Pertea G M,Antonescu C M,et al.StringTie enables improvedreconstruction of a transcriptome from RNA-seq reads[J].Nat Biotechnol.2015,33(3):290-295.
[6]Love M I,Huber W,Anders S.Moderated estimation of fold change anddispersion for RNA-seq data with DESeq2[J].Genome Biol.2014,15(12):550.
发明内容
本发明针对现有转录组标志物发现方法的不足,提供一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法。这一方法可以发挥组学技术的优势,克服劣势。充分利用组学技术的高通量,关注组学数据中具有类似生理意义的指标变化,通过指标间的协同变化克服单个指标的测量误差。将协同变化的指标对应的生理意义抽提出来,与目标病理生理变化对照,找到与目标病理生理变化相关的变化,再追踪这些相关变化是由哪些指标反映的。通过这一方法,我们可以在组学数据中找到那些稳健可靠、同时具有病理生理意义的标志物。
本发明是通过以下技术方案来实现的:
本发明公开了一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法,通过对样本的多时序采样和高通量转录组测序技术,获得目标病理生理过程、对照病理生理过程的动态时序转录组变化;通过基因差异共表达分析,获得在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;通过基因功能关联分析,评估差异共表达集与目标病理生理过程基因集之间的关联强度,从而发现组学数据中稳健可靠、同时具有病理生理意义的标志物。
作为进一步地改进,本发明的具体步骤包括:
步骤1)通过动态时序转录组变化分析和整理模块获得目标病理生理过程的动态时序转录组变化和对照病理生理过程的动态时序转录组变化;
步骤2)通过差异共表达分析模块找到在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;
步骤3)通过功能分析模块判断步骤2)中找到的基因集是否在基因功能关联网络上与目标病理生理过程功能相关,筛选掉不相关的基因集;
步骤4)通过因子贡献度分析模块对于步骤3)筛选剩余的基因集中的基因,评估每个基因与目标病理生理过程在基因功能关联网络上的关联强度;
步骤5)通过差异表达分析模块,评估在需要作出判断的时间点中目标病理生理过程和对照病理生理过程的表达差异;
步骤6)通过候选标志物筛选模块对步骤4、步骤5得到的基因与目标病理生理过程的功能关联强度、基因表达差异设置显著性阈值,报告超过显著性阈值的基因作为候选标志物;
步骤7)通过候选标志物验证模块组织样本验证候选标志物,报告通过验证的候选标志物为目标病理生理过程的标志物。
作为进一步地改进,本发明所述的步骤2)具体包括以下步骤:
2.1通过层次聚类法分别计算目标、对照动态时间序列下各基因间的关系,并剪枝形成各组动态时间序列下的共表达基因类;
2.2寻找在一组动态时间序列下表达是显著关联的、而在另一组动态时序中表达是非显著关联的差异共表达基因类:
所述的步骤2.1中的层次聚类使用Pearson相关性计算样本与样本的距离;所述的步骤2.2中差异共表达基因类的获得使用随机重新采样方法计算各基因类差异显著性(P值),如使用R语言的coXpress包。
作为进一步地改进,本发明所述的步骤3)具体包括以下步骤:
3.1对步骤2)分析获得的差异共表达基因集,选择合适的匀质基因功能关联网络;
3.2在基因功能关联网络中,评估差异共表达基因集与目标病理生理过程基因集间的功能协同,找到与目标病理生理功能相关的差异共表达基因集。
作为进一步地改进,本发明所述的步骤3.1是使用人类互作资源HIR,评估差异共表达基因集与目标病理生理过程基因集间的基因相互作用密度是否显著超过随机基因集之间的基因互作密度,以及这一相互作用密度差异是否确实与基因功能关联网络的拓扑结构相关;所述的步骤3.2是使用GSLA方法来给差异共表达基因集与目标病理过程基因集间的整体功能关联打分,代表目标病理生理过程的基因集通过文献查阅。
作为进一步地改进,本发明所述的步骤4)中,依据步骤3)得到的与目标病理生理功能相关的差异共表达基因集,基于匀质基因功能关联网络(如人类互作资源HIR),追踪哪些基因对目标病理生理功能基因集的整体功能关联贡献最大,通过简单地计数一个差异共表达基因与所有目标病理生理过程基因集中的基因之间的相互作用数目,来代表这个差异共表达基因对基因集之间整体功能关联的贡献。
作为进一步地改进,本发明所述的步骤6)具体包括以下步骤:
6.1依据步骤4)得到的差异共表达基因与目标病理生理功能关联强度结果,设置显著性阈值;
6.2依据步骤5)得到基因表达差异结果,设置显著性阈值;
6.3仅报告均通过两个显著性阈值的基因作为候选标志物。
作为进一步地改进,本发明所述的步骤7)中,组织样本验证候选标志物,通过组织新的样本进行差异和表型关联的验证或使用细胞、动物等模型进行机制验证。
本发明还公开了一种实现基于转录组时序动态变化和基因功能关联发现分子标志物的方法的系统,具体包括如下模块:
动态时序转录组变化分析和整理模块,用于获取mRNA表达数据,包括目标病理生理过程和对照病理生理过程;
差异共表达分析模块,用于对mRNA数据进行差异共表达分析,寻找在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;
功能分析模块,用于寻找在基因功能关联网络上与目标病理生理过程功能相关的基因集;
因子贡献度分析模块,用于评估每个基因与目标病理生理过程功能在基因功能关联网络上的关联强度;
差异表达分析模块,用于在需要做出判断的时间点,评估目标病理生理过程和对照病理生理过程的表达差异;
候选标志物筛选模块,用于根据因子贡献度分析模块和差异表达分析模块,设置显著性阈值,报告超过显著性阈值的基因作为候选标志物;
候选标志物验证模块,用于组织样本验证分析所得的候选标志物的有效性和可靠性。
本发明还公开了一种基于转录组时序动态变化和基因功能关联发现分子标志物的系统发现的分子标志物的应用,通过系统发现天然产物提取物LS0236选择性杀伤肺腺癌细胞的分子标志物,明确定义LS0236适用的肺腺癌患者亚群。
本发明具有以下创新点:
1、基于组学指标变化的连续性,本方法及系统利用时间序列数据与共表达分析连接了标志物要预测的事件与事件发生前的组学指标变化,将事件发生的信息从晚期带入早期,使标志物的筛选具有了导向性;
2、使用基因功能关联网络筛选掉与目标病理生理变化无关的基因,可以排除由于样本组成差异导致的干扰因素,在小样本下也可以发现具有病理生理意义的稳定生物标志物,且验证成功率高;
3、发现的标志物直接参与病理生理机制或干预机制的概率高,经近一步研究验证,可能发展为干预靶标;
4、充分利用时序关联和组学指标的内在关联,标志物发现的过程经过不同类型信息的相互支撑和验证,鲁棒性高,因此对数据的测量精度和方法的参数优化要求较低,没有丰富实验经验和分析经验的人员也能有效使用。
附图说明
图1为本发明技术方案系统的流程示意图;
图2是实施例中GeneX在天然产物提取物LS0236处理后的H460细胞中的表达情况图,图2A是GeneX mRNA表达情况图,图2B是GeneX蛋白表达情况图。
图3LS0236处理后的H460细胞中GeneX在不同时间点、不同细胞组分的表达及分布情况图:线粒体、细胞质、细胞核;
图4是GeneX SiRNA处理后GeneX蛋白表达情况图;
图5是GeneX SiRNA处理后凋亡细胞结果图。
具体实施方式
本发明公开了一种基于转录组时序动态变化和基因功能关联发现分子标志物的系统及方法及应用,本发明提供的方法可以通过以下系统实现:
动态时序转录组变化分析和整理模块,用于获取mRNA表达数据,包括目标病理生理过程和对照病理生理过程;
差异共表达分析模块,用于对mRNA数据进行差异共表达分析,寻找在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;
功能分析模块,用于寻找在基因功能关联网络上与目标病理生理过程功能相关的基因集;
因子贡献度分析模块,用于评估每个基因与目标病理生理过程功能在基因功能关联网络上的关联强度;
差异表达分析模块,用于在需要做出判断的时间点,评估目标病理生理过程和对照病理生理过程的表达差异;
候选标志物筛选模块,用于根据因子贡献度分析模块和差异表达分析模块,设置显著性阈值,报告超过显著性阈值的基因作为候选标志物;
候选标志物验证模块,用于组织样本验证分析所得的候选标志物的有效性和可靠性。
图1为本发明技术方案系统的流程示意图,本发明提供的方法可以通过以下步骤实现:
步骤1)获得目标病理生理过程的动态时序转录组变化和对照病理生理过程的动态时序转录组变化;
步骤2)通过差异共表达分析找到在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;
步骤3)判断步骤3中找到的基因集是否在基因功能关联网络上与目标病理生理过程功能相关,筛选掉不相关的基因集;
步骤4)对于步骤4筛选剩余的基因集中的基因,评估每个基因与目标病理生理过程在基因功能关联网络上的关联强度;
步骤5)通过差异表达分析,评估在需要作出判断的时间点中目标病理生理过程和对照病理生理过程的表达差异;
步骤6)对步骤4、步骤5得到的基因与目标病理生理过程的功能关联强度、基因表达差异设置显著性阈值,报告超过显著性阈值的基因作为候选标志物;
步骤7)组织样本验证候选标志物,报告通过验证的候选标志物为目标病理生理过程的标志物。
步骤1)中的目标病理生理过程和对照病理生理过程指标志物需要区分病理或生理过程,包括但不限于,同一疾病的患者最终生存或死亡(本方法应用于疾病预后);同一人群最终发生或不发生特定疾病(本方法应用于疾病预警预测);同一人群经过一定时间的训练后表现出不同的天赋(本方法应用于天赋预测);同一疾病人群对药物(含中药)或其他治疗产生了不同的响应,或细胞系、动物等疾病模型对药物(含中药)或其他治疗产生了不同的响应(本方法应用于作用机制关键基因鉴定,以及精准治疗伴随检测试剂研发);临床表现疑似特定疾病的患者经过一段时间后确诊是或不是该疾病(用于诊断试剂研发),目标病理生理过程和对照病理生理过程转录组采样的时间点由生理过程变化的速度、复杂性决定,不少于3个时间点。
步骤2)具体包括以下步骤:
2.1通过层次聚类法分别计算目标、对照动态时间序列下各基因间的关系,并剪枝形成各组动态时间序列下的共表达基因类;
2.2寻找在一组动态时间序列下表达是显著关联的、而在另一组动态时序中表达是非显著关联的差异共表达基因类。
步骤2.1中的层次聚类可以使用Pearson相关性计算样本与样本的距离;步骤2.2中差异共表达基因类的获得使用随机重新采样方法计算各基因类差异显著性(P值),例如使用R语言的coXpress包,也可以使用其他能够在一个变量集中寻找在一组时间序列内有协同变化特性在另一组时间序列内没有协同变化特性的变量子集的其他方法。
步骤3)具体包括以下步骤:
3.1对步骤2)分析获得的差异共表达基因集,选择合适的匀质基因功能关联网络;
3.2在基因功能关联网络中,评估差异共表达基因集与目标病理生理过程基因集间的功能协同,找到与目标病理生理功能相关的差异共表达基因集。
步骤3.1可以使用人类互作资源HIR,评估差异共表达基因集与目标病理生理过程基因集间的基因相互作用密度是否显著超过随机基因集之间的基因互作密度,以及这一相互作用密度差异是否确实与基因功能关联网络的拓扑结构相关;步骤3.2可以使用GSLA方法来给差异共表达基因集与目标病理过程基因集间的整体功能关联打分。代表目标病理生理过程的基因集可以通过文献查阅或其他方法确定。
步骤4)中,依据步骤3)得到的与目标病理生理功能相关的差异共表达基因集,基于匀质基因功能关联网络(如人类互作资源HIR),追踪哪些基因对目标病理生理功能基因集的整体功能关联贡献最大,可以通过简单地计数一个差异共表达基因与所有目标病理生理过程基因集中的基因之间的相互作用数目,来代表这个差异共表达基因对基因集之间整体功能关联的贡献。
步骤5)中需要作出判断的时间点指标志物的应用场景对应的样本时间点。通常是第一个时间点,即标志物用于预测将要发生的事件。也可以不是第一个时间点,例如,对于用于疾病预后或精准治疗伴随检测的标志物,可以在治疗开始后,利用患者对治疗的早期响应进行结局和副作用的预测。此类标志物的应用场景对应的样本时间点可能是应用干预后的某个时间点,而不是第一个时间点。表达差异的评估可以是基于统计显著性的方法,使用基于统计显著性的方法时,多重统计检验校正时检验次数应按所有差异共表达基因集中的非冗余基因数确定;也可以使用不基于统计显著性的方法,例如使用表达差异倍数来评估。
步骤6)具体包括以下步骤:
6.1依据步骤4)得到的差异共表达基因与目标病理生理功能关联强度结果,设置显著性阈值;
6.2依据步骤5)得到基因表达差异结果,设置显著性阈值;
6.3仅报告均通过两个显著性阈值的基因作为候选标志物。
步骤6.1应考虑使最终得到的候选标志物在基因功能关联网络上与所有目标病理生理过程基因集间均有功能关联,以确保候选标志物集与整体病理生理的相关性。
步骤7的验证一方面可以通过组织新的样本进行差异和表型关联的验证,也可以使用细胞、动物等模型进行机制验证,例如敲除或过表达候选标志物基因是否对表型或治疗效果的产生影响。如标志物确与表型或治疗效果的产生机制相关,则不但其作为标志物可靠性更高,还可能进一步发展成为干预的靶标。
本方法可以拓展应用于在一组指标中寻找与特定事件发生相关的指标的一般方法。在生物医学领域,可用于代谢组、蛋白组等其他组学标志物的发现。在生物医学以外,可用于金融、工业生产等多个领域。目前这一方法在其他领域广泛开展应用的瓶颈在于缺乏步骤3使用的匀质基因功能关联网络,即指标全集中每个指标之间的意义关联网络。当这一网络被构建出来之后,相关领域的普通技术人员可以从本发明公开的内容直接导出或联想到该领域的应用变形。目前,人工智能领域的知识图谱工作正在建立众多领域的指标意义关联网络。可以预见本方法在其他领域应用的基础将会很快出现。
在生物医学的代谢组、蛋白组等其他组学标志物发现领域,除了建立代谢物功能关联网络、蛋白功能关联网络等其他组学指标功能关联网络之外,还可以将代谢物、蛋白质等其他分子映射到相关的基因,利用已经存在的基因功能关联网络完成本方法的分析。通过此方法发现天然产物提取物LS0236选择性杀伤肺腺癌细胞的分子标志物,可以明确定义LS0236适用的肺腺癌患者亚群,对于基于LS0236的后续新药研发工作具有重要意义。
以下以天然产物提取物LS0236选择性杀伤肺腺癌细胞标志物的发现过程为例,说明本发明的具体实施方式。
肺腺癌是肺癌的主要亚型,是一项重大的健康威胁。天然产物提取物LS0236具有选择性杀伤肺腺癌细胞的活性。因此,发现天然产物提取物LS0236选择性杀伤肺腺癌细胞的分子标志物,可以明确定义LS0236适用的肺腺癌患者亚群,对于基于LS0236的后续新药研发工作具有重要意义。
1、获得目标病理生理过程的动态时序转录组变化。
已知肺腺癌细胞系H460能够被天然产物提取物LS0236杀伤。以LS0236处理H460细胞系杀伤H460细胞的过程为目标病理生理过程。LS0236处理前的时间点为0点,取0小时,2小时,6小时,9小时,18小时,24小时5个时间点的细胞样本,进行转录组测序。采用双端测序,读长150bp,每个样本的目标测序深度为6千万双端读段。
2、获得对照病理生理过程的动态时序转录组变化。
已知肺成纤维细胞系IMR-90不能被天然产物提取物LS0236杀伤。以LS0236处理IMR-90细胞系不能杀伤IMR-90细胞的过程为目标病理生理过程。取同样6个时间点的细胞样本,进行转录组测序。采用双端测序,读长150bp,每个样本的目标测序深度为6千万双端读段。
3、通过基因差异共表达分析找到在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集。
采用HISAT、StringTie转录组分析流程[3]得到每个基因的表达量。其中使用HISATv1.3.0软件[4]默认参数将干净的序列读段比对到人类基因组(GRCh 38.87)上;使用StringTie v1.3.0软件[5]默认参数完成每个样本的基因和转录本的组装和量化。使用R语言的coXpress包分析,获得在H460样本中共表达,IMR-90样本中没有表达相关性的基因集合36个(表1)。
表1 36个在H460样本中共表达、在IMR90样本中没有表达相关性的基因集
Ta:T统计量;prb:随机概率分布;mean.corrc:相关性均值
4、判断每个步骤3中找到的基因集是否在基因功能关联网络上与目标病理生理过程功能相关,筛选掉不相关的基因集。
因当前应用的目标是发现天然产物提取物LS0236选择性杀伤H460细胞的标志物,因此,选择GO途径中包含Apoptosis、Necrosis、Autophagy这三个关键字或其形容词形式(例如Apoptotic)的途径作为目标病理生理过程。考虑到应用方便,可直接使用GSLA的在线工具分析步骤3中获得的差异共表达基因集与所有GO生物途径基因集之间的显著功能关联,再筛选其中的目标病理生理过程。结果发现,2个差异共表达基因集与17个目标病理生理过程具有显著的功能关联(表2、3)。
表2 基因集(ID23)中与目标病理生理过程显著关联的功能情况
表3 基因集(ID15)中与目标病理生理过程显著关联的功能情况
5、对于步骤4筛选剩余的基因集中的基因,评估每个基因与目标病理生理过程在基因功能关联网络上的关联强度。
对每个步骤4中发现的差异共表达基因集与目标病理生理过程基因集,查询GSLA输出结果,统计差异共表达基因集中的基因与目标病理生理过程基因集的基因的功能关联数目,得表4。
表4 基因与目标病理生理过程基因集的基因的功能关联强度
6、在需要作出判断的时间点,对于步骤4筛选剩余的基因集中的基因,评估其在目标病理生理过程和对照病理生理过程的表达差异。
因当前应用的目标是发现天然产物提取物LS0236选择性杀伤H460细胞的标志物,标志物应在用药前可以判断H460细胞系是否能被LS0236杀伤,即需要标志物作出判断的时间点为0小时。使用R语言DESeq2包[6],步骤3中发现的差异共表达基因集中的基因,在时间点0得到H460与IMR-90的显著差异表达基因共357个(P-value<0.05);差异倍数(FoldChange)>4的基因有2946个。
7、对步骤5、步骤6得到的基因与目标病理生理过程的功能关联强度,基因表达差异设置显著性阈值,报告超过显著性阈值的基因作为候选标志物。
以步骤5发现差异共表达基因具有至少31个与目标病理生理过程基因集中的基因的相互作用、步骤6发现差异共表达基因在时间点0存在至少4倍表达量差异为阈值,可以得到6个候选标志物(表5)。
表5 6个候选标志物
8、组织样本验证候选标志物,报告通过验证的候选标志物为目标病理生理过程的标志物。
步骤7中发现的候选标志物,排名第一的DDR2在时间点0已经表现出表达差异的显著性,其功能有大量文献证明可能与肿瘤死亡相关。排名第二的GeneX时间点0虽有4倍以上的表达量差异,但并未通过统计显著性检验。因此,本应用实例重点验证了GeneX在天然产物提取物LS0236选择性杀伤H460细胞过程中的作用。
首先,通过实时荧光定量PCR分析GeneX在经过天然产物提取物LS0236处理后的H460细胞中的表达情况。图2是实施例中GeneX在天然产物提取物LS0236处理后的H460细胞中的表达情况图,结果显示,在LS0236处理6小时以后,GeneX mRNA的表达量显著增加;在9小时表达量下降,但仍维持较高水平(图2A是GeneX mRNA表达情况图)。同时,基于Westernblot发现GeneX的表达量在LS0236处理后的H460细胞中也存在上调,与GeneX mRNA表达量趋势相类似(图2B是GeneX蛋白表达情况图)。
图2是GeneX在天然产物提取物LS0236处理后的H460细胞中的表达情况,A.实时定量PCR结果图显示LS0236处理后H460细胞中GeneX mRNA的表达量增加;B.Western blot检测结果图显示LS0236处理后H460细胞中GeneX蛋白表达量增加;
已知GeneX从细胞核到线粒体或ER的转位在自噬和细胞凋亡等相关生物学过程中发挥关键作用。因此在后续实验过程中,需要确定在天然产物提取物LS0236处理后的H460细胞中是否存在GeneX转移到线粒体现象。通过时间过程分析(图3是LS0236处理后的H460细胞中GeneX在不同时间点、不同细胞组分的表达及分布情况图),结果显示,在LS0236处理后的1小时,GeneX位于线粒体富集的HM组分(图3,左侧,是在线粒体GeneX的表达及分布情况图)。该结果表明,在处理后的1小时,GeneX已开始转移至线粒体。同时,结果表明GeneX转移到线粒体的数量随着LS0236距离处理的时间点越久而增多。相应的,位于细胞核的GeneX随着时间的推移而减少(图3,右侧,是在细胞核GeneX的表达及分布情况图)。同时在细胞质中也发现存在部分GeneX(图3,中侧,是在细胞质GeneX的表达及分布情况图)。此外,使用SiRNA方法确定GeneX是否介导/参与天然产物提取物LS0236对H460细胞的凋亡或自噬等相关的生物学过程。结果显示,在GeneX SiRNA瞬时转染H460细胞后,天然产物提取物LS0236诱导的GeneX表达量降低(图4是GeneX SiRNA处理后GeneX蛋白表达情况图),同时DAPI细胞核染色结果显示凋亡细胞数减少(图5是Gene SiRNA处理后凋亡细胞结果图)。基于上述结果可得,GeneX是天然产物提取物LS0236诱导H460细胞凋亡的关键基因。
因此,GeneX确在天然产物提取物LS0236选择性杀伤H460的过程中起机制性作用,是优秀的标志物,也是天然产物提取物LS0236选择性杀伤肺腺癌的机制的关键基因之一。
以上所述仅是本发明的优选实施方式,应当指出,对于本技术领域的普通技术人员,在不脱离本发明原理的前提下,还可以做出若干改进和补充,这些改进和补充也应视为本发明的保护范围。
Claims (10)
1.一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,通过对样本的多时序采样和高通量转录组测序技术,获得目标病理生理过程、对照病理生理过程的动态时序转录组变化;通过基因差异共表达分析,获得在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;通过基因功能关联分析,评估差异共表达集与目标病理生理过程基因集之间的关联强度,从而发现组学数据中稳健可靠、同时具有病理生理意义的标志物。
2.根据权利要求1所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,本发明的具体步骤包括:
步骤1)通过动态时序转录组变化分析和整理模块获得目标病理生理过程的动态时序转录组变化和对照病理生理过程的动态时序转录组变化;
步骤2)通过差异共表达分析模块找到在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;
步骤3)通过功能分析模块判断步骤2)中找到的基因集是否在基因功能关联网络上与目标病理生理过程功能相关,筛选掉不相关的基因集;
步骤4)通过因子贡献度分析模块对于步骤3)筛选剩余的基因集中的基因,评估每个基因与目标病理生理过程在基因功能关联网络上的关联强度;
步骤5)通过差异表达分析模块,评估在需要作出判断的时间点中目标病理生理过程和对照病理生理过程的表达差异;
步骤6)通过候选标志物筛选模块对步骤4、步骤5得到的基因与目标病理生理过程的功能关联强度、基因表达差异设置显著性阈值,报告超过显著性阈值的基因作为候选标志物;
步骤7)通过候选标志物验证模块组织样本验证候选标志物,报告通过验证的候选标志物为目标病理生理过程的标志物。
3.根据权利要求2所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,所述的步骤2)具体包括以下步骤:
2.1通过层次聚类法分别计算目标、对照动态时间序列下各基因间的关系,并剪枝形成各组动态时间序列下的共表达基因类;
2.2寻找在一组动态时间序列下表达是显著关联的、而在另一组动态时序中表达是非显著关联的差异共表达基因类:
所述的步骤2.1中的层次聚类使用Pearson相关性计算样本与样本的距离;所述的步骤2.2中差异共表达基因类的获得使用随机重新采样方法计算各基因类差异显著性(P值),如使用R语言的coXpress包。
4.根据权利要求2或3所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,所述的步骤3)具体包括以下步骤:
3.1对步骤2)分析获得的差异共表达基因集,选择合适的匀质基因功能关联网络;
3.2在基因功能关联网络中,评估差异共表达基因集与目标病理生理过程基因集间的功能协同,找到与目标病理生理功能相关的差异共表达基因集。
5.根据权利要求4所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,所述的步骤3.1是使用人类互作资源HIR,评估差异共表达基因集与目标病理生理过程基因集间的基因相互作用密度是否显著超过随机基因集之间的基因互作密度,以及这一相互作用密度差异是否确实与基因功能关联网络的拓扑结构相关;所述的步骤3.2是使用GSLA方法来给差异共表达基因集与目标病理过程基因集间的整体功能关联打分,代表目标病理生理过程的基因集通过文献查阅。
6.根据权利要求2所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,所述的步骤4)中,依据步骤3)得到的与目标病理生理功能相关的差异共表达基因集,基于匀质基因功能关联网络(如人类互作资源HIR),追踪哪些基因对目标病理生理功能基因集的整体功能关联贡献最大,通过简单地计数一个差异共表达基因与所有目标病理生理过程基因集中的基因之间的相互作用数目,来代表这个差异共表达基因对基因集之间整体功能关联的贡献。
7.根据权利要求2所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,所述的步骤6)具体包括以下步骤:
6.1依据步骤4)得到的差异共表达基因与目标病理生理功能关联强度结果,设置显著性阈值;
6.2依据步骤5)得到基因表达差异结果,设置显著性阈值;
6.3仅报告均通过两个显著性阈值的基因作为候选标志物。
8.根据权利要求2或3或5或6或7所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法,其特征在于,所述的步骤7)中,组织样本验证候选标志物,通过组织新的样本进行差异和表型关联的验证或使用细胞、动物等模型进行机制验证。
9.一种实现如权利要求1或2或3或5或6或7所述的基于转录组时序动态变化和基因功能关联发现分子标志物的方法的系统,其特征在于,具体包括如下模块:
动态时序转录组变化分析和整理模块,用于获取mRNA表达数据,包括目标病理生理过程和对照病理生理过程;
差异共表达分析模块,用于对mRNA数据进行差异共表达分析,寻找在目标病理生理过程中共表达、在对照病理生理过程中没有表达相关性的基因集;
功能分析模块,用于寻找在基因功能关联网络上与目标病理生理过程功能相关的基因集;
因子贡献度分析模块,用于评估每个基因与目标病理生理过程功能在基因功能关联网络上的关联强度;
差异表达分析模块,用于在需要做出判断的时间点,评估目标病理生理过程和对照病理生理过程的表达差异;
候选标志物筛选模块,用于根据因子贡献度分析模块和差异表达分析模块,设置显著性阈值,报告超过显著性阈值的基因作为候选标志物;
候选标志物验证模块,用于组织样本验证分析所得的候选标志物的有效性和可靠性。
10.一种如权利要求9所述的基于转录组时序动态变化和基因功能关联发现分子标志物的系统发现的分子标志物的应用,其特征在于,通过系统发现天然产物提取物LS0236选择性杀伤肺腺癌细胞的分子标志物,明确定义LS0236适用的肺腺癌患者亚群。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910742777.8A CN110379459B (zh) | 2019-08-13 | 2019-08-13 | 一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910742777.8A CN110379459B (zh) | 2019-08-13 | 2019-08-13 | 一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110379459A true CN110379459A (zh) | 2019-10-25 |
CN110379459B CN110379459B (zh) | 2021-06-29 |
Family
ID=68259011
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910742777.8A Active CN110379459B (zh) | 2019-08-13 | 2019-08-13 | 一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110379459B (zh) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111261243A (zh) * | 2020-01-10 | 2020-06-09 | 华南理工大学 | 一种基于相对熵指标检测复杂生物系统相变临界点的方法 |
CN114882955A (zh) * | 2022-04-08 | 2022-08-09 | 广州国家实验室 | 转录组图像生成装置、方法和应用 |
CN115828093A (zh) * | 2022-11-02 | 2023-03-21 | 四川帕诺米克生物科技有限公司 | 组学样本的分析方法、装置、电子设备及存储介质 |
CN116453594A (zh) * | 2023-06-15 | 2023-07-18 | 北京望石智慧科技有限公司 | 基因共表达状态的量化分析方法及装置、设备和介质 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20150038353A1 (en) * | 2013-08-01 | 2015-02-05 | National Taiwan University | Molecular markers for detecting nosema infection |
CN106126893A (zh) * | 2016-06-17 | 2016-11-16 | 浙江大学 | 一种基于基因功能关联网络发现慢性病机制及其预警干预策略的方法 |
CN107463796A (zh) * | 2017-07-12 | 2017-12-12 | 北京航空航天大学 | 基于基因共表达网络传播分析的早期致病因子探测方法 |
CN108732349A (zh) * | 2017-04-20 | 2018-11-02 | 中国科学院上海生命科学研究院 | Lta4h作为指示肝内结节及早期预警肝癌的生物标志物 |
CN109841281A (zh) * | 2017-11-29 | 2019-06-04 | 郑州大学第一附属医院 | 基于共表达相似性识别肺腺癌早期诊断标识及风险预测模型的构建方法 |
-
2019
- 2019-08-13 CN CN201910742777.8A patent/CN110379459B/zh active Active
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20150038353A1 (en) * | 2013-08-01 | 2015-02-05 | National Taiwan University | Molecular markers for detecting nosema infection |
CN106126893A (zh) * | 2016-06-17 | 2016-11-16 | 浙江大学 | 一种基于基因功能关联网络发现慢性病机制及其预警干预策略的方法 |
CN108732349A (zh) * | 2017-04-20 | 2018-11-02 | 中国科学院上海生命科学研究院 | Lta4h作为指示肝内结节及早期预警肝癌的生物标志物 |
CN107463796A (zh) * | 2017-07-12 | 2017-12-12 | 北京航空航天大学 | 基于基因共表达网络传播分析的早期致病因子探测方法 |
CN109841281A (zh) * | 2017-11-29 | 2019-06-04 | 郑州大学第一附属医院 | 基于共表达相似性识别肺腺癌早期诊断标识及风险预测模型的构建方法 |
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111261243A (zh) * | 2020-01-10 | 2020-06-09 | 华南理工大学 | 一种基于相对熵指标检测复杂生物系统相变临界点的方法 |
CN111261243B (zh) * | 2020-01-10 | 2023-04-21 | 华南理工大学 | 一种基于相对熵指标检测复杂生物系统相变临界点的方法 |
CN114882955A (zh) * | 2022-04-08 | 2022-08-09 | 广州国家实验室 | 转录组图像生成装置、方法和应用 |
CN115828093A (zh) * | 2022-11-02 | 2023-03-21 | 四川帕诺米克生物科技有限公司 | 组学样本的分析方法、装置、电子设备及存储介质 |
CN115828093B (zh) * | 2022-11-02 | 2024-04-05 | 四川帕诺米克生物科技有限公司 | 组学样本的分析方法、装置、电子设备及存储介质 |
CN116453594A (zh) * | 2023-06-15 | 2023-07-18 | 北京望石智慧科技有限公司 | 基因共表达状态的量化分析方法及装置、设备和介质 |
CN116453594B (zh) * | 2023-06-15 | 2023-11-21 | 北京望石智慧科技有限公司 | 基因共表达状态的量化分析方法及装置、设备和介质 |
Also Published As
Publication number | Publication date |
---|---|
CN110379459B (zh) | 2021-06-29 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110379459A (zh) | 一种基于转录组时序动态变化和基因功能关联发现分子标志物的方法及系统 | |
Domínguez-Salas et al. | Impact of general cognition and executive function deficits on addiction treatment outcomes: Systematic review and discussion of neurocognitive pathways | |
Liu et al. | Single-cell RNA-seq analysis of the brainstem of mutant SOD1 mice reveals perturbed cell types and pathways of amyotrophic lateral sclerosis | |
CN103415624B (zh) | 胰腺癌生物标记及其用途 | |
CN106599616B (zh) | 基于duplex-seq的超低频突变位点检测分析方法 | |
CN109872776A (zh) | 一种基于加权基因共表达网络分析对胃癌潜在生物标志物的筛选方法及其应用 | |
CN103958662A (zh) | 用于疾病生物标记鉴定和样品质量评价的优选样品处理和加工方案的选择 | |
CN109686439A (zh) | 遗传病基因检测的数据分析方法、系统及存储介质 | |
CN113053535B (zh) | 一种医疗信息预测系统及医疗信息预测方法 | |
JPWO2014080447A1 (ja) | データ解析装置、データ解析プログラム | |
JP2009505231A (ja) | 複数のサンプルから得られる代謝物のデータを、コンピュータシステムのデータベースを用いて比較および編集するためのシステム、方法、ならびにコンピュータプログラム | |
Barnhart‐Magen et al. | Differential diagnostics of thalassemia minor by artificial neural networks model | |
CN110021346A (zh) | 基于RNAseq数据的基因融合与突变检测方法及系统 | |
CN108038352A (zh) | 结合差异化分析和关联规则挖掘全基因组关键基因的方法 | |
CN112669960A (zh) | 一种基于机器学习方法的肝脏纤维化预测模型的构建方法、预测系统、设备和存储介质 | |
van Albada et al. | Bringing anatomical information into neuronal network models | |
CN116597916A (zh) | 一种基于器官芯片和深度学习的抗肿瘤化合物预后药效的预测方法 | |
EP2701579A2 (en) | Stratifying patient populations through characterization of disease-driving signaling | |
CN110322963A (zh) | 一种新生儿遗传代谢病检测分析方法、装置及系统 | |
CN110111890A (zh) | 一种基于基因测序技术的个体精准养生方法 | |
CN114822690A (zh) | 应用于全基因组表达谱数据的多类别多功能智能分类方法 | |
CN110396538A (zh) | 偏头痛生物标志物及其用途 | |
Birk et al. | Large-scale characterization of cell niches in spatial atlases using bio-inspired graph learning | |
CN105243294B (zh) | 一种用于预测癌症病人预后相关的蛋白质对的方法 | |
CN110070942A (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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |