基因融合在许多癌症中都会作为致癌驱动因素发挥重要作用，所以识别基因融合具有重要意义。而今天小编要和大家分享的是一篇今年三月发表在Genome research（IF:9.043）杂志上关于检测基因融合的文章。

Accurate and efficient detection of gene fusions from RNA sequencing data

从RNA测序数据中准确高效地检测基因融合

一． 研究背景

从RNA测序数据中识别基因融合对癌症研究和肿瘤精准治疗至关重要。而基因融合对患者的治疗有直接的影响是由于许多基因融合可以用靶向药物治疗。然而，尽管如今有许多计算工具可用，但融合检测仍然具有挑战性。现有的许多方法预测精度差，计算要求高。此外，虽然多年来，各种检测基因融合的计算工具已经发展起来，但仍然没有金标准。由于文库制备和序列比对过程中引入了大量的人工制品，从短链RNA-seq中可靠地预测基因融合是十分困难的。因此今天介绍的文章开发了一种基因融合检测算法Arriba。

二． 数据及方法

1. Arriba流程：许多融合检测算法试图通过精密的比对方法来提高灵敏度。尽管这些技术提高了融合的发现，但运行时间较长。相比之下，Arriba工作流程是线性的，只有一个单一比对步骤及过滤步骤，大大提高了运行效率，其主要工作流程如图1所示。

2. 嵌合reads的提取：Arriba建立在超快STAR RNA-seq比对的基础上，当运行参数chimSegmentMin时，STAR搜索两类嵌合比对：split reads及discordant mates。嵌合比对被收集到一个名为chimeric.out.sam的单独输出文件或当指定参数chimOutType WithinBAM时的主输出文件Aligned.out.bam中。 Arriba从这些文件中提取嵌合序列，并整合它们来识别基因融合。只有当一个不连续的片段在一个合理的距离内没有比对到下游的外显子时，STAR才报告比对为嵌合。否则，它假设在比对中的间隙代表内含子跳过并创建一个间隙比对。一些众所周知的致癌融合来自focal缺失，它将上游基因的5端和下游基因的3端拉在一起。不同于创造的嵌合比对，STAR比对支持这种融合，如同融合基因是通过剪接连接起来的，这是由于STAR仅仅根据间隙的大小而不是基因注释来定义比对类型。除了提取嵌合比对外，Arriba还对跨越注释基因边界的比对进行筛选，以避免focal缺失导致的融合缺失。与许多其他融合检测pipelines不同，Arriba可以重复利用现有的STAR比对，而不是为了进行基因融合而要求reads进行专门比对。在所有提出的方法中，Arriba是唯一一个提供无缝集成到标准RNA-seq比对工作流的方法。

3. 过滤人工噪声：在收集候选比对后，Arriba应用一个过滤集来删除人工影响并提供高可信度。作者采用正向及负向过滤器，负向过滤器移除了人为因素产生的候选物，如过量不匹配的候选物等。此外，应用位置特异的黑名单来移除重复出现的人为因素及良性组织中的转录本。如果有强有力的证据表明某个候选项被错误丢弃，比如在注释的剪接位点上有断点的候选项，用户定义的已知或高度重复融合的白名单，或者通过WGS检测到的相关结构变体等，正向选择过滤器则可以补救被负向选择过滤器丢弃的候选项。正向选择滤器和根据支持读数筛选候选的统计模型是实现Arriba高灵敏度的关键特征。Arriba假设支持读取数和背景噪声水平之间的多项式关系。只有具有比估计的背景噪声水平更多的支持读数的候选会被报告(图2A)。此外，该模型还包含了几个与背景噪声水平相关的变量，包括测序深度、断点距离(图2B)、文库制备方案(图2C)和断点位置(内含子与外显子与剪接位点)。基于支持候选的读数，使用公式e-value = base_level_bg_noise * depth_penalty * distance_penalty * inv_to_dup_ratio * intron_to_exon_ratio计算背景噪声的预期水平(值)。其中每个变量有对应计算方程，这些方程中的小写分量表示动态计算变量；大写的是经验确定的常量，这些常量经过NCT及DKTK MASTER队列的RNA-seq样本训练，证明在不同的数据集之间相当稳定。

