【论文标题】Protein engineering via Bayesian optimization-guided evolutionary algorithm and robotic experiments

【作者团队】Ruyun Hu, Lihao Fu, Yongcan Chen, Junyu Chen, Yu Qiao, Tong Si

【发表时间】2022/08/12

【机 构】中科院

【论文链接】https://doi.org/10.1101/2022.08.11.503535

蛋白质工程的目的是在一个巨大的设计空间中找到更好的功能序列。对于这样一个昂贵的 "黑箱 "功能优化问题，贝叶斯优化是一种效率采样方法，它由目标函数的代用模型指导。不幸的是，贝叶斯优化在巨大的搜索空间中是难以计算的。更糟糕的是，它是按顺序提出序列的，这使得它与成批的湿实验室测量不兼容。在本文报告了一种可扩展的、分批进行的方法，即贝叶斯优化引导的进化（BO-EVO）算法，在预训练模型的蛋白质适应度景观基础上探索组合突变库，并指导多轮机器人实验。本文首先研究了基于蛋白质GB1基准的景观的各种设计，然后将BO-EVO成功地推广到另一个大肠杆菌激酶PhoQ上。这种方法随后被应用于指导机器人库的创建和筛选，以设计RhlA的酶特异性，RhlA是鼠李糖脂生物表面活性剂的一种关键生物合成酶。经过4次迭代，本算法在检查了所有可能的突变体中不到1%的突变体后，在生产目标鼠李糖脂同系物方面取得了4.8倍的改进。总的来说，BO-EVO被证明是一种高效和通用的方法，可以在没有先验知识的情况下指导蛋白工程。

BO-EVO整合了具有互补优势的BO和进化算法，在评估整个蛋白质组合空间（20^N）时，BO的计算时间随着目标残基的数量（N）呈指数级增长。为了提高算法的可扩展性，BO-EVO将搜索限制在由父序列的随机突变（进化算法）产生的子空间。父序列或者通过基于适应度的抽样从测量的序列中抽出，或者根据代用模型的置信度设置为一个新提出的序列。另一方面，与贪婪的进化算法相比，BO允许加强探索。本文使用GPR（高斯过程回归）作为代用模型来量化不确定性。GPR是在由蛋白质序列表征法，如预训练模型ESM-1v编码的序列上训练的。为了确定实验中突变体序列的优先次序，本文利用UCB作为获取函数来平衡探索和利用。一旦达到每轮设定的预算（即24、96、384个序列），就会在实验中创建一个突变体批次并进行分析。由此产生的序列功能数据被用来完善代用模型，以指导进一步的下回合的蛋白工程。

上图展示了在benchmark GB1上的测试，展示了所有推出的变体在100次模拟中的最大（红色）和平均（绿色）适应度。

第一个实验（左图）将迭代数固定为4，这样蛋白质工程可以在一个月左右完成。作者观察到，最大适应度随着批次大小的增加而增加，当批次大小等于或大于384，或总样本预算等于或大于1,536时，可以持续达到全局最优。有趣的是，所有推出的变体的平均适应度先是增加，然后随着batch的增大而减少，这可能是由于序列空间中高适应度的变体数量有限（只有2.31%的变体表现出高于WT的适应度），以及BO-EVO在batch较大的情况下加强了对不适应度区域的探索。

对于第二个实验（右图），本文将总的样本预算固定为1,536（占整个设计空间160,000的0.096%），这是在BO-EVO四次迭代中持续达到全局最优所需的最小样本数。可以注意到，随着batch大小的增加和迭代次数的减少，1,536个建议的突变体的最大和平均适应度都在恶化。这些结果表明，算法-实验反馈对于提高BO-EVO在序列设计上的性能是必要的。总的来说，本文由此决定继续使用适度的迭代数（4-5）和batch大小（384）进行BO-EVO实验。

上图展示了各种序列编码方法的评价。

首先，a图展示了编码策略对代用模型质量的影响，以Pearson相关系数为评价标准。本文首先从GB1景观中抽取了384个样本，然后用剩余的数据建立了一个新的景观。在这个新的实验景观上进行了模拟，并在保留的测试集上评估了代用模型。本文发现，对于所有的编码策略来说，随着更多轮的迭代，模型的性能不断提高，这表明模型的改进是由与景观数据的相互作用造成的。另一方面，内容丰富的编码并不能保证有好的模型。与Onehot编码相比，使用UniRep、ESM-1v和ProtT5-XL-U50等学习型编码是可行的，而使用TAPE或Georgiev编码是不可行的。

