A Rigorous Framework for Calculating Protein-Protein Binding Affinities in Membranes.

Calculating the binding free energy of integral transmembrane (TM) proteins is crucial for understanding the mechanisms by which they recognize one another and reversibly associate. The glycophorin A (GpA) homodimer, composed of two α-helical segments, has long served as a model system for studying TM protein reversible association. The present work establishes a methodological framework for calculating the binding affinity of the GpA homodimer in the heterogeneous environment of a membrane. Our investigation carefully considered a variety of protocols, including the appropriate choice of the force field, rigorous standardization reflecting the experimental conditions, sampling algorithm, anisotropic environment, and collective variables, to accurately describe GpA dimerization via molecular dynamics-based approaches. Specifically, two strategies were explored: (i) an unrestrained potential mean force (PMF) calculation, which merely enhances sampling along the separation of the two binding partners without any restraint, and (ii) a so-called "geometrical route", whereby the α-helices are progressively separated with imposed restraints on their orientational, positional, and conformational degrees of freedom to accelerate convergence. Our simulations reveal that the simplified, unrestrained PMF approach is inadequate for the description of GpA dimerization. Instead, the geometrical route, tailored specifically to GpA in a membrane environment, yields excellent agreement with experimental data within a reasonable computational time. A dimerization free energy of -10.7 kcal/mol is obtained, in fairly good agreement with available experimental data. The geometrical route further helps elucidate how environmental forces drive association before helical interactions stabilize it. Our simulations also brought to light a distinct, long-lived spatial arrangement that potentially serves as an intermediate state during dimer formation. The methodological advances in the generalized geometrical route provide a powerful tool for accurate and efficient binding-affinity calculations of intricate TM protein complexes in inhomogeneous environments.

Improving Speed and Affordability without Compromising Accuracy: Standard Binding Free-Energy Calculations Using an Enhanced Sampling Algorithm, Multiple-Time Stepping, and Hydrogen Mass Repartitioning.

Accurate evaluation of protein-ligand binding free energies in silico is of paramount importance for understanding the mechanisms of biological regulation and providing a theoretical basis for drug design and discovery. Based on a series of atomistic molecular dynamics simulations in an explicit solvent, using well-tempered metadynamics extended adaptive biasing force (WTM-eABF) as an enhanced sampling algorithm, the so-called "geometrical route" offers a rigorous theoretical framework for binding affinity calculations that match experimental values. However, although robust, this strategy remains expensive, requiring substantial computational time to achieve convergence of the simulations. Improving the efficiency of the geometrical route, while preserving its reliability through improved ergodic sampling, is, therefore, highly desirable. In this contribution, having identified the computational bottleneck of the geometrical route, to accelerate the calculations we combine (i) a longer time step for the integration of the equations of motion with hydrogen-mass repartitioning (HMR), and (ii) multiple time-stepping (MTS) for collective-variable and biasing-force evaluation. Altogether, we performed 50 independent WTM-eABF simulations in triplicate for the "physical" separation of the Abl kinase-SH3 domain:p41 complex, following different HMR and MTS schemes, while tuning, in distinct protocols, the parameters of the enhanced-sampling algorithm. To demonstrate the consistency and reliability of the results obtained with the best-performing setups, we carried out quintuple simulations. Furthermore, we demonstrated the transferability of our method to other complexes by triplicating a 200 ns separation simulation of nine chosen protocols for the MDM2-p53:NVP-CGM097 complex. [Holzer et al. J. Med. Chem. 2015, 58, 6348-6358.] Our results, based on an aggregate simulation time of 14.4 µs, allowed an optimal set of parameters to be identified, able to accelerate convergence by a factor of three without any noticeable loss of accuracy.

When the Dust Has Settled: Calculation of Binding Affinities from First Principles for SARS-CoV-2 Variants with Quantitative Accuracy.

