来源论文: https://arxiv.org/abs/2603.26411v1 生成时间: Mar 31, 2026 03:48

0. 执行摘要

预测蛋白质突变,特别是远离结合位点的远端突变如何影响药物结合,是药物发现领域长期面临的重大挑战。传统的分子动力学(MD)模拟虽然能提供构象采样的动态信息,但由于缺乏对电子的显式建模,其在揭示药物-蛋白质相互作用的精细电子结构变化方面存在局限性。为弥补这一不足,本研究提出了一种创新的耦合模拟工作流,旨在结合长时程分子动力学(MD)与高通量量子力学(QM)分析,以揭示 HIV-1 蛋白酶中远端突变导致的达芦那韦(Darunavir)耐药性的电子结构特征。

该工作流充分利用了现代异构超级计算机的优势:通过图形处理器(GPU)加速的 MD 模拟生成丰富的构象集合,并在此基础上,在中央处理器(CPU)节点上并行执行基于线性标度密度泛函理论(DFT)的量子力学计算。计算结果通过量子力学复杂度降低(QM-CR)框架进行后处理,从而量化并可视化药物片段与蛋白质残基之间的静电相互作用能(Eel)和片段键序(FBO)。

研究以 HIV-1 蛋白酶及其三种变体(野生型、2 个活性位点突变、11 个包含远端突变的复杂变体)与达芦那韦的结合为案例。核心发现是,在包含 11 个突变的变体中,远端突变导致了药物-蛋白质相互作用的系统性减弱。具体表现为总静电相互作用能的增加(吸引力减弱或排斥力增强)和总片段键序的减少(电子耦合减弱),这与实验观察到的结合亲和力大幅下降高度一致。分析进一步揭示,药物骨架(APC 片段)与催化区域的相互作用减弱,以及苯胺(ANL)片段与特定残基的短程电子耦合降低是主要的贡献因素。这些结果突出了构象采样和量子洞察在理解远端突变效应中的关键作用。

这项研究不仅提供了一种可扩展的计算策略,用于大规模、原子级地研究药物耐药性的复杂生物物理机制,而且为设计能够有效应对突变诱导失稳的新型抑制剂奠定了基础。通过将药物发现范式从静态结构分析转向动态电子耦合映射,本工作为开发更具鲁棒性的药物提供了新的方向。

1. 核心科学问题、理论基础、技术难点与方法细节

1.1 核心科学问题

在药物发现和设计领域,预测蛋白质突变如何影响药物结合亲和力,尤其是当这些突变远离药物结合位点时,是一个长期存在的重大挑战。传统的药物设计往往侧重于结合位点的局部结构互补性,但越来越多的证据表明,远端突变能够通过变构效应或构象重塑等间接机制,显著影响药物的结合力,甚至导致耐药性。这些远端效应通常不涉及结合位点的直接结构改变或明显的空间位阻,因此难以用静态结构模型或基于几何形状的经典分子动力学(MD)模拟准确预测。

以 HIV-1 蛋白酶(HIV-1 PR)为例,抗逆转录病毒药物达芦那韦(Darunavir, DRV)旨在通过与保守残基结合来降低病毒产生耐药性的几率。然而,临床观察发现,HIV-1 病毒可以通过积累多个远端突变来产生对 DRV 的高水平耐药性。Schiffer 团队等人的研究 [6] 发现,某些 HIV-1 PR 变体即使只有少量活性位点突变,其抑制效力也仅中度降低,但当远端突变与更广泛的遗传背景共同发生时,药物结合亲和力会急剧下降,跨越五个数量级。这些远端突变并未直接破坏抑制剂-蛋白质的接触,而是通过改变蛋白酶的构象集合来间接调节抑制剂结合。理解这些长程、非加性(上位性)机制的分子基础,特别是它们如何通过电子结构的变化影响药物结合,是开发更有效、更具耐药性鲁棒性的抗 HIV-1 药物的关键。本研究的核心科学问题正是:如何通过结合先进的 MD 和 QM 方法,在原子和电子层面揭示远端突变对 HIV-1 蛋白酶-达芦那韦复合物电子结构和相互作用模式的精确影响,从而提供设计新型抑制剂所需的深入洞察。

1.2 理论基础

本研究工作流的理论基础主要建立在分子动力学(MD)模拟和密度泛函理论(DFT)的基础上,并通过一个定制的量子力学复杂度降低(QM-CR)框架进行结果分析。

1.2.1 分子动力学 (MD)

MD 是一种基于经典力场的计算机模拟方法,用于模拟原子和分子系统的物理运动。它通过求解牛顿运动方程,在给定温度和压力的条件下,追踪系统中每个原子的运动轨迹。MD 在本研究中扮演了关键角色,用于:

  • 构象采样: 蛋白质-配体结合是一个高度动态的过程,涉及柔性基团在有限温度下的运动。MD 能够有效地探索蛋白质-配体复合物的构象集合,捕获其动态行为和构象变化。这些构象集合是后续 QM 分析的基础,确保分析结果反映真实的热力学涨落。
  • 环境描述: MD 模拟提供了蛋白质-配体复合物在水溶液环境中的动态描述,包括显式溶剂分子和离子的影响,这对于准确模拟生物大分子系统至关重要。

然而,经典的 MD 模拟依赖于经验力场,这些力场将原子建模为“球棒”结构,缺乏对电子的显式描述。这意味着 MD 无法直接提供关于化学键形成、电荷转移或极化等电子结构效应的详细信息,而这些效应对于理解药物结合的本质至关重要。

1.2.2 密度泛函理论 (DFT)

DFT 是一种广泛应用的量子力学方法,能够精确地计算多电子体系的电子结构和能量。与传统基于波函数的方法不同,DFT 通过电子密度来描述体系,显著降低了计算复杂度,使其能够应用于相对较大的分子系统。在本研究中,DFT 的应用至关重要,因为它能够:

  • 显式电子建模: DFT 能够明确地计算电子波函数和密度,从而提供详细的电子结构信息,如原子电荷、偶极矩、分子轨道等。这对于理解药物结合界面上的精确电子相互作用至关重要,尤其是在经典力场可能失效的情况下。
  • 相互作用能量计算: DFT 可以提供比经典力场更准确的分子间相互作用能量,包括静电、交换-相关和色散相互作用。这使得研究人员能够量化突变对结合亲和力的电子层面影响。
  • 克服 MD 局限: 通过对 MD 轨迹中的关键构象进行 DFT 计算,可以弥补经典力场在描述电子结构相关现象(如电荷重新分布、化学键合强度变化)方面的不足。

