PRODIGY: personalized prioritization of driver genes

摘要

动机

癌症的进化是由少数破坏细胞进程，导致异常增殖和肿瘤发展的体细胞突变驱动的，而大多数体细胞突变对进展没有影响。 区分那些驱动患者肿瘤发生的突变基因是癌症治疗的主要目标：对这些基因及其作用途径的了解可以阐明疾病的机制，并指明潜在的治疗方法和药物靶标。 当前的研究主要集中在队列水平的驱动基因识别上，但是患者特异性驱动基因的识别仍然是一个挑战。

方法

我们开发了一种针对患者基因的驱动基因特异性排序的新算法。 PRODIGY算法可分析患者的表达和突变情况，以及有关已知通路和蛋白质-蛋白质相互作用的数据。 PRODIGY利用获奖的斯坦纳树模型来量化每一个突变基因对每一个非控制通路的影响。 突变基因按其对所有失控途径的总体影响进行排序。

结果

在对超过2500名患者的5个TCGA癌症队列进行测试并与经过验证的驱动基因进行比较后，Prodigy胜过了现有方法和基于网络中心性测度的排名。 我们的结果指出了驱动基因的多效性，并表明Prodigy能够识别甚至非常罕见的驱动。 因此，Prodigy向个性化医学和治疗迈进了一步。

1介绍

癌症是一个进化过程，正常细胞会累积各种基因组和表观基因组改变，包括单核苷酸变异（SNVs）和染色体畸变。 这些改变中的一些赋予突变的细胞生长和阳性选择优势，从而引起密集的增殖和肿瘤（Stratton等，2009）。 这种改变可以通过种系突变而遗传，如乳腺癌中BRCA1和BRCA2的情况（King等，2003）或体细胞发生（Stratton等，2009）。 虽然体细胞突变也发生在正常细胞中，但它们是中性的或引起凋亡，但不会导致转化为癌细胞。

1.1驱动突变

使细胞获得选择性生长优势从而将其“驱动”成肿瘤的突变事件称为驱动突变(或驱动事件)，发生这些突变的基因称为驱动基因。相比之下，与驱动突变在同一基因组中发生但不会影响适应性的其他突变称为乘客突变（Vogelstein等，2013）。

观察到的突变的总数在肿瘤组织之间变化。 Kim＆Kim（Kim and Kim，2018）分析了TCGA的数十名癌症患者队列（Weinstein等人，2013），发现在某些癌症亚型中，每个肿瘤的平均体细胞突变数可以达到数百个。关于每个肿瘤中观察到的突变中的驱动子突变的数量，存在着非常广泛的争论（Stratton等，2009； Vogelstein等，2013； Anna C. Schinzel，2008），但共识是这个数字非常低 。显然，有许多因素会导致驱动子数量的变化，包括肿瘤的进展阶段（Vogelstein和Kinzler，2015），肿瘤的起源组织（Kandoth等，2013），环境特性（例如吸烟） （Govindan等人，2012）和其他因素（例如年龄）（Xie等人，2014）。 Tomasetti等。 （Tomasetti等人，2015）表明，仅三个驱动基因突变就足以发展为肺癌和大肠癌。 诺丁（Nordling，1953）和阿米蒂奇（Armitage and Doll，1954）提出了六到七种典型的驱动子数量。

因此，从乘客突变中区分驱动子是一个挑战。 出于以下几个原因，在癌症研究中，尤其是在个性化癌症医学中，这样做的优先级很高：1）对驱动程序及其操作机制的了解可以建议潜在的治疗方法和药物目标。 2）根据分子特征而不是疾病器官进行癌症治疗，可以为尚未接受针对其特定类型癌症的治疗方案的患者提供治疗机会。 例如，许多“篮子”临床试验均根据特定的生物标记物为患有多种癌症的患者提供了特定的药物，这些试验表明，如果检测到正确的突变，则相同的药物在不同类型的药物中都具有很高的效率（Hyman等人 。，2017）。

1.2 大型人群的驱动基因鉴定

