Basic Information 英文标题:Nucleotide Transformer: building and evaluating robust foundation models for human genomics 中文标题:核苷酸 Transformer:构建和评估用于人类基因组学的稳健基础模型 文章作者:Hugo Dalla-Torre | Thomas Pierrot 文章链接:https://www.nature.com/articles/s41592-024-02523-z Abstract Para_01 从DNA序列预测分子表型仍然是基因组学中的一个长期挑战,通常由有限的注释数据和无法在任务之间转移学习所驱动。 在这里,我们对预训练在DNA序列上的基础模型进行了广泛的研究,这些模型被称为核苷酸转换器,参数范围从5000万到25亿,并整合了来自3,202个人类基因组和850个不同物种基因组的信息。 这些转换器模型生成了特定上下文的核苷酸序列表示,即使在数据量较少的情况下也能进行准确预测。 我们展示了开发的模型可以以低成本进行微调,以解决各种基因组学应用问题。 尽管没有监督,这些模型学会了关注关键基因组元素,并可用于改进遗传变异的优先排序。 基因组学中基础模型的训练和应用提供了一种广泛适用的方法,可以从DNA序列中准确预测分子表型。 Main Para_01 人工智能(AI)中的基础模型以其大规模特性为特征,包含在大量数据集上训练的数百万个参数。这些模型可以适应各种后续预测任务,并深刻改变了AI领域。 自然语言处理(NLP)中的著名例子包括所谓的语言模型(LMs)BERT和GPT。近年来,由于它们能够在未标记的数据上进行训练,创建能够解决下游任务的通用表示,语言模型变得非常流行。 它们通过解决数十亿个完形填空测试来实现对语言的全面理解,在这些测试中,它们预测给定句子中空白处应填入的正确单词。这种方法被称为掩码语言建模。 早期将这一目标应用于生物学的基础模型实例包括在蛋白质序列上训练语言模型,其中它们的任务是预测大型蛋白质序列数据集中被掩码的氨基酸。 当这些蛋白质语言模型使用迁移学习应用于下游任务时,展示了与先前方法竞争甚至超越的能力,例如在数据稀缺的情况下预测蛋白质结构和功能等任务。 Para_02 除了蛋白质序列之外,DNA序列中编码的依赖模式在理解基因组过程方面起着基础性作用,从表征调控区域到评估单个变异在其单倍型背景下的影响。 在此背景下,专门的深度学习(DL)模型已被训练用于揭示DNA中的有意义模式。 例如,深度学习模型已被用于根据DNA序列预测基因表达,最近的进步结合了卷积神经网络和变压器架构,能够对位于上游多达100千碱基(kb)的调控元件进行编码。 现代基因组学研究产生的大量数据既是一个机会也是一个挑战。 一方面,跨物种和种群的自然变异复杂模式是现成可用的;另一方面,处理大规模数据的强大深度学习方法对于从未标记数据集中准确提取信号是必要的。 基于核苷酸序列训练的大规模基础模型似乎是一种值得探索的方法,以应对这一挑战。 Para_03 在这里,我们构建了稳健的基础模型来编码基因组序列,命名为核苷酸变换器(NT),并进行了系统的研究和基准测试以评估其性能。 我们通过构建四个不同大小的语言模型开始了我们的研究,这些模型的参数范围从5亿到25亿不等。
这些模型在三个不同的数据集上进行了预训练,包括人类参考基因组、3,202个不同人类基因组的集合以及来自各种物种的850个基因组。 训练后,我们以两种方式利用每个模型的表示(嵌入)。 为了评估NT在适应各种任务时性能的稳定性,我们在18个基因组精选预测任务上训练了每个模型,并使用系统的十倍交叉验证程序将它们与三种替代的DNA基础模型以及一种最先进的非基础模型进行比较。 此外,为了扩大我们的评估范围,我们将表现最佳的模型与三种经过优化的特定任务最先进监督基线模型进行了比较。 为了破译预训练期间学到的序列特征,我们探索了模型的关注图和困惑度,并对其嵌入进行了数据降维。 此外,我们通过零样本得分评估了嵌入的能力,以建模对人类功能重要遗传变异的影响。 在扩展初始实验结果的基础上,我们开发了第二组四个语言模型,参数大小从5亿减少到5千万,以研究此类模型的缩放规律。 我们成功构建了一个模型,该模型仅用十分之一的参数数量和加倍的感知场大小就达到了我们之前最佳模型的性能。 Results The Nucleotide Transformer model accurately predicts genomics tasks 核苷酸转换模型准确预测基因组学任务
Para_01 我们开发了一组基于Transformer的DNA语言模型NT,这些模型通过从6千碱基对未注释的基因组数据中学习得到了通用的核苷酸序列表示(图1a和方法)。 受自然语言处理领域趋势的启发,在该领域中,更大的训练数据集和模型尺寸已显示出性能的提升,我们构建了具有不同参数大小和数据集的Transformer模型: (1) 一个使用从人类参考基因组提取的序列进行训练的5亿参数模型(‘Human ref 500M’); (2) 一个使用来自3,202个遗传多样的人类基因组的序列进行训练的5亿参数模型(‘1000G 500M’),以及 (3) 一个同样使用这些序列进行训练的25亿参数模型(‘1000G 2.5B’); (4) 另外还有一个25亿参数模型,涵盖了850种来自不同门的物种(‘Multispecies 2.5B’),包括11种模式生物(图1c和补充表1-4)。 Fig. 1: Nucleotide Transformer: effective methodology to pre-train, fine-tune, analyze and compare foundational models for genomics.
- 图片说明 - a,b,NT训练概述(a)及其通过微调在下游基因组预测任务中的应用(b)。通过探测进行的下游任务预测类似,但没有NT中的重新缩放权重。c,NT模型与其他基础基因组学模型在感知场大小、参数数量和我们基准中包含的18个精选下游任务性能方面的比较。d,为下游任务考虑的基因组特征的图形表示(改编自其他地方48)。 - ,
Para_02 为了评估这些模型在预测多种分子表型方面的有效性,我们从公开资源中整理了18个基因组数据集,包括剪接位点预测任务(GENCODE29)、启动子任务(真核启动子数据库30)和组蛋白修饰及增强子任务(ENCODE31,32,33),每个任务设计为合理大小以实现快速且严格的交叉验证程序(图1d,补充表5和方法)。 虽然有更大规模的监督模型数据集可用,但这18个基因组数据集的编纂提供了一个多样且稳健的选择,用于统计严格地审查模型在不同任务中的适应性,并与其他DNA自监督基础模型进行比较。 这些数据集被处理成标准化格式,以促进实验并确保大规模语言模型性能评估的可重复性(方法)。 我们通过两种不同的技术评估了我们的变压器模型:探测和微调(图1b)。 探测指的是使用学习到的DNA序列LM嵌入作为输入特征,以预测基因组标签的简单模型。 具体来说,我们使用逻辑回归或最多两层隐藏层的小型多层感知器(MLP)探测了LM的十个任意选择的层。 在微调的情况下,LM头被分类或回归头取代,并采用参数高效的重新训练技术(方法)。 为了确保不同模型之间的公平和准确比较,我们实施了十折交叉验证策略。 Para_03 为了将我们的预训练基础模型方案与该领域中的标准监督方法进行比较,我们在每个18个任务上从头开始训练了BPNet卷积架构的不同变体(方法)。 BPNet架构在基因组学中被广泛使用,并代表了一种非常强大的默认架构,通过监督学习从头开始对小型数据集进行建模。 我们观察到原始BPNet模型在各项任务上的表现非常强劲(平均马修斯相关系数(MCC)为0.665),通过将其大小增加到2800万个参数,我们可以进一步提高其性能(平均MCC为0.683),这证实了直接监督的卷积架构在基因组学任务中表现出色(图2a,b)。 接下来,我们评估了NT模型的探测和微调与这些监督基线模型在我们的基准数据集上的比较情况。 如果结果的两个标准偏差要么重叠要么优于报告的基线值,我们将这些模型视为等同于或优于其他模型。 Fig. 2: The Nucleotide Transformer model accurately predicts diverse genomics tasks after fine-tuning.
- 图片说明 - a,基于MCC对微调后的NT模型以及HyenaDNA、DNABERT和Enformer预训练模型在下游任务中的性能结果进行了比较。我们还从头训练了BPNet模型以进行对比(原始模型有121,000个参数;大型模型有2800万个参数)。数据以十折交叉验证程序的平均MCC ± 2 × 标准差表示(每个点n = 10)。 - b,所有语言模型在微调后按类别划分的下游任务中MCC性能的归一化均值。 - c,Multispecies 2.5B模型在不同人类细胞和组织中对DNase I超敏感位点(DHSs)、组蛋白标记(HMs)和转录因子位点预测的性能与基线DeepSEA模型的比较。每个点代表不同基因组图谱的受试者工作特征(ROC)曲线下面积(AUC)。每个模型的平均AUC已标注。 - d,Multispecies 2.5B模型在预测人类基因组剪接位点方面的性能与SpliceAI和其他剪接模型的比较。 - e,Multispecies 2.5B模型在预测果蝇S2细胞发育和管家增强子活性方面的性能与基线DeepSTARR模型的比较。
Para_04
使用这一标准,NT 模型在 5 项任务中与基线 BPNet 模型相当,并且通过单独探针在 18 项任务中的 8 项中超过了它们(补充图 1 和补充表 6),并且显著优于从原始标记进行探针的结果。 与最近的工作一致,我们观察到最佳性能既依赖于模型也依赖于层(补充表 8)。 我们还注意到,最高模型性能从未通过使用来自最终层的嵌入来实现,如早期工作所示。 例如,在增强子类型预测任务中,我们观察到最高和最低性能层之间的相对差异高达 38%,表明各层之间学习表示存在显著变化(补充图 3)。 与我们的探针策略相比,我们的微调模型要么匹配(n = 6)要么超越(n = 12)了 18 个基线模型中的大多数(图 2a、b 和补充表 7 和 9)。 值得注意的是,微调后的 NT 模型表现优于探针模型,并且更大更多样化的模型始终优于其较小的同类模型。 这些结果支持了将 NT 基础模型微调至特定任务以实现优越性能的必要性。 我们的结果还表明,在多样化数据集上训练,由 Multispecies 2.5B 模型代表,优于或匹配 1000G 2.5B 模型在几个基于人类实验的任务上的表现(图 2a、b)。 这表明增加序列多样性而不是仅仅增加模型大小的策略可能会带来更好的预测性能,特别是在计算资源有限的情况下。 Para_05 微调在之前的研究中并未得到广泛探索,可能是因为其计算需求较高。我们通过采用一种最近的参数高效微调技术克服了这一限制,该技术仅需要总模型参数的0.1%(图1b和方法)。这种方法使得在单个GPU上更快地进行微调,相比所有微调参数减少了1000倍的存储需求,同时仍然提供了可比的性能。 在实践中,我们观察到严格的探测比微调更慢且计算强度更高,尽管使用简单的下游模型在嵌入上的表观简单性。这种差异源于诸如层选择、下游模型选择和超参数等因素对性能的显著影响。 此外,微调表现出较小的性能方差,增强了方法的稳健性。总体而言,这种通用方法具有多功能性和适应性,可以适用于各种任务,而无需调整模型架构或超参数。 这与监督模型形成了对比,后者通常具有不同的架构,并且需要从头开始为每个任务进行训练。 Para_06 最后,我们旨在评估大型语言DNA模型在使用大量数据集和优化架构进行监督训练的稳健基线中竞争的潜力。为此,我们将Multispecies 2.5B模型应用于三个额外的基因组预测任务,这些任务包括对来自多种人类细胞和组织的919个染色质谱型进行分类、预测整个人类基因组中的典型剪接受体和供体位点以及从果蝇S2细胞中预测发育和管家增强子活性(方法)。值得注意的是,尽管没有对其原始微调架构进行任何额外更改或优化,Multispecies 2.5B模型的表现水平与专用深度学习模型非常接近。 例如,在染色质特征谱型分类的情况下,我们获得的曲线下面积(AUC)值平均仅比DeepSEA低约1%左右(图2c)。关于预测每个前mRNA转录本位置是否为剪接受体、剪接供体或两者都不是的问题,我们调整了NT模型以提供核苷酸级别的剪接位点预测,并实现了95%的top-k准确率和0.98的精确-召回AUC(图2d)。 值得注意的是,我们的2.5B 6-kb上下文模型在性能上与最先进的SpliceAI-10k36相当,后者是在15-kb输入序列上进行训练的,此外还优于其他剪接基线;并且在测试6-kb输入序列时超过了SpliceAI。最后,在管家和发育增强子预测的情况下,我们的模型分别略微超越了1%并获得了较低4%的相关值(图2e),相比于DeepSTARR12而言。 在这三个不同的任务中,我们还比较了参数高效的微调和全模型微调(训练整个模型的所有参数以优化其在特定任务或数据集上的表现)。值得注意的是,我们在染色质和剪接预测方面未观察到显著改进,仅在增强子活性预测方面有适度的3%提升(补充图2),这支持了我们高效微调方法的使用。总而言之,我们的广泛基准和结果证明了NT作为解决许多不同基因组任务的通用方法具有灵活性和高性能。 Benchmark of genomics foundational models 基因组基础模型的基准测试
Para_01 我们比较了NT模型与其他基因组基础模型:DNABERT-2(参考文献23)、HyenaDNA(1-kb和32-kb上下文长度;参考文献25)以及Enformer(将其用作预训练模型的替代架构,参考文献19;图2a,b和方法)。 我们将DNABERT-1从此次比较中排除,因为它只能处理最大输入长度为512 bp的数据,因此无法用于大多数任务(参考文献20)。 为了确保公平比较,所有模型在18个下游任务中都使用相同的协议进行微调和评估(方法)。 与DNABERT-2、HyenaDNA-32-kb和Enformer相比,我们的多物种2.5B模型在各项任务中总体表现最佳(图2a,b和表9)。 尽管如此,Enformer在增强子预测和一些染色质任务上表现出色,证明它是一个强大的DNA基础模型。 我们的模型在所有启动子和剪接任务上的表现优于其他每个模型。 值得注意的是,尽管HyenaDNA是在人类参考基因组上预训练的,我们的多物种2.5B模型在所有18项任务中与其持平(n = 7)或超越了它(n = 11),突显了在多样化基因组序列上预训练的优势。 我们建立了一个交互式排行榜,包含所有模型在每项任务上的结果,以便于比较(https://huggingface.co/spaces/InstaDeepAI/nucleotide_transformer_benchmark)。 这代表了广泛的基因组基础模型基准测试,并应作为进一步开发基因组学语言模型的参考(图1c)。 Detection of known genomic elements in an unsupervised manner 以无监督的方式检测已知的基因组元素
Para_01 为了深入了解可解释性并理解NT模型在进行预测时使用的序列元素类型,我们探索了其架构的不同方面。 首先,我们评估了嵌入能够在多大程度上捕获与五个基因组元素相关的序列信息(补充表7和10)。 我们观察到,NT模型在没有任何监督的情况下学会了区分被唯一注释为基因间区、内含子区、编码区和非翻译区(UTRs)的基因组序列,尽管在不同层面上的熟练程度有所不同(图3a,补充图7和方法)。 特别是,500-M大小的模型和那些在较少多样性的序列上训练的模型,在基因组区域之间的分离度较低,这强化了最大模型在自我监督训练期间捕捉相关基因组模式的能力。 在多物种2.5B模型的情况下,在第1层观察到基因间区和基因区之间最强的分离,其次是第5层上的5′ UTR区,并且在第21层上大多数区域之间的分离(图3a)。 3′ UTR区与其他元素之间的有限分离表明该模型尚未完全学会区分这种类型的元素,或者如先前所建议的那样,这些区域中的许多可能被错误标注。 与这些观察结果一致的是,我们的有监督探测策略展示了对这些元素的高分类性能,准确值超过0.78,特别是在较深层次(图3b)。 这表明NT模型已经学会了以无监督的方式在其嵌入中检测已知的基因组元素,这可以用于高效的下游基因组任务预测。 Fig. 3: The Nucleotide Transformer models acquired knowledge about genomic elements.
- 图片说明 - a,基于Multispecies 2.5B模型,第1层、第5层和第21层中五个基因组元素的嵌入的t-SNE投影。 - b,基于探测分类五个基因组元素在各层中的准确度估计。 - c,描述评估给定基因组元素上注意力水平的示意图。 - d,在Multispecies 2.5B模型中,对5′ UTR、外显子、增强子和启动子区域计算每个头和每层的注意力百分比。每个瓷砖图右侧的条形图显示了给定层中所有头的最大注意力百分比。
Para_02 接下来,我们通过注意力分析模型,以了解哪些序列区域被注意力层捕获和利用。我们计算了包含九种不同类型的基因组元素(与基因结构和调控特征相关)的序列在每个模型头和层上的注意力百分比(图3c)。 从形式上讲,当某个注意力头的注意力百分比显著超过该元素在预训练数据集中自然出现的频率时,就认为它能够识别特定元素(方法)。 例如,50%的百分比意味着,在整个人类基因组中,该特定头部的注意力平均有50%是针对感兴趣的元素类型。通过这种方法,我们在大约10,000个不同的6-kb序列中对每种类型的元素进行应用,其中元素可以位于各种位置,并占序列的2-11%(补充表10),我们发现注意力在其多样化的头部和层中明显集中在特定类型的基因组元素上(图3d和补充图8-16)。 跨层的重要注意力头部数量在不同模型之间差异显著,对于内含子(117个中的640个头部),外显子(n = 72)和转录因子(TF)结合位点(n = 74),多物种2.5B模型观察到最多的重要注意力头部(补充图8和9以及补充表12),尽管含有外显子和TF基序的序列比例相对较小。 关于增强子,最大注意力百分比最高的为最大的模型,例如1000G 2.5B模型几乎达到了100%的注意力(补充图15)。 对于其他基因组元素,如3′ UTR、启动子和TF结合位点,也观察到了类似的模式,其中1000G 2.5B模型展示了高度专业化的头部,具有高注意力,特别是在前几层(补充图8-16)。 Para_03 为了在更高分辨率下深入了解预训练的NT多物种2.5B模型(重点关注更局部的序列特征),我们检查了不同类型的基因组元素中的令牌概率,作为序列约束和模型学习到的重要性指标。具体来说,我们计算了第22号染色体中每个6千碱基窗口内的六聚体令牌概率(基于每次屏蔽一个令牌)。我们的发现表明,除了重复元件(这些元件如预期那样被模型很好地重建)之外,预训练模型还学习了多种基因结构和调控元件。 这些包括受体和供体剪接位点、polyA信号、CTCF结合位点和其他基因组元件(补充图17a-d)。此外,我们将我们的令牌预测与MST1R基因第11个外显子的实验饱和诱变剪接分析进行了比较(数据来自Braun等39)。该分析揭示了实验突变效应与Multispecies 2.5B模型的令牌预测之间存在显著相关性(皮尔逊相关系数(PCC) = 0.44;补充图17e)。 该模型不仅捕捉到了不同剪接连接处的约束条件,还识别出了第二个内含子中间的一个区域,这个区域对这个外显子的剪接至关重要。这些结果为NT模型在无监督预训练期间获得的生物学知识提供了强有力的验证。 Para_04 最后,对在DeepSTARR增强子活性数据上进行了全面微调的Multispecies 2.5B模型进行了检查,以确定其是否专门学习了转录因子(TF)基序及其相对重要性,特别是针对增强子活性。我们使用了一个包含数百个增强子序列中五种不同TF基序类型的数百个单独实例的实验突变数据集,并评估了该模型在预测这些突变效应方面的准确性。 与最先进的增强子活性DeepSTARR模型相比,我们的模型在四种TF基序上取得了相似的表现,并且在Dref基序上表现出色(补充图18)。 总体而言,这些结果展示了NT模型如何获得了恢复基因结构和基因组序列功能特性的能力,并将其直接整合到其注意力机制中。 这种编码的信息应在评估遗传变异的重要性时提供帮助。 The pre-trained embeddings predict the impact of mutations 预训练的嵌入预测了突变的影响
Para_01 此外,我们评估了NT模型评估各种遗传变异严重程度并优先考虑具有功能意义的变异的能力。 我们首先研究了零样本得分的使用,这些得分用于预测模型在训练期间未见过的类别。 具体来说,我们计算了使用嵌入空间中不同向量距离方面的零样本得分,以及从损失函数派生的得分,并比较了它们在十种不同类型、严重程度各异的遗传变异中的分布(图4a和方法)。 令人鼓舞的是,这些零样本得分中有几个在模型中显示出与严重程度的适度相关性(补充图19)。 这说明了仅通过无监督训练就捕获了与遗传突变潜在严重性相关的信息,并强调了评估不同评分方法的效用。 分数之间的相关性的高变异性也表明,嵌入空间的不同方面可能更有效地捕捉与严重性相关的信息。 在这些分数中,余弦相似度在模型中显示出最高的与严重性的相关性,r2值范围从−0.35到−0.3(P < 6.55 × 10−186;补充图19)。 在模型中,最低的余弦相似度分数被分配给影响蛋白质功能的遗传变异,如终止获得变异,以及同义和错义变异(图4b)。 相反,我们注意到较高的分数被分配给可能功能重要性较低的变异,如基因间变异,突出了其潜在用途来捕捉遗传变异严重性的影响。 Fig. 4: Prioritizing functional genetic variants.
- 图片说明 - a,零样本预测在NT应用中的概述。 - b,基于余弦相似度度量的模型中,变体后果术语在十分位数中的比例。根据Ensembl估计,后果术语按严重程度(从较轻到较重)排序。 - c,基于不同距离度量的零样本预测在优先考虑功能变异方面的比较。 - d,基于GRASP eQTLs和meQTLs、ClinVar和HGMD注释突变的微调模型和可用方法在优先考虑功能变异方面的比较。模型性能通过ROC AUC进行测量。展示了三个表现最佳模型的AUC。
Para_02 接下来,我们还探索了零样本得分在优先考虑功能变异以及具有致病效应的变异方面的潜力。 具体来说,我们评估了模型对影响基因表达调控(表达数量性状位点(eQTLs))、与DNA甲基化变异相关的遗传变异(甲基化数量性状位点(meQTLs))、在ClinVar数据库中标注为致病的遗传变异和在人类基因突变数据库(HGMD)中报告的遗传变异进行分类的能力。 零样本得分展示了高分类性能,在四个任务中的最高AUC值范围从0.7到0.8(图4c)。 对于ClinVar变异获得的最高性能(Multispecies 2.5B模型的AUC为0.80),表明至少对于高度致病的变异,零样本得分可能可以直接应用。 最后,为了更正式地评估这些模型的有效性,我们还基于微调模型进行了预测,并将其性能与其他几种方法进行了比较。 这些方法包括测量基因组保守水平的方法,以及从基于功能特征训练的模型中获得的分数。
值得注意的是,微调模型的表现要么略微优于其他模型,要么与其相当(图4d和补充图20)。 用于优先考虑分子表型(eQTLs和meQTLs)的最佳模型是那些基于人类序列训练的模型,而用于优先考虑致病变异的最佳模型则是基于多物种序列的模型。 鉴于最严重的致病变异往往由于氨基酸变化影响基因功能,可能是多物种模型利用跨物种的序列变异来学习不同位点的保守程度。 我们的结果还表明,通过更好地学习来自增加的人类遗传变异性的序列变异,可以实现对非编码变异如eQTLs和meQTLs的更高预测能力。 此外,与零样本得分相比,点积对eQTLs和meQTLs分别产生了0.73和0.71的AUC值,略高于或匹配由微调模型获得的结果。 鉴于大多数这些遗传变异倾向于位于调控区域,可能是模型在没有任何监督的情况下学会了区分与基因表达和甲基化变异相关的相关调控基因组特征。 这与观察到的各层和头的注意力水平一致,特别是对于相关的调控序列如增强子和启动子,它们已被证明富含meQTLs和eQTLs。 总体而言,这些结果说明了基于DNA的转换模型如何有助于揭示并理解与分子表型和疾病相关的变异的潜在生物学意义。 Model optimization for cost-effective predictions in genomics 基因组学中用于成本效益预测的模型优化
Para_01 最后,我们探讨了通过引入当代架构进步和延长训练时间来优化我们表现最佳的模型的潜力。我们开发了四个新的NT模型(NT-v2),参数数量从5000万到5亿不等,并引入了一系列架构改进(补充表1和方法)。这些改进包括引入旋转嵌入、实施swiGLU激活函数以及消除MLP偏置和dropout机制,符合最新的研究结果。 此外,我们将上下文长度扩展到12千碱基以适应更长的序列并捕捉更远距离的基因组相互作用。我们将2.5亿和5亿参数模型的训练时间延长到涵盖1万亿个标记,与文献中的最新建议一致(图5a)。 在相同的多物种数据集上进行预训练后,所有四个NT-v2模型都进行了微调并在相同的18个下游任务中进行了评估,其结果与最初的四个NT模型的结果进行了比较(图5b和补充表7和9)。 Fig. 5: Efficient model architecture allows matching performance while strongly reducing the number of model parameters.
- 图片说明 - a,NTr-v2模型在训练过程中将损失值的变化表示为迄今为止看到的标记数量的函数。 - b,归一化后的NT-v2模型MCC性能均值表示为预训练期间看到的标记数量的函数。 - c,所有NT(灰色)和-v2(蓝色阴影)NT模型在微调后在下游任务(按类别划分)上的归一化MCC性能均值。黑色虚线表示NT 2.5B多物种模型的表现。 - d,NT在核苷酸水平剪接位点预测任务上微调的概述。预训练权重和从头开始训练的权重被突出显示。 - e,Multispecies v2 500M模型在预测人类基因组剪接位点方面的表现,与其2.5B对应模型和SpliceAI进行比较。
Para_02 我们观察到,具有50M参数的NT-v2模型在性能上与我们的两个拥有5亿参数的NT模型以及在1000基因组数据集上训练的25亿参数模型相当。这表明,优越的预训练数据集与训练技术和架构的进步相结合,可以在显著减少模型参数的同时提高性能(图5c)。 事实上,拥有2.5亿和5亿参数的NT-v2模型在保持明显更少的参数数量并使感知范围加倍的情况下,实现了优于25亿参数的多物种模型的性能。 特别值得注意的是,拥有2.5亿参数的NT-v2模型在我们的基准测试中取得了最佳性能(平均MCC为0.769),同时其规模仅为25亿参数模型的十分之一(图5c)。 Para_03 为了进一步了解需要更长时间的预训练,我们对NT-v2模型在预训练期间看到的令牌数量方面的性能进行了系统评估(图5b)。 这表明,只有在训练了9000亿个令牌之后,2.5亿参数的模型才比5亿参数的模型有小幅改进。 总之,配备12千字节上下文长度的NT-v2模型由于其紧凑的尺寸,适合部署在经济高效的加速器上。 因此,它们为寻求在其下游应用中利用前沿基础模型的用户提供了一种经济可行且实用的替代方案。 Para_04 作为一个相关的用例,我们使用一个已建立的高性能模型,通过在SpliceAI剪接任务上评估500M模型来评估NT-v2模型更长上下文长度的优势。 我们评估了500M模型,因为该模型在与剪接相关的任务中表现出最高的性能(补充表7和9)。 我们调整了分类头以预测核苷酸水平的概率,即成为剪接供体、受体或无(图5d和方法部分)。 与6-kb NT模型相比,我们的12-kb NT-v2 500M模型将性能提高了1%,达到了96%的top-k准确率和0.98的精确召回AUC(图5e)。 这一性能超过了最先进的SpliceAI-10k36,后者是在15-kb输入序列上训练的。 值得注意的是,我们并未专门针对剪接预测任务优化我们的模型架构;相反,我们采用了与其他下游任务类似的微调方法,并对分类头进行了调整以生成核苷酸级别的预测。 进一步针对特定任务如剪接进行的架构改进可能会提高性能。 总之,这些结果肯定了NT v1和v2模型在广泛的基因组学任务中的实用性和有效性,只需少量修改和计算能力即可实现高精度。 Discussion Para_01 本研究旨在探讨用于预训练同等规模的变换器模型的不同数据集对DNA序列的影响。我们的结果基于不同的基因组预测任务,表明无论是同物种(当在单一物种的多个基因组上进行训练时)还是跨物种(在不同物种的基因组上)的变异性都显著影响了任务的准确性(图1c和2a)。在大多数考虑的人类预测任务中,使用不同物种基因组训练的模型优于仅使用人类序列训练的模型。这表明,使用多种物种训练的变换器模型已经学会了捕捉可能具有跨物种功能重要性的基因组特征,从而在各种基于人类的预测任务中实现了更好的泛化。
基于这一发现,我们预计未来的研究可能会受益于利用跨物种的遗传变异性,包括确定采样这种变异性的最佳方法。另一个有趣的方向是探索编码同物种变异性的不同方法。我们混合所有个体基因组序列的方法只带来了有限的改进,因此表明当大多数基因组共享时,利用不同个体的基因组可能并不那么简单。 Para_02 本研究中训练的变压器模型参数范围从5000万到25亿,这比DNABERT-2(参考文献23)大20倍,比Enformer19骨干模型大十倍。 正如在自然语言处理研究中先前所证明的那样(参考文献27),我们在18个基因组预测任务上的结果证实,增加模型大小可以持续提高性能。 为了训练具有最大参数规模的模型,我们利用了16个计算节点上的总共128个GPU进行了28天的训练。 我们在工程高效的训练例程方面进行了大量投资,以充分利用基础设施,强调了专用基础设施和专用软件解决方案的重要性。 然而,一旦训练完成,这些模型可以以相对较低的成本用于推理,并且我们提供了笔记本以将这些模型应用于任何感兴趣的下游任务,从而促进进一步的研究。 Para_03 之前的工作基于主要使用蛋白质序列的生物数据训练的语言模型,通过探测最后一层变压器来评估下游性能。这一选择可能是由于其易于使用、相对较好的性能和较低的计算复杂性。 在本研究中,我们的目标是通过对不同变压器层、下游模型和超参数扫描进行计算密集型和彻底的探测来评估下游准确性。 我们观察到,最佳探测性能是由中间变压器层实现的(补充图 1),这与计算生物学中的最新工作和自然语言处理中的常见实践一致。 通过仅探测,Multispecies 2.5B 模型在 18 项任务中的 8 项上超过了 BPNet 基线模型。 此外,我们探索了一种最近的下游微调技术,该技术在变压器中引入了少量可训练权重。 这种方法提供了相对较快且资源高效的微调程序,与全模型微调相比差异较小(IA3)。 值得注意的是,这种微调方法只需要总参数数量的 0.1%,即使是我们最大的模型也可以在单个 GPU 上不到 15 分钟内完成微调。 与广泛的探测练习相比,该技术在使用较少计算资源的情况下产生了更好的结果,证实了下游模型工程可以带来性能改进。 使用这种技术使微调在操作层面与探测具有竞争力,无论是在训练还是推理方面,同时实现了性能提升。 Para_04 虽然经过优化架构并在大型数据集上进行广泛训练的监督模型继续展示出卓越的性能,但DNA语言方法在包括组蛋白修饰、剪接位点和调控元件表征等多样任务中仍然与这些模型具有竞争力。 此外,NT模型代表了一种多功能的方法,可以无缝地针对人类和非人类物种的各种任务进行调整。 我们的无监督预训练方法的价值在处理较小数据集时尤为明显,在这种情况下,从头开始训练监督模型通常只能产生有限的结果。 认识到基础基因组学模型在该领域中的关键作用,我们进行了广泛的比较和基准测试研究,评估了我们的模型与四种不同架构的预训练模型:DNABERT-2、HyenaDNA-1 kb、HyenaDNA-32 kb和Enformer19。 这一稳健的基准将成为未来基因组学语言模型发展的参考点。 Para_05 通过对与变压器架构相关的各种分析,我们已经证明这些模型不仅获得了重建遗传变异的能力(补充说明),还能够识别关键的调控基因组元素。 这种特性在注意力图、嵌入空间、标记重建和概率分布的分析中得到了体现。 控制基因表达的关键调控元件,如增强子和启动子,被所有模型在多个头部和层中一致检测到。 此外,我们观察到每个模型至少包含一个层,该层生成的嵌入能够清晰地区分所分析的五种基因组元素。 鉴于自我监督训练促进了这些元素的检测,我们预计这种方法可以在未来的研究中用于表征或发现新的基因组元素。 Para_06 我们已经证明,变压器模型可以匹配甚至超越其他方法来预测变异效应和有害性。除了开发监督变压器模型外,我们还展示了基于零样本得分的效用,特别是用于预测非编码变异效应。 鉴于这些基于零样本的得分仅可以从基因组序列中得出,我们鼓励将其应用于非人类生物体,尤其是那些功能注释有限的生物体。 Para_07 最后,我们展示了增强模型架构的潜力,以实现减少模型大小和提高性能的双重好处。我们的后续一系列NT-v2模型具有12 kb的上下文长度。 为了更好地理解这一点,这个长度分别是DNABERT-1(500 bp)和DNABERT-2(3 kb)各自平均上下文长度的24倍和4倍。 这些先进模型不仅表现出更好的下游性能,还具备适合在经济硬件上执行和微调的优势。 这一优势源于它们的小巧尺寸,并结合了我们引入的高效微调技术。 Para_08 虽然我们的模型的注意力跨度仍然有限,但最近开发的Enformer模型表明,将感知场增加到200 kb对于捕捉人类基因组中的长程依赖关系是必要的。作者认为,这些依赖关系对于准确预测基因表达是必需的,而基因表达是由通常位于距离转录起始位点超过10-20 kb的远端调控元件控制的。 由于自注意力机制在序列长度上的二次缩放,使用标准变压器架构处理如此大的输入是不可行的。Enformer模型通过在到达变压器层之前通过卷积层来减少输入维度,从而解决了这个问题。 然而,这种选择削弱了其在语言建模中的有效性。另一方面,最近的HyenaDNA模型在训练时使用了高达1百万bp的感知场,但我们在下游任务中的基准分析表明,当训练中使用的感知场增加时,这些模型的性能迅速下降。 根据我们的结果,我们建议开发能够在处理长输入的同时保持高短输入性能的变压器模型是一个有前景的方向。在一个多组学数据迅速扩展的时代,我们最终期望这里提出的方法,以及可用的基准和代码,将促进基因组学中大型基础语言模型的采用、发展和改进。 Methods
Models 模型
Para_01 LMs 主要是在自然语言处理领域内开发的,用于模拟口语1,2。 一个 LM 是一种概率分布,它对任意序列的标记(通常是单词)给出概率;一个 LM 会返回该句子存在的概率。 由于 LMs 能够利用大量的未标注数据集生成通用表示,从而在很少有监督数据的情况下也能解决下游任务,它们因此变得流行起来49。 一种训练 LM 模型的技术是预测序列中被掩盖位置上最有可能的标记,这通常被称为掩码语言建模(MLM)。 受到在蛋白质研究领域中使用 MLM 所取得结果的启发4,5,在这个领域中蛋白质被视为句子,氨基酸被视为单词,我们将 MLM 应用于训练基因组学中的 LMs 变压器,将核苷酸序列视为句子,k-mer(k = 6)视为单词。 变压器是一类深度学习模型,在机器学习领域取得了突破性进展,包括自然语言处理和计算机视觉。 它们由一个初始嵌入层组成,该层将输入序列中的位置转换为嵌入向量,然后是一系列自我注意层,这些层依次细化这些嵌入。 使用 MLM 训练 LM 变压器的主要技术称为来自变压器的双向编码器表示(BERT)1。 在 BERT 中,序列中的所有位置都可以相互关注,允许信息在两个方向上流动,这在 DNA 序列的上下文中至关重要。 在训练过程中,网络的最终嵌入被馈送到语言模型头中,该头将其转换为输入序列上的概率分布。 Architecture 建筑
Para_01 我们所有的模型都遵循一种仅编码器的变压器架构。一个嵌入层将标记序列转换为嵌入序列。 然后,将位置编码添加到序列中的每个嵌入中,以向模型提供位置信息。 我们使用了一个可学习的位置编码层,该层最多接受1,000个标记。 我们使用了六聚体标记作为序列长度(最多6 kb)和嵌入大小之间的折衷,并且与其他标记长度相比,它实现了最高的性能。 每个变压器层通过一个层归一化层和一个多头自注意力层对其输入进行变换。 此操作的结果然后通过一个新的层归一化层和一个具有GELU激活函数的两层感知器传递。 每个模型的头数、嵌入维度、感知器隐藏层中的神经元数量以及总层数可以在补充表1中找到。 在自监督训练期间,由堆栈的最后一层返回的嵌入通过语言模型头部转换为序列中每个位置上现有标记的概率分布。 Para_02 我们的第二个版本 NT-v2 模型包括一系列被证明更有效的架构变化:我们没有使用学习到的位置嵌入,而是使用在每个注意力层中使用的旋转嵌入;我们使用没有偏置的带有 Swish 激活函数的门控线性单元,使自然语言处理更加高效。 这些改进的模型还接受最多 2,048 个标记的序列,从而导致长达 12 千字节的上下文窗口。 Training 培训
Para_01 这些模型是按照BERT的方法训练的。在每个训练步骤中,会采样一批标记化的序列。批量大小根据可用硬件和模型大小进行调整。我们在A100 GPU集群上进行了所有实验,并分别以14和2个序列的批量大小来训练‘500M’和‘2.5B’参数模型。 在序列中的15%的标记子集中,80%被替换为特殊的掩码(‘MASK’)标记。对于在人类参考基因组和多物种数据集上的训练运行,额外有10%的15%标记子集被替换为随机选择的标准标记(任何不同于类(‘CLS’)、填充(‘PAD’)或MASK标记的标记),这与BERT中的操作一致。 对于在1000G数据集上的训练运行,我们跳过了这个额外的数据增强,因为添加的噪声大于人类基因组中存在的自然突变频率。对于每个批次,损失函数被计算为预测概率与真实标记之间的交叉熵损失之和,在每个选定位置上。 梯度被累积以达到每批次100万个标记的有效批量大小。我们使用了Adam优化器,学习率调度和标准值用于指数衰减率和epsilon常数,β1 = 0.9, β2 = 0.999 和 ϵ = 1 × 10−8。在初始预热期间,学习率在16,000步内从5 × 10−5线性增加到1 × 10−4,然后按照平方根衰减直到训练结束。 Para_02
我们对 NT-v2 模型的超参数进行了轻微修改:优化器和学习率计划保持不变;然而,我们将批量大小增加到了 512(每批次 100 万个标记)。受到 Chinchilla 缩放定律的启发,我们还比其他深度学习模型训练了更长时间的 NT-v2 模型。具体来说,我们预训练了具有 5000 万和 2.5 亿参数的 NT-v2 模型,使用了 3000 亿个标记,而我们的‘2.5 亿’和‘5 亿’参数模型则训练了多达 1 万亿个标记,以了解其中的缩放定律。 相比之下,NT-v1 的 25 亿参数模型训练了 3000 亿个标记,其 5 亿参数的对应模型则训练了 500 亿个标记。最终,我们为 NT-v2 模型使用了以下模型检查点:‘5000 万’和‘1 亿’模型使用了 3000 亿个标记的检查点,‘2.5 亿’模型使用了 8000 亿个标记的检查点,‘5 亿’模型使用了 9000 亿个标记的检查点。 Probing 探测
Para_01 我们参考了通过探测模型嵌入的质量来解决下游任务的评估。训练完成后,对于每个任务,我们探测模型的每一层,并比较几种下游方法,以深入评估模型的表示能力。 换句话说,给定一个用于下游任务的核苷酸序列数据集,我们计算并存储由模型的十层返回的嵌入。 然后,将每一层的嵌入作为输入,我们训练了几个下游模型来解决下游任务。 我们测试了scikit-learn默认超参数的逻辑回归和一个多层感知器(MLP)。 正如我们观察到的那样,诸如学习率、激活函数和每层隐藏层的数量等超参数的选择会影响最终性能,因此我们也对每个下游模型进行了超参数扫描。 我们使用了一个十折验证方案,其中训练数据集被分成十个训练和验证集,这些集包含不同的洗牌,分别占初始集的90%和10%。 对于一组给定的超参数,在十个分割上训练了十个模型,并对其验证性能进行了平均。 这个过程用树结构的Parzen估计器求解器引导超参数空间搜索,运行了100次,然后在测试集上评估表现最佳的模型组。 因此,对于每个下游任务,对于每个预训练模型的十层,记录了在超参数搜索结束时的测试集性能。 在预训练模型及其各层中表现最佳的探针的超参数报告在补充表8中。 这种探测策略导致训练了760,000个下游模型,这提供了关于训练和使用语言模型的各种方面的详细分析,例如不同层对下游任务性能的影响。 Para_02 作为基线,我们评估了逻辑回归模型的性能,该模型将标记化的序列作为输入,然后通过变压器层传递这些标记。使用原始标记化序列作为输入的表现远远优于使用一个向量的情况,其中标记ID被独热编码并通过池化层(在序列长度轴上求和或平均)。 Fine-tuning 微调
错误!!! - 待补充 错误!!! - 待补充
Para_03 在实践中,通过重新缩放权重和分类/回归头权重引入的额外参数数量大约占模型总权重数量的0.1%。这加快了微调速度,因为只需要更新一小部分参数。 同样,这也减轻了存储需求,对于每个使用传统微调的下游任务,只需为5亿到25亿参数中的0.1%的新参数腾出空间。 例如,对于拥有25亿参数的模型,其权重占用9.5GB的空间。考虑到18个下游任务,经典微调需要9.5 × 18 = 171 GB的存储空间,而参数高效的微调仅需171 MB。 Para_04 与探测方案类似,训练数据集被分成十次,分为训练集和验证集,其中包含初始集的90%和10%的不同随机排列。 对于每次划分,模型经过了10,000步的微调,并使用在验证集上产生最高得分的参数来评估测试集上的模型。 我们使用了批量大小为八和学习率为3×10^-3的Adam优化器。 每个任务中的每个模型都进行了10,000步的微调。 这些超参数是根据它们在自然语言处理领域的有希望的结果选择的。 在我们的实验中,不同的超参数选择并没有带来显著的改进。 我们还比较了这种方法与从随机初始化检查点进行微调的方法。 Comparison with supervised BPNet baseline 与监督学习的BPNet基线进行比较
Para_01
作为基线监督模型,我们在每个18个任务上从头开始训练了不同变体的BPNet卷积架构。我们测试了原始架构(121,000个参数)、一个较大的架构(2800万)和一个特大架构(1.13亿)。对于每一个架构,我们手动调整超参数以在验证集上获得最佳性能。 我们实施了一种十折交叉验证策略来测量每种架构的18个模型的性能。 Comparison with published pre-trained genomics models 与已发布的预训练基因组学模型的比较
Para_01 我们比较了三种不同的预训练模型在18个下游任务上微调后的NT模型性能:DNABERT-2(参考文献23)、HyenaDNA(上下文长度为1-kb和32-kb;文献25)和Enformer19。由于DNABERT-1只能处理最大输入长度为512 bp的数据,因此无法用于大多数任务,所以我们在对比中排除了DNABERT-1。 我们将每个模型的架构和训练权重移植到我们的代码框架中,并按照上述方法对每个模型的transformer部分进行了参数高效的微调,使用相同的交叉验证方案以确保公平比较。 所有结果都可以在一个交互式的排行榜上可视化(https://huggingface.co/spaces/InstaDeepAI/nucleotide_transformer_benchmark)。 仅对于HyenaDNA,我们进行了完整的微调,因为我们的参数高效微调方法与该模型架构不兼容。 Para_02 请注意,Enformer 最初是通过监督学习来解决染色质和基因表达任务的。为了进行基准测试,我们重新使用了提供的模型躯干作为预训练模型,这并不是原始论文所推荐的用途;然而,我们认为这种比较有助于突出自我监督学习与监督学习在预训练方面的差异,并观察到即使对于与基因表达不同的任务,Enformer 仍然是一个非常有竞争力的基线。尽管如此,我们注意到 Enformer 最初是用不同的数据分割进行预训练的,因此它在我们的基准评估中的表现可能会因潜在的数据泄漏而被夸大。 Pre-training datasets 预训练数据集
The Human reference genome dataset 人类参考基因组数据集
Para_01 人类参考数据集是通过考虑参考组装GRCh38/hg38(https://www.ncbi.nlm.nih.gov/assembly/GCF_000001405.26)中的所有常染色体和性染色体序列构建的,总共有32亿个核苷酸。 The 1000G dataset 1000G 数据集
Para_01 为了向模型提供有关人类自然发生的遗传多样性信息,我们构建了一个训练数据集,其中包括来自不同人类群体的遗传变异。具体来说,我们从1000基因组(1000G)项目中下载了变体调用格式(VCF)文件(http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/data_collections/1000G_2504_high_coverage/working/20201028_3202_phased/),该项目旨在记录在人类群体中频率至少为1%的遗传变异。 该数据集包含3,202个高覆盖率的人类基因组,这些基因组源自27个地理结构化的人群,包括非洲、美洲、东亚和欧洲血统,详细信息见补充表2,总共构成20.5万亿个核苷酸。 为了允许从VCF文件中以FASTA格式进行单倍型重建,我们考虑了数据的相位版本,这对应于总共1.25亿个突变,其中1.11亿和1400万分别是单核苷酸多态性(SNPs)和插入缺失。 The Multispecies dataset 多物种数据集
Para_01 多物种数据集的基因组选择主要基于两个因素:(1)可用参考基因组的质量和(2)所用物种之间的多样性。 该数据集所选的基因组来自NCBI的参考序列(RefSeq)集合(https://www.ncbi.nlm.nih.gov/)中可用的基因组。 为了确保基因组的多样性,我们从RefSeq中每个主要分类群(古菌、真菌、脊椎动物哺乳类、其他脊椎动物等)中随机选取了一个属水平的基因组。 然而,由于可用细菌基因组数量庞大,我们选择仅包含其中的一个随机子集。 植物和病毒基因组未被考虑在内,因为它们的调控元件与本研究感兴趣的不同。 最终收集的基因组被降采样到总共850个物种,其基因组总共有1740亿个核苷酸。 每种类别的最终贡献,在核苷酸数量方面,与从NCBI解析的原始集合中的总数相同,显示在补充表3中。 最后,我们通过选择几个在文献中广泛研究的基因组来丰富这个数据集(补充表4)。 Data preparation 数据准备
Para_01 一旦收集了每个基因组/个体的FASTA文件,它们就被组装成一个唯一的FASTA文件,每个数据集在训练前都会进行预处理。在这个数据处理阶段,所有非A、T、C、G的核苷酸都被替换为N。 使用了一个分词器将字母字符串转换为标记序列。这个分词器使用的字母表是由A、T、C、G组合得到的46 = 4,096种可能的六聚体组合,以及五个额外的标记来表示独立的A、T、C、G和N。 它还包括三个特殊标记,即填充(PAD)、掩码(MASK)和序列开始(也称为类(CLS))标记。这使得词汇量达到了4,104个标记。 为了对输入序列进行分词,分词器会从一个类标记开始,然后从左到右转换序列,在可能的情况下匹配六聚体标记,或者在需要时退回到独立标记(例如,当出现字母N或序列长度不是6的倍数时)。 Para_02 对于多物种和人类参考数据集,基因组被分割成重叠的6,100个核苷酸片段,每个片段分别与前一个和后一个片段共享第一个和最后一个50个核苷酸。 作为数据增强练习,对于每个周期和片段,随机采样起始核苷酸索引在0到100之间,并从该核苷酸开始进行标记化,直到达到1,000个标记。 周期的数量根据数据集确定,以便模型在训练过程中总共处理3000亿个标记。 在每一步中,随机采样的序列批次在周期集合内被输入到模型中。 对于1000G数据集,按照上述方法准备的人类参考基因组的序列批次在每一步中被采样。 然后,对于每个采样的片段,从1000G数据集中随机选择一个个体,如果该个体在该片段对应的染色体和位置上携带突变,则将这些突变引入序列并替换相应的标记。 这种数据处理技术确保了在训练期间对基因组和个体的均匀采样,并使我们能够高效地仅存储每个个体的突变,而不是完整的基因组。 Hardware 硬件
Para_01 所有模型都在剑桥-1 Nvidia 超级计算机系统上进行训练,使用了16个节点,每个节点配备了八个 A100 GPU,总共使用了128个 A100 GPU。 在训练过程中,模型权重在每个 GPU 上复制,而批次则分布在各个 GPU 之间。 梯度在每个分片上计算,并在平均后在整个设备上反向传播之前累积。 我们依赖于 Jax 库(https://jax.readthedocs.io/en/latest/_autosummary/jax.pmap.html),该库依赖于 NCCL(https://developer.nvidia.com/nccl)协议来处理节点和设备之间的通信,并观察到随着可用 GPU 数量的增加,训练时间几乎呈线性减少。 具有5亿参数的模型在一个节点上训练了一天,而具有25亿参数的模型则需要整个集群训练28天。 具有不同参数数量(从5千万到5亿)的 NT-v2 模型同样在一个节点上训练了一天。 所有的微调运行都是在一个配备八个 A100 GPU 的节点上进行的。 对于训练运行,模型权重被复制,批次分布在各个 GPU 之间。 由于我们在微调时使用了批量大小为八,每个 GPU 在平均梯度并应用它们之前处理一个样本。 平均而言,具有5亿参数的模型的微调运行持续了20分钟,而具有25亿参数的模型则持续了50分钟。 Para_02 对于探测实验,所有嵌入(针对所有下游任务中的所有序列,以及每个模型的选定层)都在一个配备了八个A100 GPU的节点上计算并存储,计算过程需要2天。然后,在一个包含3,000个CPU的集群上拟合了760,000个下游模型,这需要2.5天。 Nucleotide Transformer downstream tasks
核苷酸转换器下游任务
Epigenetic marks prediction 表观遗传标记预测
Para_01 从ENCODE31(https://www.encodeproject.org/)获得了K562人类细胞系中十个组蛋白标记的组蛋白ChIP-seq数据。我们下载了带有以下标识符的bed narrowPeak文件:H3K4me3(ENCFF706WUF)、H3K27ac(ENCFF544LXB)、H3K27me3(ENCFF323WOT)、H3K4me1(ENCFF135ZLM)、H3K36me3(ENCFF561OUZ)、H3K9me3(ENCFF963GZJ)、H3K9ac(ENCFF891CHI)、H3K4me2(ENCFF749KLQ)、H4K20me1(ENCFF909RKY)和H2AFZ(ENCFF213OTI)。 对于每个数据集,我们选择了包含峰的1-kb基因组序列作为正例,并选择了所有不与峰重叠的1-kb序列作为负例。 Promoter sequence prediction 启动子序列预测
Para_01 我们构建了一个启动子序列数据集,以评估模型识别启动子基序的能力。我们从真核启动子数据库(https://epd.expasy.org/epd/)下载了所有人源启动子,涵盖了转录起始位点上游49个碱基对和下游10个碱基对(https://epd.expasy.org/ftp/epdnew/H_sapiens/006/Hs_EPDnew_006_hg38.bed)。这产生了29,598个启动子区域,其中3,065个是TATA盒启动子(使用基序注释https://epd.expasy.org/ftp/epdnew/H_sapiens/006/db/promoter_motifs.txt)。 我们选择了包含启动子的300个碱基对基因组序列作为正例,并选择了所有不与启动子重叠的300个碱基对序列作为负例。 这些正例和负例被用来创建三个不同的二分类任务:存在任何启动子元件(所有启动子)、具有TATA盒基序的启动子(TATA启动子)或没有TATA盒基序的启动子(非TATA启动子)。 Enhancer sequence prediction 增强子序列预测
Para_01 人类增强子元件是从ENCODE的SCREEN数据库(https://screen.wenglab.org/)中检索出来的。远端和近端增强子被合并在一起。 根据Meuleman等人(https://www.meuleman.org/research/dhsindex/)的词汇表,增强子被分为组织特异性和组织不变性。 与分类为组织不变性的区域重叠的增强子被定义为组织不变性增强子,而所有其他增强子被定义为组织特异性增强子。 我们选择了包含增强子的400个碱基对的基因组序列作为正例,并选择所有不与增强子重叠的400个碱基对序列作为负例。 我们创建了一个二元分类任务,用于判断序列中是否存在增强子元件(增强子),并创建了一个多标签预测任务,标签为组织特异性增强子、组织不变性增强子或无(增强子类型)。 Splice site prediction 剪接位点预测
Para_01 我们从GENCODE29 V44基因注释中获得了所有人工标注的剪接位点。注释经过过滤,排除了三级转录本(自动注释),以确保所有训练数据均由人工标注。 我们使用了HISAT2(参考文献56)中的extract_splice_sites.py脚本(https://github.com/DaehwanKimLab/hisat2/blob/master/hisat2_extract_splice_sites.py)来提取相应的剪接位点注释。 我们选择了包含中心剪接受体或供体位点的600个碱基对的基因组序列作为正例,并将所有不与剪接位点重叠的600个碱基对序列作为负例。 我们使用这些序列创建了三个不同的任务来评估剪接预测:一个多标签预测任务,标签为受体、供体或无(所有剪接位点);一个用于预测剪接受体的二分类任务(剪接受体);以及一个用于预测剪接供体的二分类任务(剪接供体)。 Dataset splits and performance evaluation
数据集划分和性能评估
Para_01 模型训练和性能评估是在人类基因组的不同染色体集上进行的。具体来说,来自20号和21号染色体的序列用于测试,其余的用于训练不同的模型。 我们通过将负样本子采样到与正样本相同数量来平衡每个数据集。 为了获得可以快速用于基准测试任何新设计选择或模型的小型数据集,我们进一步随机子采样示例,最多为30,000个用于训练,3,000个用于验证和测试(详见补充表5)。 Para_02 我们使用了十折交叉验证方案来评估每个模型,其中训练集被分成十个部分,每次我们用其中的九部分进行训练,并保留十分之一用于验证和选择最终检查点。 我们对每个模型重复这个过程十次,每次保留不同的十分之一用于验证,并在相同的保留测试集上评估最终性能。 我们使用十个部分中的中位数性能作为每个模型在给定任务上的最终性能指标。 Para_03 为了在不同任务中进行一致的性能评估,我们使用了二元或多元的MCC作为指标,因为它对不平衡的标签比例具有鲁棒性。 为了最终比较不同模型,我们计算了每个模型在三种不同任务类别(染色质谱型、调控元件和剪接)中的平均MCC,在每个类别中我们使用了跨任务的中位数MCC。 Additional downstream tasks 额外的下游任务
Chromatin profile prediction 染色质谱型预测
Para_01 我们使用了由Zhou等人编译的DeepSEA数据集(http://deepsea.princeton.edu/media/code/deepsea_train_bundle.v0.9.tar.gz)来进行染色质谱预测。该数据集包含240万个序列,每个序列长度为1000个核苷酸,并与919个染色质特征相关联。 这些特征包括690个转录因子、125个DNase和104个组蛋白特征。 与原始出版物一样,我们的模型同时在919个分类任务上进行训练,具有919个独立的分类头,损失是交叉熵损失的平均值。 由于每个标签高度不平衡且主要由负样本组成,因此正样本相关的损失被放大了8倍。 与DeepSEA方法不同,后者独立训练了两个模型,一个用于正向序列,一个用于相应的反向互补序列,并评估它们预测结果的平均值,本文提出的模型仅在正向序列上进行了训练。 SpliceAI benchmark SpliceAI 基准测试
Para_01 我们使用了Illumina Basespace平台(https://basespace.illumina.com/projects/66029966/)上可用的脚本来重现SpliceAI36中展示的训练数据集。简而言之,这个训练数据集是使用GENCODE V24lift37注释和来自GTEx队列的RNA-seq数据构建的,仅关注主转录本的剪接位点注释。 训练数据集包括位于染色体2、4、6、8和10-22以及X和Y染色体上的基因注释,而其余染色体上的基因注释(非旁系同源基因)构成了测试数据集。 每个通过该流程生成的序列长度为15,000 bp,其中中央5,000 bp包含需要预测的位点。 有关训练数据集构建的更多细节可以在原始出版物中找到。 我们调整了这个原始的SpliceAI数据集以运行我们的模型,将序列长度减少到6,000 bp(用于NT-v1模型)和12,000 bp(用于NT-v2模型),减少了侧翼上下文但保留了中央5,000 bp。 在与SpliceAI的第一个数据集进行比较时,我们称之为SpliceAI-6k,我们附加了9,000个‘N’核苷酸作为侧翼序列,因为SpliceAI基于一个具有15,000-bp输入的模型。 在与SpliceAI的第二个数据集进行比较时,我们报告了原始出版物中的性能表现,与该数据集相比,其中包括15,000 bp的序列长度而不是12,000 bp。
Para_02 这项任务是对每个输入序列的核苷酸进行多标签分类,类似于SpliceAI36(图5d)。从变压器模型输出的每个嵌入中,一个头部预测由令牌嵌入表示的六个核苷酸中的每一个的三个标签概率:剪接受体、剪供体或无。 这个头部是一个简单的分类层,它预测18个类别(每个核苷酸有三个标签)。 为了确保每个嵌入都与一个六聚体相关联,序列被截断,使其长度可以被6整除。 此外,训练集和测试集中的所有包含Ns的序列都被移除,这些序列在数据中所占比例微不足道。 请注意,如果我们要使用字节对编码标记器,例如DNABERT-2(参考文献23),每个嵌入所表示的核苷酸数量会有所不同,这将使核苷酸级别的预测任务显著变得更加复杂。 Enhancer activity prediction 增强子活性预测
错误!!! - 待补充
Additional performance analysis 额外的性能分析
t-SNE projections of embeddings 嵌入的t-SNE投影
Para_01 t-SNE(t-分布随机邻域嵌入)被用来将NT内部嵌入降维为二维向量,以可视化不同基因组元素的分离情况。 每个基因组元素的NT嵌入在多个变换层输出中计算,并计算出与该元素对应的序列位置上的平均嵌入。 然后,这些平均嵌入作为输入传递给sklearn Python包中的默认参数t-SNE降维对象。 Reconstruction accuracy and perplexity 重建准确性和困惑度
错误!!! - 待补充
Para_02 其中({{\mathcal{P"}}}_{{\rm{masked"}}}) 是被掩码位置的集合,({\mathcal{V"}}) 是词汇表(所有现有标记的集合)。困惑度通常在自回归生成模型的上下文中定义。在这里,我们依赖于Rives4中使用的另一种定义,并将其定义为在被掩码位置上计算的损失函数的指数: Para_03 困惑度衡量模型重构掩码位置的能力,它比准确率更精细,因为它还考虑了幅度。 与准确率相反,较低的困惑度表明更好的重构能力,因此性能更好。 Reconstruction of tokens in different genomics elements 在不同基因组元素中重构标记
Para_01
我们还在整个22号染色体上以6千碱基的窗口执行了这种标记重建方法。我们只保留了没有N的窗口。对于每个窗口,每次掩蔽一个标记,恢复序列中原标记的预测概率。我们将这些分数显示为WashU表观基因组浏览器会话(https://shorturl.at/jov28)。 为了获得不同类型基因组元素的平均标记概率,我们从Ensembl检索了基因注释区域(外显子、内含子、剪接受体和供体、5′ UTR、3′ UTR、转录起始位点和转录终止位点)(https://www.ensembl.org/info/data/ftp/index.html),从GENCODE检索了polyA信号位点(https://www.gencodegenes.org/human/),并从ENCODE检索了调控元件(来自SCREEN数据库的增强子、启动子和CTCF结合位点)(https://api.wenglab.org/screen_v13/fdownloads/GRCh38-ccREs.bed)。 Functional variant prioritization 功能变异优先排序
Para_01 为了获得具有不同相关严重程度的遗传变异,我们使用了Variant Effect Prediction(VEP)软件,并对人类基因组中的序列进行了注释。具体来说,我们在整个基因组中随机采样序列,并保留那些被注释为以下类别的遗传变异:‘内含子变异’、‘基因间变异’、‘调控区变异’、‘错义变异’、‘3'非翻译区变异’、‘同义变异’、‘转录因子结合位点变异’、‘5'非翻译区变异’、‘剪接区域变异’和‘终止增益变异’。 在仅保留单核苷酸多态性(SNPs)并过滤掉被注释为多个后果的变异(例如同时被注释为终止增益和剪接变异的变异)之后,我们最终获得了每种类别包含920个遗传变异的数据集。 Para_02 作为功能性遗传变异的正集,我们从四个不同的资源中编译了SNPs。我们使用了与基因表达(eQTLs)相关的SNPs和来自GRASP数据库中P值<10^-12的meQTLs,ClinVar中标注为‘可能致病’的SNPs以及在HGMD(公共版v.2020.4)中报告的SNPs。 经过这些过滤后,我们分别为eQTLs、meQTLs、ClinVar和HGMD SNP数据集保留了总计80,590、11,734、70,224和14,626个遗传变异。 对于这四个数据集中的每一个,然后基于来自1000G项目的次要等位基因频率>5%且不与测试数据集中报告的任何变异重叠并且位于相关变异100 kb范围内的SNPs构建了一组负变异,从而产生了四个平衡的数据集。 Para_03 为了计算给定兴趣位点的零样本得分,我们进行了以下操作。对于每个SNP,我们基于人类参考基因组获得了以该SNP为中心的6,000个碱基对序列。 然后,我们创建了两个序列,一个携带参考等位基因,另一个在SNP位置携带替代等位基因。 接着,我们计算了多个零样本得分,这些得分捕捉了这两个序列在嵌入空间中向量距离的不同方面,即:(1) L1距离(曼哈顿距离);(2) L2距离(欧几里得距离);(3) 余弦相似度;以及(4) 点积(未归一化的余弦相似度)。 我们还计算了替代等位基因的损失以及携带替代和参考等位基因的序列之间损失的差异,作为另外两个零样本得分。 在功能性变异的情况下,除了零样本得分外,我们还微调了变换器模型来分类正变体和负变体。 我们采用了之前描述的类似策略,主要区别在于训练集和测试集按染色体划分,并严格保持不重叠。 具体来说,我们将22条染色体分为五个集合,并依次将每个集合用作测试集,其余四个集合用作训练集。 通过微调训练集,我们可以推导出测试集中每个序列成为正变体的概率。 Para_04 为了将这些预测与其他方法进行比较,我们从四个数据集中各随机抽取了10,000个正负SNP。 然后,我们使用Combined Annotation Dependent Depletion工具(v.GRCh38-v1.6)计算了CADD、GERP、phastCons和phyloP分数。 DeepSEA分数是使用在https://hb.flatironinstitute.org/sei/上可用的Beluga模型计算的。 Attention maps analysis 注意力图分析
Para_01 我们分析了从预训练模型中收集的注意力图如何捕捉关键的基因组元素。我们遵循了之前工作中提出的方法论。 对于一个基因组元素,我们定义了一个在标记上的指示函数 f(i),如果标记 i 中的一个或多个核苷酸属于该元素,则 f(i) 等于 1,否则为 0。 我们在一个注意力头中计算了关注该基因组元素的平均比例,通过在一个核苷酸序列数据集 X 上进行聚合计算如下: 错误!!! - 待补充
Para_03 我们计算了所有模型的所有头部和所有层的pα(f)值,并考虑了九个元素(‘5′ UTR’, ‘3′ UTR’, ‘外显子’, ‘内含子’, ‘增强子’, ‘启动子’, ‘CTCF结合位点’, ‘开放染色质’ 和 ‘转录因子结合位点’)。 我们在由90,000个序列组成的数据集上进行了这些分析,每个特征有10,000个序列,长度为6 kb,从人类参考基因组中提取。 对于每个序列,在数据集创建期间,特征在序列中的位置被均匀采样。 正如先前的工作所建议的那样,我们为所有实验选择了一个置信度阈值μ = 0.3。 Para_04 我们认为,如果特征的量pα(f)显著大于数据集中该特征的自然出现频率(补充表10),则该特征被一个注意力头捕获。 为了验证这一点,我们进行了一个两比例z检验,零假设为特征的自然频率,备择假设为pα(f)。 每个模型的总头数被用作对显著性水平(α)为0.05的Bonferroni校正。 我们计算了每个模型中每个基因组元素的每个头的z分数和相应的P值,如下所示: Data availability Para_01 NT 预训练序列来自公开资源。1000G 序列来自 http://ftp.1000genomes.ebi.ac.uk/vol1/ftp/data_collections/1000G_2504_high_coverage,人类和多物种参考基因组来自 https://ftp.ncbi.nlm.nih.gov/genomes/refseq/。 基因注释来自 GENCODE 数据库(https://www.gencodegenes.org/)和 Ensembl 数据库(https://www.ensembl.org)。 变异效应预测使用 Ensembl 的 VEP API 获得(https://www.ensembl.org/info/docs/tools/vep/index.html)。 致病性和调控变异从 ClinVar(https://ftp.ncbi.nlm.nih.gov/pub/clinvar/vcf_GRCh38/)、GRASP(https://grasp.nhlbi.nih.gov)和 HGMD(公共 v.2020.4)通过 Ensembl Biomart 提取。 NT 下游任务整理自以下公开资源:ENCODE31,32,33(https://www.encodeproject.org/, https://screen.wenglab.org/ 和 https://www.meuleman.org/research/dhsindex/),真核启动子数据库30(https://epd.expasy.org/epd/)和 GENCODE29(补充表 5)。 我们提供了我们的多物种和人类预训练数据集的 HuggingFace 版本(https://huggingface.co/collections/InstaDeepAI/nucleotide-transformer-65099cdde13ff96230f2e592)以及我们的下游任务基准(https://huggingface.co/datasets/InstaDeepAI/nucleotide_transformer_downstream_tasks_revised)。 我们建立了一个交互式排行榜(https://huggingface.co/spaces/InstaDeepAI/nucleotide_transformer_benchmark),其中包含每个任务中所有模型的结果,以便于比较。 我们还在 WashU 表观基因组浏览器上创建了一个交互式浏览会话,展示了整个 22 号染色体上的预训练模型标记概率,网址为 https://shorturl.at/jov28。 Code availability Para_01 预训练变压器模型的模型代码和权重以及Jax中的推理代码可以在https://github.com/instadeepai/nucleotide-transformer上用于研究目的。PyTorch版本的模型可以在https://huggingface.co/collections/InstaDeepAI/nucleotide-transformer-65099cdde13ff96230f2e592找到。示例笔记本可以在HuggingFace上找到,网址为https://huggingface.co/docs/transformers/notebooks#pytorch-bio。