1.2.3 量子力学复杂度降低 (QM-CR) 框架

为了从复杂的 QM 计算结果中提取有意义的化学洞察力,本研究采用了其课题组提出的 QM-CR 框架 [23, 25, 26]。QM-CR 框架通过对密度矩阵的分析,将大的蛋白质-配体复合物分解为多个化学上可理解的片段,并量化这些片段之间的相互作用。其核心概念包括:

  • 片段纯度指标 (|Π|) [Appendix, Eq. 4]: |Π| 量化了给定原子集合作为独立电子子系统的程度。它基于片段投影密度矩阵的非幂等性来衡量,本质上是片段电子布居方差的负归一化值 [Appendix, Eq. 26]。较小的 |Π| 值表示片段行为更独立,可以可靠地视为一个单元,从而确保片段化是化学上有意义的。本研究利用该指标对达芦那韦进行可靠的片段化,并将其与氨基酸序列分解的蛋白质片段相结合,实现了统一的解析度。
  • 片段键序 (FBO) [Appendix, Eq. 5]: FBO 量化了两个片段之间共享的密度矩阵部分,直接衡量了它们之间的电子耦合强度。它推广了 Mayer 键序概念到片段层面,提供了一种短程化学键合强度的度量。FBO 的减少通常意味着电子耦合的削弱。
  • 静电相互作用能 (Eel) [Appendix, Eq. 28, 45, 46]: Eel 量化了两个片段之间基于其电荷分布的多极展开的静电相互作用。它描述了长程相互作用,可以根据电荷中心和相对取向而具有吸引性或排斥性。Eel 的变化可以指示电荷分布的重新分配或静电匹配的变化。
  • 相互作用图: QM-CR 框架能够生成药物片段与蛋白质氨基酸之间的相互作用图,将相互作用分类为短程化学键合和长程静电相互作用,从而直观地展示结合界面的相互作用网络。

1.3 技术难点

将长时程 MD 模拟与高精度 QM 计算相结合,特别是应用于大规模生物分子体系,带来了多方面的技术挑战:

  • 计算成本与可扩展性: 高精度 QM 方法(如 DFT)对计算资源需求巨大,即使是中等大小的蛋白质-配体复合物(例如本研究中约 5200 个原子)也需要高性能计算。线性标度 DFT 方法 [19] 虽然能有效降低计算复杂度,但要处理 MD 模拟产生的大量构象(本研究中每种体系 40 帧),仍然需要大规模并行计算资源。如何有效地利用这些资源,并在有限的时间内完成大量计算是关键。
  • 异构超算架构的利用: 现代超级计算机越来越多地采用异构架构,即融合了中央处理器(CPU)和图形处理器(GPU)。MD 模拟通常受益于 GPU 的高度并行计算能力,而许多 QM 代码(包括本研究中使用的 BigDFT)则在 CPU 上表现更佳。如何在一个统一的工作流中无缝地协调和调度在不同架构上运行的软件,最大化两种资源的利用率,同时避免数据传输瓶颈和调度复杂性,是一个技术挑战。本研究的解决方案是在 GPU 分区运行 MD,在 CPU 分区运行 QM。
  • 工作流的自动化与编排: 传统上,MD 和 QM 计算是分步进行的,MD 运行完成后才能进行 QM 分析,中间涉及大量数据传输、格式转换和手动脚本。这种“事后”分析方式效率低下,且难以扩展到大规模研究。为了实现高通量、连续性的分析,需要开发一个高度自动化的工作流管理系统,能够在 MD 模拟进行中“即时”(in-operando)或“动态”(on-the-fly)地提取快照,触发 QM 计算,并进行后处理,从而实现“数据流导向”的分析范式。
  • 数据量与数据管理: 长时程 MD 轨迹和大量 QM 计算结果会产生庞大的数据集。如何有效地存储、管理和分析这些数据,特别是从 QM-CR 框架中提取有意义的相互作用数据,并将其转化为可解释的化学洞察,也需要精巧的设计。
  • 准确性与效率的平衡: 在选择 QM 方法和计算参数时,需要在计算成本和结果准确性之间取得平衡。例如,PBE 泛函 [20] 在描述电子密度相关相互作用时表现良好,但可能不如更昂贵的泛函在结合自由能的绝对值上精确。如何确保所选方法的组合足以回答科学问题,同时保持计算可行性,是方法论上的考量。

1.4 方法细节

本研究采用了一个精心设计的耦合计算工作流,将 MD 模拟、QM 计算和 QM-CR 分析无缝结合起来,并在异构超算平台上高效执行。

1.4.1 HIV-1 PR 数据集准备

研究对象是 HIV-1 蛋白酶与达芦那韦(DRV)的复合物,共分为三种变体 [Page 3]:

  • 野生型(Wild Type, WT): 基于 PDB 结构 6OPS。
  • 2 个突变体(6OPT): 包含活性位点突变 V82F 和 I84V。基于 PDB 结构 6OPT。
  • 11 个突变体(6OPZ): 包含活性位点(V32I, V82F, I84V)、近端/翻盖(L33F, K45I, M46I, I54L, L76V)和远端(I13V, G16E, A71V)突变。基于 PDB 结构 6OPZ。论文为了与 Henes 等人研究 [6] 保持一致,将所有不在活性位点内的突变都定义为远端突变,不再细分“近端”类别。

所有结构均使用 CHARMM-GUI 溶液构建工具 [8-10] 进行 MD 模拟准备:

  • 水分子与溶剂环境: 排除晶体水,体系在显式 0.15 M NaCl 溶液中进行溶剂化。
  • 模拟盒子: 盒子尺寸设置为 12 Å 的边缘距离,以匹配之前的研究,导致 WT、2 突变和 11 突变变体的盒子尺寸分别为 79 Å、81 Å 和 80 Å。
  • 蛋白酶形式: 模拟 HIV-1 PR 的二聚体形式。
  • 质子化状态: 催化三联体中一个链的 D25 残基建模为去质子化状态 [11],与中子散射研究 [12] 的发现一致。

1.4.2 MD 模拟

MD 模拟部分使用 GENESIS 软件包 [13] 进行 [Page 3]:

  • 力场: 蛋白质使用 CHARMM36m 力场 [14],配体使用 CHARMM 通用力场(CGenFF)[15, 16],水分子使用 TIP3P 模型 [17]。
  • 模拟阶段:
    • 优化: 40,000 步最陡下降能量最小化,以移除体系中的结构应力。
    • 平衡: 250,000 步(总计 250 ps)在 303.15 K 的 NVT 恒温恒容系综下进行,采用 1 fs 的时间步长和 Velocity Verlet 算法。
    • 生产: 50,000,000 步(总计 80 ns)在 303.15 K 的 NPT 恒温恒压系综下进行,采用 2 fs 的时间步长和 RESPA 算法。
  • 算法细节: SHAKE 算法用于约束涉及氢原子的键长,截断距离为 12.0 Å。
  • 数据输出: 在生产运行期间,每 1,000,000 步(即每 1 ns)写入一次中间轨迹快照,共生成 40 个快照用于 QM 后处理。

