Nat. Biomed. Eng. | 传统分子动力学模拟和新型深度生成模型加速抗菌肽的发现

方法概览

在此，研究者提供了一种有效的经电脑模拟的筛选方法，该方法使用了，经高通量物理驱动的分子模拟增强的深度学习分类器(图1)。这种用于从头抗菌设计的计算方法，明确说明了广谱效价和低毒性，并且研究者通过实验，证实了这些特性。

研究者通过筛选的20个候选序列(从近90000个生成的序列池中)的合成，从而发现了两个新的短肽(YI12和FK13)，它们对不同的病原体有较强的抗菌活性，包括一株难以治疗的耐多药革兰氏阴性肺炎克雷伯菌。众多实验表明，YI12和FK13都是有希望的治疗候选者，值得进一步研究。

图1. 方法概览

图片来源于Nat. Biomed. Eng.

肽的自编码

为了对肽潜在空间进行建模，研究者使用了基于由编码器和解码器，两个神经网络组成的深度自编码器的生成模型。用ϕ参数化的编码器qϕ(z|x)学习，将输入x映射到一个变分分布，用θ参数化的解码器pθ(x|z)的目标，是根据所学的分布，给定潜在向量z，重构输入x，如图2a所示。变分自编码器(VAE)是这一类中最流行的模型，它假设潜变量z~p(z)，遵循简单的先验分布(例如高斯分布)。然后解码器，在连续表示z的序列上生成一个分布。生成过程指定为，从而对潜在变量进行积分。另一种标准的VAE(被认为是改进的变体)是Wasserstein自动编码器(WAE)。在VAE/WAE框架内，肽生成被表述为密度建模问题，即估计p(x)，其中x是短的可变长度氨基酸串。密度估计程序，必须指定已知肽的高可能性。因此，模型泛化意味着可信的新肽，可以从具有高概率密度的区域在模型下生成。肽序列，是由20个天然氨基酸组成的文本字符串。只考虑长度≤25的序列，进行模型训练和生成，因为仅需较短的AMPs。 

图2. 属性控制肽序列生成的阶段

图片来源于Nat. Biomed. Eng.

在此，研究者在已知的肽序列上，训练了一个密度模型。事实上，无监督表示学习，是通过在一个大的文本库上进行预训练，最近在文本和语音的下游任务以及蛋白质生物学中，都有令人印象深刻的结果，从而启发了这种方法。此外，与蛋白质序列生成的类似模型相比，研究者并不局限于学习，单个蛋白质家族或特定三维折叠的密度。相反，研究者训练了一个，在所有已知的短肽序列表达在不同的有机体上的全局模型。这种全局方法，应该能够实现跨多个科的有意义的密度建模，它们之间的插值，更好地学习看似合理的肽的“语法”，并超越已知的抗菌模板探索。

在肽序列上训练WAE(而不是简单的VAE)的优势，从众多的评价指标中，得到了明显的体现。此外，还观察到，当WAE训练所有的肽序列，而不是只训练AMP序列时，产生的序列具有很高的重建精度和多样性。受自然语言过程的启发，研究者进一步分析了肽WAE的信息含量。使用“探测”的方法，已经证明编码的句子，可以保留很多语言信息。同样，研究者分析了序列之间的相似性是否可以通过潜伏z空间的编码来捕获，因为已知这些信息可以指定肽序列的生物学功能和折叠度。图3a显示WAE模型的序列相似性，与z空间的欧氏距离呈负相关关系(Pearson相关系数=−0.63)，这表明WAE本质上，捕获了肽空间内的序列关系。VAE潜在空间，未能捕捉到这种关系。 

图3. 生成式自编码器潜在空间的特征

图片来源于Nat. Biomed. Eng.

为了有条件生成新肽序列的最终目标，确保在z空间中学习到的编码，保留原始序列功能属性的可识别信息，是至关重要的。为此，研究者分析了空间是否线性可分成不同的属性，从而使用序列的z编码训练用于二元(yes/no)函数属性预测的线性分类器(图2b)。