Accurate determination of binding free energy is pivotal for the study of many biological processes and has been applied in a number of theoretical investigations to compare the affinity of severe acute respiratory syndrome coronavirus 2 variants toward the host cell. Diversity of these variants challenges the development of effective general therapies, their transmissibility relying either on an increased affinity toward their dedicated human receptor, the angiotensin-converting enzyme 2 (ACE2), or on escaping the immune response. Now that robust structural data are available, we have determined with utmost accuracy the standard binding free energy of the receptor-binding domain to the most widespread variants, namely, Alpha, Beta, Delta, and Omicron BA.2, as well as the wild type (WT) in complex either with ACE2 or with antibodies, namely, S2E12 and H11-D4, using a rigorous theoretical framework that combines molecular dynamics and potential-of-mean-force calculations. Our results show that an appropriate starting structure is crucial to ensure appropriate reproduction of the binding affinity, allowing the variants to be compared. They also emphasize the necessity to apply the relevant methodology, bereft of any shortcut, to account for all the contributions to the standard binding free energy. Our estimates of the binding affinities support the view that while the Alpha and Beta variants lean on an increased affinity toward the host cell, the Delta and Omicron BA.2 variants choose immune escape. Moreover, the S2E12 antibody, already known to be active against the WT (Starr et al., 2021; Mlcochova et al., 2021), proved to be equally effective against the Delta variant. In stark contrast, H11-D4 retains a low affinity toward the WT compared to that of ACE2 for the latter. Assuming robust structural information, the methodology employed herein successfully addresses the challenging protein-protein binding problem in the context of coronavirus disease 2019 while offering promising perspectives for predictive studies of ever-emerging variants.

Hazardous Shortcuts in Standard Binding Free Energy Calculations.

Calculating the standard binding free energies of protein-protein and protein-ligand complexes from atomistic molecular dynamics simulations in explicit solvent is a problem of central importance in computational biophysics. A rigorous strategy for carrying out such calculations is the so-called "geometrical route". In this method, two molecular objects are progressively separated from one another in the presence of orientational and conformational restraints serving to control the change in configurational entropy that accompanies the dissociation process, thereby allowing the computations to converge within simulations of affordable length. Although the geometrical route provides a rigorous theoretical framework, a tantalizing computational shortcut consists of simply leaving out such orientational and conformational degrees of freedom during the separation process. Here the accuracy and convergence of the two approaches are critically compared in the case of two protein-ligand complexes (Abl kinase-SH3:p41 and MDM2-p53:NVP-CGM097) and three protein-protein complexes (pig insulin dimer, SARS-CoV-2 spike RBD:ACE2, and CheA kinase-P2:CheY). The results of the simulations that strictly follow the geometrical route match the experimental standard binding free energies within chemical accuracy. In contrast, simulations bereft of geometrical restraints converge more poorly, yielding inconsistent results that are at variance with the experimental measurements. Furthermore, the orientational and positional time correlation functions of the protein in the unrestrained simulations decay over several microseconds, a time scale that is far longer than the typical simulation times of the geometrical route, which explains why those simulations fail to sample the relevant degrees of freedom during the separation process of the complexes.

Accurate determination of protein:ligand standard binding free energies from molecular dynamics simulations.

Designing a reliable computational methodology to calculate protein:ligand standard binding free energies is extremely challenging. The large change in configurational enthalpy and entropy that accompanies the association of ligand and protein is notoriously difficult to capture in naive brute-force simulations. Addressing this issue, the present protocol rests upon a rigorous statistical mechanical framework for the determination of protein:ligand binding affinities together with the comprehensive Binding Free-Energy Estimator 2 (BFEE2) application software. With the knowledge of the bound state, available from experiments or docking, application of the BFEE2 protocol with a reliable force field supplies in a matter of days standard binding free energies within chemical accuracy, for a broad range of protein:ligand complexes. Limiting undesirable human intervention, BFEE2 assists the end user in preparing all the necessary input files and performing the post-treatment of the simulations towards the final estimate of the binding affinity.