1.4.3 QM 计算与复杂度降低 (QM-CR)

QM 计算部分使用 BigDFT 软件包进行 [Page 4]:

  • 软件: BigDFT [18],一个基于 Daubechies 小波基组的 DFT 代码。
  • 体系提取: 从 MD 快照中提取蛋白质-配体复合物,并包含复合物原子 3 Å 范围内的所有水分子。这导致了大约 5200 个原子的系统规模,可以通过 BigDFT 的线性标度模式 [19] 进行高效处理。
  • 电荷处理: QM 计算中未明确包含中和离子,体系的净电荷根据质子化状态设置。
  • 泛函与赝势: 采用 PBE 泛函 [20],网格间距为 0.45 Bohr,HGH 赝势 [21](包含非线性核修正项)。
  • QM-CR 后处理: QM 计算完成后,数据通过 QM-CR 框架 [23, 25, 26] 进行后处理 [Page 4]。
    • 药物分子片段化: 首先将目标药物分子分解为一系列化学上有意义的片段。这基于片段纯度指标 (|Π|) 进行优化,以确保片段化是可靠的。
    • 相互作用图生成: 框架能够生成药物片段与蛋白质氨基酸之间的相互作用图,并将相互作用分为短程化学键合(通过片段键序 FBO 量化)和长程静电相互作用(通过静电相互作用能 Eel 量化)。
    • 数学细节: Appendix 中详细介绍了 QM-CR 方案的数学细节,包括片段投影密度矩阵、纯度指标、片段键序及其代数关系,以及静电相互作用和电子密度分解方法。这些描述符允许在原子分辨率上分析突变诱导的相互作用模式。

1.4.4 多异构超算平台与工作流编排

本研究的工作流在一个异构超算平台——日本东京大学的 Wisteria/BDEC-0 超算上实现 [Page 4]:

  • 平台组成: Wisteria 包含两个主要分区——Aquarius(配备 NVIDIA A100 GPU)和 Odyssey(配备 A64FX CPU),两者通过共享文件系统连接,并支持 MPI 跨分区通信。
  • 资源分配:
    • MD 计算: 使用 Aquarius 分区的一个完整节点(8 块 A100 GPU)运行 GENESIS MD 模拟。
    • QM 计算: 每个 BigDFT 计算在 Odyssey 分区的 64 个 CPU 节点上运行。
    • 后处理: QM-CR 后处理在 Odyssey 分区的一个 CPU 节点上进行。
  • 工作流编排: 整个工作流通过 remotemanager [37] 软件包进行编排 [Page 5, Figure 2]。
    • 管理进程: 一个管理进程在用户的本地计算机上运行,负责协调 MD、QM 和 QM-CR 的耦合。
    • MD 任务提交: GENESIS MD 任务被提交,并以固定频率写入新的轨迹文件。
    • 实时快照提取与 QM 提交: 一个子例程周期性地使用 MDAnalysis [38] 检查新的 MD 帧。一旦有新的位置数据可用,就会生成并提交 BigDFT QM 计算到 CPU 分区。这种“即时”处理避免了等待整个 MD 模拟完成后再进行 QM 分析的传统模式。
    • QM 结果序列化: BigDFT 计算完成后,一个序列化步骤被提交,使用 QM-CR 框架将数据后处理成便于分析的格式。
    • 异步与并行: remotemanager 支持任务的异步提交,一旦单个任务完成,其结果即可用于后续分析。这使得 MD、QM 和 QM-CR 任务可以高度并行地进行,特别是多个 BigDFT 计算可以同时运行,显著提高了整体吞吐量。

这种设计允许了持续性的分析,而无需中断生产运行或在集群之间传输大量数据集,从而实现了数据流导向的工作流,特别适合于大规模和复杂的生物分子模拟研究。

2. 关键 benchmark 体系、计算所得数据与性能数据

2.1 基准体系

本研究的核心基准体系是人免疫缺陷病毒 1 型(HIV-1)蛋白酶与其抑制剂达芦那韦(Darunavir, DRV)的复合物。HIV-1 蛋白酶是病毒复制的关键酶,其活性位点突变是导致药物耐药性的常见原因。然而,本研究特别关注远离活性位点的“远端”突变如何间接影响药物结合。

研究中考虑了 HIV-1 蛋白酶的三个变体 [Page 3]:

  1. 野生型(Wild Type, WT): 这是未经突变的蛋白酶,作为所有比较的基准。其晶体结构信息来源于 PDB 数据库条目 6OPS。
  2. 2 个活性位点突变体(6OPT): 该变体包含两个位于活性位点附近的突变:V82F 和 I84V。这些突变通常被认为是直接影响药物结合的突变。其晶体结构信息来源于 PDB 数据库条目 6OPT。
  3. 11 个突变体(6OPZ): 这是最复杂的变体,包含 11 个突变,旨在模拟高度耐药的临床相关 HIV-1 毒株。这些突变分为三类:
    • 活性位点突变: V32I, V82F, I84V。这三个突变直接改变了酶的活性口袋的形状和特性。
    • 近端/翻盖突变: L33F, K45I, M46I, I54L, L76V。这些突变靠近活性位点或位于构成酶活性口袋的“翻盖”区域,可能通过构象变化间接影响结合。
    • 远端突变: I13V, G16E, A71V。这些突变远离活性位点,其作用机制被认为是长程变构效应,而非直接空间阻碍。本研究的重点之一是理解这些远端突变如何通过电子结构变化引起结合力下降。

所有蛋白酶都以其二聚体形式进行模拟,这反映了其在生物学上的活性状态。体系的准备过程包括使用 CHARMM-GUI 工具进行溶剂化(在 0.15 M NaCl 溶液中),并调整盒子尺寸以保持 12 Å 的边缘距离,以确保模拟环境的充分性和一致性。催化三联体中的 D25 残基被建模为去质子化状态,这与中子散射研究的发现相符 [11, 12]。

