Overcoming Free-Energy Barriers with a Seamless Combination of a Biasing Force and a Collective Variable-Independent Boost Potential

Lambda-ABF: Simplified, Portable, Accurate, and Cost-Effective Alchemical Free-Energy Computation

We introduce the lambda-Adaptive Biasing Force (lambda-ABF) method for the computation of alchemical free-energy differences. We propose a software implementation and showcase it on biomolecular systems. The method arises from coupling multiple-walker adaptive biasing force with λ-dynamics. The sampling of the alchemical variable is continuous and converges toward a uniform distribution, making manual optimization of the λ schedule unnecessary. Contrary to most other approaches, alchemical free-energy estimates are obtained immediately without any postprocessing. Free diffusion of λ improves orthogonal relaxation compared to fixed-λ thermodynamic integration or free-energy perturbation. Furthermore, multiple walkers provide generic orthogonal space coverage with minimal user input and negligible computational overhead. We show that our high-performance implementations coupling the Colvars library with NAMD and Tinker-HP can address real-world cases including ligand-receptor binding with both fixed-charge and polarizable models, with a demonstrably richer sampling than fixed-λ methods. The implementation is fully open-source, publicly available, and readily usable by practitioners of current alchemical methods. Thanks to the portable Colvars library, lambda-ABF presents a unified user interface regardless of the back-end (NAMD, Tinker-HP, or any software to be interfaced in the future), sparing users the effort of learning multiple interfaces. Finally, the Colvars Dashboard extension of the visual molecular dynamics (VMD) software provides an interactive monitoring and diagnostic tool for lambda-ABF simulations.

Understanding the Reversible Binding of a Multichain Protein-Protein Complex through Free-Energy Calculations

We demonstrate that the binding affinity of a multichain protein-protein complex, insulin dimer, can be accurately predicted using a streamlined route of standard binding free-energy calculations. We find that chains A and C, which do not interact directly during binding, stabilize the insulin monomer structures and reduce the binding affinity of the two monomers, therefore enabling their reversible association. Notably, we confirm that although classical methods can estimate the binding affinity of the insulin dimer, conventional molecular dynamics, enhanced sampling algorithms, and classical geometrical routes of binding free-energy calculations may not fully capture certain aspects of the role played by the noninteracting chains in the binding dynamics. Therefore, this study not only elucidates the role of noninteracting chains in the reversible binding of the insulin dimer but also offers a methodological guide for investigating the reversible binding of multichain protein-protein complexes utilizing streamlined free-energy calculations.

Collective Variable-Based Enhanced Sampling: From Human Learning to Machine Learning

Enhanced-sampling algorithms relying on collective variables (CVs) are extensively employed to study complex (bio)chemical processes that are not amenable to brute-force molecular simulations. The selection of appropriate CVs characterizing the slow movement modes is of paramount importance for reliable and efficient enhanced-sampling simulations. In this Perspective, we first review the application and limitations of CVs obtained from chemical and geometrical intuition. We also introduce path-sampling algorithms, which can identify path-like CVs in a high-dimensional free-energy space. Machine-learning algorithms offer a viable approach to finding suitable CVs by analyzing trajectories from preliminary simulations. We discuss both the performance of machine-learning-derived CVs in enhanced-sampling simulations of experimental models and the challenges involved in applying these CVs to realistic, complex molecular assemblies. Moreover, we provide a prospective view of the potential advancements of machine-learning algorithms for the development of CVs in the field of enhanced-sampling simulations.

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.

Standard Binding Free-Energy Calculations: How Far Are We from Automation?

Recent success stories suggest that in silico protein-ligand binding free-energy calculations are approaching chemical accuracy. However, their widespread application remains limited by the extensive human intervention required, posing challenges for the neophyte. As such, it is critical to develop automated workflows for estimating protein-ligand binding affinities with minimum personal involvement. Key human efforts include setting up and tuning enhanced-sampling or alchemical-transformation algorithms as a preamble to computational binding free-energy estimations. Additionally, preparing input files, bookkeeping, and postprocessing represent nontrivial tasks. In this Perspective, we discuss recent progress in automating standard binding free-energy calculations, featuring the development of adaptive or parameter-free algorithms, standardization of binding free-energy calculation workflows, and the implementation of user-friendly software. We also assess the current state of automated standard binding free-energy calculations and evaluate the limitations of existing methods. Last, we outline the requirements for future algorithms and workflows to facilitate automated free-energy calculations for diverse protein-ligand complexes.

Deep-Learning-Assisted Enhanced Sampling for Exploring Molecular Conformational Changes

