Multidimensional Potential Of Mean Force Research Articles

Molecular transport through membranes: Accurate permeability coefficients from multidimensional potentials of mean force and local diffusion constants.

Estimating the permeability coefficient of small molecules through lipid bilayer membranes plays an important role in the development of effective drug candidates. In silico simulations can produce acceptable relative permeability coefficients for a series of small molecules; however, the absolute permeability coefficients from simulations are usually off by orders of magnitude. In addition to differences between the lipid bilayers used in vitro and in silico, the poor convergence of permeation free energy profiles and over-simplified diffusion models have contributed to these discrepancies. In this paper, we present a multidimensional inhomogeneous solubility-diffusion model to study the permeability of a small molecule drug (trimethoprim) passing through a POPC (1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine) lipid bilayer. Our approach improves the permeation model in three ways: First, the free energy profile (potential of mean force, PMF) is two-dimensional in two key coordinates rather than simply one-dimensional along the direction normal to the bilayer. Second, the 2-D PMF calculation has improved convergence due to application of the recently developed transition-tempered metadynamics with randomly initialized replicas, while third, the local diffusivity coefficient was calculated along the direction of the minimum free energy path on the two-dimensional PMF. The permeability is then calculated as a line integral along the minimum free energy path of the PMF. With this approach, we report a considerably more accurate permeability coefficient (only 2-5 times larger than the experimental value). We also compare our approach with the common practice of computing permeability coefficients based only on the translation of the center of mass of the drug molecule. Our paper concludes with a discussion of approaches for minimizing the computational cost for the purpose of more rapidly screening a large number of drug candidate molecules.

ChemInform Abstract: Quantum and Molecular Mechanical Monte Carlo Techniques for Modeling Condensed‐phase Reactions

AbstractReview: 98 refs.

Quantum and Molecular Mechanical (QM/MM) Monte Carlo Techniques for Modeling Condensed-Phase Reactions.

A recent review (Acc. Chem. Res. 2010, 43:142-151) examined our use and development of a combined quantum and molecular mechanical (QM/MM) technique for modelling organic and enzymatic reactions. Advances included the PDDG/PM3 semiempirical QM (SQM) method, computation of multi-dimensional potentials of mean force (PMF), incorporation of on-the-fly QM in Monte Carlo simulations, and a polynomial quadrature method for rapidly treating proton-transfer reactions. The current article serves as a follow up on our progress. Highlights include new reactions, alternative SQM methods, a polarizable OPLS force field, and novel solvent environments, e.g., "on water" and room temperature ionic liquids. The methodology is strikingly accurate across a wide range of condensed-phase and antibody-catalyzed reactions including substitution, decarboxylation, elimination, isomerization, and pericyclic classes. Comparisons are made to systems treated with continuum-based solvents and ab initio or density functional theory (DFT) methods. Overall, the QM/MM methodology provides detailed characterization of reaction paths, proper configurational sampling, several advantages over implicit solvent models, and a reasonable computational cost.

Quantum and Classical Simulations of Orotidine Monophosphate Decarboxylase: Support for a Direct Decarboxylation Mechanism

Orotidine 5'-monophosphate (OMP) decarboxylase (ODCase) catalyzes the decarboxylation of OMP to uridine 5'-monophosphate (UMP). Numerous studies of this reaction have suggested a plethora of mechanisms including covalent addition, ylide or carbene formation, and concerted or stepwise protonation. Recent experiments and simulations present strong evidence for a direct decarboxylation mechanism, although direct comparison between experiment and theory is still lacking. In the current work we present hybrid quantum mechanics-molecular mechanics simulations that address the detailed decarboxylation mechanisms for OMP and 5-fluoro-OMP by ODCase. Multidimensional potentials of mean force are computed as functions of structural progress coordinates for the Methanobacterium thermoautotrophicum ODCase reaction: the decarboxylation reaction coordinate, an orbital rehybridization coordinate, and the proton transfer coordinate between Lys72 and the substrate. The computed free energy profiles are in accord with the available experimental data. To facilitate further direct comparison with experiment, we compute the kinetic isotope effects (KIEs) for the enzyme-catalyzed reactions using a mass-perturbation-based path-integral method. The computed KIE provide further support for a direct decarboxylation mechanism. In agreement with experiment, the data suggest a role for Lys72 in stabilizing the transition state in the catalysis of OMP and, to a somewhat lesser extent, in 5-fluoro-OMP.