达芦那韦(DRV)分子在 QM-CR 框架下被分解为七个化学片段,以便进行精细的相互作用分析 [Figure 4, Page 7]:

  • APC (Core Scaffold): 核心骨架,是 DRV 的主要结构部分。
  • THF (bis-THF): 双四氢呋喃环。
  • TOL (Benzene): 苯环。
  • IBA (Isobutyl): 异丁基。
  • SO2 (Sulfonyl Connector): 磺酰基连接器。
  • ANL (Aniline): 苯胺。

这些片段化使得研究人员能够识别 DRV 的哪些特定化学基团受到突变诱导的结合力损失的影响,从而为抑制剂设计提供靶点。

2.2 计算所得数据

2.2.1 DRV 片段化与纯度指标

QM-CR 框架首先对达芦那韦(DRV)分子进行了片段化处理,以识别化学上独立的模块元素。基于纯度指标(|Π|)的优化算法,DRV 被分解为七个核心片段:APC(核心骨架)、THF(双四氢呋喃)、TOL(苯)、IBA(异丁基)、SO2(磺酰连接器)和 ANL(苯胺) [Figure 4, Page 7]。这些片段的 |Π| 值均低于 0.075 的阈值,表明这种片段化是可靠且具有化学意义的。图 5 [Page 8] 显示了在 MD 运行过程中,DRV 片段、蛋白质氨基酸和水分子的 |Π| 值分布。|Π| 值通常高于气相结果,这归因于与蛋白酶的相互作用。SO2 和 ANL 片段的变异性最大,表明它们在构象集合中与环境的相互作用具有更高的动态性。

2.2.2 突变诱导的电子相互作用全局修饰

为了量化突变对 DRV 和 PR 之间相互作用的影响,研究分析了通过 QM-CR 框架获得的静电相互作用能(Eel)和片段键序(FBO)的变化 [Figure 6, Page 9]。

  • 6OPT (2 个突变体): 该体系的 Eel 和 FBO 分布与野生型相比,表现出较小的偏差。这表明活性位点附近的两个突变对 DRV-PR 电子相互作用网络的整体影响相对有限,网络结构基本保持。
  • 6OPZ (11 个突变体): 相比之下,11 个突变体(6OPZ)体系显示出显著的系统性变化:
    • Eel 变化: 静电相互作用项的分布整体向更正值移动。这意味着 DRV 与 PR 之间的吸引力减弱,或者排斥力增加,导致总体的静电结合力减弱。
    • FBO 变化: 片段键序的分布整体减小。由于 FBO 衡量的是片段之间共享的电子密度部分,其减少表明 DRV 与周围残基之间的电子耦合减弱。这与实验观察到的 6OPZ 体系对 DRV 结合亲和力的大幅下降高度一致。

这些结果清晰地表明,多重突变(尤其是远端突变)通过重塑电子相互作用网络,系统性地削弱了药物的结合力,而非仅仅改变局部的几何形状。

2.2.3 相互作用指纹与残基解析

为了识别对 DRV-PR 相互作用贡献最大的蛋白酶区域,研究分析了野生型复合物中残基分辨的 Eel 和 FBO [Figure 7, Page 10]。

  • 野生型:
    • Eel 分布: 最强的静电相互作用主要来自抑制剂的中央片段(特别是 APC 和 SO2 基团),它们直接与蛋白酶催化区域(D25, T26, G27 附近)和翻盖铰链区域(~45 位点)的残基相互作用。THF 和 ANL 片段的贡献中等,而 IBA 和 TOL 片段的静电相互作用较小。
    • FBO 分布: APC 片段表现出最大的 FBO 贡献,表明它是抑制剂在结合口袋内的主要电子锚点。THF 和 ANL 片段也对电子耦合有贡献,但程度较小,而 IBA、SO2 和 TOL 片段的键序相对较小。
  • 突变体比较 [Figure 8, Page 12]:
    • 2 个突变体(6OPT): 在突变位点 V82 和 I84 处,FBO 相互作用减弱,翻盖区域(I47, G48, G49)的 FBO 也有所降低。然而,整体变化不足以预测 DRV 抑制作用的显著减弱。
    • 11 个突变体(6OPZ): 在催化三联体区域(特别是 D25, D29, D30),FBO 显著减弱。值得注意的是,V32 处的相互作用增强,这可能是一种补偿机制,旨在减轻其他地方的结合损失。

2.2.4 抑制剂片段水平变化 (11 个突变体系统)

针对 11 个突变体(6OPZ)体系的详细分析 [Figure 9, Page 13] 揭示了结合力损失在 DRV 不同片段间的分布:

  • 静电相互作用: Eel 的变化主要通过抑制剂的中央片段(SO2 和 APC)进行调节。这些片段是配体的主要静电相互作用中心,将突变诱导的扰动传递到整个结合口袋。例如,APC 与 G16E(突变点)以及对催化区域的 G27 和翻盖区域的 G48, G49 的间接相互作用均表现出减弱。
  • 片段键序: FBO 的减少主要集中在 DRV 的核心骨架(APC 片段)。在催化区域的残基中,APC 片段的键序下降最为明显。此外,ANL 片段与蛋白酶的相互作用也有显著减弱。这意味着针对 ANL 基团的修饰可能成为克服耐药性的一种潜在途径。

2.2.5 相互作用图谱

FBO 相互作用网络被可视化为图谱,其中节点代表 DRV 片段或蛋白质残基,边表示它们之间的电子耦合强度 [Figure 10, Page 14]。

  • 野生型(6OPS): 展示了 DRV 片段与核心活性位点残基及翻盖区域的紧密电子耦合。V82 和 I84 与 DRV 的两个片段直接连接。远端突变 L33、I54、L76 围绕 V32 形成星形拓扑。K45 与 D30 有强相互作用。
  • 突变体(6OPT 和 6OPZ): 图谱显示,随着突变的引入,特别是 6OPZ 体系中,许多相互作用的强度(由边粗细和节点颜色表示)显著减弱或消失。这直观地展示了远端突变如何重塑结合口袋内的电子耦合网络。

2.2.6 相互作用网络的稳定性

时间分辨的相互作用曲线(Figure 11, Page 15)显示,在 MD 模拟过程中,不同配体片段的相对贡献保持稳定。尽管热力学涨落导致瞬时相互作用值有中等程度的变化,但片段贡献的层次结构基本保持不变。晶体结构与 MD 构象集合的相互作用值进行比较后发现,许多片段的相互作用值接近 MD 分布的中位数,表明实验晶体构象能够合理地代表主要的相互作用模式。然而,APC 和 THF 片段与蛋白酶的 FBO 相互作用在晶体结构中略有减弱,这表明相互作用网络应被视为一个动态集合,而非单一的静态构象。MD 轨迹提供了关于结合口袋内抑制剂稳定性的相互作用强度范围的重要信息。