有关驱动基因的计算研究首先侧重于区分一群患者（通常是同一来源组织）中的驱动突变与乘客突变：MuSiC（Dees等人，2012年）使用了高于预期突变率的统计显着性，并与途径突变率以及与临床特征的相关性以检测驱动因素。 MutSigCV（Lawrence et al。，2013）估算每个基因的背景突变率，并鉴定出明显偏离该基因率的突变。 MEMo（Ciriello等人，2012）试图找到属于同一途径并表现出内部互斥模式的基因的小型子网络。 HotNet2（Leiserson等人，2015）结合了蛋白质-蛋白质相互作用（PPI）网络的知识，使用热扩散过程找到了经常突变的基因的小型子网络。TieDie（Paull等人，2013）还结合了PPI和mRNA表达数据，利用热扩散来发现具有高度突变和表达值的重叠子网络。 DriverNet（Bashashati et al。，2012）试图找到一组简约的突变基因，这些基因与给定PPI网络中mRNA表达失控的基因相关。 Paradigm-Shift（Ng et al。，2012）利用SNV，拷贝数变异，表达和已知途径来推断单例患者中突变基因的功能增强或丧失。 Cheng等人研究了队列中更多的驱动基因检测方法。 （2016）和Tokheim等。 （2016年）

以上方法侧重于一般驱动基因的检测，但并不旨在提供个性化的诊断或治疗手段：个别患者可能具有不同的突变驱动基因组成（补充图S1）。 此外，这些方法依赖于大型人群的统计能力，因此不可避免地低估了仅在少数患者中发生的稀有驱动基因的重要性[也称为“长尾巴现象”（Garraway and Lander，2013）。 ）]，并且仅对他们重要。 在这里，我们专注于特定于患者的驱动子优先级。

尽管已通过实验验证了许多驱动程序突变（Futreal等人，2004），但出于以下几个原因，仍需要对驱动程序进行个性化设置优先级：（i）一些患者在数十种已知驱动程序中进行了突变（补充图S1），了解哪些是病人真正的驱动基因是很重要的。 （ii）有些患者在任何已知的驱动程序中都没有突变（补充图S1），因此必须从头寻找假定的突变。 （iii）即使患者在已知驱动程序中只有很少的突变，并且假设它们都处于活动状态，我们仍然需要在内部对其进行排名，因为由于毒性和不良事件，可以同时给予个体的治疗数量非常少（Kroschinsky等人，2017; Park等人，2013）

1.3 个性化的驱动基因图谱

为了满足个性化驱动基因识别和优先级排序的需求，必须开发一种可以对单个患者的数据进行操作的方法。 为此，已经进行了一些尝试：DawnRank（Hou和Ma，2014年）使用Google的PageRank变体，根据个体对大型定向PPI网络中下游基因表达失调的影响，对个体的突变基因谱进行排名。 它通过使用扩散过程量化每个基因对差异表达基因（DEG）的影响来对基因进行排名。 SCS（Guo et al。，2018）试图在大型定向PPI网络中找到与下游DEG相关联的简约突变基因集。 这些方法将患者的假定驱动基因排序。相比之下，Hitn’DRIVE （Shrestha et al。，2014，2017）输出的一组候选驱动基因没有内部排名。它试图找到一组简短的突变，这些突变的预期路径长度短于一组DEG。 从治疗的角度来看，缺乏排名是一个缺点，特别是当预测基因的数量很大时。

1.4 这项研究

在这里，我们开发了一种用于对个体驱动基因进行排名的新算法。 该算法称为PRODIGY（驱动基因分析的个性化排名），它通过对多个已知途径解除管制的影响来评估突变。 与上述方法不同，Prodigy从同一肿瘤的许多局部视图而不是一个全局视图收集多个信号。 这些局部视图基于选定的路径，每个视图都反映了肿瘤失调状态的不同方面。 因此，一个给定的突变解释多途径失控的程度可代表该突变确实是驱动因素之一的可能性。 我们的算法假设驱动程序突变会影响受影响路径中其他基因的失调。 特别是，真正的驱动程序将与这些路径具有良好的连通性，并且我们的方法旨在使用奖赏性斯坦纳树（PCST）问题的变体来正确对此类连接进行评分。 通过汇总一个人的所有变异的许多局部视图，可以获得全局图片，并且可以对驱动基因的个性化景观进行组装和排名。

