Implementation of adaptive integration method for free energy calculations in molecular systems.
Ontology highlight
ABSTRACT:
SUBMITTER: Mirabzadeh CA
PROVIDER: S-EPMC7808261 | biostudies-literature |
REPOSITORIES: biostudies-literature
Similar Datasets
Project description:Relative free energy calculations are fast becoming a critical part of early stage pharmaceutical design, making it important to know how to obtain the best performance with these calculations in applications that could span hundreds of calculations and molecules. In this work, we compared two different treatments of long-range electrostatics, Particle Mesh Ewald (PME) and Reaction Field (RF), in relative binding free energy calculations using a nonequilibrium switching protocol. We found simulations using RF achieve comparable results to those using PME but gain more efficiency when using CPU and similar performance using GPU. The results from this work encourage more use of RF in molecular simulations.
| S-EPMC8208087 | biostudies-literature
Project description:There has been a resurgence of interest in free energy methods motivated by the performance enhancements offered by molecular dynamics (MD) software written for specialized hardware, such as graphics processing units (GPUs). In this work, we exploit the properties of a parameter-interpolated thermodynamic integration (PI-TI) method to connect states by their molecular mechanical (MM) parameter values. This pathway is shown to be better behaved for Mg2+ → Ca2+ transformations than traditional linear alchemical pathways (with and without soft-core potentials). The PI-TI method has the practical advantage that no modification of the MD code is required to propagate the dynamics, and unlike with linear alchemical mixing, only one electrostatic evaluation is needed (e.g., single call to particle-mesh Ewald) leading to better performance. In the case of AMBER, this enables all the performance benefits of GPU-acceleration to be realized, in addition to unlocking the full spectrum of features available within the MD software, such as Hamiltonian replica exchange (HREM). The TI derivative evaluation can be accomplished efficiently in a post-processing step by reanalyzing the statistically independent trajectory frames in parallel for high throughput. We also show how one can evaluate the particle mesh Ewald contribution to the TI derivative evaluation without needing to perform two reciprocal space calculations. We apply the PI-TI method with HREM on GPUs in AMBER to predict p Ka values in double stranded RNA molecules and make comparison with experiments. Convergence to under 0.25 units for these systems required 100 ns or more of sampling per window and coupling of windows with HREM. We find that MM charges derived from ab initio QM/MM fragment calculations improve the agreement between calculation and experimental results.
| S-EPMC5849537 | biostudies-literature
Project description:Free energy simulations are an established computational tool in modelling chemical change in the condensed phase. However, sampling of kinetically distinct substates remains a challenge to these approaches. As a route to addressing this, we link the methods of thermodynamic integration (TI) and swarm-enhanced sampling molecular dynamics (sesMD), where simulation replicas interact cooperatively to aid transitions over energy barriers. We illustrate the approach by using alchemical alkane transformations in solution, comparing them with the multiple independent trajectory TI (IT-TI) method. Free energy changes for transitions computed by using IT-TI grew increasingly inaccurate as the intramolecular barrier was heightened. By contrast, swarm-enhanced sampling TI (sesTI) calculations showed clear improvements in sampling efficiency, leading to more accurate computed free energy differences, even in the case of the highest barrier height. The sesTI approach, therefore, has potential in addressing chemical change in systems where conformations exist in slow exchange.
| S-EPMC4676921 | biostudies-literature
Project description:Accelerated molecular dynamics (aMD) has been shown to enhance conformational space sampling relative to classical molecular dynamics; however, the exponential reweighting of aMD trajectories, which is necessary for the calculation of free energies relating to the classical system, is oftentimes problematic, especially for systems larger than small poly peptides. Here, we propose a method of accelerating only the degrees of freedom most pertinent to sampling, thereby reducing the total acceleration added to the system and improving the convergence of calculated ensemble averages, which we term selective aMD. Its application is highlighted in two biomolecular cases. First, the model system alanine dipeptide is simulated with classical MD, all-dihedral aMD, and selective aMD, and these results are compared to the infinite sampling limit as calculated with metadynamics. We show that both forms of aMD enhance the convergence of the underlying free energy landscape by 5-fold relative to classical MD; however, selective aMD can produce improved statistics over all-dihedral aMD due to the improved reweighting. Then we focus on the pharmaceutically relevant case of computing the free energy of the decoupling of oseltamivir in the active site of neuraminidase. Results show that selective aMD greatly reduces the cost of this alchemical free energy transformation, whereas all-dihedral aMD produces unreliable free energy estimates.
| S-EPMC2976571 | biostudies-literature
Project description:Free Energy Perturbation with Replica Exchange Molecular Dynamics (FEP/REMD) offers a powerful strategy to improve the convergence of free energy computations. In particular, it has been shown previously that a FEP/REMD scheme allowing random moves within an extended replica ensemble of thermodynamic coupling parameters "lambda" can improve the statistical convergence in calculations of absolute binding free energy of ligands to proteins [J. Chem. Theory Comput. 2009, 5, 2583]. In the present study, FEP/REMD is extended and combined with an accelerated MD simulations method based on Hamiltonian replica-exchange MD (H-REMD) to overcome the additional problems arising from the existence of kinetically trapped conformations within the protein receptor. In the combined strategy, each system with a given thermodynamic coupling factor lambda in the extended ensemble is further coupled with a set of replicas evolving on a biased energy surface with boosting potentials used to accelerate the inter-conversion among different rotameric states of the side chains in the neighborhood of the binding site. Exchanges are allowed to occur alternatively along the axes corresponding to the thermodynamic coupling parameter lambda and the boosting potential, in an extended dual array of coupled lambda- and H-REMD simulations. The method is implemented on the basis of new extensions to the REPDSTR module of the biomolecular simulation program CHARMM. As an illustrative example, the absolute binding free energy of p-xylene to the nonpolar cavity of the L99A mutant of T4 lysozyme was calculated. The tests demonstrate that the dual lambda-REMD and H-REMD simulation scheme greatly accelerates the configurational sampling of the rotameric states of the side chains around the binding pocket, thereby improving the convergence of the FEP computations.
| S-EPMC3157700 | biostudies-literature
Project description:Conformational changes that occur upon ligand binding may be too slow to observe on the time scales routinely accessible using molecular dynamics simulations. The adaptive integration method (AIM) leverages the notion that when a ligand is either fully coupled or decoupled, according to λ, barrier heights may change, making some conformational transitions more accessible at certain λ values. AIM adaptively changes the value of λ in a single simulation so that conformations sampled at one value of λ seed the conformational space sampled at another λ value. Adapting the value of λ throughout a simulation, however, does not resolve issues in sampling when barriers remain high regardless of the λ value. In this work, we introduce a new method, called Accelerated AIM (AcclAIM), in which the potential energy function is flattened at intermediate values of λ, promoting the exploration of conformational space as the ligand is decoupled from its receptor. We show, with both a simple model system (Bromocyclohexane) and the more complex biomolecule Thrombin, that AcclAIM is a promising approach to overcome high barriers in the calculation of free energies, without the need for any statistical reweighting or additional processors.
| S-EPMC4025579 | biostudies-literature
Project description:Galectin-1, a member of the conserved family of carbohydrate-binding proteins with affinity for ?-galactosides, is a key modulator of diverse cell functions such as immune response and regulation. The binding affinity and specificity of galectin-1 for eight different ?-galactosyl terminal disaccharides was studied using molecular-dynamics simulations in which the ligand was pulled away from the binding site using a mechanical force. We present what we believe to be a novel procedure, based on combinations of multistep trajectories, that was used to estimate the binding free energy (?G) of each disaccharide. The computed binding free energy differences show excellent correlation with experimental values determined previously. The small differences in affinity among the disaccharides are the result of an exquisite balance between the strengths of the galectin-sugar H-bonds and the H-bonds the protein and the disaccharides make with the solvent. Analysis of the free energies along the reaction coordinate shows that disaccharide unbinding/binding presents no energetic barrier and, therefore, is diffusion-limited. In addition, the calculations revealed that as the ligand is undocked from the binding site, breaking of protein-disaccharide H-bonds takes place in stages with intermediate states in which the interactions are bridged by water molecules.
| S-EPMC3149269 | biostudies-literature
Project description:The histone acetylation-deacetylation at lysine regulates the functions of many cellular proteins. An increased expression of HDAC6 can cause an increased amount of deacetylated histones, which leads to an inhibition of gene expression and has been associated with cancer cell proliferation. The present study screened the ZINC database to find novel HDAC6 inhibitors using virtual high-throughput screening techniques. The docking score, free energy, and binding pattern of the complexes were used to select a best ligand for further study. Molecular dynamic simulations, binding interactions, and the stability of docked conformations were investigated. Several parameters that determine protein-ligand interactions, such as root-mean-square deviation (RMSD), root-mean-square fluctuation (RMSF), radius of gyration (Rg), and binding pattern, were observed. Hydrogen bonds were observed at His 573 and Gly 582 after a 150 ns simulation with identified compound ZINC000002845205, and they were similar to known inhibitor Panobinostat. The molecular mechanics with generalised Born and surface area solvation (MM/GBSA) free energy was comparable to known inhibitor Panobinostat. ZINC000002845205 qualifies drug-likeness according to Lipinski's rule-of-five, rule-of-three, and the World Drug Index (WDI)-like rule, but there is one violation in the lead-like rule.
| S-EPMC9178742 | biostudies-literature
Project description:The absolute (standard) binding free energy of eight FK506-related ligands to FKBP12 is calculated using free energy perturbation molecular dynamics (FEP/MD) simulations with explicit solvent. A number of features are implemented to improve the accuracy and enhance the convergence of the calculations. First, the absolute binding free energy is decomposed into sequential steps during which the ligand-surrounding interactions as well as various biasing potentials restraining the translation, orientation, and conformation of the ligand are turned "on" and "off." Second, sampling of the ligand conformation is enforced by a restraining potential based on the root mean-square deviation relative to the bound state conformation. The effect of all the restraining potentials is rigorously unbiased, and it is shown explicitly that the final results are independent of all artificial restraints. Third, the repulsive and dispersive free energy contribution arising from the Lennard-Jones interactions of the ligand with its surrounding (protein and solvent) is calculated using the Weeks-Chandler-Andersen separation. This separation also improves convergence of the FEP/MD calculations. Fourth, to decrease the computational cost, only a small number of atoms in the vicinity of the binding site are simulated explicitly, while all the influence of the remaining atoms is incorporated implicitly using the generalized solvent boundary potential (GSBP) method. With GSBP, the size of the simulated FKBP12/ligand systems is significantly reduced, from approximately 25,000 to 2500. The computations are very efficient and the statistical error is small ( approximately 1 kcal/mol). The calculated binding free energies are generally in good agreement with available experimental data and previous calculations (within approximately 2 kcal/mol). The present results indicate that a strategy based on FEP/MD simulations of a reduced GSBP atomic model sampled with conformational, translational, and orientational restraining potentials can be computationally inexpensive and accurate.
| S-EPMC1578458 | biostudies-literature
Project description:We have calculated relative binding affinities for eight tetrafluorophenyl-triazole-thiogalactoside inhibitors of galectin-3 with the alchemical free-energy perturbation approach. We obtain a mean absolute deviation from experimental estimates of only 2-3 kJ/mol and a correlation coefficient (R2) of 0.5-0.8 for seven relative affinities spanning a range of up to 11 kJ/mol. We also studied the effect of using different methods to calculate the charges of the inhibitor and different sizes of the perturbed group (the atoms that are described by soft-core potentials and are allowed to have differing coordinates). However, the various approaches gave rather similar results and it is not possible to point out one approach as consistently and significantly better than the others. Instead, we suggest that such small and reasonable variations in the computational method can be used to check how stable the calculated results are and to obtain a more accurate estimate of the uncertainty than if performing only one calculation with a single computational setup.
| S-EPMC5889414 | biostudies-literature