2.3 性能数据

本研究的工作流设计旨在最大化计算资源的利用率和吞吐量,尤其是在异构超算环境中 [Page 6, Figure 3]。

  • 工作流效率时间线: Figure 3 展示了野生型变体第一个生产阶段(40 ns)的计算时间线。它清晰地描绘了 GENESIS(MD)计算的运行、快照的生成、以及随后 BigDFT(QM)和序列化(QM-CR)任务的提交和执行。
  • 即时处理(In-operando processing): 实验结果显示,BigDFT 和序列化任务在 MD 快照可用后几乎立即开始执行。这种“即时”处理模式显著减少了传统工作流中存在的等待时间和数据传输延迟。例如,不需要等到整个 MD 运行完成后才开始 QM 分析,而是可以在 MD 生产的同时并行进行 QM 计算。
  • 计算重叠与并发性: 时间线显示,工作流中存在多个 BigDFT 计算和多个序列化计算同时进行的重叠。在低系统利用率时期,这种并发性使得 BigDFT 和序列化任务能够几乎在提交后立即启动。这种并发处理能力对于提高整体吞吐量至关重要。
  • 资源利用率: 论文指出,如果所有 20 个 BigDFT 计算同时进行,大约需要 Odyssey CPU 分区六分之一的资源。通过“即时”提交计算,并允许它们在资源可用时尽快运行,工作流在实际机器负载下表现出更高的吞吐量,因为一次性预留如此大的计算资源通常不太可能。
  • MDAnalysis 延迟: 值得注意的是,由于 MDAnalysis 处理不完整轨迹的方式,最后两个快照的计算直到整个 MD 运行完成后才启动。这表明在 MD 软件与分析工具之间的接口优化仍有改进空间。

总体而言,这些性能数据验证了所提出工作流在利用异构超算资源方面是高效且有效的,尤其是在需要将长时程 MD 与高通量 QM 分析相结合的场景中。

3. 代码实现细节、复现指南、所用的软件包及开源 repo link

本研究的强大之处在于其精心设计的计算工作流,它无缝整合了多个高性能计算软件包,并通过一个智能的编排工具实现自动化。以下是其代码实现细节、复现指南以及所使用的主要软件包和开源资源。

3.1 工作流编排细节

整个工作流的自动化和耦合是通过 remotemanager [37] Python 软件包实现的 [Page 5]。remotemanager 能够将任意 Python 函数远程执行到高性能计算集群的作业调度系统上。其核心功能包括:

  • 计算机定义: 用户可以定义远程计算机,包括作业脚本模板、主机 URL 和提交器(例如 pjsub)。这使得工作流能够适应不同的超算环境和调度系统。例如,Listing 1 [Page 5] 展示了如何在 Wisteria/BDEC-0 机器上定义一个计算机,包括 PJM 调度器的资源请求(如节点数、OMP 线程、运行时间等)。
  • 远程函数执行: 用户可以编写一个普通的 Python 函数,remotemanager 会将其转换为可在远程计算机上执行的作业。这极大地简化了远程执行的复杂性。
  • 数据集管理: remotemanager 允许用户创建“数据集”,将多个远程函数调用组织在一起。每个函数调用可以有不同的参数,从而实现高通量计算。
  • 异步提交与结果获取: 任务可以异步提交,并且一旦单个任务完成,其结果即可立即获取。这种设计对于本研究的“即时”(in-operando)工作流至关重要,因为它允许在 MD 模拟进行中,一旦快照可用,就立即启动 QM 计算,而无需等待整个 MD 运行完成。

整体耦合策略 [Figure 2, Page 6]:

  1. 管理进程: 一个 Python 管理进程在用户本地计算机上运行,负责协调整个 MD + QM + QM-CR 工作流。
  2. MD 任务: 首先提交一个 GENESIS 任务,它在 Aquarius GPU 分区上运行分子动力学模拟,并以固定频率(每 1 ns)写入新的轨迹文件。
  3. 快照提取与 QM 提交: 一个独立的 Python 子例程周期性地使用 MDAnalysis [38] 检查 GENESIS 轨迹文件,以提取新的快照。一旦检测到新的快照,remotemanager 就会生成并提交相应的 BigDFT QM 计算到 Odyssey CPU 分区。
  4. QM 结果序列化: 当 BigDFT QM 计算完成后,remotemanager 会自动提交一个序列化任务。这个任务使用 QM-CR 框架对 QM 结果进行后处理,将其转换为易于分析的格式。
  5. 阶段性工作流: 论文为 GENESIS 的不同计算阶段(优化、平衡、生产)生成了单独的工作流,增强了灵活性。

这种设计使得工作流能够实现高度自动化、并行化和“即时”处理,显著提高了计算效率和数据分析的及时性。

3.2 所用的主要软件包

