ActivePathways整合多维组学通路分析

Integrative pathway enrichment analysis of multivariate omics data

2020 Feb 5发表于Nature Communications

导语

GUIDE ╲

通路富集分析是利用现有基因和生物过程知识解释高通量（组学）数据的一个重要步骤。常见的应用是对大量的候选基因进行分子通路、生物过程和其他功能注释的统计学富集。基因组学、转录组学、蛋白质组学和表观基因组学实验能够作为基础生物学的一个互补领域，如很多大型研究计划都是这种整合分析，比如The Cancer Genome Atlas (TCGA)、Clinical Proteome Tumor Analysis Consortium (CPTAC)、International Cancer Genome Consortium (ICGC)和Genotype-Tissue Expression (GTEx)等。越来越需要同时分析多个候选基因列表的特征通路。

背景介绍

有许多方法可以用来解释单基因列表。如GSEA算法可以识别基因表达数据集中上调和下调的通路。基于网络的方法，如Panther、ToppCluster和g:Profiler，可以在排序或未排序的基因列表中识别显著富集的通路，通常适用于基因和蛋白质的各种分析。一些方法允许对多个输入基因列表进行分析，但这些方法主要依赖于可视化而不是数据整合来评估不同基因列表的贡献。最后，目前尚无方法对编码和非编码突变的全基因组测序（WGS）数据的统一通路分析，或者将其与其他类型的DNA畸变（如拷贝数改变和基因组平衡重排）整合起来。

癌症基因组的特征是多种类型的突变，包括单核苷酸变异(SNVs)、小插入缺失(indels)、拷贝数改变和易位。编码和非编码driver突变的挖掘是大型癌症全基因组测序工作的主要目标。PCAWG汇总了来自ICGC和TCGA项目产生的38种肿瘤类型的2658个癌症样本的全基因组测序数据，识别生殖系变异和体细胞获得性突变，PCAWG肿瘤变异的一致分析，生成了一个高可信度的蛋白编码driver基因(CDS)、5 '和3 '非编码非翻译元件（UTRs）、启动子和增强子中的driver基因突变目录。使用已知分子通路和基因互作网络，对PCAWG driver突变的通路和网络分析，能够进一步发现包括非编码基因组在内的罕见候选driver突变。

本工作使用数据融合技术，目的是解决多组数据的综合通路分析，研发ActivePathways方法。它在多个数据集中检测到显著富集的通路，包括那些在单个数据集中不明显的通路。然后提出了几个分析来证明这种方法。首先，将癌症driver基因与PCAWG数据集预测的编码和非编码突变整合在一起，揭示了大量的编码和非编码突变的过程和附加基因。其次，将患者的临床信息和METABRIC据中乳腺癌的转录组和拷贝数的改变相结合，发现乳腺癌亚型的预后通路和过程。第三，将GTEx project正常组织的转录组数据与ChIP-seq数据进行整合，推断出Hippo通路下游的基因调控网络和组织生长和再生的生物学过程。因此，ActivePathways是一个可以结合各种多组数据集的多用途的方法。

数据介绍

1. 从g:Profiler web server下载Gene Ontology(GO)的生物学过程和来自Reactome数据库的分子通路数据。

2. ICGC-TCGA PCAWG project的47套2583个样本的全基因数据，包括蛋白质编码序列，启动子，增强子和非翻译区域等。

3. 已知癌症基因：COSMIC Cancer Gene Census (CGC) database

4. METABRIC 数据库的四种亚型(HER2-enriched, basal-like, luminal-A, luminal-B)的1780个乳腺癌样本的mRNA表达和拷贝数变异数据。

5. GTEx v7 data portal 的RNA-seq数据，包含来自53种组织的11,688个样本中的21,518个蛋白编码基因。

6. KEGG 数据库的308个 Hippo通路相关基因。

7. ENCODE project转录因子结合位点 (TFBS) 对应基因ChIP-seq数据。

结果解析

01

