Treating Entropy and Conformational Changes in Implicit Solvent Simulations of Small Molecules

Implicit solvent models are increasingly popular for estimating aqueous solvation (hydration) free energies in molecular simulations and other applications. In many cases, parameters for these models are derived to reproduce experimental values for small molecule hydration free energies. Often, these hydration free energies are computed for a single solute conformation, neglecting solute conformational changes upon solvation. Here, we incorporate these effects using alchemical free energy methods. We find significant errors when hydration free energies are estimated using only a single solute conformation, even for relatively small, simple, rigid solutes. For example, we find conformational entropy (TΔS) changes of up to 2.3 kcal/mol upon hydration. Interestingly, these changes in conformational entropy correlate poorly (R2 = 0.03) with the number of rotatable bonds. The present study illustrates that implicit solvent modeling can be improved by eliminating the approximation that solutes are rigid

Confine-and-Release Method: Obtaining Correct Binding Free Energies in the Presence of Protein Conformational Change

Free energy calculations are increasingly being used to estimate absolute and relative binding free energies of ligands to proteins. However, computed free energies often appear to depend on the initial protein conformation, indicating incomplete sampling. This is especially true when proteins can change conformation on ligand binding, as free energies associated with these conformational changes are either ignored or assumed to be included by virtue of the sampling performed in the calculation. Here, we show that, in a model protein system (a designed binding site in T4 lysozyme), conformational changes can make a difference of several kcal/mol in computed binding free energies and that they are neglected in computed binding free energies if the system remains kinetically trapped in a particular metastable state on simulation timescales. We introduce a general “confine-and-release” framework for free energy calculations that accounts for these free energies of conformational change. We illustrate its use in this model system by demonstrating that an umbrella sampling protocol can obtain converged binding free energies that are independent of the starting protein structure and include these conformational change free energies

Treating Entropy and Conformational Changes in Implicit Solvent Simulations of Small Molecules

Implicit solvent models are increasingly popular for estimating aqueous solvation (hydration) free energies in molecular simulations and other applications. In many cases, parameters for these models are derived to reproduce experimental values for small molecule hydration free energies. Often, these hydration free energies are computed for a single solute conformation, neglecting solute conformational changes upon solvation. Here, we incorporate these effects using alchemical free energy methods. We find significant errors when hydration free energies are estimated using only a single solute conformation, even for relatively small, simple, rigid solutes. For example, we find conformational entropy (TΔS) changes of up to 2.3 kcal/mol upon hydration. Interestingly, these changes in conformational entropy correlate poorly (R2 = 0.03) with the number of rotatable bonds. The present study illustrates that implicit solvent modeling can be improved by eliminating the approximation that solutes are rigid

Absolute Binding Free Energy Calculations for Buried Water Molecules

Water often plays a key role in mediating protein–ligand interactions. Understanding contributions from active-site water molecules to binding thermodynamics of a ligand is important in predicting binding free energies for ligand optimization. In this work, we tested a non-equilibrium switching method for absolute binding free energy calculations on water molecules in binding sites of 13 systems. We discuss the lessons we learned about identified issues that affected our calculations and ways to address them. This work fits with our larger focus on how to do accurate ligand binding free energy calculations when water rearrangements are very slow, such as rearrangements due to ligand modification (as in relative free energy calculations) or ligand binding (as in absolute free energy calculations). The method studied in this work can potentially be used to account for limited water sampling via providing endpoint corrections to free energy calculations using our calculated binding free energy of water

Small Molecule Solvation Free Energy: Enhanced Conformational Sampling Using Expanded Ensemble Molecular Dynamics Simulation

We present an efficient expanded ensemble molecular dynamics method to calculate the solvation free energy (or residual chemical potential) of small molecules with complex topologies. The methodology is validated by computing the solvation free energy of ibuprofen in water, methanol, and ethanol at 300 K and 1 bar and comparing to reference simulation results using Bennett’s acceptance ratio method. Difficulties with ibuprofen using conventional molecular dynamics methods stem from an inadequate sampling of the carboxylic acid functional group, which, for the present study, is subject to free energy barriers of rotation of 14–20 kBT. While several advances have been made to overcome such weaknesses, we demonstrate how this shortcoming is easily overcome by using an expanded ensemble methodology to facilitate conformational sampling. Not only does the method enhance conformational sampling but it also boosts the rate of exploration of the configurational phase space and requires only a single simulation to calculate the solvation free energy. Agreement between the expanded ensemble and the reference calculations is good for all three solvents, with the reported uncertainties of the expanded ensemble being comparable to the uncertainties of the reference calculations, while requiring less simulation time; the reduced simulation time demonstrates the improved performance gained from the expanded ensemble method

To Design Scalable Free Energy Perturbation Networks, Optimal Is Not Enough

Drug discovery is accelerated with computational methods such as alchemical simulations to estimate ligand affinities. In particular, relative binding free energy (RBFE) simulations are beneficial for lead optimization. To use RBFE simulations to compare prospective ligands in silico, researchers first plan the simulation experiment, using graphs where nodes represent ligands and graph edges represent alchemical transformations between ligands. Recent work demonstrated that optimizing the statistical architecture of these perturbation graphs improves the accuracy of the predicted changes in the free energy of ligand binding. Therefore, to improve the success rate of computational drug discovery, we present the open-source software package High Information Mapper (HiMap)a new take on its predecessor, Lead Optimization Mapper (LOMAP). HiMap removes heuristics decisions from design selection and instead finds statistically optimal graphs over ligands clustered with machine learning. Beyond optimal design generation, we present theoretical insights for designing alchemical perturbation maps. Some of these results include that for n number of nodes, the precision of perturbation maps is stable at n·ln­(n) edges. This result indicates that even an “optimal” graph can result in unexpectedly high errors if a plan includes too few alchemical transformations for the given number of ligands and edges. And, as a study compares more ligands, the performance of even optimal graphs will deteriorate with linear scaling of the edge count. In this sense, ensuring an A- or D-optimal topology is not enough to produce robust errors. We additionally find that optimal designs will converge more rapidly than radial and LOMAP designs. Moreover, we derive bounds for how clustering reduces cost for designs with a constant expected relative error per cluster, invariant of the size of the design. These results inform how to best design perturbation maps for computational drug discovery and have broader implications for experimental design