探索由WAE所模拟的z空间发现，空间确实是线性可分的不同功能属性。使用WAE z分类器和序列级分类器，对测试数据进行属性AMP的分类预测准确率分别为87.4%0和88.0%。研究者的序列级LSTM模型，显示了88%的精确度。与文献报道的或内部训练的能够访问原始序列的分类器相比，这一比较表明了z级分类器的性能非常接近。研究者的目标是在潜在特征上，训练一个预测器，从而产生可比较的精度，它可以通过使用CLaSS(条件潜在(属性)空间采样)，直接从潜在空间进行条件采样，来自动生成新的AMP候选对象。值得注意的是，比较不同的AMP预测模型并不简单，因为训练AMP数据集大小、序列长度、AMP和非AMP的不同定义，以及其他数据管理标准，不同的方法差异很大。此外，与许多现有的预测工具相比，这里使用的z级和序列级AMP分类器，都不需要任何手动定义的特征集。

对于毒性分类任务，与报告准确率高达90%的相似序列级深度分类器(图3)相比，使用基于潜在特征训练的模型的准确率要低得多。这些结果表明，一些属性，如毒性，从学习的潜在肽表示，对其预测更具挑战性；一个可能的原因是，在训练数据中存在较高的类别不平衡。此外，研究者还通过分析，在两个遥远的训练序列之间，沿z空间的线性插值向量生成的序列，来研究潜在空间的平滑性(图3b,c)。序列相似性、功能属性(AMP和毒性类概率)，以及一些物理化学性质，包括芳香性、电荷和疏水矩(表明螺旋的两亲性)，在插值过程中平稳变化。这些结果令人欣喜，因为在更大数量的未标记数据上训练的WAE潜在空间，似乎在功能、物理化学和序列相似性方面，具有重要的结构。图3c还表明，在线性插值过程中，可以识别出与两个端点序列明显不同的序列(s)，这表明经过训练的潜在空间，具有生成新的、唯一的序列的潜力。

CLaSS的控制序列生成

对于控制生成，研究者的目标是控制一组感兴趣的二元(yes/no)属性，如抗菌功能和/或毒性。为此，研究者建议使用CLaSS。CLaSS利用在肽z空间上直接训练的属性分类器，因为它们可以捕获重要的属性信息(图3)。目标，是对指定的目标属性组合at，有条件地采样p(x|at)。利用CLaSS(图2c)对该任务进行了求解，假设属性条件密度因子为：。研究者近似地从潜在z空间的模型中采样p(z|at)，使用贝叶斯规则，并使用属性分类器建模p(at|z)(图2b, c)。由于CLaSS只使用简单的属性预测模型和来自z空间模型的拒绝采样，因此它是一种简单有效的纯正向筛选方法。与现有的可控生成方法(如BO、RL或SS生成模型)相比，它不需要对潜在空间进行复杂的优化。CLaSS易于再使用，同时具有高度的并行性，并且不需要在潜在空间中定义起点。

由于在潜在特征上训练的毒性分类器似乎较弱(图3)，抗菌功能(yes/no)被用作控制，从潜在肽空间采样的唯一条件。生成的抗菌素候选菌，在后代筛选过程中使用序列级分类器筛选毒性。值得注意的是，CLaSS并不对潜在空间执行启发式搜索，相反，它依赖于基于概率拒绝抽样的方案，来进行属性条件生成。CLaSS也不同于基于局部密度的抽样方法，因为这些方法依赖于，对标记数据进行聚类，然后通过相似度，搜索找到测试样本的聚类分配；因此，它适合于向前设计任务。与此相反，对于反向设计问题，研究者提出了CLaSS，并通过对潜在空间进行属性条件采样，然后进行解码，实现了有针对性的生成。

CLaSS生成放大器的特点

为了检查CLaSS生成的AMP序列与训练数据的同源性，研究者进行了BLAST序列相似性搜索。通过使用实际序列的比对评分期望值(E值)，来评估序列的同源性，同时使用其他比对指标，如原始比对评分、百分比身份和正匹配、比对中的差距和覆盖率，来获得序列相似性的整体感觉。E值，表示查询与来自特定大小数据库的序列之间，匹配的统计意义(也称为生物学意义)。E值越大，表示目标与查询之间的相似度仅仅是巧合的可能性就越大——也就是说，查询与目标没有同源性或相关性。

此外，研究者分析了匹配得分最高的匹配值。通常，在查询包含约2.2亿个序列的UniProt非冗余数据库时，使用≤0.001的E值来推断同源性，这里，研究者采用≤10−6的E值来表示同源性。于是，在考虑与最高比对得分的匹配时，约有14%的生成序列E值≥10，另有36%的生成序列E值为>1，表明生成序列与训练序列的相似性不显著。如果只考虑与score >20的比对，则发现E的平均值为>2，进一步说明生成的序列具有非同源性。