作为基本斯坦纳树问题的推广，PCST问题是NP难题，但是可以使用有效的解决方法（Bailly-Bechet等，2011； Ljubi c等，2006）。 该问题的变体以前曾用于生物信息学，尤其是E. Fraenkel的研究小组：Huang和Fraenkel（2009）将PCST模型应用于转录组，磷酸化蛋白质组学和基因筛选数据，以检测酵母中调控和信号通路的变化。 Bailly-Bechet等。 使用PCST分析与酵母信息素反应相关的转录数据。 Tuncbag等。 （2013年）使用转录组和蛋白质组学数据，使用奖赏性斯坦纳森林（PCSF）配方来发现酵母中信息素反应诱导的多种途径改变。 吉特等。 （2014年）对PCSF问题进行了概括，以找到跨多个患者共享的“共识网络”。 在这里，由于Akhmedov等人，我们使用了一种算法。 （2017a）在R（Akhmedov et al。，2017b）中实现以解决PCST问题。

在对超过2500名患者的5个TCGA癌症队列进行测试，并与经过验证的驱动基因进行比较后，PRODIGY胜过了现存的方法和基于网络中心性测度的排名。 我们的结果强调了驱动基因的多效性，并表明PRODIGY能够识别甚至非常罕见的驱动。 因此，肿瘤学家可以根据预测的驱动因素及其影响的路径，在个性化治疗的决策中使用PRODIGY。

1.4.1 注意事项

请注意，尽管我们偶尔会谈论驱动程序突变，但我们所有的分析都是在基因水平上进行的，如在SCS和DawnRank中一样，同一基因中的不同突变无法区分。 由于患者中每个突变基因的突变数通常是一个（补充表S1），因此这种区分对个性化排名的重要性不如对队列分析的重要性。 而且，正如我们将看到的，我们通常会每位患者识别并排列10个或更多的基因，因此本研究中的驱动程序概念比文献中常见的宽松一些。 但是，我们的结果表明，实际上有大量的预测驱动因素对性能有贡献。

材料与方法

给定一组突变基因的集合和个体的表达谱，我们希望对个体中的突变基因进行排名。 我们的假设是驱动基因的影响是沿着途径传播的，并由DEGs表现出来。 通过汇总来自突变基因的多种途径的证据，我们对它解释途径失控的程度进行评分。 该分数代表该基因是患者体内驱动因子的可能性的代表。在数学上，我们使用无向PCST模型对突变对解除调控途径的影响进行评分。

图1. Prodigy的方法概述。 （A）对突变的基因g在途径p上的影响进行评分：与之或距离全局网络中的突变基因g距离1的途径和基因，以及它们之间的全局边缘，构成了网络Gp， 分析（请参阅第2节）。 这是此处显示的网络。 节点奖赏（正值）反映了DEG在p中的差异表达程度，节点罚分反映了其他节点的度（由指数a校准）。 边罚反映了交互的信心。 目的是在以突变基因g为根的网络中找到最大权重的子树。 它的权重是PCST解决方案的分数。 在此示例中，用橙色虚线标记的子树是PCST解决方案，得分为9-3a。 该对的影响分数（p，g）是PCST解决方案的分数，除以属于p的DEG值的总和（此处为10.2）。 （B）在计算了所有对的影响力分数（p，g）之后，我们从得分矩阵中筛选出一些途径和基因（请参阅第2节）。 最终输出是通过剩余途径上剩余基因的总得分对它们进行排序。PRODIGY 1833从https://academic.oup.com/bioin

2.1 PCST模型

在这个问题中（图1A），目标是在加权图中找到一个子树，该子树最大化节点权重之和减去其中边的成本。 输入是无向图G=（V，E，W，P）。W：E！Rþ是边上的正权函数，P：V！ R是节点上的权重函数。 在我们的上下文中，边缘权重是反映交互可靠性的惩罚，正节点权重是给予DEG的奖励，而其他可以充当树中间节点的节点（Steiner节点）被赋予非正值作为惩罚。 给定节点g 2V; 目的是找到包含g并最大化的G的子树T：分数
换句话说，T的分数是预定义奖品的总利润减去使用中间边缘和节点的惩罚。 在我们的案例中，我们还要求T中包含特定的突变基因g。这样，T就可以解释g中的突变如何导致网络中的失调：如果g是驱动基因，那么应该进行很多失调 被树解释。

2.2 数据和相关工作