Self-Learning Adaptive Umbrella Sampling Method for the Determination of Free Energy Landscapes in Multiple Dimensions.

The potential of mean force describing conformational changes of biomolecules is a central quantity that determines the function of biomolecular systems. Calculating an energy landscape of a process that depends on three or more reaction coordinates might require a lot of computational power, making some of multidimensional calculations practically impossible. Here, we present an efficient automatized umbrella sampling strategy for calculating multidimensional potential of mean force. The method progressively learns by itself, through a feedback mechanism, which regions of a multidimensional space are worth exploring and automatically generates a set of umbrella sampling windows that is adapted to the system. The self-learning adaptive umbrella sampling method is first explained with illustrative examples based on simplified reduced model systems, and then applied to two non-trivial situations: the conformational equilibrium of the pentapeptide Met-enkephalin in solution and ion permeation in the KcsA potassium channel. With this method, it is demonstrated that a significant smaller number of umbrella windows needs to be employed to characterize the free energy landscape over the most relevant regions without any loss in accuracy.

Studies of Mg2+/Ca2+ complexes of naturally occurring dinucleotides: potentiometric titrations, NMR, and molecular dynamics

Dinucleotides (Np(n)N'; N and N' are A, U, G, or C, n = 2-7) are naturally occurring physiologically active compounds. Despite the interest in dinucleotides, the composition of their complexes with metal ions as well as their conformations and species distribution in living systems are understudied. Therefore, we investigated a series of Mg(2+) and Ca(2+) complexes of Np(n)N's. Potentiometric titrations indicated that a longer dinucleotide polyphosphate (N is A or G, n = 3-5) linker yields more stable complexes (e.g., log K of 2.70, 3.27, and 3.73 for Ap(n)A-Mg(2+), n = 3, 4, 5, respectively). The base (A or G) or ion (Mg(2+) or Ca(2+)) has a minor effect on K(M)(ML) values. In a physiological medium, the longer Ap(n)As (n = 4, 5) are predicted to occur mostly as the Mg(2+)/Ca(2+) complexes. (31)P NMR monitored titrations of Np(n)N's with Mg(2+)/Ca(2+) ions showed that the middle phosphates of the dinucleotides coordinate with Mg(2+)/Ca(2+). Multidimensional potential of mean force (PMF) molecular dynamics (MD) simulations suggest that Ap(2)A and Ap(4)A coordinate Mg(2+) and Ca(2+) ions in both inner-sphere and outer-sphere modes. The PMF MD simulations additionally provide a detailed picture of the possible coordination sites, as well as the cation binding process. Moreover, both NMR and MD simulations showed that the conformation of the nucleoside moieties in Np(n)N'-Mg(2+)/Ca(2+) complexes remains the same as that of free mononucleotides.

Automated Umbrella Sampling Simulations for the Calculation of Multidimensional Potentials of Mean Force

The potential of mean force (PMF) describing the conformational changes of biomolecules is a central quantity that determines the function of biomolecular systems. Calculating a multidimensional PMF (free energy landscape) is a time consuming process. Each additional reaction coordinate drastically increases the required simulation time for conformation sampling, making calculations in four or more dimensions practically impossible. However, in most cases, only a small fraction of such multidimensional space is energetically relevant. PMF calculations in high dimensionality could thus be achieved if one could effectively focus the simulation effort on those region of high interest. We have developed a method based on umbrella sampling (US) that determines, using a feedback mechanism, which regions of the multidimensional space is worth exploring. The first aim of such application is to manage the creation and analysis of sampling windows for the calculation of PMFs involving up to four reaction coordinates. While this approach could in some cases be used to find the minimum free energy pathway underlying large conformational changes, it is not as general as other approaches that were specifically developed in that purpose. However, contrary to most other approaches, our method allows for the simultaneous characterization of several pathways, and not only the most probable one. The current implementation consists in a C++ application with a python interface that manages MD simulations performed with the biomolecular simulation program CHARMM. The feedback mechanism and final PMF calculation involves the unbiasing of simulations using the weighted-histogram analysis methods (WHAM).

Free energy surfaces for the interaction of D-glucose with planar aromatic groups in aqueous solution.