ActivePathways多组学通路富集分析

ActivePathways方法主要分为三步(Fig. 1)，是对该研究团队之前的工作一个扩展，它需要两个输入数据集。第一个输入是一个p值表(Fig. 1a)，其中行是基因，列是数据集。列可以代表差异基因表达、基因必要性、突变或拷贝数变异负荷等p值。第二个输入是一组基因集，是对应相应的生物学过程的通路基因，通常用来自 gene ontology(GO)的生物学过程和来自Reactome数据库的分子通路，也可用许多其他类型的基因集，如转录因子或microRNAs的靶点基因集。

在ActivePathways的第一步（Fig. 1b），使用Brown, M. B.建立的方法（是Fisher’s联合概率P值方法的扩展）对每个基因计算不同的组学数据集基因的联合P值，Brown的方法考虑了数据集之间的相关性，从而为多个相似的组学数据集支持的基因提供了更保守的重要性估计。整合后输入基因列表按重要性（Brown P-values）递减排序，并且使用P < 0.1过滤筛选出重要基因。第二步，用排序超几何检验对候选基因的排序表中的通路进行统计富集。该测试每次只考虑一个通路的基因集（生物过程、分子途径和其他基因注释），并分析从排名最高的基因列表中不断增加的输入基因子集。排序超几何检验（下图）的目的是获取与少数 top-ranking基因（P值显著）紧密相关的通路，以及与较大的输入基因子集相关的更广泛的通路。然后使用HB方法校正P值，选择整合后的基因列表显著富集的通路(Qpathway < 0.05)。

在第三步（Fig. 1c），分别对单个组学数据集的基因列表进行类似的分析，以确定能够支持第二步中确定的整合通路分析结果的组学证据。第三步还强调了只有通过整合数据集才能识别的通路，而单独在任何单个组数据集中都无显著性。该方法提供了用于 Cytoscape的EnrichmentMap app的输入文件，用于可视化结果通路。

02

整合2658个样本癌症基因组的编码driver和非编码driver

接下来，基于体细胞SNVs和indels数据对PCAWG项目中预测的driver基因进行了整合通路分析，这里涉及多种肿瘤类型的47套样本数据。

ActivePathways在89%的患者中发现了至少一个显著富集的生物学过程或通路(42/47,Qpathway <0.05，排序超几何检验) (Fig. 2a)。大多数样本群组在基因蛋白编码突变支持的通路中显示富集(37/47或79%)，而大多数目前已知的癌症driver基因有频繁的蛋白质编码突变。

基因的非编码突变富集分析显示，24/47(51%)的样本群显示出显著的富集通路，若只分析与UTRs、启动子或增强子对应的非编码driver因子得分时，这些通路很显著。一些频繁突变的通路在大多数PCAWG的肿瘤样本群(41/47或87%)中，这些通路在整合编码和非编码突变时很明显，但在分别考虑编码或非编码突变时却未被检测到，这强调了本工作整合方法的价值，较大的数据集能够更好地区分涉及生物学通路和过程的罕见突变基因。

接下来研究了16种肿瘤类型的1773个样本的腺癌样本群，整合通路分析发现有432个基因在526条通路中显著富集（Qpathway <0.05）（Fig. 2b）。大多数通路仅由具有频繁突变的编码基因（328/526或62％）富集，但是编码和非编码基因突变均在另外的101条通路（19％）中富集，仅在编码和非编码突变的综合分析中具有显著性的是72条通路（14％），而25条（5 ％）仅富集非编码突变。

03

基于通路的罕见突变癌基因优先排序

通路分析可以识别一些在基于基因的分析中无法检测到的基因。ActivePathways通过跨多个数据集更加增强了这种基因挖掘。在PCAWG编码和非编码突变的通路分析中，本工作重点关注了一组通过编码和非编码突变整合检测到的肾脏发育过程(Fig. 2 c, d)。ActivePathways发现了18个参与这些过程的基因，其中只有5个在PCAWG项目的分析中被预测为driver基因。ActivePathways可以利用功能基因注释和多重组学信号来发现更多的候选基因。

