从零基础到AI蛋白设计需要多久？

奋斗的小鸟

从零基础到AI蛋白设计需要多久？
原文地址：https://www.zhihu.com/question/517473317

感恩由您

2024年科学界翘首以盼的诺贝尔化学奖盛大揭晓，这一奖项毫无悬念地聚焦于AI蛋白设计领域，瞬间成为全球瞩目的焦点。美国科学家大卫・贝克凭借在蛋白质结构设计领域的开创性贡献，揽获一半奖金，其对生命构建模块的精准把控，创造出全新蛋白质的壮举，让世人惊叹。另一半奖金则由英国科学家德米斯・哈萨比斯和美国科学家约翰・江珀共同收入囊中，他们利用人工智能技术，成功预测几乎所有已知蛋白质的结构，为蛋白质研究开辟了全新路径。


这一成果的诞生，是人工智能与科学研究深度融合的里程碑。近年来，人工智能以令人惊叹的速度渗透到科研的各个角落，从繁杂的数据收集整理，到精巧的实验设计，都有AI的身影。此次诺贝尔化学奖的颁发，不仅是对科学家们卓越成就的认可，凸显了AI在蛋白质研究领域的巨大潜力，更标志着AI在蛋白质研究领域的重大突破，为未来生物技术的发展指明了方向。


AI蛋白设计的核心原理
▷蛋白质，作为生命活动的主要承担者，其功能与结构紧密相关。传统上，科学家们通过漫长而复杂的实验方法来研究蛋白质的结构与功能，探索新的蛋白质设计。而AI技术的介入，为这一领域带来了革命性的变革。


▷AI蛋白设计的核心在于利用深度学习算法，对海量的蛋白质数据进行分析与学习。这些数据包括蛋白质的氨基酸序列、三维结构以及它们在不同环境下的功能表现等。通过对这些数据的深度挖掘，AI模型能够理解氨基酸序列与蛋白质结构之间的复杂关系，从而实现对蛋白质结构的精准预测与设计。
▷以深度学习算法中的神经网络为例，它通过构建多层的神经元网络结构，模拟人类大脑的学习过程。在蛋白质结构预测方面，神经网络可以将已知的蛋白质氨基酸序列作为输入，经过层层计算与分析，最终输出对应的蛋白质三维结构预测结果。这种方法大大提高了预测的准确性和效率，使得科学家们能够在短时间内获得大量蛋白质的结构信息。
▷在蛋白质设计方面，AI则基于对蛋白质结构与功能关系的理解，通过优化算法，设计出具有特定功能的全新蛋白质序列。例如，科学家可以设定目标蛋白质的功能需求，如催化特定化学反应、结合特定分子等，然后AI模型会根据这些需求，在大量的可能序列中进行搜索与优化，生成最有可能满足功能要求的蛋白质序列。


2024诺贝尔化学奖
体现AI蛋白设计核心
▷2024年诺贝尔化学奖得主的研究成果，正是AI蛋白设计原理的杰出体现。美国科学家大卫・贝克在蛋白质结构设计领域取得了开创性进展。他的团队通过开发先进的计算方法，能够精确地设计出具有特定结构和功能的全新蛋白质。他们利用AI模型，对蛋白质的折叠过程进行模拟和预测，从而实现了对蛋白质结构的精准控制。这一成果不仅加深了我们对生命基本构建模块的理解，更为开发新型生物材料、药物等提供了强大的工具。


▷另一半奖金获得者，英国科学家德米斯・哈萨比斯和美国科学家约翰・江珀，则在蛋白质结构预测方面取得了重大突破。他们带领团队开发的人工智能程序，能够利用深度学习算法，根据蛋白质的氨基酸序列快速而准确地预测其三维结构。这一成果打破了传统实验方法在蛋白质结构解析上的瓶颈，极大地推动了蛋白质研究的进展。


美尚洁：AI蛋白设计产业化的先行者
▷在全球因诺贝尔化学奖为AI蛋白设计欢呼时，中国的美尚洁生物科技有限公司早已在这一领域默默耕耘多年，成为AI蛋白设计产业化的先驱者。早在2019-2021年期间，美尚洁凭借在AI技术在氨基酸序列研究的深厚积累，成功取得国内首创 “重组纤连人源化胶原蛋白” 专利技术。他们利用AI强大的数据分析与模拟能力，深入探究蛋白质结构与功能的内在联系，从分子层面精准设计出全新的重组纤连人源化胶原蛋白。


▷2023年美尚洁生物科技就已成为国内唯一一家将AI技术用于氨基酸序列优化测序的重组蛋白企业。2024年在诺贝尔化学奖揭晓的次日，美尚洁生物科技有限公司与途深智合公司强强联合，聚焦 “AI 技术用于重组蛋白氨基酸序列研究领域”，借助途深智合的TourSynbio™大模型，进一步为美尚洁的重组融合蛋白研发注入强大动力。
▷同年，美尚洁在张家港建立了自己的GMP管线，这一生产线严格遵循最严格的药品生产质量管理规范，从原材料采购到产品出厂，每一个环节都有精准的AI技术监控与质量把控。这条管线将公司前期研发的各类重组蛋白产品，实现了从实验室到规模化生产的无缝衔接，确保高品质、稳定供应的产品源源不断地走向市场，满足医疗美容、生物医药等多个领域对于高端重组蛋白产品的迫切需求。


从诺奖到产业
美尚洁开启产业化新篇
▷随着美尚洁生物科技的不断发展壮大，其影响力也将持续辐射。在医疗美容领域，基于AI设计的重组蛋白产品将为肌肤修复、抗衰老等提供更安全、更高效的解决方案；在生物医药领域，有望借助AI优化的重组蛋白药物靶点设计，攻克一系列疑难病症的药物研发难题，加速新药上市进程，为全球患者带来新的希望之光。


▷美尚洁生物科技有限公司始终专注于通过基因编辑以及合成生物技术，矢志不渝地实现重组蛋白在各应用领域对传统动物源胶原蛋白的全面替代。在这一宏伟愿景的驱动下，美尚洁将继续砥砺前行，为人类健康与科技进步贡献更多来自中国的智慧与力量。

卡卡

1.背景介绍
深度学习领域的快速发展对蛋白质设计产生了重大影响。深度学习方法最近在蛋白质结构预测方面取得了突破，形成了可用于数百万种蛋白质的高质量模型。随着生成建模和序列分析的方法的发展，深度学习在过去几年提高了识别新蛋白质序列和结构的准确性和能力，从而极大地促进了蛋白质设计发展。深度神经网络现在可以学习和提取蛋白质结构的基本特征，预测它们如何与其他生物分子相互作用，并能够创造新的治疗疾病的潜在有效药物。随着它们在蛋白质设计中的适用性的迅速增长，作者回顾了深度学习方法的最新发展和技术，并提供了它们在生成新的功能蛋白方面的例子。
在蛋白质设计中，人们通过修改一个自然序列（定向进化）或从头开始设计一个新的序列（从头设计）来探索蛋白质空间。在定向进化中，遗传多样性被引入自然蛋白质，继续通过迭代高通量或虚拟筛选/突变，以实现模拟自然选择的理想功能。相比之下，从头设计—本综述的主要重点—并不依赖于自然蛋白质序列，而是使用计算方法来设计新的蛋白质。考虑到蛋白质的序列-结构-功能关系，三种从头设计类型如图1所示。这三个主要任务包括从结构设计序列、从功能设计序列和从功能设计结构。通常从头设计方法的目的是为一个给定的蛋白质结构构建一个新的蛋白质序列。这个问题的设置被称为“逆蛋白质折叠”。然而，了解目标蛋白的期望功能，人们可以探索其结构或直接探索序列空间来设计功能。所有这些方法的困难在于需要探索可能的蛋白质序列的巨大空间，这使得蛋白质设计具有计算密集型的特质，因此需要有效的数据驱动方法。这篇综述不仅旨在概述深度学习驱动的蛋白质设计的当前能力，而且还为以应用为重点的研究人员在深入学习中一些最具影响力的概念提供更深入的见解，帮助他们理解和评估这个快速发展的领域的新方法。该综述的结构如下：首先描述了支持大多数最先进的蛋白质设计工具的相关深度学习方法。这部分是具有相当技术性的，并且应该是对深度学习方法的一个高级的概述。然后，总结了该工具集在蛋白质设计方面的最新应用，并围绕序列-结构-功能关系的不同方式进行构建（见图1）。