类似的标准，也被用于检测设计的短抗菌剂的新锐性。CLaSS生成的AMPs也更加多样化，因为唯一的(即在序列集合中只发现一次)k-mers (k = 3-6)，比训练序列或其片段更丰富。这些结果突出了，当前方法生成短长度AMP序列的能力，平均而言，相对于训练数据而言，它们是独特的，而且它们本身也是多样化的。

图4a-d中比较了训练和生成的AMPs的关键分子特征分布，如氨基酸组成、电荷、疏水性(H)和疏水性矩(μH)。CLaSS生成的AMP序列表现出明显的特征：与训练抗菌序列相比，这些序列的Arg、Leu、Ser、Gln和Cys含量更丰富，而Ala、Gly、Asp、His、Asn和Trp含量减少(图4a)。此外，研究者给出了最丰富的k-mers (k = 3,4)，分析表明在生成的和训练的AMPs中，最常见的3-和4-mers是富含Lys和leu的，在生成的序列中频率更高。

生成的AMPs的特征是：全局净正电荷和芳香性介于未标记和AMP-标记的训练序列之间，而疏水力矩与已知AMPs相当(图4b-d)。这些趋势表明，产生的抗菌素仍然是阳离子的，可以形成假定的两亲性α-螺旋，类似于大多数已知的抗菌素。有趣的是，与训练序列相比，它们也表现出较高的疏水性和脂肪族指数。这些观察突出了CLaSS生成AMP序列的独特物理化学性质，这是研究者学习范式的SS本质的结果，这可能有助于它们的治疗应用。例如，已知较低的芳香性和较高的脂肪族指数，会诱导短肽具有较好的氧化敏感性和较高的热稳定性，而较低的疏水性则与较低的毒性有关。

图4. 理化性质比较

图片来源于Nat. Biomed. Eng.

经电脑模拟的后代筛选

为了筛选大约90,000个CLaSS生成的AMP序列，研究者首先使用了一组独立的二进制(yes/no)序列水平的深层神经网络分类器，用于筛选抗菌功能、广谱有效性、二级结构和毒性(图1)。通过这次筛选的163个候选分子，随后经过处理，进行了肽-膜相互作用的粗粒度分子动力学模拟。这些模拟的计算效率，使它们成为选择高通量和物理启发肽序列过滤的一个有吸引力的选择。

由于没有使用分子模拟筛选抗菌候选药物的标准方案，因此，研究者对已知序列进行了一组具有或不具有抗菌活性的控制模拟。从这些对比运行中，发现阳性残留物和膜脂之间接触数量的差异，是活性抗菌的预测。具体来说，接触方差区分高效AMPs和非抗菌序列的敏感性为88%，特异性为63。根本上，这一特征可以解释为，测量肽序列与模型膜的强大结合趋势。因此，研究者使用接触方差截止值2，来进一步过滤，通过分类器筛选的163个产生的AMPs。

湿实验室验证

最后一组20个CLaSS生成的AMP序列，通过了上述接触方差筛选，以及它们的模拟和物理化学特性。这些序列在湿实验室中，进行了抗菌活性测试，使用最低抑菌浓度(MIC；越低越好)对革兰氏阳性金黄色葡萄球菌和革兰氏阴性大肠杆菌，进行了测量。同时，11条产生的非AMP序列，也进行了抗菌活性筛选。没有一个设计的非AMP序列，显示出足够低MIC值，从而可认为是抗菌剂，这意味着研究者的方法，不容易出现假阴性预测。研究者还推测，AMP测试集和从已标记和未标记的肽序列上训练的潜在空间生成的序列之间的域移位，可能导致假阳性预测(18/20)。