对于每个患者，Prodigy使用两种类型的基因组数据：突变的基因列表，即在编码区域中具有SNV或小插入/缺失的所有基因，以及mRNA表达谱。 来自健康组织样品的mRNA表达谱也用于差异表达分析。 Prodigy还使用两种类型的无向交互网络：（i）摘自STRING v10.5的全球PPI网络（Szklarczyk等，2015）。 在这里，我们仅使用经过实验验证的物理交互，以及来自其他精选数据库的，可信度> 0.7的交互，因此仅包括高度可靠的交互。 生成的网络具有11302个节点和273210个边。 （ii）各种途径。 在这里，我们使用了Reactome（Joshi-Tope等，2005），NCI PID（Schaefer等，2009）或KEGG（Ogata等，1999）。 有关途径数据库的信息在补充表S2中给出。

2.3 天才算法

该算法的示意图如图1所示。该算法的工作原理如下

2.3.1 预处理

给定患者的mRNA表达谱（以读数计数），使用DeSEQ2进行差异表达分析（Love等人，2014），方法是将该谱图与来自同一来源组织的健康样品的背景表达分布进行比较。 将具有在loglog-fold绝对变化（FC）> b且差异表达在FDR <c时具有统计学意义的所有基因鉴定为DEG。 使用超几何分数测试了每种途径的基因组中DEG的富集，将显着富集的途径（FDR <d）称为去调控。 我们使用了b¼2和c¼d¼0:05：

2.3.2 驾驶员—道路分数

我们使用全球互动网络G¼ðV； E;WÞ其中W是边缘置信度得分。 对于失调的路径p，我们也有其网络Gp¼ðVp;EpÞ。 这两个网络都是无向的。 突变基因g对途径p的影响得分计算如下：
1.我们构造一个新的网络Gp; g¼ðVp; g; Ep; g; Wp; g; Pp;gÞ衍生自G; Gp和g，如下所示：网络的节点是解除管制的路径g的节点，和NðVp[gÞ-它们在G中的距离为1的邻居：

它的边缘是解除管制的路径的边缘加上全局网络的所有边缘，两端均在Vp; g中：

从p开始的边成本为0.1。 对于源自全局网络G的其他边缘，其成本取决于其在该网络中的置信度得分，而较高置信度的边缘成本较低。

由于路径数据库没有提供交互作用的置信度得分，因此将路径的边沿始终分配为0.1的罚分，但是这些路径是高度策划的。 相反，全局网络边缘的置信度得分的上限为0.8，因此其在Gp; g中的成本至少为0.2。 这样做的理由是，我们希望引导算法优先考虑原始路径边缘，同时允许进行一些更改。 最后，取决于该路径的FC，每个属于该路径的DEG均具有正（奖）分，取决于其在Gp中的程度，每个其他节点v均具有负（惩罚）分； g如下

注意，Vp； g Vp中的DEG为负值。 PCST问题旨在收集尽可能多的奖励节点值，同时为中间边缘和节点付出最少的代价。 具有较高程度（“集线器”）的中间节点会打开更多的连接选项，因此根据其程度会受到更高的惩罚。 参数a控制该惩罚。 g的总体影响力得分为ing g，其中DP是患者解除管制的路径集。

路径过滤

我们计算所有突变基因和所有失控途径的驾驶员途径影响力得分。 对于最终得分，我们排除了超过一半的基因得分为阳性的途径。 这些主要是非常大的途径，在全球网络中具有很高的连通性，因此，某些基因可能偶然获得积极的影响力得分。

基因过滤

在许多途径中获得阳性评分的基因比仅在少数途径中获得阳性评分的基因更有可能代表对肿瘤的真正影响，这可能是由于网络的拓扑结构所致。 在某些患者中，当绘制所有突变基因g上InflðgÞ分数的分布图（经过过滤途径后）时，我们观察到了双峰分布（参见补充图S2）。 通常，一种分布包含从许多途径收集的高分基因，而另一种分布则包含从一些途径收集的低分基因。 我们将此分布建模为两个高斯的混合，并使用期望最大化来计算其最大似然参数（McLachlen和Peel，2000年）。 然后，我们排除了具有较高后验概率的所有基因来自均值较低的分布（补充图S2）。 如果拟合的双峰分布的可能性低于拟合的单峰正态分布，则我们不过滤任何基因。

最后排名

在过滤步骤之后，根据基因的总体影响力得分对其进行排名。

2.4 与其他个性化方法的比较