We present a novel strategy to explore conformational changes and identify stable states of molecular objects, eliminating the need for a priori knowledge. The approach applies a deep learning method to extract information about the movement modes of the molecular object from a short, high-dimensional, and parameter-free preliminary enhanced-sampling simulation. The gathered information is described by a small set of deep-learning-based collective variables (dCVs), which steer the production-enhanced-sampling simulation. Considering the challenge of adequately exploring the configurational space using the low-dimensional, suboptimal dCVs, we incorporate a method designed for ergodic sampling, namely, Gaussian-accelerated molecular dynamics (MD), into the framework of CV-based enhanced sampling. MD simulations on both toy models and nontrivial examples demonstrate the remarkable computational efficiency of the strategy in capturing the conformational changes of molecular objects without a priori knowledge. Specifically, we achieved the blind folding of two fast folders, chignolin and villin, within a time scale of hundreds of nanoseconds and successfully reconstructed the free-energy landscapes that characterize their reversible folding. All in all, the presented strategy holds significant promise for investigating conformational changes in macromolecules, and it is anticipated to find extensive applications in the fields of chemistry and biology.

MemCross: Accelerated Weight Histogram method to assess membrane permeability

Passive permeation events across biological membranes are determining steps in the pharmacokinetics of xenobiotics. To reach an accurate and rapid prediction of membrane permeation coefficients of drugs is a complex challenge, which can efficiently support drug discovery. Such predictions are indeed highly valuable as they may guide the selection of potential leads with optimum bioavailabilities prior to synthesis. Theoretical models exist to predict these coefficients. Many of them are based on molecular dynamics (MD) simulations, which allow calculation of permeation coefficients through the evaluation of both the potential of mean force (PMF) and the diffusivity profiles. However, these simulations still require intensive computational efforts, and novel methodologies should be developed and benchmarked. Free energy perturbation (FEP) method was recently shown to estimate PMF with a significantly reduced computational cost compared to the adaptive biasing force method. This benchmarking was achieved with small molecules, namely short-chain alcohols. Here, we show that to estimate the PMF of bulkier, drug-like xenobiotics, conformational sampling is a critical issue. To reach a sufficient sampling with FEP calculations requires a relatively long time-scale, which can lower the benefits related to the computational gain. In the present work, the Accelerated Weight Histogram (AWH) method was employed for the first time in all-atom membrane models. The AWH-based protocol, named MemCross, appears affordable to estimate PMF profiles of a series of drug-like xenobiotics, compared to other enhanced sampling methods. The continuous exploration of the crossing pathway by MemCross also allows modeling subdiffusion by computing fractional diffusivity profiles. The method is also versatile as its input parameters are largely insensitive to the molecule properties. It also ensures a detailed description of the molecule orientations along the permeation pathway, picturing all intermolecular interactions at an atomic resolution. Here, MemCross was applied on a series of 12 xenobiotics, including four weak acids, and a coherent structure-activity relationship was established.

Enhancing sampling with free-energy calculations

In recent years, considerable progress has been made to enhance sampling and help address biological questions, including, but not limited to conformational transitions in biomolecules and protein-ligand reversible binding, hitherto intractable by brute-force computer simulations. Many of these advances result from the development of a palette of methods aimed at exploring rare events through reliable free-energy calculations. The advent of new, often conceptually related methods has also rendered difficult the choice of the best suited option for a given problem. Here, we focus on geometrical transformations and algorithms designed to enhance sampling along adequately chosen progress variables, tracing their theoretical foundations, and showing how they are connected and can be blended together for improved performance.

Accurate Description of Solvent-Exposed Salt Bridges with a Non-polarizable Force Field Incorporating Solvent Effects

The strength of salt bridges resulting from the interaction of cations and anions is modulated by their environment. However, polarization of the solvent molecules by the charged moieties makes the accurate description of cation-anion interactions in an aqueous solution by means of a pairwise additive potential energy function and classical combination rules particularly challenging. In this contribution, aiming at improving the representation of solvent-exposed salt-bridge interactions with an all-atom non-polarizable force field, we put forth here a parametrization strategy. First, the interaction of a cation and an anion is characterized by hybrid quantum mechanical/molecular mechanics (QM/MM) potential of mean force (PMF) calculations, whereby constantly exchanging solvent molecules around the ions are treated at the quantum mechanical level. The Lennard-Jones (LJ) parameters describing the salt-bridge ion pairs are then optimized to match the reference QM/MM PMFs through the so-called nonbonded FIX, or NBFIX, feature of the CHARMM force field. We apply the new set of parameters, coined CHARMM36m-SBFIX, to the calculation of association constants for the ammonium-acetate and guanidinium-acetate complexes, the osmotic pressures for glycine zwitterions, guanidinium, and acetate ions, and to the simulation of both folded and intrinsically disordered proteins. Our findings indicate that CHARMM36m-SBFIX improves the description of solvent-exposed salt-bridge interactions, both structurally and thermodynamically. However, application of this force field to the standard binding free-energy calculation of a protein-ligand complex featuring solvent-excluded salt-bridge interactions leads to a poor reproduction of the experimental value, suggesting that the parameters optimized in an aqueous solution cannot be readily transferred to describe solvent-excluded salt-bridge interactions. Put together, owing to their sensitivity to the environment, modeling salt-bridge interactions by means of a single, universal set of LJ parameters remains a daunting theoretical challenge.