Multidimensional potentials of mean force for the interactions in aqueous solution of both anomers of D-glucopyranose with two planar aromatic molecules, indole and para-methyl-phenol, have been calculated using molecular dynamics simulations with umbrella sampling and were subsequently used to estimate binding free energies. Indole and para-methyl-phenol serve as models for the side chains of the amino acids tryptophan and tyrosine, respectively. In all cases, a weak affinity between the glucose molecules and the flat aromatic surfaces was found. The global minimum for these interactions was found to be for the case when the pseudoplanar face of β-D-glucopyranose is stacked against the planar surfaces of the aromatic residues. The calculated binding free energies are in good agreement with both experiment and previous simulations. The multidimensional free energy maps suggest a mechanism that could lend kinetic stability to the complexes formed by sugars bound to sugar-binding proteins.

Advances in Quantum and Molecular Mechanical (QM/MM) Simulations for Organic and Enzymatic Reactions

Application of combined quantum and molecular mechanical (QM/MM) methods focuses on predicting activation barriers and the structures of stationary points for organic and enzymatic reactions. Characterization of the factors that stabilize transition structures in solution and in enzyme active sites provides a basis for design and optimization of catalysts. Continued technological advances allowed for expansion from prototypical cases to mechanistic studies featuring detailed enzyme and condensed-phase environments with full integration of the QM calculations and configurational sampling. This required improved algorithms featuring fast QM methods, advances in computing changes in free energies including free-energy perturbation (FEP) calculations, and enhanced configurational sampling. In particular, the present Account highlights development of the PDDG/PM3 semi-empirical QM method, computation of multi-dimensional potentials of mean force (PMF), incorporation of on-the-fly QM in Monte Carlo (MC) simulations, and a polynomial quadrature method for efficient modeling of proton-transfer reactions. The utility of this QM/MM/MC/FEP methodology is illustrated for a variety of organic reactions including substitution, decarboxylation, elimination, and pericyclic reactions. A comparison to experimental kinetic results on medium effects has verified the accuracy of the QM/MM approach in the full range of solvents from hydrocarbons to water to ionic liquids. Corresponding results from ab initio and density functional theory (DFT) methods with continuum-based treatments of solvation reveal deficiencies, particularly for protic solvents. Also summarized in this Account are three specific QM/MM applications to biomolecular systems: (1) a recent study that clarified the mechanism for the reaction of 2-pyrone derivatives catalyzed by macrophomate synthase as a tandem Michael-aldol sequence rather than a Diels-Alder reaction, (2) elucidation of the mechanism of action of fatty acid amide hydrolase (FAAH), an unusual Ser-Ser-Lys proteolytic enzyme, and (3) the construction of enzymes for Kemp elimination of 5-nitrobenzisoxazole that highlights the utility of QM/MM in the design of artificial enzymes.

Investigation of the Regulatory Mechanism of the ZAP-70 Immunological Signaling Enzyme

ZAP-70 is a critical Syk-family protein tyrosine kinase (PTK) that functions in the initial step of T-cell receptor (TCR) signaling. Its importance is highlighted by mutations in ZAP-70 that cause severe combined immunodeficiency (SCID). Moreover, elevated levels of ZAP-70 are associated with T-cell proliferative diseases and chronic lymphocytic leukemia (CLL). Under normal conditions, T-cell activation by foreign antigen results in co-localization with ZAP-70, due to binding of its associated tandem SH2 domains (tSH2) to phosphorylated ITAM (Immunoreceptor Tyrosine Activation Motif) sequences present on the TCR ζ subunits. This leads to up-regulation of kinase activity and downstream initiation of an immune response. Experimental structures reveal that upon ITAM binding, the tSH2 must undergo a large conformational change. Of interest, therefore, is how ITAM-binding controls the tSH2 conformational equilibrium, and the subsequent proposed disassembly of the kinase domain. To this end, we have carried out 1.5 μs of molecular dynamics (MD) simulations of various ZAP-70 systems, composed of 50,000 to 70,000 atoms. These include the isolated tSH2, and variants of the ITAM-bound state. Collectively, our results suggest that one phosphotyrosine site regulates the tSH2 conformational equilibrium, whilst the other controls localization of the tSH2 to membrane-proximal ITAM targets. Using these results, along with those from targeted MD simulations, we have identified a set of order parameters that describe the tSH2 conformational transition. We have subsequently performed multidimensional potential of mean force (PMF) calculations, via umbrella-sampling, to describe the transition; a single such PMF for the tSH2 “switch” constitutes over 3 μs of simulation time. Thus, we have achieved an atomic level characterization of the free-energy landscape of the tSH2 in the presence and absence of the catalytic domain, with important consequences for the mechanism of signaling and autoinhibition in Syk-family kinases.

