针对SARS-CoV-2刺突糖蛋白激活变构抑制的分子数据库筛选集成计算流程

《Journal of Biomolecular Structure and Dynamics》：Integrated computational protocol for sampling molecular databases towards allosteric inhibition of SARS-CoV-2 spike activation

　　

SARS-CoV-2刺突（S）糖蛋白受体结合结构域（RBD）的开放是病毒进入所需的关键构象变化。本工作中，研究人员采用集成计算方法来评估小分子如何变构调控该过程。从初始约60,000个ZINC化合物的库中，通过物理化学、吸收和空间标准筛选出739个候选分子，经聚类后选取9个代表物，研究其对蛋白质激活（开放）的变构效应。经典分子动力学（MD）模拟揭示了配体稳定性的显著差异，而牵引（非平衡）模拟量化了各化合物对RBD开放运动的影响。通过伞形采样构建的自由能剖面，以及基于主成分分析（PCA）推导的构象群体密度图表明：小尺寸配体抑制效果有限，中等尺寸配体产生最一致的稳定效应，而两个最重配体则通过稳定闭合或中间构象显著重塑了激活路径。总体而言，所提出的方法为分析小分子对病毒融合蛋白大尺度构象变化的影响提供了实用的筛选框架，并允许在纳秒时间尺度上评估反应能垒。

一、研究背景

新型冠状病毒SARS-CoV-2自2019年底出现以来引发了全球COVID-19大流行。病毒通过其表面的刺突（S）糖蛋白与宿主细胞受体ACE2结合，介导膜融合，这是感染的关键步骤。S蛋白以三聚体形式存在，其受体结合结构域（RBD）可在“开放”（受体可及）和“闭合”（受体不可及）构象之间动态转换，RBD的开放是病毒进入的必要条件。因此，靶向S蛋白构象变化的抑制剂，特别是通过变构机制稳定其闭合构象，成为抗病毒药物研发的重要策略。虽然已有大量计算和实验研究探索S蛋白的结构与功能，但系统评估小分子变构调控其激活过程的计算筛选流程仍不完善，尤其是在纳秒尺度量化能垒变化方面存在挑战。

二、研究目的与意义

本研究旨在开发一个集成计算协议，用于从大规模分子数据库中筛选能变构抑制SARS-CoV-2 S蛋白RBD开放的小分子。通过结合分子对接、多步筛选、分子动力学模拟及自由能计算，研究人员希望评估不同尺寸和性质的配体对激活路径和能量景观的影响，为设计新型变构抑制剂提供理论基础和方法框架。论文发表在《Journal of Biomolecular Structure and Dynamics》上，其方法学具有通用性，可推广至其他病毒融合蛋白的构象动力学研究。

三、关键技术方法

研究人员采用以下主要计算技术：

1. 1.
   数据集准备：从ZINC15数据库获取60,411个“in vivo”标记化合物，进行结构校验、电荷校正（使用Gasteiger、MMFF94、Mulliken等方法），并使用Open Babel和Gaussian 16软件处理。
2. 2.
   分子对接与筛选：针对S蛋白闭合构象（PDB ID: 6vxx），使用AutoDock 4.2.6进行半柔性分子对接。对接后应用Lipinski五规则、BOILED-Egg模型（用于类药性评估）和空间邻近性筛选（配体几何中心需在S1/S2切割位点Asn⁵³⁶和Glu⁶¹⁹15 ?内，且至少有一个原子在5 ?内），得到739个候选物。进一步通过结构聚类（ChemmineR库，相似性阈值0.5）和分子量分箱（A: ≤250 Da, B: 250-500 Da, C: 500-750 Da），从每箱中选取对接得分最优、最差和中位数的代表物各一，最终得到9个化合物（A1-A3, B1-B3, C1-C3）进行后续模拟。
3. 3.
   分子动力学模拟：使用GROMACS 2023.2和CHARMM36力场进行全原子MD模拟。包括：
   • •
     平衡MD模拟：对S蛋白开/闭构象（6vyb/6vxx）进行5×50 ns模拟以评估构象稳定性；对9个配体-蛋白复合物进行3×30 ns模拟，评估配体在结合位点的稳定性（通过配体重原子RMSD、氢键分析和短程非键相互作用能）。
   • •
     牵引分子动力学（SMD）模拟：对apo（无配体）和配体结合系统进行闭合→开放的5 ns非平衡牵引模拟，施加谐和势（牵引速率10 nm/ns，力常数50 kJ/mol/nm²），以诱导RBD开放。
   • •
     伞形采样（US）：沿RBD开放反应坐标，在窗口间隔0.1 nm处进行谐和约束模拟（力常数50 kJ/mol/nm²），每个窗口进行3 ns NPT模拟（含0.3 ns平衡）。
4. 4.
   构象与能量分析：
   • •
     主成分分析（PCA）：对SMD轨迹进行PCA，以前两个主成分（PC1: RBD开放，PC2: RBD倾斜）构建构象群体密度图，分析配体对激活路径的影响。
   • •
     平均力势（PMF）计算：使用加权直方图分析方法（WHAM）处理伞形采样数据，计算沿反应坐标的PMF，量化配体对开放能垒的影响。

四、研究结果

分子对接和数据集筛选

分子对接显示，对接得分（结合亲和力预测）随分子量增加而改善，在C箱达到峰值。应用Lipinski规则、BOILED-Egg模型和空间距离筛选后，从60,407个初始化合物中选出739个（1.2%），并聚类为10个簇（I-IX为显著簇，X为小簇）。最终选取的9个代表物化学结构多样，涵盖不同分子量范围和对接得分。