接下来评估了数据整合方法的效果，分析在腺癌样本群中检测到的333个与通路相关的候选基因(Fig. 2e)，结果发现这些基因包括了相当一部分已知的癌症基因，60/64个显著突变的基因也在PCAWG分析中被鉴定出来，47/333个在COSMIC Cancer Gene Census中注释为癌基因。

ActivePathways能够检测这些额外的基因有几个原因。首先，对整合的基因列表使用宽松的统计域值进行筛选。其次通过数据融合的方法，将基因的编码区、启动子、UTRs、增强子等对应的多个较弱的p值结合起来，得到单个较强的p值，从而对某些基因进行了更新。因此，与单个突变数据集中的原始排序相比，整合程序优先选择特定的通路相关基因。ActivePathways在逐个基因分析中发现了其他未被检测到的候选基因，是由于它们在多重组学数据中的通路信号中，所以能够突出表现。

04

整合乳腺癌预后的CNA和mRNA信号

接下来研究了与乳腺癌患者预后相关的通路和过程。利用METABRIC 数据库的四种亚型的1780个人乳腺癌样本的mRNA表达和拷贝数变异数据。根据基因的表达值中值将样本分为两组（高表达组和低表达组），然后进行生存分析，识别与生存相关的基因。根据基因的拷贝数状态将样本分为两组（扩增和其他；缺失和其他），然后进行生存分析。

ActivePathways在乳腺癌4个亚型中发现了192个显著富集的 GO生物学过程和Reactome通路，其中9个通路在多种癌症亚型中富集，33个通路仅在整合通路分析中发现。乳腺癌患者预后相关基因主要富集的通路有免疫反应、细胞凋亡、核糖体生物发生和染色体分离(Fig. 3a)。

basal-like 和HER2-enriched乳腺癌中，免疫活性与预后基因相关。在肿瘤细胞或肿瘤邻近细胞中，高表达与改善患者预后相关的基因，HER2-enriched中有50/61个基因，basal型中有78/113个基因(Fig. 3b)。这些基因中只有少数(10个)在两种乳腺癌亚型中都是显著的，这表明了不同亚型中免疫活动的不同模式。基于通路的研究结果说明，乳腺肿瘤细胞和周围微环境的免疫活动对肿瘤的进展和预后有负面影响。

HER2-enriched和luminal-A乳腺癌富集到GO的凋亡相关的过程，如“凋亡过程负调控”和“程序性细胞死亡”(Fig. 3c)。抗凋亡通路仅在整合分析中检测到，而在基因组和转录组基因签名中未检测到。在负调控凋亡的基因中，DUSP1在HER2-enriched的乳腺癌中提供了最强的预后信号，这个显著性在基因表达和拷贝数扩增上都有体现(Fig. 3d)。ActivePathways可用于整合临床数据和分子改变的多组学信息，这些分析可以为功能研究和生物标志物的开发提供线索。

05

Hippo通路的TFs的共表达和DNA结合靶点

为了证明ActivePathways在基因调控研究中的作用，接下来分析了GTEx proiect的非癌组织的转录组数据。重点研究了Hippo信号通路在器官大小控制、组织稳态和癌变中的作用，并研究了两种转录因子(TFs) YAP和TAZ（由YAP1和WWTR1编码）的下游调控网络。在哺乳动物中，YAP和TAZ是进化上保守的Hippo信号的主调节因子。

对Hippo通路的两个主调节因子的转录组学和表观基因组学数据进行了整合通路富集。首先，共表达分析和Robust rank aggregation(RRA)方法预测了YAP和TAZ的转录靶基因(分别为1898和1319个基因)。其次，研究了来自ReMap数据库中重新分析的YAP ChIP-seq study的YAP的YAP个靶基因。将三个基因列表以及相应的显著性值输入到ActivePathways中进行整合分析。综合分析YAP/TAZ的转录和DNA结合的靶基因，发现有225个显著富集的GO生物学过程和Reactome通路(Fig.4a)。富集结果包括Hippo信号通路，以及发育和形态发生、细胞运动、肌动蛋白骨架和细胞间连接的组织、EGFR、Wnt、Robo、TGF-beta、rho等信号转导通路。使用扩展Hippo通路基因列表（106个），有36/308基因在通路富集(Fig.4b)。