图1：从头设计的三种类型。左边的方框列出了在现代计算蛋白质设计中发挥关键作用的深度学习概念（更多细节请参见“深度学习的最新创新”部分）。右边的图表显示了这些想法是如何被应用于各种蛋白质设计任务中的。主要的结构到序列的方法通过图表示来利用结构信息，并结合基于几何和物理的归纳偏差来预测相关的序列。功能到结构和功能到顺序的方法，在不同的部分中讨论，大多基于生成建模，并利用最适合相应的目标数据域的数据表示和操作。
2.深度学习
在讨论蛋白质设计师使用深度学习的不同方式之前，作者总结了他们可以使用的一套工具。本文涵盖了专门为处理图结构数据而开发的神经网络体系结构、具有精心设计的输出约束的数据高效模型，序列数据的Transfomer和生成模型。
图已经被应用于不同抽象层次的生物学问题，从表示分子结构到细胞中的大型相互作用网络。在蛋白质设计中，图通常将节点分配给原子或氨基酸，这些氨基酸由表示化学键或空间关系的边缘连接，如图2A所示。这种为节点和边缘分配生物单位的灵活性允许基于图的学习算法根据特定任务的需求定制局部性的概念。这一转变背后的主要驱动力之一是人们对图形机学习和消息传递神经网络的兴趣迅速增长，这导致了其在许多蛋白质设计任务上取得了重大进步和卓越的性能。




图2：蛋白质残基图(A)的例子从蛋白质结构构建图（PDB：1FSD）。节点用残差标识进行注释。边缘的定义是特定于用例的，例如，可以基于肽序列（蓝色）或空间接近度（洋红色）。(B)消息传递。计算所有邻居的消息，并聚合以更新节点i的特征。
归纳偏差是对作为学习算法设计基础的问题的先验信念或知识。用于蛋白质设计的神经网络的诱导偏差是由对输入表示（序列、结构、表面等）的选择来定义的，以及数据是如何传递给算法的。例如，序列信息可以编码为一个单个氨基酸串、一个多序列比对（MSA）或一个特定位置的评分矩阵（PSSM），而结构信息的选项包括成对距离图、离散体素网格、具有原子分辨率的图或粗粒度图。另一个例子是MaSIF框架所采用的表面表示，以解释蛋白质之间的相互作用主要涉及其表面这一事实。通过提取蛋白质核心的细节，MaSIF在分子识别任务上取得了最先进的结果，如蛋白质-蛋白质相互作用预测。数据表示的选择强烈地影响算法将产生的预测，必须为给定的任务仔细选择。




图3:蛋白质语言模型(A)蛋白质语言模型是根据目前超过2亿个序列的序列数据进行训练的。(B)对具有掩蔽标记预测目标的语言模型进行自我监督训练。(C)蛋白质语言模型可以应用于逆折叠问题、蛋白质结构预测，以及其他各种任务，包括评估蛋白质嵌入（TAPE）的任务。此外，预先训练的语言模型产生的序列嵌入可以帮助探索“蛋白质宇宙”。
在过去的几十年里，新的测序技术使蛋白质序列数据库呈指数级增长，每年序列的总数每年翻一番，如图3A所示。随着时间的推移，测序成本变得越来越便宜，蛋白质序列数据库的增长速度甚至超过了计算能力。因此，利用序列信息来研究蛋白质的结构和功能特性已经成为一个很有前途的选择，大量可用的蛋白质序列使得通过训练蛋白质语言模型（PLMs）捕获蛋白质组中的统计模式成为可能-学习观察到的氨基酸序列的概率分布。蛋白质语言模型在大量的蛋白质序列语料库上进行训练后，能够学习蛋白质序列和进化的规则，这些规则通常被称为生命的“语言”。
深度生成模型深度生成模型是一类近似概率分布的神经网络，通过大数据集进行经验定义。将它们看作是概率框架而不是神经网络架构是很有用的。在下面概述的每一种情况下，不同的神经网络主干（如cnn、gnn、Transformer）可以根据手头的任务使用。根据生成方法的类型，数据分布可以显式或隐式地建模。然而，所有的生成算法的设计方式都可以从学习到的分布中采样新的数据点。在蛋白质设计的背景下，生成模型很有潜力，因为它们不依赖于能量函数的准确性，因此可以被视为对传统设计方法的补充方法。随着更大的数据集的使用，这种方法在基于结构的设计中越来越受欢迎。目前，最流行的深度生成模型类型包括了变分自编码器（VAEs）、生成对抗网络（GANs）和扩散模型。




图4：深度生成模型(A)变分自编码器通过联合学习潜在变量z的分布和将这些潜在向量转换为数据点的概率解码器来优化数据似然的下界。为了生成新数据，从正态分布中采样随机向量，并由解码器模块进行处理。(B)生成对抗网络由两个部分组成：生成器试图模拟数据分布，而鉴别器被训练将样本分类为“真实”或“假”。经过训练后，生成器可以用来将随机噪声样本转换为与训练数据集的统计数据相似的真实数据。(C)扩散模型学习逐渐去噪训练数据的扰动版本。通常，原始信号x被一个固定的扩散过程破坏，该过程将数据转化为噪声，并训练模型通过具有参数化转移概率的马尔可夫链恢复数据样本。
3.基于深度学习的蛋白质设计
3.1结构到序列
Rosetta的早期版本是蛋白质结构预测和设计的领先生物信息学工具之一，它能够使用具有相似局部序列的结构片段组合和模拟退火方案准确预测小的蛋白质结构域。不久之后，该框架迅速配备了蒙特卡罗优化程序，用于寻找自然折叠和新结构的稳定序列，从而将其转化为蛋白质设计工具。最近，AlphaFold前所未有的成功引发了另一波这类浪潮。
一类流行的序列设计方法试图通过序列空间的优化来反转结构预测模型-逆折叠任务。一种常见的方法是从一个随机的氨基酸序列开始，这将导致一个随机的结构，并使用一个目标折叠的置信度和结构损失来优化它。通常，结构相似性是用旋转和平移不变的方式来测量的。例如，一个蛋白质可以被表示为每个残基对的距离、角度和二面体的二维矩阵，这些残基间特征之间的差异可以用来定义一个优化目标。结构和稳定性损失已成功地应用于各种蛋白质设计任务，如自由幻觉蛋白设计。如图5所示。