To Design Scalable Free Energy Perturbation Networks, Optimal Is Not Enough

Drug discovery is accelerated with computational methods such as alchemical simulations to estimate ligand affinities. In particular, relative binding free energy (RBFE) simulations are beneficial for lead optimization. To use RBFE simulations to compare prospective ligands in silico, researchers first plan the simulation experiment, using graphs where nodes represent ligands and graph edges represent alchemical transformations between ligands. Recent work demonstrated that optimizing the statistical architecture of these perturbation graphs improves the accuracy of the predicted changes in the free energy of ligand binding. Therefore, to improve the success rate of computational drug discovery, we present the open-source software package High Information Mapper (HiMap)a new take on its predecessor, Lead Optimization Mapper (LOMAP). HiMap removes heuristics decisions from design selection and instead finds statistically optimal graphs over ligands clustered with machine learning. Beyond optimal design generation, we present theoretical insights for designing alchemical perturbation maps. Some of these results include that for n number of nodes, the precision of perturbation maps is stable at n·ln­(n) edges. This result indicates that even an “optimal” graph can result in unexpectedly high errors if a plan includes too few alchemical transformations for the given number of ligands and edges. And, as a study compares more ligands, the performance of even optimal graphs will deteriorate with linear scaling of the edge count. In this sense, ensuring an A- or D-optimal topology is not enough to produce robust errors. We additionally find that optimal designs will converge more rapidly than radial and LOMAP designs. Moreover, we derive bounds for how clustering reduces cost for designs with a constant expected relative error per cluster, invariant of the size of the design. These results inform how to best design perturbation maps for computational drug discovery and have broader implications for experimental design

Fragment Pose Prediction Using Non-equilibrium Candidate Monte Carlo and Molecular Dynamics Simulations

Part of early stage drug discovery involves determining how molecules may bind to the target protein. Through understanding where and how molecules bind, chemists can begin to build ideas on how to design improvements to increase binding affinities. In this retrospective study, we compare how computational approaches like docking, molecular dynamics (MD) simulations, and a non-equilibrium candidate Monte Carlo (NCMC)-based method (NCMC + MD) perform in predicting binding modes for a set of 12 fragment-like molecules, which bind to soluble epoxide hydrolase. We evaluate each method’s effectiveness in identifying the dominant binding mode and finding additional binding modes (if any). Then, we compare our predicted binding modes to experimentally obtained X-ray crystal structures. We dock each of the 12 small molecules into the apo-protein crystal structure and then run simulations up to 1 μs each. Small and fragment-like molecules likely have smaller energy barriers separating different binding modes by virtue of relatively fewer and weaker interactions relative to drug-like molecules and thus likely undergo more rapid binding mode transitions. We expect, thus, to see more rapid transitions between binding modes in our study. Following this, we build Markov State Models to define our stable ligand binding modes. We investigate if adequate sampling of ligand binding modes and transitions between them can occur at the microsecond timescale using traditional MD or a hybrid NCMC+MD simulation approach. Our findings suggest that even with small fragment-like molecules, we fail to sample all the crystallographic binding modes using microsecond MD simulations, but using NCMC+MD, we have better success in sampling the crystal structure while obtaining the correct populations

Fragment Pose Prediction Using Non-equilibrium Candidate Monte Carlo and Molecular Dynamics Simulations

Part of early stage drug discovery involves determining how molecules may bind to the target protein. Through understanding where and how molecules bind, chemists can begin to build ideas on how to design improvements to increase binding affinities. In this retrospective study, we compare how computational approaches like docking, molecular dynamics (MD) simulations, and a non-equilibrium candidate Monte Carlo (NCMC)-based method (NCMC + MD) perform in predicting binding modes for a set of 12 fragment-like molecules, which bind to soluble epoxide hydrolase. We evaluate each method’s effectiveness in identifying the dominant binding mode and finding additional binding modes (if any). Then, we compare our predicted binding modes to experimentally obtained X-ray crystal structures. We dock each of the 12 small molecules into the apo-protein crystal structure and then run simulations up to 1 μs each. Small and fragment-like molecules likely have smaller energy barriers separating different binding modes by virtue of relatively fewer and weaker interactions relative to drug-like molecules and thus likely undergo more rapid binding mode transitions. We expect, thus, to see more rapid transitions between binding modes in our study. Following this, we build Markov State Models to define our stable ligand binding modes. We investigate if adequate sampling of ligand binding modes and transitions between them can occur at the microsecond timescale using traditional MD or a hybrid NCMC+MD simulation approach. Our findings suggest that even with small fragment-like molecules, we fail to sample all the crystallographic binding modes using microsecond MD simulations, but using NCMC+MD, we have better success in sampling the crystal structure while obtaining the correct populations

Comparison of Charge Models for Fixed-Charge Force Fields: Small Molecule Hydration Free Energies in Explicit Solvent

Comparison of Charge Models for Fixed-Charge Force Fields: Small Molecule Hydration Free Energies in Explicit Solven