在人工智能设计的20个候选AMP中，两个序列YLRLIRYMAKMI-CONH2 (YI12, 12个氨基酸)和FPLTWLKWWKWKK-CONH2 (FK13, 13个氨基酸)的MIC值最低(表1)。这两种多肽都带正电荷，并具有非零疏水力矩，这表明它们的阳离子和两亲性，与已知的抗菌剂是一致的。研究者利用这些多肽，对更难治疗的革兰氏阴性铜绿假单胞菌和鲍曼不动杆菌，以及耐多药的革兰氏阴性肺炎K.肺炎杆菌，进行了进一步的评估。如图5所示，YI12和FK13均表现出较强的广谱抗菌活性，MIC值相当。已知，LLKKLLKKLLKK是一种现有的α-螺旋形成AMP，具有良好的抗菌活性和选择性；因此，研究者进一步比较了YI12和FK13与LLKKLLKKLLKK的MIC值。LLKKLLKKLLKK对金黄色葡萄球菌、大肠杆菌和铜绿假单胞菌的MIC值，分别为>500、>500和63。FK13和FY12对铜绿假单胞菌的MIC值与LLKKLLKKLLKK相当。然而，FK13和YI12对金黄色葡萄球菌和大肠杆菌的MIC值，明显低于LLKKLLKKLLKK，说明本研究发现的AMPs具有更高的效率。此外，研究者还报告了几种现有AMP预测方法对YI12和FK13的结果(表1)。iAMPPred和CAMP-RF预测YI12和FK13都是AMPs。另外两种方法都存在误分类，如基于局部相似度搜索的DBAASP-SP误分类FK13。

与此同时，研究者进一步进行了体外和体内毒性试验。根据在50%溶血(HC50)时测定的活性和致死剂量(LD50)毒性值(表1)，这两种多肽似乎是生物相容的(因为HC50和LD50值远高于MIC值)；FK13比YI12更具有生物相容性。重要的是，这两种多肽的LD50值优于多粘菌素B，多粘菌素B是临床上用于治疗耐药革兰氏阴性细菌感染的抗菌药物。

表1. YI12和FK13的抗菌活性和毒副作用

表格来源于Nat. Biomed. Eng.

图5. 原子模拟和位阻捕获研究

图片来源于Nat. Biomed. Eng.

序列相似性分析

为了研究YI12和FK13在训练序列方面的相似性，研究者详细分析了BLAST同源搜索返回的比对分数。评分指标包括，原始比对评分，E值，比对覆盖率的百分比，身份的百分比，正匹配或相似度的百分比，以及比对间隙的百分比。对训练数据库进行e值阈值为10的BLAST搜索，没有发现任何YI12的匹配，说明YI12不存在有统计学意义的匹配。因此，研究者在更大的UniProt数据库中，搜索YI12的相关序列，该数据库由近2.235亿个非冗余序列组成，而研究者的模型训练只包含了其中的一小部分。YI12与其最接近的匹配E值为2.9，是细菌含有EAL-结构域的蛋白的11个残基片段。

以上结果表明，即使考虑到UniProt中的所有蛋白质序列，YI12的相似性也很低。研究者还对包含约6550万专利肽的PATSEQ数据库进行了YI12的BLAST搜索，仍然获得了最小E值1.66。在PATSEQ序列中与YI12最接近的序列，是一个由79个氨基酸组成的8个氨基酸片段，相似性为87.5%，覆盖率仅为66.7%，进一步证实了YI12与已知序列的较低相似性。

与训练数据库中最接近的匹配相比，FK13的特征小于75%，在匹配上存在差距，覆盖率为85%，这意味着FK13与训练序列的相似度也很低。然而，YI12比FK13更独特。在训练数据库中与FK13最接近的匹配，是小麦胚乳中嘌呤吲哚碱A蛋白的一个13个氨基酸杀菌域(PuroA: FPVTWRWWKWWKG)的合成变体。FK13的抗菌活性和溶血活性与PuroA接近。然而，FK13与PuroA有本质上的不同；FK13富含Lys，色氨酸含量低，导致较低的整体亲水性评分(−0.854 vs −0.962)、较高的脂肪指数(60.0 vs 22.3)和较低的不稳定性指数(15.45 vs 58.30)，这些共同表明了较高的肽稳定性。

事实上，在湿法实验中，较低的色氨酸含量有利于FK13的稳定，因为色氨酸在空气中容易被氧化。较低的色氨酸含量，也涉及到改善体内肽的稳定性。综上所述，这些结果表明，WAE模型的潜在肽空间CLaSS，能够通过有效地学习肽中复杂的序列-功能关系，并将这些知识用于控制探索，从而生成独特和最佳的抗菌序列。当与后续的电脑模拟筛选相结合时，新的和最佳的候选物，具有实验证实的高广谱有效性和选择性，成功率为10%。整个周期(从数据库管理到湿实验室确认)，在一次迭代中总共花费了48天(图1)。