图5：受深度学习结构预测算法启发的设计管道在基于MCMC的生成中，在初始序列中引入随机突变，并使用损失函数对预测进行评分。当损失减少或满足接受标准时，接受突变序列。这些步骤被重复进行多次迭代，从而得到一个满足损失函数的序列。基于梯度下降的生成（红色箭头）计算输入在损失之上的梯度，从而产生一个误差梯度矩阵。这个矩阵用于优化PSSM，从中可以提取最可能的序列进行下一次迭代。几轮优化结果是PSSM收敛，可以采样可能的蛋白质序列。
解决“逆折叠”问题的另个方法是训练模型可以明确地预测给定一个结构的序列，例如，ProteinMPNN是一个基于图的神经网络。作者和其他快速采用该工具的研究人员报告了在其各种蛋白质设计任务上的高实验成功率，如单体蛋白质、对称重复结构、纳米颗粒和目标蛋白质结合物。
3.2功能到结构
在许多蛋白质设计应用中，目标主干是未知的，必须从头生成。考虑到定制的功能，这里的主要目标是通过优化选定的结构来设计新的蛋白质，这些结构可以执行特定的任务或表现出特定的特性，如酶活性、与目标分子的结合亲和力、溶解度或稳定性。迄今为止，许多工作将这项任务定义为一个分布学习问题，并利用了“深度生成模型”部分中讨论的概率框架之一。如图6所示，生成模型（VAEs、GANs和扩散模型）的应用跨越了整个蛋白质设计子任务范围，如从头骨架，小分子或蛋白质结合物、对称低聚物、表位特异性抗体和基序支持支架。人们通常非常关注的是允许用户引导生成过程，从而根据自己的规范设计蛋白质。到目前为止，基于扩散的蛋白质设计方法通过只设计蛋白质的部分，或者通过操作更粗粒度的蛋白质表示来克服了这个问题。例如，ProtDiff引入了一种基于粒子滤波的新采样方案，该方案可用于从无条件训练的扩散模型中生成近似条件样本，例如，用于基于给定功能主题的支架采样。这在实践中很有用，因为它减少了为条件生成任务创建新的、专门的数据集的负担。Chroma以低温采样方法有原则地交换样本多样性和质量，并使用分类器指导让用户控制蛋白质的形状、对称性或类，甚至支持以文本输入为条件的蛋白质设计。这一类的一个具体例子是RFdiffusion，它得到了广泛的实验验证的支持。RFdiffusion是上述讨论的许多想法与预先训练的蛋白质结构预测模型（RoseTTAfold）的结合。它与ProteinMPNN一起形成了一个成熟的相关方法工具箱，用于各种蛋白质设计应用，如稳定单体设计、对称低聚物设计、功能基序或酶活性位点的支架设计，以及蛋白质结合物设计。




图6：功能到结构/序列设计的一般工作流程可考虑多种蛋白质功能为起点，包括抗体设计、基于蛋白质的生物发光（或荧光素酶)、酶活性、结构稳定性和折叠、设计特定低聚态和蛋白质结合物。流行的深度学习方法的选择是VAEs、GANs、扩散模型和蛋白质语言模型。
3.3功能到序列
与以结构为中心的蛋白质设计不同，功能到序列设计从所需的功能或活性开始，目的是识别可以执行该功能的氨基酸序列，而无需明确建模蛋白质的折叠状态。
生成的VAE和GAN模型在这个领域中发现了许多应用程序（图6）。然而，最先进的蛋白质序列建模，是大型的蛋白质语言模型。一系列的工作已经表明，蛋白质语言模型是有效探索蛋白质序列空间的良好候选模型。如ProGen，给定一系列输入的蛋白质属性，如生物功能或相关家族，它以进行下一个标记预测任务进行训练。在推理时，它将一系列属性作为输入，并从头开始迭代地生成蛋白质序列。
4.讨论
深度学习对蛋白质设计景观产生了变革性的影响，可以说是当今研究最积极的计算蛋白质设计方法家族。虽然严格的数学往往不能描述生物系统的复杂性，因此不适合工程化，但深度学习已经证明了它在复杂的高维空间中寻找模式的有效性。在本文中，作者描述了最新的深度学习技术的工具箱计算蛋白质设计者和他们目前如何应用于创建(1)序列为给定的骨干，(2)新结构为特定目的，和(3)序列基于函数规范没有显式建模的结构。虽然有重要的证据表明，深度学习方法将成为这个具有挑战性且不断发展的蛋白质设计领域的正确工具箱，但仍有几个局限性需要考虑。从技术的角度来看，第一个挑战是数据的可用性。深度学习模型需要大量的标记数据来进行训练。然而，获得这些具有良好质量的蛋白质数据是具有挑战性和昂贵的。另一个尚未被该领域最终回答的相关问题是，大多数机器学习模型对蛋白质结构的大量依赖是否总是一个合理的假设。事实上，在结构水平上的工作似乎提供了一种归纳倾向，即针对更有物理意义的解决方案，但同时也是生物现实的重要部分（溶剂相互作用、翻译后修饰、动力学等）。在本文综述中讨论的方法中仍然缺乏理解，或者至少没有明确地建模，这对这种模型偏差的有效性提出了质疑。该领域的一个更广泛的问题是机器学习方法的泛化和分布外外推。由于监督学习算法从根本上将复杂的非线性函数拟合到所提供的数据点，随着我们远离训练例子，它们的不确定性会大幅增加。作为最后的技术挑战，作者提到深度学习模型的可解释性。主要的“黑盒”方法不能为其预测提供任何解释，因此阻碍了人类有针对性的干预。由于可解释性的缺乏，也限制了我们对蛋白质功能的新见解的程度，并阻止我们在运行耗时和昂贵的体外实验之前失败的情况。其中一些问题通过置信估计技术来解决，这些技术在越来越多的数据驱动的蛋白质设计工具中实现，但需要更多的工作来提高计算蛋白质设计的鲁棒性。
蛋白质设计中深度学习方法面临的第二个主要挑战是验证策略。研究人员利用不同的设计蛋白来评估他们的模型的性能。这些目标通常是根据它们与实际用例的相关性或它们所构成的挑战来选择的。由于缺乏严格的基准测试，使不同的方法和蛋白质设计策略之间不能够进行公平的比较。在少数几篇通过实验验证了其基于深度学习的蛋白质设计的论文中，对实验成功率（如成功表达、折叠或靶点结合）的定义相当模糊。目前的计算分析价值有限，因为评估指标和基准数据集都没有标准化，这意味着新方法不能轻易地与现有方法进行比较。
-------------------------------------------
欢迎点赞收藏转发！
下次见！

检验医师

【爆肝万字】基于Rosetta的配体结合位点的设计（Rosetta and the Design of Ligand Binding Sites ）

1. 前言

Rosetta是蛋白质设计中的重要平台，在蛋白质设计领域中占有重要地位。本文以2016年的文献Rosetta and the Design of Ligand Binding Sites进行实践，探索用Rosetta进行配体结合位点的设计。经过三个月实践成功，得到本文档，在此分享给读者，希望能够帮到有需要的读者完成实验工作，早发文章！
觉得有用还请点个赞呀！ヾ(◍°∇°◍)ﾉﾞ
文献信息：Methods Mol Biol. 2016 ; 1414: 47–62. doi:10.1007/978-1-4939-3569-7_4.
Rosetta and the Design of Ligand Binding Sites - PubMed


总体设计流程如下图所示。


2. 材料

（1）运行类Unix操作系统（如Linux或MacOS）的计算机。尽管测试运行和小型生产运行可以在传统的笔记本电脑和台式机系统上执行，但建议在生产运行中使用多处理器计算群集。
（2）Rosetta。Rosetta建模包可以从RosettaCommons网站（https://www.rosettacommons.org/software/license-and-download）获得。Rosetta许可证免费提供给学术用户。Rosetta以源代码形式提供，必须在使用前编译。有关如何编译Rosetta的说明，请参阅Rosetta文档（https://www.rosettacommons.org/docs/latest/）。本文中的方案已使用Rosetta每周发布版本2021.16.61629进行了测试。
特别注意：版本不同，使用时会有所差别，比如原文章中使用的是2016版本的Rosetta，具体到2021年的版本是有些差异的。
（3）处理小分子的程序。OpenBabel是一个自由软件包，允许操作许多小分子文件格式。有关下载和安装信息，请访问http://openbabel.org获取下载和安装信息。本文中的协议已在OpenBabel版本2.3.1中进行了测试。也可以使用其它小分子操作程序。
（4）配体构象生成程序。我们推荐BCL，可从http://meilerlab.org/index.php/bclcommons免费获得，用于学术用途，但需要剑桥结构数据库的额外许可证才能生成构象异构体。本文中的协议已使用BCL 3.2版进行了测试。也可以使用其他构象异构体生成程序，如Omega、MOE或RDKit。
（5）标准格式的靶小分子结构，如SDF或SMILES。
虽然Rosetta可以忽略远离结合位点的链断裂和缺失环，但蛋白质的结构在配体结合区域应该是完整的。如果结合口袋缺失残基，使用起始结构作为模板，用比较建模方案重塑这些残基。
（6）待重新设计的蛋白质结构，PDB格式。
可接受的格式取决于小分子处理程序的能力。OpenBabel可用于将大多数小分子表示（包括SMILES和InChI）转换为Rosetta所需的sdf格式。
在与密切相关的配体的络合物中测定的高分辨率实验结构是最理想的，但不是必需的。可以使用未配体化蛋白质的实验结构甚至同源性模型。
3. 方法