本研究集成了多个高性能计算和数据分析软件包:

  • 分子动力学 (MD) 软件:
    • GENESIS [13]: (General-purpose, ENabled by an E-platform for biomolecular Simulations) — 这是一个高性能的分子动力学模拟软件,特别针对 GPU 加速进行了优化。它用于生成蛋白质-配体复合物的长时间构象轨迹。 (非开源,但其论文是公开的)
  • 量子力学 (QM) 软件:
    • BigDFT [18]: — 这是一个基于 Daubechies 小波基组的密度泛函理论代码。其显著特点是支持线性标度算法 [19],使得它能够高效处理包含数千个原子的大型体系,这是本研究中 QM 计算能够实现的关键。BigDFT 在 CPU 节点上运行。 (开源,项目主页:https://bigdft.org/)
  • 力场 (Force Fields):
    • CHARMM36m [14]: — 用于蛋白质的最新 CHARMM 泛用力场版本,包含对蛋白质折叠和构象采样的改进。 (商业力场,但相关参数文件和使用指南广泛可用)
    • CHARMM General Forcefield (CGenFF) [15, 16]: — 用于配体(达芦那韦)的力场,兼容 CHARMM 蛋白质力场。 (商业力场)
    • TIP3P [17]: — 一种广泛使用的三点水模型,用于描述溶剂水分子。 (广泛应用于 MD 软件中)
  • 数据分析与处理工具:
    • MDAnalysis [38]: — 一个开源的 Python 库,用于读取、操作和分析分子动力学轨迹。在本研究中,它用于从 GENESIS 轨迹中提取快照。 (开源,项目主页:https://www.mdanalysis.org/)
    • QM-CR 框架 (Complexity Reduction): — 这是作者课题组开发的一个定制的后处理框架,用于分析 BigDFT 结果中的片段间相互作用。它并非一个独立的通用软件包,而是作为其研究工作流的一部分。 (相关理论和实现细节在 [23, 25, 26] 论文和本文附录中描述)
  • 工作流管理:
    • remotemanager [37]: — 一个开源的 Python 库,用于在远程 HPC 集群上提交、管理和监控计算任务,实现了本研究中 MD 和 QM 任务的无缝耦合。 (开源,项目主页:https://remotemanager.readthedocs.io/)

3.3 复现指南 (概念性)

为了复现本研究的工作流,使用者需要具备访问 Wisteria/BDEC-0 或类似异构超算平台(具备 GPU 和 CPU 分区,并支持跨分区通信和作业调度)的权限。以下是大致的复现步骤:

  1. 获取代码与数据:

    • 开源仓库: 所有的计算设置,包括输入参数,都可在 Gitlab 仓库中获取:https://gitlab.com/wddawson/genesis-bigdft-coupling [Page 16]。下载此仓库以获取所需脚本和配置。
    • 初始结构: 获取 HIV-1 PR-DRV 复合物的 PDB 结构文件 (6OPS, 6OPT, 6OPZ)。
  2. 环境配置:

    • 超算登录: 登录到目标超算平台(例如 Wisteria/BDEC-0)。
    • 软件模块加载: 加载必要的环境模块,例如 Python 环境、MPI 库、编译器、以及 GENESIS 和 BigDFT 运行时所需的特定库。
    • Python 包安装: 确保 Python 环境中安装了 remotemanagerMDAnalysis 等所需的 Python 库。可能需要创建虚拟环境。
    • GENESIS 和 BigDFT 可执行文件: 确保 GENESIS 和 BigDFT 的可执行文件已在系统上正确编译并配置,路径可被 remotemanager 作业脚本调用。
  3. 准备 MD 模拟体系:

    • 使用 CHARMM-GUI: 按照论文 [8-10] 中描述的方法,使用 CHARMM-GUI 在线工具准备 MD 模拟体系。这包括蛋白酶和配体的质子化、溶剂化(添加水和 NaCl 离子)、以及力场参数的生成。
    • 生成输入文件: 根据 CHARMM36m、CGenFF 和 TIP3P 力场,生成 GENESIS 所需的拓扑和参数文件。
  4. 配置工作流脚本:

    • 修改 remotemanager 配置: 根据您所使用的超算环境,调整 Gitlab 仓库中的 remotemanager 配置文件和作业脚本模板(例如 Listing 1)。这包括:
      • 作业调度器参数(PJM、Slurm 或 LSF 等)。
      • 资源请求(GPU/CPU 节点数量、核心数、内存、运行时间)。
      • 软件路径和模块加载命令。
    • 定义 QM-CR 逻辑: QM-CR 后处理的细节可能需要根据具体的分析需求进行配置。
  5. 执行工作流:

    • 启动管理进程: 在本地或一个稳定的远程会话中,运行 remotemanager 的主控制脚本。这个脚本将启动整个工作流。
    • MD 阶段: GENESIS 任务将首先提交到 GPU 分区进行能量最小化、平衡和生产模拟。
    • QM + QM-CR 阶段: 管理进程将监控 MD 轨迹的生成。一旦有新的快照可用,它将自动提取快照,生成 BigDFT 输入,并提交 QM 计算到 CPU 分区。QM 计算完成后,将自动触发 QM-CR 后处理步骤。
    • 数据收集: remotemanager 将异步收集所有完成任务的结果,并在本地进行进一步分析和可视化。

注意事项:

  • 超算资源: 确保您有足够的计算配额和访问权限,尤其是在 GPU 和大规模 CPU 分区上。
  • 队列等待: BigDFT 任务在 64 个 CPU 节点上运行,可能需要较长的队列等待时间。
  • 故障排除: 监控作业日志是解决任何计算问题(如资源不足、软件错误、路径问题)的关键。

遵循上述概念性指南,并参考 Gitlab 仓库中的详细脚本和说明,理论上可以复现本研究的计算结果。

4. 关键引用文献与对这项工作局限性的评论

4.1 关键引用文献

本研究综合了分子动力学、量子力学和高性能计算的多个前沿领域,因此引用了大量基础性及前沿性文献,这些文献共同构成了该工作的理论和方法支撑:

方法论与理论基础:

  • MD 模拟指南: [1] Brázdová, V.; Bowler, D. R. (2013) 的《Atomistic computer simulations: a practical guide》提供了原子模拟的实践指导,是 MD 领域的重要参考。
  • QM 在计算分子科学中的地位: [2] Dawson, W. et al. (2022) 在《WIREs Comput. Mol. Sci.》上发表的文章强调了 QM 在计算分子科学中提供电子结构洞察的独特价值。
  • DFT 理论的奠基: [3] Hohenberg, P.; Kohn, W. (1964) 和 [4] Kohn, W.; Sham, L. J. (1965) 在《Phys. Rev.》上发表的经典论文,奠定了密度泛函理论(DFT)的理论基础。
  • BigDFT 软件与线性标度: [18] Ratcliff, L. E. et al. (2020) 和 [19] Mohr, S. et al. (2015) 的工作详细介绍了 BigDFT 代码及其在处理大体系时的线性标度算法,这是本研究 QM 计算可行性的关键。
  • QM-CR 框架: [23] Dawson, W. et al. (2023), [25] Mohr, S. et al. (2017) 和 [26] Dawson, W. et al. (2020) 的系列论文详细阐述了量子力学复杂度降低(QM-CR)框架的理论和应用,包括片段纯度指标、片段键序和静电相互作用的计算。

MD 软件与力场:

  • GENESIS 软件: [13] Jung, J. et al. (2024) 在《J. Phys. Chem. B》上发表的文章介绍了 GENESIS 软件包,本研究利用其进行 GPU 加速的 MD 模拟。
  • CHARMM 力场: [14] Huang, J. et al. (2017) 在《Nat. Methods》上发表的文章介绍了 CHARMM36m 蛋白质力场。此外,[15, 16] Vanommeslaeghe, K. et al. (2012) 的两篇论文详细描述了 CHARMM 通用力场(CGenFF)的开发,用于配体参数化。
  • 水模型: [17] Jorgensen, W. L. et al. (1983) 在《J. Chem. Phys.》上发表的文章介绍了 TIP3P 水模型,这是 MD 模拟中广泛使用的显式水模型。

案例研究背景与 HIV-1 蛋白酶:

  • 达芦那韦抗性研究的灵感来源: [6] Henes, M. et al. (2019) 在《ACS Chem. Biol.》上发表的综合实验和计算研究,分析了 HIV-1 PR 对达芦那韦的抗性,特别是远端突变的作用,为本研究提供了重要的案例背景。
  • HIV-1 PR 相关研究: [11] Smith, R. et al. (1996) 涉及 HIV-1 PR 的结构。 [12] Gerlits, O. et al. (2016) 和 [39] L. N. M. N. et al. (2010) 的研究讨论了 HIV-1 PR 的突变与抑制剂结合机制。
  • 药物设计与抗性: [40] Ghosh, A. K. et al. (2008) 在《Acc. Chem. Res.》中讨论了 HIV-1 PR 抑制剂的设计,强调了其结构与功能的关键方面。

超算平台与工作流管理:

  • Wisteria/BDEC-0 超算: [33] 东京大学的官方介绍提供了 Wisteria/BDEC-0 超算平台的信息。 [34] Sumimoto, S. et al. (2022) 讨论了在异构计算环境下多 MPI 程序耦合的系统级通信。
  • 异构计算与多物理场: [35] Caviedes-Voullième, D. et al. (2025) 和 [36] Hamamura, I. et al. (2026) 的文章探讨了在 Wisteria 上的多物理场应用和量子经典耦合,体现了异构超算的潜力。
  • remotemanager 软件包: [37] Dawson, W. et al. (2024) 介绍了 remotemanager 软件,它是本研究实现跨异构平台工作流编排的关键工具。

高级采样方法:

  • [41] Yang, Y. I. et al. (2019) 的文章讨论了化学物理中的高级采样方法,这提示了未来研究中可以进一步整合的潜在方向。

4.2 对这项工作局限性的评论

尽管本研究在结合 MD 和 QM 揭示远端突变效应方面取得了显著进展,并提供了一个高效的异构超算工作流,但仍存在一些潜在的局限性,值得在未来的工作中加以解决或进一步探讨:

  1. DFT 泛函的选择与精度:

    • PBE 泛函的通用性: 本研究使用了 PBE 泛函 [20]。PBE 是一种广受认可的广义梯度近似(GGA)泛函,在描述电子密度和某些相互作用(如氢键)方面表现良好。对于本研究关注的基于密度的描述符(如 Eel 和 FBO)以及相互作用模式的相对变化,PBE 泛函可能是足够的。
    • 结合能的局限性: 然而,PBE 泛函通常不足以精确预测结合自由能的绝对值,尤其是在包含色散相互作用的生物分子体系中。更先进的泛函,如包含色散校正(DFT-D)的泛函,或专门为结合能校正设计的泛函 [22],可能会提供更准确的绝对结合能数据。因此,本研究的结论更侧重于相互作用模式的改变和相对强度,而非对结合自由能的精确量化。
  2. QM-CR 片段化参数:

    • 纯度指标阈值: QM-CR 框架通过纯度指标 (|Π|) 确保片段划分的化学意义和可靠性。论文设定了一个阈值(0.075),低于该阈值的片段被认为是独立的。虽然这种方法是系统性的,但特定阈值的选择可能会影响片段的细节和分析的颗粒度。论文讨论了更宽松阈值(0.1)可能带来的额外片段化 [Figure 4],但未深入探讨不同阈值对最终结果解释的敏感性。
  3. 计算成本与采样广度:

    • MD 采样时间: 尽管 80 ns 的 MD 模拟对于传统的蛋白质-配体体系来说是“长时程”的,但对于某些涉及罕见构象转变或需要更全面探索自由能景观的复杂变构过程,这个采样时间可能仍然有限 [41]。更长的 MD 模拟或结合增强采样技术(如伞形采样、副本交换 MD 等)可能会捕获更广泛的构象,从而提供更全面的电子结构洞察。
    • QM 计算开销: 即使采用了线性标度 DFT 和异构超算策略,每个 BigDFT 计算仍需 64 个 CPU 节点,并且对 20 帧快照进行计算,总计算资源消耗仍然巨大。这限制了在极大规模(例如,对数百个突变体或数千个 MD 构象进行 QM 分析)上的应用。
  4. 因果关系与关联性:

    • 机制解释: 本研究成功地将远端突变导致的电子结构变化(Eel 和 FBO 的改变)与实验观察到的药物结合力下降关联起来。然而,从这些关联性推断严格的因果关系,即精确量化某种电子结构变化对结合自由能贡献的定量因果,仍需更深入的理论推导和计算方法(例如,结合自由能分解分析)。
  5. 水分子在 QM 体系中的处理:

    • 截断距离: QM 计算中只包含了蛋白质-配体复合物 3 Å 范围内的水分子。虽然这是在计算成本和环境效应之间取得平衡的常见做法,但更远处的水分子,特别是那些在结合口袋或通道中形成稳定水桥的水分子,其长程静电效应和动态重排可能未被完全捕获。这可能对 Eel 的计算精度产生一定影响。
  6. 缺乏直接的自由能计算:

    • 本研究主要关注电子结构描述符(Eel, FBO)的变化,这些变化与药物结合亲和力的损失高度相关。然而,它并未直接计算结合自由能(ΔG)。虽然电子结构洞察对于理解机制至关重要,但直接的自由能计算将提供更全面的热力学描述,并与实验测量的结合亲和力直接对比。

总而言之,本研究提出了一种强大且创新的工作流,为理解远端突变效应提供了前所未有的电子结构视角。未来的工作可以在上述局限性方面进行拓展,以进一步提升方法的精度、广度和应用范围。

5. 其他必要的补充

5.1 远端突变与变构效应的深远意义

本研究的核心价值之一在于其对远端突变(distal mutations)效应的深入探讨。长久以来,药物设计主要聚焦于结合位点及其附近的相互作用。然而,本研究强调了远端突变在药物耐药性中的关键作用,这代表了对传统药物设计范式的深刻挑战。这些突变通过非加性机制(如上位性效应)和集体重塑构象集合来影响药物结合,而并非直接改变结合位点的几何形状。这种长程的变构(allostery)效应使得基于静态结构或简单几何匹配的药物设计难以应对。

通过耦合 MD 与 QM,本研究提供了一种前所未有的能力,能够直接观察和量化远端突变如何通过电子结构层面(如电荷重新分布、键合强度变化)进行信息传递和影响药物结合。这不仅揭示了 HIV-1 蛋白酶中达芦那韦耐药性的复杂机制,也为理解其他生物系统中的变构调节和远端效应提供了通用框架。这种对电子耦合网络进行显式映射的方法,有望将药物发现从纯粹的“结构驱动”转向“电子驱动”,从而设计出更具鲁棒性和广谱性的抑制剂。

5.2 异构超算与数据流导向工作流的强大优势

本研究成功地展示了异构超级计算机(结合 GPU 和 CPU 架构)在解决复杂科学问题中的巨大潜力。Wisteria/BDEC-0 超算上的 Aquarius(GPU)分区和 Odyssey(CPU)分区分别针对 MD 和 QM 计算进行了优化,这种异构性反映了当前高性能计算硬件发展的趋势。

关键在于 remotemanager [37] 软件包所实现的“数据流导向”(data-flow-oriented)工作流和“即时”(in-operando)分析模式。传统的工作流通常是串行的:MD 模拟完成后,所有数据传输到另一个位置,再进行 QM 计算和后处理。这种模式效率低下,且数据传输量巨大。本研究的创新之处在于:

  • 无缝集成与自动化: remotemanager 作为管理层,自动化了 MD 模拟、快照提取、QM 计算提交和 QM-CR 后处理的整个过程。用户只需配置一次,系统即可自行调度和执行。
  • 实时分析能力: 允许在 MD 模拟仍在进行中时,一旦有新的快照可用,就立即启动 QM 计算。这极大地减少了整体处理时间,实现了“持续分析”,避免了长时间的等待。这对于快速迭代药物设计周期或在有限的超算资源下实现高吞吐量至关重要。
  • 最大化资源利用: 通过将 GPU 用于 MD(例如 GENESIS)和 CPU 用于 QM(例如 BigDFT),工作流能够充分利用异构超算中不同硬件的优势,同时处理不同的计算任务。 Figure 3 [Page 7] 清晰地展示了 MD、QM 和 QM-CR 任务是如何重叠执行的,从而提高了整体吞吐量和资源利用效率。
  • 可扩展性: 这种模块化和异步提交的设计使得工作流具有良好的可扩展性,能够轻松适应更大规模的模拟(例如更多突变体、更长的模拟时间)和更复杂的分析任务。

这种异构超算上数据流导向的工作流,不仅是计算化学领域的一个技术突破,也为其他需要大规模、多尺度模拟的科学领域提供了借鉴,例如材料科学、气候模拟等。

5.3 药物发现范式的潜在变革

本研究不仅解决了 HIV-1 达芦那韦耐药性的具体问题,更重要的是,它预示着药物发现范式的一次潜在变革。

  • 从静态结构到动态电子耦合: 传统的药物设计多依赖于晶体结构等静态信息,侧重于几何匹配和形状互补性。本研究通过对动态构象集合进行电子结构分析,将重点从“静态结构”转移到“动态电子耦合”的映射。这种转变允许研究者更深入地理解药物-靶点相互作用的本质,包括电荷转移、极化效应、以及共价键形成等电子层面的微观细节。
  • 识别“可成药”的电子景观: 通过 QM-CR 框架识别药物片段上相互作用损失集中的特定化学基团,为先导化合物的优化提供了“处方模型”。例如,发现 ANL 片段的相互作用减弱,可能提示针对该区域进行化学修饰以恢复或增强结合力。这种基于电子结构洞察的策略,远比简单的结构修改更具指导性。
  • 设计抗性鲁棒性抑制剂: 理解系统性、突变诱导的失稳机制是设计能够维持结合稳定性,对抗病毒变异的抑制剂的关键。通过揭示远端突变如何通过长程电子相互作用影响结合口袋,研究者可以开发出更具鲁棒性的药物,降低耐药性产生的几率。
  • 通用性与未来拓展: 该方法不仅限于 HIV-1,还可推广到其他生物系统,研究由远端突变或变构效应引起的结构-活性关系。随着便携式和可扩展量子方法的发展,这种工作流将成为计算药物发现日益核心的工具。

5.4 可视化与数据解释

本研究利用了一系列直观且信息丰富的图表来展示其复杂的计算结果,这对于理解和解释科学发现至关重要。这些可视化工具有效地将高维的电子结构数据转化为易于理解的化学洞察:

  • 工作流时间线 (Figure 3): 清晰地展示了 MD、QM 和 QM-CR 任务的提交和执行顺序,以及它们之间的重叠,直观地证明了工作流的效率和并行化能力。
  • DRV 片段化 (Figure 4): 通过化学结构图和表格,明确展示了达芦那韦被 QM-CR 框架分解的七个化学片段,并报告了每个片段的纯度指标,为后续分析奠定基础。
  • 纯度指标分布 (Figure 5): 使用小提琴图或盒形图,展示了 DRV 片段、氨基酸和水分子在 MD 模拟过程中纯度指标的变异性,突出了动态环境中电子结构的柔性。
  • 全局相互作用修饰 (Figure 6): 采用盒形图比较了不同突变体中总静电相互作用能(Eel)和片段键序(FBO)的分布变化,直观地揭示了突变对整体药物结合力损失的影响方向和程度。
  • 残基分辨相互作用热力图 (Figure 7, S1, S2): 这些热力图清晰地展示了 DRV 片段与蛋白质每个氨基酸之间的平均 Eel 和 FBO,通过颜色深浅和分布图直观地识别出结合界面的关键相互作用“热点”。这使得研究者能够精确识别哪些蛋白质残基和药物片段对结合贡献最大,以及它们如何因突变而变化。
  • 氨基酸水平相互作用分布 (Figure 8, S4, S5): 详细的盒形图比较了不同突变体中关键氨基酸与 DRV 整体或特定片段之间的 Eel 和 FBO 变化,特别突出了统计学上显著的差异,帮助定位突变影响的具体作用点。
  • 抑制剂片段分解 (Figure 9): 针对 11 个突变体体系,该图表进一步分解了特定蛋白质残基与各个 DRV 片段之间的 Eel 和 FBO 变化,从而更细致地揭示了相互作用损失在药物分子上的分布。
  • FBO 相互作用网络 (Figure 10): 以节点和边的形式表示了 DRV 片段与蛋白质残基之间的 FBO 相互作用网络。节点颜色和边粗细代表相互作用强度,突变位点用不同颜色标记,直观地展示了突变如何重塑电子耦合网络,以及远端突变如何通过多跳连接影响活性位点。
  • 相互作用时间稳定性 (Figure 11): 累积平均图展示了不同 DRV 片段的 Eel 和 FBO 随 MD 模拟时间的变化,评估了相互作用模式的稳定性,并与晶体结构值进行比较,强调了动态模拟的重要性。

这些可视化工具的结合,使得本研究不仅能生成大量复杂的计算数据,还能以清晰、易懂的方式呈现这些数据,从而促进了对 HIV-1 蛋白酶中药物耐药性机制的深入理解。