4. 基准：研究中所有的融合检测工具都使用默认参数运行，但有以下例外:PRADA的参数junL没有默认值，被设置为开发人员推荐的读数长度的80%。为了检测基因间断点融合的基准，nFusion使用参数--allow-intronic，--allow-intergenic，和--allow-non-coding执行。默认情况下，FusionCatcher使用内部已知的致癌融合列表来提高敏感性。为了更好地反映FusionCatcher对新融合发现的敏感度，作者使用skip-known-fusion参数调用了FusionCatcher，从而禁用了这个列表。另外，脚本extract_fusion_genes.py的参数allowed-labels的默认值必须被清空，参数skip-known-fusion才会生效。而clock时间、CPU时间和内存消耗是由GNU时间实用程序度量的。作者认为，如果融合partners与一组已验证的融合匹配，或者断点距离匹配的WGS样本中检测到的结构变异在100 kb以内，则预测为真阳性。由于FusionCatcher没有报告这一信息，所以并不要求基因组断点的方向与转录组断点的方向匹配。用GENCODE v19基因模型对所有工具的预测断点进行重新注释，以匹配基因名称。如果一个工具报告了涉及同一对基因的多个选择性剪接，由于产生于相同的基因组重排，只有一个转录本被计数。类似地，如果一对断点与多个基因重叠，并且以不同的基因名报告了不止一次，则只计算其中一个实例。PRADA和SOAPfuse不会根据置信对输出进行排序。因此，这些工具的预测根据支持读取的数量按降序排列。

5. MCF-7细胞系的融合预测验证：对于每种融合检测方法，如果预测没有在之前的研究中得到验证或结构变异证实，作者用Sanger测序对MCF-7细胞系的预测进行了实验验证。

6. 样本收集：样本收集过程和伦理批准可在本研究中分析的样本的相关出版物中找到。除了已发表的样本，作者还纳入了一名KRAS野生型胰腺癌患者的样本，该患者是NCT/DKTK MASTER队列中招募的。

7. 胰腺癌样本融合的识别：作者运行2.5.3版本的STAR来比对RNA-seq读数。基因融合工具的运行参数与基准测试相同。对于某些工具，需要人工干预才能在一小部分样本上成功完成。作者使用BCFtools的mpileup调用和过滤模块，以及Annovar来识别KRAS突变。此外，通过检测IGV中支持读数，人工筛选除11、12、13、61外的密码子突变。当RNA-seq数据中没有发现KRAS错义突变时，作者从各自的研究中获取KRAS的突变状态。 为了在收集的队列中确定重复，作者比较了所有样本的1000个常见SNP位点的基因型。使用基于欧式距离的层次聚类将样本聚在一起的被认为是重复的。作者从特征的组合推断基因融合是否应该被认为是一个驱动因素，包括表达水平，Arriba的置信评分，阅读框的保存，致癌活性的基本域的保留，以及这些基因之前是否曾被描述参与胰腺癌或其他实体的致癌融合。使用R包PBase，从蛋白坐标到基因组坐标绘制Pfam蛋白结构域。跨膜结构域的基因组坐标来自UniPro。在IGV中可视化潜在的融合候选体，以识别潜在的比对噪音。PACACA队列患者PCSI_0326进行了TRIM24-BRAF融合。Arriba只报告了一份融合记录。进一步观察BRAF的软剪接读数发现，一些读数将TRIM24的第9外显子连接到BRAF的第8外显子，这是IGV内置的BLAT功能所揭示的。据推测，STAR未能对这些读数进行比对是因为它们包含了BRAF内含子中的20个碱基，这些碱基不能单独定位于人类基因组。接着在WebGestalt 的帮助下，对通路中过多代表的基因进行了分析。作者将所有人类蛋白编码基因作为KEGG数据库的背景和通路作为基因集进行过度表达测试。

8. 实验验证:慢病毒转导、定量RT-PCR、克隆形成实验、蛋白免疫印记。

9. 药物敏感性：对于剂量反应曲线，将 2000 个 MCF10A 细胞接种在 96 孔板中，在不含EGF的培养基中用特定浓度的 MAPK (ERK) 抑制剂 FR180204 或 RAF1 抑制剂sorafenib进行处理，48小时后用CellTiter 96 Aqueous One Solution Cell Proliferation Assay (Promega) MTS法检测细胞活力。

三． 研究的主要内容及结果

1. 精度基准

在文章的第一部分作者为了展示Arriba方法在不同类型输入数据上的鲁棒性，在四种类型的基准数据集上评估了其精确性。作者首先使用的数据集是模拟的150个融合转录本，并将它们整合到良性组织的RNA-seq样本中，作为背景表达。为了测量作为融合转录本表达水平函数的方法的灵敏度，模拟了9个不同的表达水平，从5倍到200倍不等。作者使用的第二个数据集是在黑色素瘤细胞系COLO-829的RNA文库中加入合成的RNA分子。这种合成的RNA分子模仿了九种不同癌症类型中发现的致癌融合的转录序列。第三个数据集是来自4个细胞系的8个样本。最后一个数据集是来自ICGC早发前列腺癌队列和TCGA弥漫巨B细胞淋巴瘤队列的患者数据。结果图3A中ROC等曲线表明Arriba显示了良好的准确性。在所有四种基准数据集中，Arriba的敏感性最高(图4)。在特异性方面，Arriba也可以与最先进的方法相媲美。Arriba为其预测指定了三个置信等级:低、中、高。使用者可以通过选择高于某个置信级别的事件，在敏感性和特异性之间选择他们喜欢的平衡。