Statistically optimal analysis of the extended-system adaptive biasing force (eABF) method

The extended-system adaptive biasing force (eABF) method and its newer variants offer rapid exploration of the configuration space of chemical systems. Instead of directly applying the ABF bias to collective variables, they are harmonically coupled to fictitious particles, which separates the problem of enhanced sampling from that of free energy estimation. The prevalent analysis method to obtain the potential of mean force (PMF) from eABF is thermodynamic integration. However, besides the PMF, most information is lost as the unbiased probability of visited configurations is never recovered. In this contribution, we show how statistical weights of individual frames can be computed using the Multistate Bennett's Acceptance Ratio (MBAR), putting the post-processing of eABF on one level with other frequently used sampling methods. In addition, we apply this formalism to the prediction of nuclear magnetic resonance shieldings, which are very sensitive to molecular geometries and often require extensive sampling. The results show that the combination of enhanced sampling by means of extended-system dynamics with the MBAR estimator is a highly useful tool for the calculation of ensemble properties. Furthermore, the extension of the presented scheme to the recently published Gaussian-accelerated molecular dynamics eABF hybrid is straightforward and approximation free.

Conformational Change from U- to I-Shape of Ion Transporters Facilitates K+ Transport across Lipid Bilayers

We have investigated, at the atomic level, the ion-fishing mechanism underlying the ion transport across membranes mediated by an artificial ion transporter composed of a hydroxyl-rich cholesterol group, a flexible alkyl chain, and a crown ether. Our results show that the transporter can spontaneously insert into the membrane and switch between the folded (U-shaped) and extended (I-shaped) conformations. The free-energy profile associated with the conformational transition indicates that compared with the U-shaped conformation of the transporter, the I-shaped one is thermodynamically more favorable. Furthermore, the free-energy profiles describing the ion translocation reveal that the transporter capturing the ion in U-shape on one side of the membrane and releasing it in I-shape on the other side constitutes a key way for the highly efficient transport of K⁺ ions. We present herewith a rigorous and rational framework to decipher the detailed ion-fishing mechanism of transmembrane ion transport with exceptionally high activity.

Computer-aided design of molecular machines: techniques, paradigms and difficulties

With their development in the past decade, molecular machines, which achieve specific tasks by responding to external stimuli, have gradually come to be regarded as powerful tools for a wide range of applications, rather than interesting molecular toys. This conceptual change in turn motivates scientists to design molecular machines with complex architectures. Due to the lack of general principles bridging the functions and the chemical structures of molecular machines, experience-based design becomes difficult with the increase of size and complexity of the architectures. Computer-aided molecular-machine design, therefore, has attracted widespread attention on account of its ability to model and investigate complex molecular architectures without too much time and expense required for synthetic experiments. Using leading-edge numerical-simulation techniques, the mechanisms underlying achieving tasks through response to external stimuli of a large number of existing molecular machines have been successfully explored. Based on the experience of studying existing molecular machines, generalized methodologies of predicting the properties and working principles of molecular candidates have been established, paving the way for de novo computer-aided design of molecular machines. In this perspective, we introduce cutting-edge techniques that have been applied for investigating and designing molecular machines. We show paradigms of computer-aided design of molecular machines, which can serve as guidelines for the investigation of new supramolecular architectures. Moreover, we discuss the limitations and possible future developments of current techniques and methodologies in the field of computer-aided design of molecular machines.

Gaussian-Accelerated Molecular Dynamics with the Weighted Ensemble Method: A Hybrid Method Improves Thermodynamic and Kinetic Sampling

Gaussian-accelerated molecular dynamics (GaMD) is a well-established enhanced sampling method for molecular dynamics simulations that effectively samples the potential energy landscape of the system by adding a boost potential, which smoothens the surface and lowers the energy barriers between states. GaMD is unable to give time-dependent properties such as kinetics directly. On the other hand, the weighted ensemble (WE) method can efficiently sample transitions between states with its many weighted trajectories, which directly yield rates and pathways. However, convergence to equilibrium conditions remains a challenge for the WE method. Hence, we have developed a hybrid method that combines the two methods, wherein GaMD is first used to sample the potential energy landscape of the system and WE is subsequently used to further sample the potential energy landscape and kinetic properties of interest. We show that the hybrid method can sample both thermodynamic and kinetic properties more accurately and quickly compared to using either method alone.