ActivePathways还可以与其他类别的功能基因集(如TF靶基因)一起解释组学数据。接下来使用了ENCODE项目的161个TFs的靶基因数据集，进一步阐明在YAP / TAZ下游的下游的网络。发现了一个由17个TFs和1426个靶基因组成的调控网络，该网络丰富了YAP / TAZ调控组（Q <0.05，排名超几何检验）(Fig.4c)。调控网络显著富集了50个Hippo相关基因和6个核心Hippo基因，类似上述基于通路的分析。

然而，这两种分析揭示了不同的基因：两种分析共同发现886个基因，1180个基因仅在基于通路的分析中发现，540个基因仅在基于TF的分析中发现(Fig.4d)。因此，TF靶基因的整合富集分析为GO通路分析提供了补充信息。证明了ActivePathways在研究基因调控网络和组学数据集研究具有互补作用。

06

评估ActivePathways鲁棒性和敏感性

接下来使用PCAWG预测的癌症driver基因数据集对ActivePathways进行测试。首先，比较了ActivePathways与六种方法（(Hierarchical HotNet,SSA−ME, NBDI, induced subnetwork

analysis, CanIsoNet, hypergeometric test）的性能。这些方法使用分子相互作用网络、功能基因集和/或转录组数据来分析预测癌症driver基因的PCAWG泛癌数据集。然后对使用方法预测的蛋白编码和非编码突变的driver基因列表进行一致性分析。ActivePathways有最高的准确率：100%的编码driver基因(87/87)和85%的非编码候选基因(79/93)被检测到(Fig.5a)。因此在使用通路和网络环境对已知和候选癌症driver基因进行优先排序时，ActivePathways建议使用其他几种方法的集成。接下来将ActivePathways的性能与考虑单个统计过滤基因列表的标准通路富集分析方法(使用排秩超几何)进行了比较(Fig.5b)。与标准方法相比，ActivePathways显示出更高的通路富集分析敏感性，尤其是涉及非编码突变的通路，而这些非编码突变在任何单基因列表中都没有明显的表达。

接下来评估了ActivePathways对参数变化和丢失数据的鲁棒性。改变了用于富集分析的排序基因列表的参数Pgene（默认阈值Pgene < 0.1），大多数PCAWG样本群(40/47或85%)有显著富集的通路，当阈值限制非常严格(Pgene < 0.001)时，富集通路减少了67%。随机删除初始基因的scores，甚至去除50%的基因driverp值时（P < 0.001），大多数样本群(37/47或79%)至少有一个显著富集通路，但平均少66%。然后评估了ActivePathways的预期假阳性率，使用破坏通路基因的注释的方法，错误率只是略有升高。总的来说，这些测试表明，ActivePathways是一种对多组学数据的综合分析以挖掘显著丰富的通路和过程的敏感而稳健的方法。

小编总结

本工作研发ActivePathways方法，是整合了多组数据然后进行富集分析的方法，可以挖掘单数据集未识别的通路等。然后从三方面进行验证，首先用PCAWG数据集识别的编码和非编码突变富集的通路。然后将患者的临床信息与乳腺癌的转录组和拷贝数的改变相结合，发现乳腺癌亚型的预后通路和过程，最后推断Hippo通路下游的基因调控网络。

引用：

Paczkowska M, Barenboim J, Sintupisut N, et al. Integrative pathway enrichment analysis of multivariate omics data. Nat Commun. 2020;11(1):735. Published 2020 Feb 5. doi:10.1038/s41467-019-13983-9