2. 运行时间和内存消耗

在文章的第二部分，作者在运行时间及内存消耗两个方面评估Arriba方法的性能，测量了所有工具的运行时间。结果发现Arriba在运行时间和CPU时间方面都是最快的(图3B)。此外，平均而言，基于Arriba的工作流消耗了38 GB的内存，比最有效的内存工具SOAPfuse多5.8倍(图3C)。

3. 在实践中使用Arriba

在文章的第三部分，作者对Arriba的应用作了进一步介绍。可以了解到为了加快基因融合相关研究任务，Arriba提供了一些其他的功能，而不仅仅是预测融合。它提供了连接位点侧面的转录序列，有助于Sanger测序验证引物的设计。它还计算了嵌合转录本产生的肽序列，这可以作为预测融合衍生的基础。此外，Arriba提供了可视化工具，以促进基因融合的解释。而且由于STAR以SAM格式存储嵌合比对，这些比对可以加载到基因组浏览器中。Arriba提供了一个蛋白质结构域的特征轨迹，它可以与序列一起装载到IGV中，以评估融合的功能影响。除了RNA-seq数据，临床研究项目偶尔会为每个患者生成WGS数据。Arriba可以通过提供从WGS得到的结构变量列表来进一步提高预测精度。

4. 胰腺癌致癌基因融合的识别

在这一部分，由于在KRAS野生型胰腺癌中发现涉及NRG1的复发融合等融合病例的研究，作者在胰腺癌中对基因融合进行了识别。作者收集了803名捐赠者的RNA-seq样本，包括了18项已发表的胰腺癌研究。其中327个样本，有匹配的WGS数据。当Arriba从转录组数据中预测基因融合时，作者在WGS数据中检查了一个相关的结构变体，以确认预测的有效性。结果在RNA-seq数据中检测到30个潜在的驱动融合(图5)，所有这些经结构变异在WGS数据可用。同时作者也对影响的功能进行了分析，其中在纳入分析的105个样本中，没有检测到KRAS改变。在这些样本中，致癌融合显著富集，在KRAS突变肿瘤中仅发现4个融合。在79例肿瘤中，既没有检测到KRAS突变，也没有检测到驱动融合。为了排除由于Arriba的缺点而忽视融合的可能性，作者也在队列上运行了其他融合检测工具，但没有一个报告表现出驱动融合超出了Arriba的设置。事实上，Arriba显示出了最高的灵敏度，比其他方法多检测了3到11个驱动融合，从而证实了基准测试结果。其中一些已确认的胰腺基因融合已经在其他癌症类型的背景下被报道。

5. 两个新的融合基因的功能验证

在文章的最后一部分作者试图通过实验验证预测的基因融合作为致癌驱动因素。作者选择了RASGRP1-ATP1A1和RRBP1-RAF1(图6A、B)，因为RASGRP1之前没有参与致癌融合，而RRBP1是RAF1的新伴侣，与接近全长的RAF1融合，而不是更常见的第8外显子融合。这些融合通过慢病毒转导进入H6c7细胞和TP53缺陷的MCF10A细胞。与空载体对照相比，两种融合显著增强了EGF独立的集落形成(图6C)。此外，当EGF退出时，融合蛋白增加了MAP2K1/2 (MEK1/2)和MAPK1/3 (ERK2/1)的磷酸化，表明MAPK通路的结构性激活(图6D)。总之，这些实验证实了RASGRP1-ATP1A1和RRBP1-RAF1的致癌活性，进一步支持了Arriba预测致癌融合的观点。为了测试融合是否能够用于治疗，作者用两种针对MAPK信号的化合物处理细胞:RAF1抑制剂sorafenib和MAPK (ERK)抑制剂FR180204。结果发现虽然细胞培养对所有化合物都有反应，但融合阳性细胞并不比空载体对照更敏感。

到这里，这篇文章的主要内容就介绍完了。文章开发了一种从RNA-seq中预测基因融合的方法，并对方法进行了比较测试应用验证等分析，对基因融合感兴趣或者做方法的小伙伴不要错过呀。