结构和机理分析

研究者对这两个序列，进行了全原子显式水模拟，其中存在一个从α-螺旋结构开始的脂膜。如图5a所示，这两个序列观察到了不同的膜结合机制。YI12使用带正电荷的n端Arg残基嵌入膜中，而FK13既嵌入N端Phe，也嵌入C端Trp和Lys残基。这些结果为研究YI12和FK13，在膜相互作用早期不同的作用模式，提供了机理上的见解。

用CD光谱对多肽进行了进一步的实验表征。YI12和FK13在水中均呈随机卷曲状结构，而在20% SDS缓冲液中形成α-螺旋。结构分类器预测和全原子模拟与CD结果一致。从CD谱上看，YI12的α-螺旋度似乎比FK13强，这与它更强的疏水矩一致。综上所示，理化分析和CD光谱分析表明，YI12和FK13的阳离子性质和两亲性螺旋结构，是诱导其抗菌性质的潜在因素。

为了深入了解YI12和FK13的抗菌作用机制，研究者进行了琼脂平板实验，发现YI12和FK13两种多肽都具有杀菌作用。在2× MIC时，菌落减少了99.9%。

耐药性分析

最后，研究者对亚胺培南(一种静脉注射的β-内酰胺类抗生素，YI12或FK13)，在亚-MIC浓度下的大肠杆菌进行耐药性-习得研究。从图5b的结果可以看出，YI12和FK13在传代25次后均未产生耐药性，而大肠杆菌在传代6次后，已开始对抗生素亚胺培南产生耐药性。同时，研究者还研究了这些多肽对耐多粘菌素b的肺炎克雷伯菌的疗效，这是一种耐多粘菌素B-a -内酰胺抗生素的菌株。表1显示了YI12、FK13和多粘菌素B的MIC值，表明与对耐多药肺炎克雷伯菌染色(来自ATCC)的MIC相比，这两种发现的多肽中的任何一种的MIC都没有增加。相比之下，多粘菌素B对同样多药耐药的肺炎克雷伯菌染色(来自ATCC)的MIC为2 μg ml−1，但当肺炎克雷伯菌对多粘菌素B产生耐药性时，MIC为>125 μg ml−1。在耐多黏菌素B的菌株中，YI12和FK13的MIC值仍然低于多黏菌素B，说明对YI12和FK13没有表现出对多黏菌素B的耐药性。综上所示，YI12和FK13具有治疗耐药菌株的潜力，因此需进一步研究。

展望与结论

在此，研究者提供了一个完全自动化的计算框架，该框架结合了可控生成建模、深度学习和物理驱动学习，用于重新设计广谱强效和选择性AMP序列，并通过实验验证这些序列的广谱有效性和毒性。此外，发现的肽对耐最后一种抗生素的菌株，表现出了很高的疗效，同时也减轻了耐药性的发作。湿实验室的结果证实了，该方法在设计新的和优化的序列时的效率，并且只合成和测试了非常少量的候选化合物。在这个概念验证研究中，目前的设计方法，获得了10%的成功率和48天的快速周转，强调了将人工智能驱动的计算策略与实验相结合，可以获得更有效的候选药物的重要性。这里提出的生成模型方法，不仅可以用于生成新的候选药物，而且可以用于设计独特的联合疗法和抗生素佐剂，以进一步推进抗生素治疗。

由于CLaSS是一种通用方法，它适用于各种受控生成任务，可以同时处理多个控件。同时，该方法易于实现、快速、高效和可扩展性，因为它不需要对潜在空间进行任何优化。CLaSS在可重复使用性方面，有额外的优势，因为添加新的约束，只需简单的预测器训练。因此，本研究的未来方向，将探索额外的相关约束条件，如诱导耐药性、动物感染模型的有效性和细粒度菌株特异性，对使用本文提出的方法设计的AMPs的影响。未来，可将CLaSS的应用扩展到，其他受控分子设计任务中，如靶向性和选择性药物样小分子生成，正在进行中。最后，AI模型，将在一个主动学习框架中，使用模拟和/或实验的反馈，以迭代的方式进一步优化。

参考文献

Payel Das, Tom Sercu, Kahini Wadhawan, et al., Accelerated antimicrobial discovery via deep generative models and molecular dynamics simulations. Nat. Biomed. Eng. 2021, ASAP. DOI: 10.1038/s41551-021-00689-x.