然后，b图研究了序列编码如何影响BO-EVO的性能。当检查所有建议的序列在第五次迭代后的平均适应度时，观察到类似于a图的趋势：ProtT5-XL-U50 ~ ESM-1v ~ UniRep > Onehot > TAPE > Georgiev编码，表明代用模型的性能与BO-EVO的性能之间有很强的相关性。然而，当c图中评估特定迭代后的成功率时，趋势变得有些不同。在前三轮中，简单的Onehot编码的表现超过了所有其他策略。在第五轮迭代后，性能排名是ProtT5-XL-U50 > ESM-1v ~ Onehot > UniRep > TAPE > Georgiev，这又与代用模型的情况很相关。总的来说，ProtT5-XL-U50和ESM-1v策略是首要选择。在随后的实验中，本文继续使用ESM-1v编码，因为与Prot5-XL-U50相比，它的模型尺寸更小，计算效率更高。值得注意的是，更复杂的编码充其量只比简单的Onehot编码表现好一点，这与最近的其他研究一致。

上图展示了优化算法的评价。

本文将BO-EVO与纯进化算法（AdaLead）和纯BO算法进行了基准测试，以考察结合这两种探索策略的必要性；随机突变（Random）也作为基线进行了评估。正如预期的那样，随机算法的表现是四种算法中最差的。从轮次成功率（a图）和五轮后所有建议序列的最大和平均适应度（b图）来看，AdaLead通过对具有高预测适应度的变体的贪婪利用，比Random好，但这种纯进化算法比BO-EVO和BO差很多，这两种算法都利用了UCB获取函数。这些结果表明，在探索崎岖不平的适应度景观时，考虑变体适应度和模型不确定性是很重要的。另一方面，通过在每个迭代中对整个设计空间（16万个序列）的详细探索，纯BO取得了比BO-EVO更好的性能，后者每轮只评估了3072个序列（整个设计空间的1.92%）。虽然BO-EVO的性能不如纯BO，但在探索组合突变景观时，BO-EVO的计算时间几乎是恒定的，而纯BO的计算时间随着目标残基数的增加呈指数级增长，并迅速变得难以解决。

上图展示了BO-EVO指导的酶工程结果。

RhlA是一种关键的酶，负责合成鼠李糖脂（RLs）的脂质分子，这是一种重要的生物表面活性剂。RhlA的酶的特异性决定了脂质分子的化学结构（a图），这进一步影响了相应RL分子的生理化学和生物活性。然而，事实证明很难通过（半）理性设计或定向进化来改变RhlA的酶的特异性。为了将BO-EVO应用于RhlA工程，本文开发了一个机器人流程来建立和测试RhlA的4残基组合SSM库中16万个成员中的任何384个。为了测量酶的特异性，作者应用MALDI-ToF-MS来量化两个RL产品Rha-（C8-C10）和Rha-C10-C10（b图）。在微孔板培养中，携带WT型RhlA的大肠杆菌细胞产生Rha-C10-C10作为主要产物，故目标是将RhlA突变体的产物特异性转向较小的Rha-（C8-C10）。

为了应用BO-EVO，本文选择了Arg74、Ala101、Leu148和Ser173（RALS）作为组合突变的四个目标残基，因为这些配体结合残基上的许多突变都能大幅提高Rha-（C8-C10）的产量。这4个残基的单残基SSM数据被用于 "热启动"，这对于提高结果是可取的。对于BO-EVO的迭代，本文观察到累积最大适应度逐轮增加，在第4轮达到7.35（AACA序列）的最高值（c图）。最终成功地确定了一个RhlA突变体，其RL-(C8-C10)的产量相对于WT有4.8倍的提高。

 

创新点

  • 本文开发了一种可扩展的、分批进行的BO算法，即BO-EVO，用于指导多轮机器人实验探索组合蛋白的适应度景观。对于4位点组合突变，BO-EVO的总实验预算减少到1,536个突变体，或不到理论库大小16万的1%，其中执行了实验室测量、模型完善和突变体设计的4次迭代，批次大小适中，为384。在这些样本预算内，BO-EVO取得了不错的成功率（75%），在一系列模拟的和实验的蛋白质适应度景观上达到了全局最优值。
  • 本文进行了一个真正的蛋白质工程任务，修改RhlA的产物特异性。BO-EVO在一个月内迅速设计出大幅改进的突变体。而据作者所知，这是第一个通过ML模型和机器人实验之间的多轮反馈进行算法指导的自动蛋白质工程的报告。
  • 对比其他方法，BO-EVO不需要事先了解目标蛋白的结构或同源序列，也不需要通过进化算法限制序列设计来评估整个虚拟设计空间。

 

 

内容中包含的图片若涉及版权问题，请及时与我们联系删除