平衡分子动力学

• •
  无配体的S蛋白在开/闭构象下均保持稳定（碳骨架RMSD约0.4 nm）。
• •
  配体稳定性分析：A箱配体（A1-A3）表现出不稳定，A1在所有模拟中均从结合口袋解离，被排除出后续分析；A2和A3也显示较高的配体RMSD和不稳定的氢键模式，且A2因带负电荷净电荷而产生静电排斥。B箱和C箱配体表现出较好的稳定性。B1、B2、C2和C3的配体RMSD值较低（0.4-0.8 nm），并形成稳定的氢键网络（主要与S1/S2位点附近的残基如Asp⁶¹⁴、Asn⁶¹⁶、Glu⁶¹⁹等相互作用）。非键相互作用能计算表明，B1、B2、C2和C3与蛋白的相互作用最强（负值最大）。

牵引（非平衡）分子动力学

牵引过程中，A2和A3配体与目标结合位点（Asn⁵³⁶-Glu⁶¹⁹）的距离增加，表明它们在平衡阶段已发生位移。B2和C3的位移最小，B1、C1、C2等次之。

构象群体密度与自由能分析

• •
  PCA群体密度图：无配体系统的SMD轨迹在PC1-PC2空间呈现连续分布，PC1对应开放运动，PC2对应倾斜运动。在PC1中段（参考构象20-25）存在低群体区域，表明快速过渡。
• •
  配体对构象景观的影响：A2和A3对RBD开放的阻碍作用有限，系统仍能达到开放构象。B箱配体影响更为显著：B1允许系统达到晶体结构开放构象，B2和B3可达到两种开放构象。B3表现出开放与倾斜的强耦合。C箱配体中，C1倾向于占据MD开放构象而非晶体结构构象；C2展现出非典型的开放路径；C3与B2类似，稳定了另一种闭合构象。总体而言，B箱和C箱的较重配体倾向于稳定闭合或中间构象，或改变开放路径。
• •
  平均力势（PMF）能垒：无配体系统的RBD开放PMF能垒约为15±5 kJ/mol。所有测试配体均增加了开放能垒：A2和A3分别增至约30±5和40±5 kJ/mol；B1的能垒最高（46±5 kJ/mol），约为无配体系统的三倍；B2和B3约为36-37 kJ/mol；C1和C3约为30-35 kJ/mol；C2的PMF剖面最接近无配体系统，但整体上移约20 kJ。能垒增加与配体在变构结合位点的稳定性相关，而非结构细节差异。这些能垒值与文献报道的S蛋白开放自由能变化范围（5-45 kJ/mol）一致。

五、讨论与结论

讨论

本研究建立的集成计算策略为探索大病毒蛋白构象转变的变构抑制提供了通用框架。该方法成功结合了大规模虚拟筛选、分子对接、药代动力学标准、平衡/非平衡MD模拟和构象景观分析。空间距离标准是筛选过程中的一个高效参数，类药性标准可用于进一步优选。

然而，该方法存在一定局限性：单反应坐标的适用性主要限于经历一维大尺度构象变化的体系；对于更复杂的构象动力学，可能需要多维自由能描述或其他增强采样策略。此外，PMF能量剖面是沿选定反应坐标的有效自由能变化，由于牵引过程的非平衡性质，不能直接定量预测生理条件下的病毒进入动力学。但本研究发现的无配体系统开放能垒与已发表数据一致，证明了该框架在纳秒尺度上比较变构效应对激活路径和自由能能垒的可行性。

结论翻译

本研究提出了一种计算策略，通过集成大规模虚拟筛选、分子对接、药代动力学标准、平衡和非平衡MD模拟以及构象景观分析，来探索SARS-CoV-2刺突糖蛋白激活的潜在变构抑制剂。从总共60,407个配体中筛选出9个配体，详细研究了它们对RBD激活的变构抑制。筛选基于独特性、类药特性（Lipinski五规则和BOILED-Egg模型）、结合亲和力（对接得分）以及与变构抑制位点的距离。距离标准似乎是一个高效参数，而类药性是选择具有所需性质配体的附加选项。通过平衡MD模拟验证对接构象，并通过非平衡SMD模拟分析这九种选定化合物对激活能量学和构象变化的影响，以探索其变构抑制作用。

研究发现，开/闭构象下的无配体系在平衡MD模拟中是稳定的，而非平衡SMD模拟代表了从闭合到开放样构象的连续转变。闭合构象下配体体系的平衡MD模拟显示，A箱配体稳定性有限（配体A1不稳定，被排除在进一步考虑之外），C箱化合物稳定性中等，而B箱配体则表现出持续强而持久的相互作用。

通过SMD模拟和从伞形采样推导的PMF（作为激活阻力增加的度量）以及描述构象景观变化的PCA群体直方图，探索了变构抑制对病毒感染性或进入效率的潜在降低。PCA构象图的变化表明，沿开放路径的闭合或中间构象稳定性增加（A3、B1-B3），或开放路径显著改变（C1-C3）。从能量角度看，激活与高达无配体系三倍的PMF需求增加相关。由于构象转变速率对激活自由能呈指数依赖，PMF能垒的观察变化可在生理条件下转化为数量级的转变概率降低。总的来说，这些效应导致闭合构象稳定性增加、向开放状态转变的PMF能垒（有效时间）增加，并可能影响最终的开放状态。

更广泛地说，本文建立的方法学为利用计算方法探索大病毒蛋白构象转变的变构抑制（涉及能量学和构象变化）提供了一个可推广的框架。