总述：在整个方案中，${ROSETTA}表示安装Rosetta的目录。文件内容和要在终端中运行的命令以斜体显示。假定使用bash shell-其他shell的用户可能需要修改命令行的语法。
参照本部分理论，使用Rosetta的relax protocol优化蛋白结构，在存放蛋白结构文件（以DEBP.pdb为例）的目录下打开Linux命令行界面。我本地Rosetta安装路径为：
/home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle注意：运行时需要找到个人安装Rosetta的路径进行替换。
3.1 将蛋白质结构预松弛到Rosetta评分函数

具体命令如下：
/home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle/main/source/bin/relax.mpi.linuxgccrelease -ignore_unrecognized_res -ignore_zero_occupancy false -use_input_sc -flip_HNQ -no_optH false -relax:constrain_relax_to_start_coords -relax:coord_constrain_sidechains -relax:ramp_constraints false -s DEBP.pdb命令正确，开始运行，输出日志！


运行结束后输出文件为DEBP_0001.pdb。


为方便起见，重命名输出结构。可以直接在xftp中更改文件名，或用mv命令更改文件名（命令如下）。
mv DEBP_0001.pdb DEBP_relaxed.pdb3.2 准备配体

3.2.1 将小分子转化为SDF格式，包括根据需要添加氢

从RCSB数据库中下载AKG的‍sdf‍小分子构象文件AKG_ideal.sdf。
使用‍OpenBabel软件为‍AKG补加氢原子，在存放配体构象文件的目录下打开‍Win10‍命令行窗口。


输入以下命令：
obabel AKG_ideal.sdf -p 7.4 -O AKG.sdf其中，-p‍指定在‍pH‍=‍7.4‍的条件下为‍AKG加氢。程序运行结束后会输出氢原子数符合体液‍pH‍条件的‍AKG构象‍AKG.sdf。
运行后显示“1 molecule converted”，即一个分子已经转化完成。


在原路径下新生成了AKG.sdf这个文件。


以下对OpenBabel软件处理前后的sdf文件做一下简单对比，可以发现处理后去掉了两个氢原子。即程序认为4个氢原子符合体液‍pH‍为7.4的条件。


分析AKG的结构可以发现一些端倪，AKG中含有两个羧基，在ph设置为7.4的弱碱性环境下，可以猜想去掉两个氢原子是合理的。


后续尝试将PH设置为14，强碱性环境下失去两个氢原子，而PH=1的强酸性环境下没有去掉氢原子，即酸性环境下AKG不会失去氢原子，这和之前的猜想是一致的。


3.2.2 生成配体构象异构体库，将构象异构体库转换为Rosetta格式的“参数文件”