Multidimensional Potentials of Mean Force from Biased Experiments along a Single Coordinate

External biasing forces are often applied to enhance sampling in regions of phase space which would otherwise be rarely observed. While the typical goal of these experiments is to calculate the potential of mean force (PMF) along the biasing coordinate, here I present a method to construct PMFs in multiple dimensions and along arbitary alternative degrees of freedom. A protocol for multidimensional PMF reconstruction from nonequilibrium single-molecule pulling experiments is introduced and tested on a series of two-dimensional potential surfaces with varying levels of correlation. Reconstruction accuracy and convergence from several methods--this new protocol, equilibrium umbrella sampling, and free diffusion--are compared, and nonequilibrium pulling is found to be the most efficient. To facilitate the use of this method, the source code for this analysis is made freely available.

Multiscale coarse-graining and structural correlations: connections to liquid-state theory.

A statistical mechanical framework elucidates the significance of structural correlations between coarse-grained (CG) sites in the multiscale coarse-graining (MS-CG) method (Izvekov, S.; Voth, G. A. J. Phys. Chem. B 2005, 109, 2469; J. Chem. Phys. 2005, 123, 134105). If no approximations are made, then the MS-CG method yields a many-body multidimensional potential of mean force describing the interactions between CG sites. However, numerical applications of the MS-CG method typically employ a set of pair potentials to describe nonbonded interactions. The analogy between coarse-graining and the inverse problem of liquid-state theory clarifies the general significance of three-particle correlations for the development of such CG pair potentials. It is demonstrated that the MS-CG methodology incorporates critical three-body correlation effects and that, for isotropic homogeneous systems evolving under a central pair potential, the MS-CG equations are a discretized representation of the well-known Yvon-Born-Green equation. Numerical calculations validate the theory and illustrate the role of these structural correlations in the MS-CG method.

Direct determination of reaction paths and stationary points on potential of mean force surfaces

A simulation approach is introduced for directly determining reaction paths and stationary points on potential of mean force (PMF) surfaces associated with molecular events that occur in complex environments. The nudged elastic band approach was employed to search for steepest descent paths on the PMF surface using the relevant PMF derivatives from a series of local simulations. The steepest descent path on the PMF surface corresponds to the minimum PMF path (i.e. the path with the lowest local PMF barrier), which contains important information about stationary points (e.g. saddle points) on the PMF surface, which in turn can provide useful insights into the thermodynamics and kinetics for the process of interest. By working with the PMF defined in a low-dimensional subspace rather than a potential energy function of full molecular dimensionality, the main features of the process under study are concisely represented and the orthogonal degrees of freedom are adequately sampled with the appropriate canonical distribution at the desired temperature (e.g. 300 K). Therefore, minimum PMF paths carry statistically meaningful mechanistic information and are complementary to reaction paths of full molecular dimensionality proposed in previous studies. The NEB based path optimization method is direct in the sense that no information regarding the global PMF surface is necessary for the determination of the local reaction path and stationary points along this path. Since only low-dimensional quantities (paths) are searched for, the PMF-path method is expected to scale better in terms of dimension of the PMF subspace than those aims to fully explore multi-dimensional PMF surfaces. Test applications on simple molecular systems, the alanine di-peptide in vacuum and in solution and a microsolvated proton-wire, indicate that reliable PMF paths can be determined for both conformational isomerization and chemical reaction processes. However, highly accurate PMF derivatives are required for determining more quantitative observables, such as the free energy profile along the minimum PMF path. Therefore, effective numerical algorithms for calculating local PMF derivatives and systematic protocols for defining the relevant subspace are the main focus in the near future. Finally, we emphasize that the minimum PMF path defined here includes thermal (e.g. entropic) effects associated with the orthogonal degrees of freedom, but finite kinetic energies associated with the PMF degrees of freedom are not included; this can be improved by adopting a different definition of the reaction path, such as the maximum flux path, on the PMF surface, or thermally sampling all degrees of freedom orthogonal to the one-dimensional path.

Implementation of an adaptive umbrella sampling method for the calculation of multidimensional potential of mean force of chemical reactions in solution.

We describe the implementation of an adaptive umbrella sampling method, making use of the weighted histogram analysis method, for computing multidimensional potential of mean force for chemical reaction in solution. The approach is illustrated by investigating the effect of aqueous solution on the free energy surface for the proton transfer reaction of [H(3)N-H-NH(3)](+) using a combined quantum mechanical and molecular mechanical AM1/TIP3P potential.