不同于蛋白大分子，配体小分子一般都具有相当大的灵活性，同一分子会以不同的构象与蛋白发生相互作用。因此，为了在设计时考虑到小分子的灵活性，需要生成小分子的构象异构体集合。
原文献中使用的是BCL的本地应用程序，而以下使用更为简便的BCL::Conf网络服务器生成配体构象异构体库。
从[文献](https://pubmed.ncbi.nlm.nih.gov/33351632/)中了解到BCL下载，以及BCL::Conf网络服务器的网址（如下），二者均可以实现生成配体构象异构体库，但实际使用时打不开。
BCL下载：http://meilerlab.org/bclcommons
BCL::Conf网络服务器：http://meilerlab.org/bclconf
之后通过检索找到可用的[BCL::Conf网络服务器](Meiler Lab)。


在 Molecules选项中上传已加氢的AKG.sdf；在Conformer Generation Mode选项中选择Standard；在Number of Conformers Desired / Molecule选项中填写100。


提交设置和任务后，提示在页面要等待大约48秒。


结果不到30秒就自动刷新页面显示了结果。


点击conformers_bcl_conf.sdf文件即可下载，该文件中含有100个加氢后的AKG分子的构象。将其重命名为AKG_conf.sdf以便进行后续操作。
接下来将该构象集合转换为 Rosetta 可读取的、用于储存小分子形状及化学性质信息的专有格式：params 文件。在存放构象异构体集合文件（AKG_conf.sdf）的目录下打开 Linux 命令行界面，输入以下命令：
${ROSETTA}/main/source/scripts/python/public/molfile_to_params.py -n LIG -p LIG -- conformers-in-one-file AKG_conf.sdf其中，-n 指定在 pdb 和 params 文件中用来表示配体名称的 3 字符缩写，这里命名为 LIG 即配体 ligand 的缩写（需要注意的是，这里不能沿用 GLY 或者其他氨基酸、金属原子的缩写，否则生成的 params 文件会在后续突变结构时，与 Rosetta 自带的氨基酸或部分金属原子的 params 文件发生冲突）；-p 指定生成文件的命名。
我实际输入命令（${ROSETTA}代表的路径不同，需替换成本地的路径）：
/home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle/main/source/scripts/python/public/molfile_to_params.py -n LIG -p LIG --conformers-in-one-file AKG_conf.sdf输出正常，很快运行完成。


在Xftp中查看输出文件，可以看到新产生了3个文件：配体形状及化学信息文件 LIG.params，Rosetta 选定的 1 个配体构象LIG.pdb，配体构象集合 LIG_conformers.pdb。


用PyMOL查看LIG.pdb和LIG_conformers.pdb的结构，分别如下图所示。




至此，配体准备工作完成。
3.3 将配体置于蛋白质中

Rosetta要求原子名称与molfile_to_params.py步骤中生成的名称匹配。即使具有正确放置了配体的起始结构，也应该将molfile_to_params.py生成的结构对齐到口袋中，以便原子命名正确。
在进行蛋白设计时会包含蛋白与配体的多次对接，因此需要先指定配体每个原子的起始坐标，这样配体才能在指定的对接盒子中移动和扭转，与结合口袋发生相互作用。
为此，只需要在PyMOL中打开蛋白结构、以及前文Rosetta选定的其中1个配体构象LIG.pdb，将配体移动至结合口袋内即可。
为了提高定位的准确性、缩小对接需要搜索的空间范围，在PyMOL中打开含原配体Glu 的蛋白结构2VHA.pdb和 LIG.pdb，将鼠标操作改为编辑模式，移动和旋转LIG.pdb，使之与2VHA.pdb中的原配体Glu在羧基处大致重合，如下图所示。


调整好LIG.pdb的位置之后，删除原配体Glu和蛋白DEBP的复合物结构2VHA.pdb，只留下已移动到结合口袋处的AKG，保存为 LIG_positioned.pdb。
这个过程相对简单，实际操作中需要注意一些细节，以下进行简单说明。
（1）如何在PyMOL中设置鼠标操作为编辑模式？
在界面右下角可以看到“Mouse Mode”，即鼠标操作模式。


只需在上图红框区域单击鼠标即可实现显示模式（3-Button Viewing）和编辑模式（3-Button Editing）的切换。


（2）在PyMOL中调整配体位置过程中的细节
因为上述操作过程中，我们只想移动配体LIG.pdb，而不想蛋白质2VHA.pdb位置发生改变，因此可以在显示模式对蛋白进行固定，固定方式为：选中蛋白结构→A(ACTION)→MOVEMENT→PROTECT。
在编辑模式下，一些鼠标操作如下：

• Shift键+鼠标左键——旋转配体；
  
• Shift键+鼠标中键——可平移配体；
  
• 长按鼠标右键滑动——缩放视野；
  
• Ctrl+Shift+鼠标中键——将鼠标所在位置设置为旋转中心。
  

（3）为蛋白质设置透明度和透明方式
在移动和旋转配体时，不透明的蛋白显示方式会造成遮挡，不便于操作。因此可以在显示模式下选中蛋白质，为其设置透明度，如果显示的方式为cartoon，将其设置成透明度设置为50%，可以按以下步骤操作：点击菜单栏中的Setting→Transparency→cartoon→50%。
更多透明度和透明方式的设置，可参见[PyMOL中文教程](入门教程 - PyMOL中文教程 2022.09 文档)。
（4）检查
调整好LIG.pdb的位置，保存为LIG_positioned.pdb之后，可以重新打开LIG_positioned.pdb和2VHA.pdb查看重合情况。如下图可以看到，LIG_positioned.pdb和原配体Glu重合情况较好。


3.4 运行Rosetta设计

3.4.1 突变位点声明

实际操作时，这个步骤很简单，就是新建一个.txt文档，在其中输入以下文字，保存后将文档后缀名由`.txt`改为`.resfile`即可。
ALLAA
AUTO
start
1 × NATAA说明：
其中，“start”以上为 Header 部分，“start”以下为 Body 部分。
Header 声明了对所有未在 Body 中声明的残基要采取何种操作：“ALLAA”即允许选定残基突变为 20 种常规氨基酸中任意一种；“AUTO”即让 Rosetta 自动检测结合界面确定待突变残基。即 Header 允许位于结合界面的残基突变为 20 种常规氨基酸中任意一种，其他位置的残基不发生突变。
Body 声明的操作作为 Header 的例外被优先确定，即当 Header 和 Body 对同一残基采取的操作发生冲突时，以 Body 为准。“1 X NATAA”即允许 X 链第 1 位残基采取 repack操作，即发生构象变化。但 1R9L 的蛋白结构只有 A 链，关于 X 链的情况需要特别说明：在 3.2.3 配体文件准备与预处理中，将构象集合 sdf 文件转换为 params 文件时，生成的LIG.pdb，LIG_conformers.pdb 都会将 GLY 放置在 X 链；LIG.params 中 MM AtomType项也设置为 X。这是因为 Rosetta 默认 X 链是配体、配体原子类型是 X。因此 X 链第 1位残基就是指目标配体 GLY。即 Body 只允许配体构象变化，不允许配体发生突变。
附：浅谈Rosetta Mutfile与Resfile用法 - 知乎
谷雨：浅谈Rosetta Mutfile与Resfile用法
补充：Commands for for controlling sequence identity（控制序列识别的命令）
参考自（下面是对内容的翻译）：Resfile syntax and conventions
每个命令都用于限制每个位置允许的氨基酸种类。如果组合多个命令，只有被每个命令单独允许的氨基酸才会被包含。这是由于对PackerTask类操作的交换性质。
注意：应该记住，resfile命令是限制性的，而不是允许性的。虽然不允许的残基不会被重新打包/重新设计，但是各种协议提供了额外的限制，可能进一步限制突变身份。特别地，许多协议会禁止重新设计二硫键半胱氨酸，即使在resfile中明确列出为可突变。
每个命令都用于限制每个位置允许的氨基酸种类。如果组合多个命令，只有被每个命令单独允许的氨基酸才会被包含。这是由于对PackerTask类操作的交换性质。 注意：应该记住，resfile命令是限制性的，而不是允许性的。虽然不允许的残基不会被重新打包/重新设计，但是各种协议提供了额外的限制，可能进一步限制突变身份。特别地，许多协议会禁止重新设计二硫键半胱氨酸，即使在resfile中明确列出为可突变。

• ALLAA … 允许所有20种氨基酸，包括半胱氨酸（与ALLAAwc相同）
  
• ALLAAwc … 允许所有20种氨基酸（默认）
  
• ALLAAxc … 允许除半胱氨酸外的所有氨基酸
  
• POLAR … 只允许极性的标准氨基酸（DEHKNQRST）
  
• APOLAR … 只允许非极性的标准氨基酸（ACFGILMPVWY）
  
• PROPERTY <property> … 不允许任何缺乏给定属性的残基类型
  
• NOTAA <list of AAs> … 只不允许指定的氨基酸（使用单字母代码，不用分隔符，如ACFYRT。对于非标准氨基酸，使用X[<full name>]。）
  
• PIKAA <list of AAs> … 只允许指定的氨基酸（使用单字母代码，不用分隔符，如ACFYRT。对于非标准氨基酸，使用X[<full name>]。）
  
• NATAA … 只允许本来的氨基酸（NATive Amino Acid）- 重新打包而不设计
  
• NATRO … 保留输入的旋转异构体（完全不打包）（NATive ROtamer）
  

3.4.2 设计protocol

准备对接和设计脚本（“design.xml”）。建议的方案基于使用RosettaScripts框架的RosettaLigand对接。它将优化配体在结合口袋中的位置（low_res_dock），重新设计周围的侧链（design_interface），并在设计的环境中优化相互作用（high_res_dock）。为了避免虚假突变，在每个位置（favor_native）给输入残基一个轻微的能量奖励。
新建一个记事本文件，将以下内容写入，并将文件另存为design.xml即可。
<ROSETTASCRIPTS>
    <SCOREFXNS>
        <ScoreFunction name="ligand_soft_rep" weights="ligand_soft_rep.wts"/>
        <ScoreFunction name="hard_rep" weights="ligandprime.wts"/>
    </SCOREFXNS>
    <TASKOPERATIONS>
        <DetectProteinLigandInterface name="design_interface" cut1="6.0" cut2="8.0" cut3="10.0" cut4="12.0" design="1" resfile="DEBP.resfile"/>
    </TASKOPERATIONS>
    <LIGAND_AREAS>
        <LigandArea name="docking_sidechain" chain="X" cutoff="6.0" add_nbr_radius="true" all_atom_mode="true" minimize_ligand="10"/>
        <LigandArea name="final_sidechain" chain="X" cutoff="6.0" add_nbr_radius="true" all_atom_mode="true"/>
        <LigandArea name="final_backbone" chain="X" cutoff="7.0" add_nbr_radius="false" all_atom_mode="true" Calpha_restraints="0.3"/>
    </LIGAND_AREAS>
    <INTERFACE_BUILDERS>
        <InterfaceBuilder name="side_chain_for_docking" ligand_areas="docking_sidechain"/>
        <InterfaceBuilder name="side_chain_for_final" ligand_areas="final_sidechain"/>
        <InterfaceBuilder name="backbone" ligand_areas="final_backbone" extension_window="3"/>
    </INTERFACE_BUILDERS>
    <MOVEMAP_BUILDERS>
        <MoveMapBuilder name="docking" sc_interface="side_chain_for_docking" minimize_water="true"/>
        <MoveMapBuilder name="final" sc_interface="side_chain_for_final" bb_interface="backbone" minimize_water="true"/>
    </MOVEMAP_BUILDERS>
    <SCORINGGRIDS ligand_chain="X" width="15">
        <ClassicGrid grid_name="vdw" weight="1.0"/>
    </SCORINGGRIDS>
    <MOVERS>
        <FavorNativeResidue name="favor_native" bonus="1.00"/>
        <Transform name="transform" chain="X" box_size="3.0" move_distance="0.1" angle="5" cycles="500" repeats="1" temperature="5" rmsd="4.0"/>
        <HighResDocker name="high_res_docker" cycles="6" repack_every_Nth="3" scorefxn="ligand_soft_rep" movemap_builder="docking"/>
        <PackRotamersMover name="design_interface" scorefxn="hard_rep" task_operations="design_interface"/>
        <FinalMinimizer name="final" scorefxn="hard_rep" movemap_builder="final"/>
        <InterfaceScoreCalculator name="add_scores" chains="X" scorefxn="hard_rep"/>
        <ParsedProtocol name="low_res_dock">
            <Add mover_name="transform"/>
        </ParsedProtocol>
        <ParsedProtocol name="high_res_dock">
            <Add mover_name="high_res_docker"/>
            <Add mover_name="final"/>
        </ParsedProtocol>
    </MOVERS>
    <PROTOCOLS>
        <Add mover_name="favor_native"/>
        <Add mover_name="low_res_dock"/>
        <Add mover_name="design_interface"/>
        <Add mover_name="high_res_dock"/>
        <Add mover_name="add_scores"/>
    </PROTOCOLS>
</ROSETTASCRIPTS>注：在XML脚本中更改resfile的名称，使其与正在使用的resfile的完整路径和文件名相匹配。切割值决定如何使用AUTO规格处理残留物。将设计在配体的切口1埃内具有C-β原子的所有AUTO残基，以及切口2内指向配体的所有残基。选择侧链的逻辑分别类似于切口3和切口4，但具有侧链灵活性而不是设计。在设计阶段将忽略剖切壳之外的任何内容，但在其他阶段可能会移动这些内容。
注：网格宽度必须足够大以容纳配体。对于较长的配体，请将该值至少增加到配体的最大延伸长度加上“变换”移动器中box_size值的两倍。
注23：box_size的值设定配体从起始位置的最大刚体置换。rmsd的值设置与起始位置的最大允许均方根偏差。如果希望使设计的配体更接近起始构象，请将这些值设置为较小的值；如果希望允许更多的移动，请将这些值设置为较大的值。这些仅是方案主动采样阶段的限值。在协议的其他阶段可能会发生额外的移动。
3.4.3 运行设计应用程序

输入以下指令运行。-nstruct 5000代表生成5000个pdb文件，可根据需要进行更改。
/home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle/main/source/bin/rosetta_scripts.mpi.linuxgccrelease -ex1 -ex2 -linmem_ig 10 -restore_pre_talaris_2013_behavior -parser:protocol design.xml -extra_res_fa LIG.params -s "DEBP_relaxed.pdb LIG_positioned.pdb" -nstruct 5000 -out:file:scorefile design_results.sc在目录下生成了所需的pdb文件。


3.5 进行两轮筛选

（1）第1次筛选
为了查看高分结构，在当前目录下打开 Linux 命令行界面，输入以下命令：
awk '{print $NF,$2}' design_results.sc | sort -nk2 > sorted_scores.txt程序运行结束后会输出20,000个结构按total_score项排列的分数记录文件，该项代表总体得分的大小，一般为负数，负得越多说明能量越小，结构越稳定。
为了获得上述各项的平均值，需要先过滤一次以查看附带输出的各项平均值信息。创建之后要写入过滤参数的空文件metric_thresholds.txt，在当前目录下打开 Linux 命令行界面，实际输入以下命令：
perl /home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle/main/source/src/apps/public/enzdes/DesignSelect.pl -d <(grep SCORE design_results.sc) -c metric_thresholds.txt -tag_column last > filtered_designs.sc命令行界面在程序一开始执行时会输出各项得分的平均值，将过滤时要用到的得分项平均值记下，删除本次输出的 filtered_designs.sc（也可不删除，因下次过滤输出的结果会覆盖此文件）。


根据过滤项平均值，打开metric_thresholds.txt写入过滤参数。写入内容如下：
req total_score value < -1238.12
req fa_rep value < 97.62
req if_X_fa_atr value < -12.93
req if_X_fa_rep value < 1.92
req ligand_is_touching_X value > 0.8
output sortmin total_score在当前目录下打开Linux命令行界面，再次输入上次过滤的命令，程序运行结束后会输出filtered_designs.sc文件，文件包含了按interface_delta_X项排序的滤过结构列表。


（2）第2次筛选
进入到要产生构象的文件目录下，输入以下命令：
awk '{print $NF ".pdb"}' filtered_designs.sc >filtered_pdbs.txtRosetta的InterfaceAnalyzer计算了许多额外的指标。不过，这些可能需要时间来评估，因此最好只在一组经过预过滤的结构上运行。在生成指标之后，可以如前面的步骤中那样过滤结构。这将生成一个分数文件（“design_interfaces.sc”），其中包含所选PDB的计算的指标值。
实际使用指令：
/home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle/main/source/bin/InterfaceAnalyzer.mpi.linuxgccrelease -interface AB_X -compute_packstat -pack_separated -score:weights ligandprime -no_nstruct_label -out:file:score_only interfaces.sc -l filtered_pdbs.txt -extra_res_fa LIG.params这步比较费时，相当于在第一步筛选结果的基础上，进一步计算新的指标并进行筛选，相当于计算一千多个构象的新的指标。运行结束之后得到interfaces.sc文件，包含第1次筛选后留下的1168个构象。


为了获得上述各项的平均值，需要先过滤一次以查看附带输出的各项平均值信息。创建之后要写入过滤参数的空文件metric_thresholds_2.txt，在当前目录下打开 Linux 命令行界面，输入以下命令：
perl /home/hdeng/software/rosetta_bin_linux_2021.16.61629_bundle/main/source/src/apps/public/enzdes/DesignSelect.pl -d <(grep SCORE interfaces.sc) -c metric_thresholds_2.txt -tag_column last > filtered_interfaces.sc运行后上述命令之后各项指标的平均值输出如下。


根据过滤项平均值，打开metric_thresholds_2.txt写入过滤参数。写入内容如下：
req packstat value > 0.75
req sc_value value > 0.75
req delta_unsatHbonds value < 4.82
req dG_separated/dSASAx100 value < -0.81
output sortmin dG_separated在当前目录下打开Linux命令行界面，再次输入第2次过滤的命令，程序运行结束后会输出filtered_interfaces.sc文件，文件包含了按dG_separated项排序的滤过结构列表和各项分数，总共筛选得到93个构象。


通过运行以下命令可以得到只包含已通过筛选的文件的名称的文件filtered_pdbs_final.txt。
awk '{print $NF ".pdb"}' filtered_interfaces.sc >filtered_pdbs_final.txt

至此，通过Rosetta进行筛选的步骤结束。
3.6 手动检查选定序列

虽然自动化程序在不断改进，并且可以在有限程度上替代，但在评估设计时，仍然没有人的专业知识可以替代。领域专家对界面的视觉检查可以捕获难以编码到自动过滤器中的系统特定需求。
3.7 重新应用设计方案，从步骤3.4开始

通过对来自前几轮设计的输出结构重复设计协议，可以获得改进的结果。设计轮次的数量取决于系统及其收敛速度，但通常需要3-5轮设计，每轮设计都从上一轮的过滤结构开始。
在后续循环中，随着配体稳定到首选结合方向，可能希望降低低分辨率采样阶段的积极性（步骤3.4.2中“变换”移动器的box_ size和rmsd值）。由于输出结构同时包含蛋白质和配体，因此传递给“-s”选项的值上不再需要引号（参见步骤3.4.3）。与之相反，你可能希望使用“-l”选项，该选项采用每行包含一个输入PDB的文本文件的名称。每个输入PDB都将生成“-nstruct”模型。减小此值，使每轮中未筛选的输出结构的总数大致相同。
3.8 设计迭代

根据计算协议进行的初始设计很少能给予精确的预期结果。通常需要执行设计和实验的迭代循环，使用从实验中获得的信息来改变设计过程。

队长是我

天然蛋白质具有临界稳定性的特征，然而临界稳定性使得蛋白质遭受胁迫压力后极易发生错误折叠并失去功能。体内蛋白质在错误折叠后产生的聚集沉淀被认为是多种疾病发生发展的原因。因此，优化蛋白质的稳定性是科学研究与工程应用领域亟待解决的关键问题。Rosetta是一种生物物理建模工具，根据蛋白质的氨基酸序列有效预测蛋白质的结构，在此基础上可以从头设计各种类型的全新蛋白质。基于 Rosetta 系列算法的蛋白设计在过去十年中在创新蛋白药物、抗体、疫苗、新型合成生物学元件及纳米药物等生物大分子研究领域中被广泛使用。
课以Rosetta软件为基础，以实例讲授和练习为主。依次讲授Rosetta软件基础、蛋白质结构viewer、Linux和Python基础、结构扰动与结构优化、蛋白质复合物预测、抗体抗原模型处理与对接、序列设计PackRotamer和FastDesign，掌握RAbD常用命令、SSD和MSD设计、CartisenDDG 突变扫描、RosettaScript开发流程、序列与结构设计、从头蛋白质设计等的操作等多个内容。




另有计算机辅助药物设计CADD专题、集成多组学数据的机器学习在生物医学中的应用专题。
了解更多欢迎私信~

计算机辅助药物设计、rosetta 蛋白抗体设计、机器学习集成多组学系列专题2023

大力水手

等下文的文章的源码公布出来，就只需要点击几下鼠标了。
引入

长期以来的研究范式告诉我们，在蛋白设计中，直接设计原子水平的蛋白结构是极其困难的：直接设计结构的话，没有那么强的结构学习网络来捕捉复杂的能量函数；设计序列时的话，误差在形成结构时会被成倍放大，使得设计成功率大幅下降。这使得长期以来只有2种研究范式：

• Inside-out策略：设计蛋白骨架->设计侧链，得益于Rosetta Commons长期以来的经验；
  
• 用强大的语言模型生成序列：得益于近年来测序数据的暴增和超大规模的transformer。
  

有没有可能实现从头的结构生成？随着AlphaFold2的公布，一些研究团队也按耐不住利用AF2来做结构生成，基于迭代优化策略的AFdesign随之诞生，那真正end-to-end的结构生成模型还会远吗？
2022.4.4，Eva-MariaStrauch团队发布了一个用3个结构不同的组合模型组合在小尺度上生成结合短肽结构的模型，给新的方法开了一个曙光；2022.5.28，一篇组合3个相似模型做结构生成的paper横空出世，作者是Namrata Anand和Tudor Achim，这篇文章不仅能直接从二级结构描述就产生有side chain的结构，还能做蛋白修补(inpainting)，简直是最近1年的里程碑之一，特此分享下阅读笔记。
这篇文章在GitHub上发布的：
Protein Structure and Sequence Generation with Equivariant Denoising Diffusion Probabilistic Models作者介绍

Namrata Anand是Po-Ssu Huang组的博士生，Harvard应用数学的本科，在Stanford读了一个AI方向的硕士后于16年进入Po-Ssu Huang组，主要从事计算蛋白设计。



Namrata Anand的LinkedIn

相信各位都听说或者读过过她在进Po-Ssu Huang组之后几篇著名的文章：《Generative modeling for protein structures》，《Fully differentiable full-atom protein backbone generation》这两篇是用GAN做骨架设计；《IG-VAE》这篇用VAE；《Protein sequence design with a learned potential》这篇则基于3D-CNN；还有《Progen》这篇基于transformer做序列生成...大有学术新星的态势。
Tudor Achim则是Stanford的一位二年级CS博士生，同时是Helm.ai的CTO，主要领域是高维概率推断，在这篇文章之前并没有蛋白设计领域相关工作。



Tudor Achim的LinkedIn

模型简介

这篇文章主要解决蛋白结构数据生成的难题：

• 与图像不同，蛋白在离散网格上没有一个自然的表示，无法直接套现有的生成模型。
  
• 与图像不同，蛋白没有自然的规范取向。因此，非旋转不变性的方法必须在模型权重中直接考虑这个变化因素，这降低了可以有益结构变化的生成能力。
  
• 在蛋白生成过程中，局部或全局结构中的重大错误会导致不可信的蛋白结构。
  

现今的结构生成模型还只能在小分子从头生成，或者如同IG-VAE那样单一拓扑结构的蛋白结构生成(还是非侧链的)。作者提出了几个解决要点：

• 一个扩散模型(去噪概率扩散模型)，扩散模型是生成模型种类的一种，知乎上有很多大佬有此类模型的介绍；
  
• 同时作者利用了AF2中的IPA模块(这个模块有很好迁移性，Jefferey J. Gray组在设计抗体时也巧妙运用了IPA，详见此处)，使得模型有SE(3)等变性；
  
• 为了处理蛋白生成所需的旋转参照系的扩散，作者使用了一个利于SO(3)差值的公式；
  
• 作者利用类似掩码语言建模的方法来生成序列；
  
• 加入足够的约束以生成最终的蛋白。
  

下面分开介绍这几点。
模型结构

蛋白数据

相信各位对蛋白结构数据已有充分了解，作者认为生成模型需要在N个残基组成的蛋白中，学习如下结构数据中变量的生成先验：

• ，蛋白骨架3D坐标；
  
• ，定义上文以 为中心的正则坐标系的全局旋转的四元数，用来恢复NCO原子的位置；
  
• ，i位置氨基酸类型
  
• ， 原子上侧链的4个 角。
  

扩散模型

扩散模型是一类潜在变量模型，将数据生成过程建模为随机先验的迭代去噪，近似后验分布的特定参数化可以解释为向固定先验分布“扩散”。从数据分布 中采样的数据点 的数据生成(反向)过程由转移核 和先验分布π递归定义，简单回顾下定义公式：


在连续情况下，近似后验，称为先验过程，将数据点 向随机先验扩散：

其中 是根据一个固定的方差计划来选择的。作者使用一个神经网络 来参数化反向转移核：
作者通过训练中最小化下面的变分边界来得到：

其中 由 对 进行噪声处理得到，FAPE见《noi19：蛋白质结构评价函数与Alphafold2中的损失函数——FAPE》。最后，采样依赖于学习的 来执行一个反向过程，将样本从先验分布映射到数据分布的样本。
扩散旋转

与坐标不同，旋转变量 和 不存在于平面几何的欧几里得流形上。因此，在训练和采样中，它们不能像坐标那样通过随机缩放和扰动它们的编码简单地扩散到先验分布。为解决这些局限性，对于先验分布 ，作者从SU(2)抽样均匀随机旋转。接下来，作者基于下表从 插入到一个样本 。作者选择球形线性插值(SLERP   ，通过因子 差值 到 )。对于转子扭转角扩散(一维旋转)，作者在 上使用均匀先验，并在单位圆上插值进行噪声和采样(“Interp”)。这些设计选择在训练和测试时都达到了预期的效果，使网络暴露在相似的随机旋转分布中。


旋转不变性

作者给去噪模型采用了一种等变transformer。该模型以中间蛋白结构 作为输入，并产生去噪后的真实数据结构估计值 。作者将transformer中的标准注意机制替换为AF2中的IPA模块。IPA将节点query和value特征划分为三维向量，并将其从目标节点的参考框架转换为全局参考框架，然后计算注意权值和注意机制的输出。注意层的输出与输入蛋白的全局方向不相关，因此在 局部坐标系中所预测的修正结果是等变的。作者使用的损失函数为FAPE。
序列扩散

作者使用类似于BERT的方法在骨架上生成序列，可以解释为带有吸收状态的扩散。在训练过程中，掩码部分作为t的函数在[0,1]中线性插值。在测试时，作者在 t = T时通过掩蔽所有残基来执行相反的过程，并从一个模型中迭代采样，该模型的输入残基以概率t/T独立掩蔽，采样过程中t从t = T逐步到t = 0。
约束条件

简洁地描述蛋白的一个方法是通过指定数量的残基N，然后一个元组的数字长度B添加到N表示二级拓扑结构块的大小，然后一块二级结构赋值 ,最后一块对称邻接矩阵在 连同并行/反并行先验对每个 配对。在实践中，作者通过在训练过程中删除块的邻接，以及默认情况下不包括任何loop块的邻接信息，来允许块邻接矩阵的不规范。这种拓扑先验的粗编码方案不会过分地约束模型，只产生一个结构；它允许出现显著变化性的结构。
训练&采样过程

生成模型的条件是对给定蛋白的紧凑规范的约束，这些约束通过一个具有三角形自我注意的transformer来embedding，以产生feature embedding，然后使用IPA来处理，以产生Cα原子局部坐标系中平移、旋转和残基的update。这些update在训练中被用来计算旋转不变的loss，在采样中逐步走向最终的结构。(采样详细过程见扩散旋转小节的图)


动图(11s)：


模型构建

作者用CATH 4.2 S95结构域数据训练模型，95%训练集，5%测试集。
约束嵌入模型

作者用5层类GPT-3 transformer和25层三角形注意力模块编码二级结构信息，然后对二级结构embed和块邻接进行向下采样。
扩散模型

扩散模型对约束网络的输出和当前结构进行了条件约束，并对最终结构构型进行了猜测。当预测 和 时，每个IPA模块产生一个中间backbone更新，在计算下一轮IPA之前将其应用于结构。
结构扩散模型
结构扩散网络有12层，4个IPA头，每个残基4个query，每个残基4个value，并使用最大256个残基的向下采样骨架进行训练。网络训练时T = 1000，batch size=160，lr=0.001。
序列扩散模型
该网络有15层，4个IPA头，每个残基4个query，每个残基4个value，并使用最大128个残基的向下采样骨架进行训练。网络训练时T = 100，batch size=160，lr=0.001。
Rotamer扩散模型
该模型将整个具有扩散侧链的全原子蛋白质结构视为输入。没有使用约束嵌入，并且在向前传递中没有完成骨架更新。该网络有6层，4个IPA头，每个残基4个query，每个残基4个value，并使用最大75个残基的向下采样骨架进行训练。网络训练时T = 500，batch size=100，lr=0.0005。
结构修补
保持整个结构固定，除了一个区域被屏蔽。对于训练中的每个数据点，以0.6的概率执行“块扩散”，以0.4的概率执行“连续扩散”。在块扩散中，每个loop以0.25的概率被屏蔽，其他块以0.05的概率被屏蔽。在连续扩散中，随机选择连续块向先验扩散，每个起始残基的概率为0.03，长度均匀分布在1到15之间。
结构修复模型由T = 1000，batch size=160，lr=0.001 ，由结构扩散模型进行微调。
联合结构&序列修补
保持预训练的结构修补网络freeze，并对结构扩散网络输出的预测结构进行微调。训练序列模型时， ，batch size=160，lr=0.0005。
训练和采样

见模型构建中的训练&采样过程小节。作者用了3个独立的模型来给结构( 和 )、序列( )、rotamer( )来做扩散：

• 结构扩散模型用15倍坐标缩放的真实结构数据训练；
  
• 序列扩散模型用真实结构训练；
  
• rotamer扩散模型用真实结构&序列数据训练。
  

作者用这套组合拳模型做了几个常见的蛋白设计任务：无序列从头生成、蛋白补全(inpainting)、序列设计、rotamer repacking。
结果

无序列从头生成

作者测试了几个典型结构：全helix、全 和 混合作为约束输入，生成结果如下图所示



看起来模型学到了不少接近原生的结构

蛋白补全

（我之前的笔记翻译为“蛋白修补”，whatever，还没有正式翻译，感觉“蛋白补全”好一点）
在这个任务中，作者训练一个附加模型 以满足现有结构的条件，在训练过程中固定结构的部分，并对固定部分的补充部分执行正向扩散过程。对于所有不向先验扩散的残基，它们的位置在训练和采样期间固定在它们的真实值位置。如下图，蛋白修补区域的键合几何分布与原生结构中的对应分布是一致的。作者还从采样中看到，该模型可以在原子水平上找到loop分布的离散模态。





蓝色为修补区域

动图(8s):


作者想看看修补过程能不能发现新拓扑结构，在采样时，作者通过简单地修改块的底层二级结构的长度来修改块的邻接条件反射信息。下图展示了如何使用模型以物理上合理的方式修改结构——生成理想化的拓扑，改变loop长度，并修改固定拓扑的次级结构长度。这些合成结构与PDB中发现的自然结构不同，这表明该模型编码了有用的物理先验。



可控生成:该模型可实现多种蛋白结构的可控生成模式。



可变长度的loop修补

序列设计和rotamer repacking

作者对50个采样序列的序列恢复率进行比较，每个序列都从原生的全原子骨架开始，没有侧链信息。结果如下图，效果还不错。


联合建模

尝试一个模型让扩散结构和序列并行。用上文描述的扩散训练和采样方法，对整个 步的结构变量 和 进行扩散。序列变量 从 扩散到 ，每步都有一个附加的网络来限制 结构分量的输出。也就是说，对给定的 ，基于来自结构网络 和 的预测对序列进行掩码预测。然后在采样的骨架和序列上运行rotamer扩散。
在下图A，作者对骨架和序列进行上下文修补，发现该模型有时几乎能够在内部修补区域的本地序列和骨架原子位置方面恢复局部解决方案。这使得模型能进行全原子环生成(下图B)，其中联合生成loop区域的loop骨架和候选序列。这种能力为有趣的工程问题开辟了道路，如免疫球蛋白(Ig)循环设计。抗体可变Ig域具有高度可变的CDR(互补决定区)环，使它们能够有选择地结合几乎任何靶点。在下图C中，作者展示了如何使用这种生成模型在一个固定的Ig骨干上共同改变CDR骨干环路和序列。



动图(8s)：


作者最后希望接下未来实现结构、序列、rotamer并行中交叉生成来提高采样质量，而不是现在按顺序生成。真是令人期待和激动，如果可行，那真是开启了AI蛋白设计新范式。
小结

作者介绍了一种蛋白结构、序列和旋转子的生成模型，该模型可以产生大量的蛋白结构域，这些结构域在物理上是合理的，并且在PDB的不同结构域类型之间有很大的差异，表明它有潜力成为一个完全end-to-end蛋白设计工具。
AF2里其实还有一些有意思的东西可以探索，作者提到，可以用本文中的扩散公式替换AlphaFold2中用于预测的监督学习“回收”过程(或者用Evoformer module的输出替换约束条件信息)。
其次还有更多的辅助能量函数可以加上来更好控制采样过程，比如那些惩罚偏离距离约束的函数。
最后期待一波源码，值得参考的地方我会在此篇笔记更新。
下篇笔记见Eva-MariaStrauch组的结合多肽生成模型。
最近Namrata做客了ML for protein engineering seminar series，带来了这篇文章的讲座视频，链接如下：
https://www.youtube.com/watch?v=i8fGzddGbU8​www.youtube.com/watch?v=i8fGzddGbU8Reference

Namrata Anand, Tudor Achim, Protein Structure and Sequence Generation with Equivariant Denoising Diffusion Probabilistic Models, https://nanand2.github.io/proteins

这篇paper也记录在了下方这个记录所有蛋白设计(protein design)和深度学习(deep learning, DL)有关文献的GitHub仓库里，各位若对AI蛋白设计感兴趣，欢迎一同维护：
https://github.com/Peldom/papers_for_protein_design_using_DL​github.com/Peldom/papers_for_protein_design_using_DL​github.com/Peldom/papers_for_protein_design_using_DL​github.com/Peldom/papers_for_protein_design_using_DL