Roadmaps through free energy landscapes calculated using the multi-dimensional vFEP approach

Unveiling Cas8 dynamics and regulation within a transposon-encoded Cascade-TniQ complex

The Vibrio cholerae Cascade-TniQ complex unveiled a new paradigm in biology, demonstrating that CRISPR-associated proteins can direct DNA transposition. Despite the tremendous potential of "knocking-in" genes at desired sites, the mechanisms underlying DNA binding and transposition remain elusive. In this system, a conformational change of the Cas8 protein is essential for DNA binding, yet how it occurs is unclear. Here, structural modeling and free energy simulations reconstruct the Cas8 helical bundle and reveal an open-closed conformational change that is key for the complex's function. We show that when Cascade-TniQ binds RNA, the Cas8 bundle changes conformation mediated by the interaction with the Cas7.1 protein. This interaction promotes the bundle's transition toward the open state, priming the complex for DNA binding. As the target DNA binds the guide RNA, the opening of the Cas8 bundle becomes more favorable, exposing positively charged residues and facilitating their interaction with DNA, which ultimately leads the DNA-binding process to completion. These outcomes provide a dynamic representation of a critical conformational change in one of the largest CRISPR systems and illustrate its role at critical steps of the Cascade-TniQ biophysical function, advancing our understanding of nucleic acid binding and transposition mechanisms.

Amber free energy tools: Interoperable software for free energy simulations using generalized quantum mechanical/molecular mechanical and machine learning potentials

We report the development and testing of new integrated cyberinfrastructure for performing free energy simulations with generalized hybrid quantum mechanical/molecular mechanical (QM/MM) and machine learning potentials (MLPs) in Amber. The Sander molecular dynamics program has been extended to leverage fast, density-functional tight-binding models implemented in the DFTB+ and xTB packages, and an interface to the DeePMD-kit software enables the use of MLPs. The software is integrated through application program interfaces that circumvent the need to perform "system calls" and enable the incorporation of long-range Ewald electrostatics into the external software's self-consistent field procedure. The infrastructure provides access to QM/MM models that may serve as the foundation for QM/MM-ΔMLP potentials, which supplement the semiempirical QM/MM model with a MLP correction trained to reproduce ab initio QM/MM energies and forces. Efficient optimization of minimum free energy pathways is enabled through a new surface-accelerated finite-temperature string method implemented in the FE-ToolKit package. Furthermore, we interfaced Sander with the i-PI software by implementing the socket communication protocol used in the i-PI client-server model. The new interface with i-PI allows for the treatment of nuclear quantum effects with semiempirical QM/MM-ΔMLP models. The modular interoperable software is demonstrated on proton transfer reactions in guanine-thymine mispairs in a B-form deoxyribonucleic acid helix. The current work represents a considerable advance in the development of modular software for performing free energy simulations of chemical reactions that are important in a wide range of applications.

Surface-Accelerated String Method for Locating Minimum Free Energy Paths

We present a surface-accelerated string method (SASM) to efficiently optimize low-dimensional reaction pathways from the sampling performed with expensive quantum mechanical/molecular mechanical (QM/MM) Hamiltonians. The SASM accelerates the convergence of the path using the aggregate sampling obtained from the current and previous string iterations, whereas approaches like the string method in collective variables (SMCV) or the modified string method in collective variables (MSMCV) update the path only from the sampling obtained from the current iteration. Furthermore, the SASM decouples the number of images used to perform sampling from the number of synthetic images used to represent the path. The path is optimized on the current best estimate of the free energy surface obtained from all available sampling, and the proposed set of new simulations is not restricted to being located along the optimized path. Instead, the umbrella potential placement is chosen to extend the range of the free energy surface and improve the quality of the free energy estimates near the path. In this manner, the SASM is shown to improve the exploration for a minimum free energy pathway in regions where the free energy surface is relatively flat. Furthermore, it improves the quality of the free energy profile when the string is discretized with too few images. We compare the SASM, SMCV, and MSMCV using 3 QM/MM applications: a ribozyme methyltransferase reaction using 2 reaction coordinates, the 2'-O-transphosphorylation reaction of Hammerhead ribozyme using 3 reaction coordinates, and a tautomeric reaction in B-DNA using 5 reaction coordinates. We show that SASM converges the paths using roughly 3 times less sampling than the SMCV and MSMCV methods. All three algorithms have been implemented in the FE-ToolKit package made freely available.

Modern Alchemical Free Energy Methods for Drug Discovery Explained

This Perspective provides a contextual explanation of the current state-of-the-art alchemical free energy methods and their role in drug discovery as well as highlights select emerging technologies. The narrative attempts to answer basic questions about what goes on "under the hood" in free energy simulations and provide general guidelines for how to run simulations and analyze the results. It is the hope that this work will provide a valuable introduction to students and scientists in the field.

AmberTools

AmberTools is a free and open-source collection of programs used to set up, run, and analyze molecular simulations. The newer features contained within AmberTools23 are briefly described in this Application note.

Polyaromatic hydrocarbon thin film layers on glass, dust, and polyurethane foam surfaces

An investigation was conducted into the dynamic behavior of two polyaromatic hydrocarbon (PAH) semi-volatile organic compound (SVOC) naphthalene (NAP) and benzo [ghi]perylene (BghiP) in air and on various surfaces including glass, dust, and polyurethane foam (PUF) to understand their interaction with different media. A confocal fluorescence microscope and an infrared microscope were employed to detect and monitor the concentration-, time-, and temperature-dependent changes of the aromatic NAP and BghiP species on the surfaces. Infrared two-dimensional mapping of the vibrational characteristic peaks was used to track the two PAHs on the surfaces. Gas chromatography-mass spectrometry (GC-MS) was employed to measure the gaseous concentrations. The sorption of NAP and BghiP on the surfaces was estimated using Arizona desert sand fine (ISO 12103-1 A2) dust and organic contaminant household (SRM 2585) dust. The surface-to-air partition coefficients of NAP and BghiP were estimated on the different surfaces of glass, dust, and PUF. Molecular dynamic simulations were performed on dust surfaces based on the Hatcher model to understand the behavior of NAP and BghiP on dust surfaces. The Weschler-Nazaroff model was introduced to predictPAH film accumulation on the surfaces, providing a better understanding of PAH interaction with different environmental media. These findings could contribute to developing effective strategies to mitigate the adverse impact of PAHs on the environment and human health.

Estimation of frequency factors for the calculation of kinetic isotope effects from classical and path integral free energy simulations

We use the modified Bigeleisen-Mayer equation to compute kinetic isotope effect values for non-enzymatic phosphoryl transfer reactions from classical and path integral molecular dynamics umbrella sampling. The modified form of the Bigeleisen-Mayer equation consists of a ratio of imaginary mode vibrational frequencies and a contribution arising from the isotopic substitution's effect on the activation free energy, which can be computed from path integral simulation. In the present study, we describe a practical method for estimating the frequency ratio correction directly from umbrella sampling in a manner that does not require normal mode analysis of many geometry optimized structures. Instead, the method relates the frequency ratio to the change in the mass weighted coordinate representation of the minimum free energy path at the transition state induced by isotopic substitution. The method is applied to the calculation of 16/18O and 32/34S primary kinetic isotope effect values for six non-enzymatic phosphoryl transfer reactions. We demonstrate that the results are consistent with the analysis of geometry optimized transition state ensembles using the traditional Bigeleisen-Mayer equation. The method thus presents a new practical tool to enable facile calculation of kinetic isotope effect values for complex chemical reactions in the condensed phase.

Altered Mechanisms for Acid-Catalyzed RNA Cleavage and Isomerization Reactions Models

RNA strand cleavage by 2'-O-transphosphorylation is catalyzed not only by numerous nucleolytic RNA enzymes (ribozymes) but also by hydroxide or hydronium ions. In experiments, both cleavage of the 5'-linked nucleoside and isomerization between 3',5'- and 2',5'-phosphodiesters occur under acidic conditions, while only the cleavage reaction is observed under basic conditions. An ab initio path-integral approach for simulating kinetic isotope effects is used to reveal the reaction mechanisms for RNA cleavage and isomerization reactions under acidic conditions. Moreover, the proposed mechanisms can also be combined through the experimental pH-rate profiles.

ACES: Optimized Alchemically Enhanced Sampling

We present an alchemical enhanced sampling (ACES) method implemented in the GPU-accelerated AMBER free energy MD engine. The methods hinges on the creation of an "enhanced sampling state" by reducing or eliminating selected potential energy terms and interactions that lead to kinetic traps and conformational barriers while maintaining those terms that curtail the need to otherwise sample large volumes of phase space. For example, the enhanced sampling state might involve transforming regions of a ligand and/or protein side chain into a noninteracting "dummy state" with internal electrostatics and torsion angle terms turned off. The enhanced sampling state is connected to a real-state end point through a Hamiltonian replica exchange (HREMD) framework that is facilitated by newly developed alchemical transformation pathways and smoothstep softcore potentials. This creates a counterdiffusion of real and enhanced-sampling states along the HREMD network. The effect of a differential response of the environment to the real and enhanced-sampling states is minimized by leveraging the dual topology framework in AMBER to construct a counterbalancing HREMD network in the opposite alchemical direction with the same (or similar) real and enhanced sampling states at inverted end points. The method has been demonstrated in a series of test cases of increasing complexity where traditional MD, and in several cases alternative REST2-like enhanced sampling methods, are shown to fail. The hydration free energy for acetic acid was shown to be independent of the starting conformation, and the values for four additional edge case molecules from the FreeSolv database were shown to have a significantly closer agreement with experiment using ACES. The method was further able to handle different rotamer states in a Cdk2 ligand identified as fractionally occupied in crystal structures. Finally, ACES was applied to T4-lysozyme and demonstrated that the side chain distribution of V111χ1 could be reliably reproduced for the apo state, bound to p-xylene, and in p-xylene→ benzene transformations. In these cases, the ACES method is shown to be highly robust and superior to a REST2-like enhanced sampling implementation alone.

Multireference Generalization of the Weighted Thermodynamic Perturbation Method

We describe the generalized weighted thermodynamic perturbation (gwTP) method for estimating the free energy surface of an expensive "high-level" potential energy function from the umbrella sampling performed with multiple inexpensive "low-level" reference potentials. The gwTP method is a generalization of the weighted thermodynamic perturbation (wTP) method developed by Li and co-workers [J. Chem. Theory Comput. 2018, 14, 5583-5596] that uses a single "low-level" reference potential. The gwTP method offers new possibilities in model design whereby the sampling generated from several low-level potentials may be combined (e.g., specific reaction parameter models that might have variable accuracy at different stages of a multistep reaction). The gwTP method is especially well suited for use with machine learning potentials (MLPs) that are trained against computationally expensive ab initio quantum mechanical/molecular mechanical (QM/MM) energies and forces using active learning procedures that naturally produce multiple distinct neural network potentials. Simulations can be performed with greater sampling using the fast MLPs and then corrected to the ab initio level using gwTP. The capabilities of the gwTP method are demonstrated by creating reference potentials based on the MNDO/d and DFTB2/MIO semiempirical models supplemented with the "range-corrected deep potential" (DPRc). The DPRc parameters are trained to ab initio QM/MM data, and the potentials are used to calculate the free energy surface of stepwise mechanisms for nonenzymatic RNA 2'-O-transesterification model reactions. The extended sampling made possible by the reference potentials allows one to identify unequilibrated portions of the simulations that are not always evident from the short time scale commonly used with ab initio QM/MM potentials. We show that the reference potential approach can yield more accurate ab initio free energy predictions than the wTP method or what can be reasonably afforded from explicit ab initio QM/MM sampling.

Who stole the proton? Suspect general base guanine found with a smoking gun in the pistol ribozyme

The pistol ribozyme (Psr) is one among the most recently discovered classes of small nucleolytic ribozymes that catalyze site-specific RNA self-cleavage through 2'-O-transphosphorylation. The Psr contains a conserved guanine (G40) that in crystal structures is in a position suggesting it plays the role of the general base to abstract a proton from the nucleophile to activate the reaction. Although some functional data is consistent with this mechanistic role, a notable exception is 2-aminopurine (2AP) substitution which has no effect on the rate, unlike similar substitutions across other so-called "G + M" and "G + A" ribozyme classes. Herein we postulate that an alternate conserved guanine, G42, is the primary general base, and provide evidence from molecular simulations that the active site of Psr can undergo local refolding into a structure that is consistent with the common "L-platform/L-scaffold" architecture identified in G + M and G + A ribozyme classes with Psr currently the notable exception. We summarize the key currently available experimental data and present new classical and combined quantum mechanical/molecular mechanical simulation results that collectively suggest a new hypothesis. We hypothesize that there are two available catalytic pathways supported by different conformational states connected by a local refolding of the active site: (1) a primary pathway with an active site architecture aligned with the L-platform/L-scaffold framework where G42 acts as a general base, and (2) a secondary pathway with the crystallographic active site architecture where G40 acts as a general base. We go on to make several experimentally testable predictions, and suggest specific experiments that would ultimately bring closure to the mystery as to "who stole the proton in the pistol ribozyme?".

Combined QM/MM, Machine Learning Path Integral Approach to Compute Free Energy Profiles and Kinetic Isotope Effects in RNA Cleavage Reactions

We present a fast, accurate, and robust approach for determination of free energy profiles and kinetic isotope effects for RNA 2'-O-transphosphorylation reactions with inclusion of nuclear quantum effects. We apply a deep potential range correction (DPRc) for combined quantum mechanical/molecular mechanical (QM/MM) simulations of reactions in the condensed phase. The method uses the second-order density-functional tight-binding method (DFTB2) as a fast, approximate base QM model. The DPRc model modifies the DFTB2 QM interactions and applies short-range corrections to the QM/MM interactions to reproduce ab initio DFT (PBE0/6-31G*) QM/MM energies and forces. The DPRc thus enables both QM and QM/MM interactions to be tuned to high accuracy, and the QM/MM corrections are designed to smoothly vanish at a specified cutoff boundary (6 Å in the present work). The computational speed-up afforded by the QM/MM+DPRc model enables free energy profiles to be calculated that include rigorous long-range QM/MM interactions under periodic boundary conditions and nuclear quantum effects through a path integral approach using a new interface between the AMBER and i-PI software. The approach is demonstrated through the calculation of free energy profiles of a native RNA cleavage model reaction and reactions involving thio-substitutions, which are important experimental probes of the mechanism. The DFTB2+DPRc QM/MM free energy surfaces agree very closely with the PBE0/6-31G* QM/MM results, and it is vastly superior to the DFTB2 QM/MM surfaces with and without weighted thermodynamic perturbation corrections. 18O and 34S primary kinetic isotope effects are compared, and the influence of nuclear quantum effects on the free energy profiles is examined.

Energetics of a protein disorder-order transition in small molecule recognition

Many proteins recognise other proteins via mechanisms that involve the folding of intrinsically disordered regions upon complex formation. Here we investigate how the selectivity of a drug-like small molecule arises from its modulation of a protein disorder-to-order transition. Binding of the compound AM-7209 has been reported to confer order upon an intrinsically disordered ‘lid’ region of the oncoprotein MDM2. Calorimetric measurements revealed that truncation of the lid region of MDM2 increases the apparent dissociation constant of AM-7209 250-fold. By contrast, lid truncation has little effect on the binding of the ligand Nutlin-3a. Insights into these differential binding energetics were obtained via a complete thermodynamic analysis that featured adaptive absolute alchemical free energy of binding calculations with enhanced-sampling molecular dynamics simulations. The simulations reveal that in apo MDM2 the ordered lid state is energetically disfavoured. AM-7209, but not Nutlin-3a, shows a significant energetic preference for ordered lid conformations, thus shifting the balance towards ordering of the lid in the AM-7209/MDM2 complex. The methodology reported herein should facilitate broader targeting of intrinsically disordered regions in medicinal chemistry.

Molecular simulations and biophysical measurements elucidate why the ligand AM-7209 orders a disordered region of the protein MDM2 on binding. This work expands strategies available to medicinal chemists for targeting disordered proteins.

Seeding the multi-dimensional nonequilibrium pulling for Hamiltonian variation: indirect nonequilibrium free energy simulations at QM levels

The combination of free energy simulations in the alchemical and configurational spaces provides a feasible route to access the thermodynamic profiles under a computationally demanding target Hamiltonian. Normally, due to the significant differences between the computational cost of ab initio quantum mechanics (QM) calculations and those of semi-empirical quantum mechanics (SQM) and molecular mechanics (MM), this indirect method could be used to obtain the QM thermodynamics by combining the SQM or MM results and the SQM-to-QM or MM-to-QM corrections. In our previous work, a multi-dimensional nonequilibrium pulling framework for Hamiltonian variations was introduced based on bidirectional pulling and bidirectional reweighting. The method performs nonequilibrium free energy simulations in the configurational space to obtain the thermodynamic profile along the conformational change pathway under a selected computationally efficient Hamiltonian, and uses the nonequilibrium alchemical method to correct or perturb the thermodynamic profile to that under the target Hamiltonian. The BAR-based method is designed to achieve the best generality and transferability and thus leads to modest (∼20 fold) speedup. In this work, we explore the possibility of further accelerating the nonequilibrium free energy simulation by employing unidirectional pulling and using the selection criterion to obtain the initial configurations used to initiate nonequilibrium trajectories following the idea of adaptive steered molecular dynamics (ASMD). A single initial condition is used to seed the whole multi-dimensional nonequilibrium free energy simulation and the sampling is performed fully in the nonequilibrium ensemble. Introducing very short ps-length equilibrium sampling to grab more initial seeds could also be helpful. The ASMD scheme estimates the free energy difference with the unidirectional exponential average (EXP), but it does not follow exactly the requirements of the EXP estimator. Another deficiency of the seeding simulation is the inherently sequential or serial pulling due to the inter-segment dependency, which triggers some problems in the parallelizability of the simulation. Numerical tests are performed to grasp some insights and guidelines for using this selection-criterion-based ASMD scheme. The presented selection-criterion-based multi-dimensional ASMD scheme follows the same perturbation network of the BAR-based method, and thus could be used in various Hamiltonian-variation cases.

Solvation Effects on the Thermal Helix Inversion of Molecular Motors from QM/MM Calculations

Molecular motors convert light and thermal energies into mechanical work, offering good opportunities to design novel molecular devices. Among them, molecular motors alternate a photoisomerization and a thermal helix inversion to achieve unidirectional rotation. The rotational speed is limited by the helix inversion step, which in turn is governed by a barrier in the electronic ground state. In this work, we systematically study the solvation effect on the thermal process of selected molecular motors, comparing reaction barriers obtained from both density functional theory (DFT) in the isolated system and umbrella sampling within a hybrid quantum mechanics/molecular mechanics (QM/MM) model in solution. We find more prominent solvation effects on those molecular motors with larger dipole moments. The results could provide insight into how to functionalize molecular motors to speed up their rotation.

Calculation of Energy for RNA/RNA and DNA/RNA Duplex Formation by Molecular Dynamics Simulation

Abstract

The development of approaches for predictive calculation of hybridization properties of various nucleic acid (NA) derivatives is the basis for the rational design of the NA-based constructs. Modern advances in computer modeling methods provide the feasibility of these calculations. We have analyzed the possibility of calculating the energy of DNA/RNA and RNA/RNA duplex formation using representative sets of complexes (65 and 75 complexes, respectively). We used the classical molecular dynamics (MD) method, the MMPBSA or MMGBSA approaches to calculate the enthalpy (ΔH°) component, and the quasi-harmonic approximation (Q-Harm) or the normal mode analysis (NMA) methods to calculate the entropy (ΔS°) contribution to the Gibbs energy ($$\Delta G_{{37}}^{^\circ }$$ ) of the NA complex formation. We have found that the MMGBSA method in the analysis of the MD trajectory of only the NA duplex and the empirical linear approximation allow calculation of the enthalpy of formation of the DNA, RNA, and hybrid duplexes of various lengths and GC content with an accuracy of 8.6%. Within each type of complex, the combination of rather efficient MMGBSA and Q-Harm approaches being applied to the trajectory of only the bimolecular complex makes it possible to calculate the $$\Delta G_{{37}}^{^\circ }$$ of the duplex formation with an error value of 10%. The high accuracy of predictive calculation for different types of natural complexes (DNA/RNA, DNA/RNA, and RNA/RNA) indicates the possibility of extending the considered approach to analogs and derivatives of nucleic acids, which gives a fundamental opportunity in the future to perform rational design of new types of NA-targeted sequence-specific compounds.

Development of Range-Corrected Deep Learning Potentials for Fast, Accurate Quantum Mechanical/Molecular Mechanical Simulations of Chemical Reactions in Solution

We develop a new deep potential─range correction (DPRc) machine learning potential for combined quantum mechanical/molecular mechanical (QM/MM) simulations of chemical reactions in the condensed phase. The new range correction enables short-ranged QM/MM interactions to be tuned for higher accuracy, and the correction smoothly vanishes within a specified cutoff. We further develop an active learning procedure for robust neural network training. We test the DPRc model and training procedure against a series of six nonenzymatic phosphoryl transfer reactions in solution that are important in mechanistic studies of RNA-cleaving enzymes. Specifically, we apply DPRc corrections to a base QM model and test its ability to reproduce free-energy profiles generated from a target QM model. We perform these comparisons using the MNDO/d and DFTB2 semiempirical models because they differ in the way they treat orbital orthogonalization and electrostatics and produce free-energy profiles which differ significantly from each other, thereby providing us a rigorous stress test for the DPRc model and training procedure. The comparisons show that accurate reproduction of the free-energy profiles requires correction of the QM/MM interactions out to 6 Å. We further find that the model's initial training benefits from generating data from temperature replica exchange simulations and including high-temperature configurations into the fitting procedure, so the resulting models are trained to properly avoid high-energy regions. A single DPRc model was trained to reproduce four different reactions and yielded good agreement with the free-energy profiles made from the target QM/MM simulations. The DPRc model was further demonstrated to be transferable to 2D free-energy surfaces and 1D free-energy profiles that were not explicitly considered in the training. Examination of the computational performance of the DPRc model showed that it was fairly slow when run on CPUs but was sped up almost 100-fold when using NVIDIA V100 GPUs, resulting in almost negligible overhead. The new DPRc model and training procedure provide a potentially powerful new tool for the creation of next-generation QM/MM potentials for a wide spectrum of free-energy applications ranging from drug discovery to enzyme design.

Extension of the Variational Free Energy Profile and Multistate Bennett Acceptance Ratio Methods for High-Dimensional Potential of Mean Force Profile Analysis

We redevelop the variational free energy profile (vFEP) method using a cardinal B-spline basis to extend the method for analyzing free energy surfaces (FESs) involving three or more reaction coordinates. We also implemented software for evaluating high-dimensional profiles based on the multistate Bennett acceptance ratio (MBAR) method which constructs an unbiased probability density from global reweighting of the observed samples. The MBAR method takes advantage of a fast algorithm for solving the unbinned weighted histogram (UWHAM)/MBAR equations which replaces the solution of simultaneous equations with a nonlinear optimization of a convex function. We make use of cardinal B-splines and multiquadric radial basis functions to obtain smooth, differentiable MBAR profiles in arbitrary high dimensions. The cardinal B-spline vFEP and MBAR methods are compared using three example systems that examine 1D, 2D, and 3D profiles. Both methods are found to be useful and produce nearly indistinguishable results. The vFEP method is found to be 150 times faster than MBAR when applied to periodic 2D profiles, but the MBAR method is 4.5 times faster than vFEP when evaluating unbounded 3D profiles. In agreement with previous comparisons, we find the vFEP method produces superior FESs when the overlap between umbrella window simulations decreases. Finally, the associative reaction mechanism of hammerhead ribozyme is characterized using 3D, 4D, and 6D profiles, and the higher-dimensional profiles are found to have smaller reaction barriers by as much as 1.5 kcal/mol. The methods presented here have been implemented into the FE-ToolKit software package along with new methods for network-wide free energy analysis in drug discovery.

Exploring the concerted mechanistic pathway for HIV-1 PR-substrate revealed by umbrella sampling simulation

HIV-1 protease (HIV-1 PR) is an essential enzyme for the replication process of its virus, and therefore considered an important target for the development of drugs against the acquired immunodeficiency syndrome (AIDS). Our previous study shows that the catalytic mechanism of subtype B/C-SA HIV-1 PR follows a one-step concerted acyclic hydrolysis reaction process using a two-layered ONIOM B3LYP/6-31++G(d,p) method. This present work is aimed at exploring the proposed mechanism of the proteolysis catalyzed by HIV-1 PR and to ensure our proposed mechanism is not an artefact of a single theoretical technique. Hence, we present umbrella sampling method that is suitable for calculating potential mean force (PMF) for non-covalent ligand/substrate-enzyme association/dissociation interactions which provide thermodynamic details for molecular recognition. The free activation energy results were computed in terms of PMF analysis within the hybrid QM(DFTB)/MM approach. The theoretical findings suggest that the proposed mechanism corresponds in principle with experimental data. Given our observations, we suggest that the QM/MM MD method can be used as a reliable computational technique to rationalize lead compounds against specific targets such as the HIV-1 protease.

Markov State Model of Lassa Virus Nucleoprotein Reveals Large Structural Changes during the Trimer to Monomer Transition

Lassa virus contains a nucleoprotein (NP) that encapsulates the viral genomic RNA forming the ribonucleoprotein (RNP). The NP forms trimers that do not bind RNA, but a structure of only the NP N-terminal domain was co-crystallized with RNA bound. These structures suggested a model in which the NP forms a trimer to keep the RNA gate closed, but then is triggered to undergo a change to a form competent for RNA binding. Here, we investigate the scenario in which the trimer is disrupted to observe whether monomeric NP undergoes significant conformational changes. From multi-microsecond molecular dynamics simulations and an adaptive sampling scheme to sample the conformational space, a Markov state model (MSM) is constructed. The MSM reveals an energetically favorable conformational change, with the most significant changes occurring at the domain interface. These results support a model in which significant structural reorganization of the NP is required for RNP formation.

Assessment of the Cruzain Cysteine Protease Reversible and Irreversible Covalent Inhibition Mechanism

Reversible and irreversible covalent ligands are advanced cysteine protease inhibitors in the drug development pipeline. K777 is an irreversible inhibitor of cruzain, a necessary enzyme for the survival of the Trypanosoma cruzi (T. cruzi) parasite, the causative agent of Chagas disease. Despite their importance, irreversible covalent inhibitors are still often avoided due to the risk of adverse effects. Herein, we replaced the K777 vinyl sulfone group with a nitrile moiety to obtain a reversible covalent inhibitor (Neq0682) of cysteine protease. Then, we used advanced experimental and computational techniques to explore details of the inhibition mechanism of cruzain by reversible and irreversible inhibitors. The isothermal titration calorimetry (ITC) analysis shows that inhibition of cruzain by an irreversible inhibitor is thermodynamically more favorable than by a reversible one. The hybrid Quantum Mechanics/Molecular Mechanics (QM/MM) and Molecular Dynamics (MD) simulations were used to explore the mechanism of the reaction inhibition of cruzain by K777 and Neq0682. The calculated free energy profiles show that the Cys25 nucleophilic attack and His162 proton transfer occur in a single step for a reversible inhibitor and two steps for an irreversible covalent inhibitor. The hybrid QM/MM calculated free energies for the inhibition reaction correspond to -26.7 and -5.9 kcal mol^-1 for K777 and Neq0682 at the MP2/MM level, respectively. These results indicate that the ΔG of the reaction is very negative for the process involving K777, consequently, the covalent adduct cannot revert to a noncovalent protein-ligand complex, and its binding tends to be irreversible. Overall, the present study provides insights into a covalent inhibition mechanism of cysteine proteases.

Reaction mechanism of nucleoside 2'-deoxyribosyltransferases: free-energy landscape supports an oxocarbenium ion as the reaction intermediate

Insight into the catalytic mechanism of Lactobacillus leichmannii nucleoside 2'-deoxyribosyltransferase (LlNDT) has been gained by calculating a quantum mechanics-molecular mechanics (QM/MM) free-energy landscape of the reaction within the enzyme active site. Our results support an oxocarbenium species as the reaction intermediate and thus an SN1 reaction mechanism in this family of bacterial enzymes. Our mechanistic proposal is validated by comparing experimental kinetic data on the impact of the single amino acid replacements Tyr7, Glu98 and Met125 with Ala, Asp and Ala/norLeu, respectively, and accounts for the specificity shown by this enzyme on a non-natural substrate. This work broadens our understanding of enzymatic C-N bond cleavage and C-N bond formation.

Genome mining and characterisation of a novel transaminase with remote stereoselectivity

Microbial enzymes from pristine niches can potentially deliver disruptive opportunities in synthetic routes to Active Pharmaceutical Ingredients and intermediates in the Pharmaceutical Industry. Advances in green chemistry technologies and the importance of stereochemical control, further underscores the application of enzyme-based solutions in chemical synthesis. The rich tapestry of microbial diversity in the oceanic ecosystem encodes a capacity for novel biotransformations arising from the chemical complexity of this largely unexplored bioactive reservoir. Here we report a novel ω-transaminase discovered in a marine sponge Pseudovibrio sp. isolate. Remote stereoselection using a transaminase has been demonstrated for the first time using this novel protein. Application to the resolution of an intermediate in the synthesis of sertraline highlights the synthetic potential of this novel biocatalyst discovered through genomic mining. Integrated chemico-genomics revealed a unique substrate profile, while molecular modelling provided structural insights into this ‘first in class’ selectivity at a remote chiral centre.

BAR-based multi-dimensional nonequilibrium pulling for indirect construction of QM/MM free energy landscapes: from semi-empirical to ab initio

The indirect method for the construction of quantum mechanics (QM)/molecular mechanics (MM) free energy landscapes provides a cheaper alternative for free energy simulations at the QM level. The indirect method features a direct calculation of the free energy profile with a computationally efficient but less accurate Hamiltonian (i.e. low-level Hamiltonian) and a low-level-to-high-level correction. In the thermodynamic cycle, the direct low-level calculation along the physically meaningful reaction coordinate is corrected via the alchemical method, which is often achieved with perturbation-based techniques. In our previous work, a multi-dimensional nonequilibrium pulling framework is proposed for the indirect construction of QM/MM free energy landscapes. Previously, we focused on obtaining semi-empirical QM (SQM) results indirectly from direct MM simulations and MM to SQM corrections. In this work, we apply this method to obtain results under ab initio QM Hamiltonians by combining direct SQM results and SQM to QM corrections. A series of SQM and QM Hamiltonians are benchmarked. It is observed that PM6 achieves the best performance among the low-level Hamiltonians. Therefore, we recommend using PM6 as the low-level theory in the indirect free energy simulation. Considering its higher similarity to the high-level Hamiltonians, PM6 corrected with the bond charge correction could be more accurate than the existing AM1-BCC model. Another central result in the current work is a basic protocol of choosing the strength of restraints and an appropriate time step in nonequilibrium free energy simulation at the stiff spring limit. We provide theoretical derivations to emphasize the importance of using a sufficiently large force constant and choosing an appropriate time step. It is worth noting that a general rule of thumb for choosing the time step, according to our derivation, is that a time step of 1 fs or smaller should be used, as long as the stiff spring approximation is employed, even in simulations with constraints on bonds involving hydrogen atoms.

Accelerated computation of free energy profile at ab initio quantum mechanical/molecular mechanical accuracy via a semi-empirical reference potential. II. Recalibrating semi-empirical parameters with force matching

An efficient and accurate reference potential simulation protocol is proposed for producing ab initio quantum mechanical/molecular mechanical (AI-QM/MM) quality free energy profiles for chemical reactions in a solvent or macromolecular environment. This protocol involves three stages: (a) using force matching to recalibrate a semi-empirical quantum mechanical (SE-QM) Hamiltonian for the specific reaction under study; (b) employing the recalibrated SE-QM Hamiltonian (in combination with molecular mechanical force fields) as the reference potential to drive umbrella samplings along the reaction pathway; and (c) computing AI-QM/MM energy values for collected configurations from the sampling and performing weighted thermodynamic perturbation to acquire an AI-QM/MM corrected reaction free energy profile. For three model reactions (identity SN2 reaction, Menshutkin reaction, and glycine proton transfer reaction) in aqueous solution and one enzyme reaction (Claisen arrangement in chorismate mutase), our simulations using recalibrated PM3 SE-QM Hamiltonians well reproduced QM/MM free energy profiles at the B3LYP/6-31G* level of theory all within 1 kcal mol-1 with a 20 to 45 fold reduction in the computer time.

Using AMBER18 for Relative Free Energy Calculations

With renewed interest in free energy methods in contemporary structure-based drug design, there is a pressing need to validate against multiple targets and force fields to assess the overall ability of these methods to accurately predict relative binding free energies. We computed relative binding free energies using graphics processing unit accelerated thermodynamic integration (GPU-TI) on a data set originally assembled by Schrödinger, Inc. Using their GPU free energy code (FEP+) and the OPLS2.1 force field combined with the REST2 enhanced sampling approach, these authors obtained an overall MUE of 0.9 kcal/mol and an overall RMSD of 1.14 kcal/mol. In our study using GPU-TI from AMBER with the AMBER14SB/GAFF1.8 force field but without enhanced sampling, we obtained an overall MUE of 1.17 kcal/mol and an overall RMSD of 1.50 kcal/mol for the 330 perturbations contained in this data set. A more detailed analysis of our results suggested that the observed differences between the two studies arise from differences in sampling protocols along with differences in the force fields employed. Future work should address the problem of establishing benchmark quality results with robust statistical error bars obtained through multiple independent runs and enhanced sampling, which is possible with the GPU-accelerated features in AMBER.

BAR-based multi-dimensional nonequilibrium pulling for indirect construction of a QM/MM free energy landscape

Construction of free energy landscapes at the quantum mechanics (QM) level is computationally demanding. As shown in previous studies, by employing an indirect scheme (i.e. constructing a thermodynamic cycle connecting QM states via an alchemical pathway), simulations are converged with much less computational burden. The indirect scheme makes QM/molecular mechanics (MM) free energy simulation orders of magnitude faster than the direct QM/MM schemes. However, the indirect QM/MM simulations were mostly equilibrium sampling based and the nonequilibrium methods were merely exploited in one-dimensional alchemical QM/MM end-state correction at two end states. In this work, we represent a multi-dimensional nonequilibrium pulling scheme for indirect QM/MM free energy simulations, where the whole free energy simulation is performed only with nonequilibrium methods. The collective variable (CV) space we explore is a combination of one alchemical CV and one physically meaningful CV. The current nonequilibrium indirect QM/MM simulation method can be seen as the generalization of equilibrium perturbation based indirect QM/MM methods. The test systems include one backbone dihedral case and one distance case. The two cases are significantly different in size, enabling us to investigate the dependence of the speedup of the indirect scheme on the size of the system. It is shown that the speedup becomes larger when the size of the system becomes larger, which is consistent with the scaling behavior of QM Hamiltonians.

Fluorinated synthetic anion carriers: experimental and computational insights into transmembrane chloride transport

A series of fluorinated tripodal tris-thioureas function as highly active anion transporters across lipid bilayers and cell membranes. Here, we investigate their mechanism of action using anion transport assays in cells and synthetic vesicles and molecular modelling of transporter-lipid interactions. When compared with non-fluorinated analogues, fluorinated compounds demonstrate a different mechanism of membrane transport because the free transporter cannot effectively diffuse through the membrane. As a result, in H+/Cl- cotransport assays, fluorinated transporters require the presence of oleic acid to form anionic oleate complexes for recycling of the transporter, whereas non-fluorinated analogues readily diffuse through the membrane as free transporters and show synergistic transport with the proton transporter gramicidin. Molecular dynamics simulations revealed markedly stronger transporter-lipid interactions for fluorinated compounds compared with non-fluorinated analogues and hence, higher energy barriers for fluorinated compounds to cross the membrane as free transporters. With use of appropriate proton transporters to ensure measurement of the correct rate-limiting steps, the transport rates determined in synthetic vesicle assays show excellent agreement with the anion transport rates determined in cell-based assays. We conclude that integration of computational and experimental methods provides a strategy to optimise transmembrane anion transporter design for biomedical applications.

Molecular features of hydration layers probed by atomic force microscopy

Structurally-ordered layers of water are universally formed on a solid surface in aqueous solution or under ambient conditions. Although such hydration layers are commonly probed via atomic force microscopy (AFM), the current understanding on how the hydration layers manifest themselves in an AFM experiment is far from complete. By using molecular dynamics simulation, we investigate the hydration layers on a hydrophilic or hydrophobic surface probed by a nanoscale tip. We study the density and molecular orientation of water, the free energy, and the force on the tip by varying the tip-surface distance. The force-distance curve oscillates due to the transition between the mono-, bi-, and tri-layers of water confined between the tip and the surface. If both the tip and the surface are hydrophobic, water confined between the tip and the surface evaporates due to the dewetting transition, giving a hydrophobic force without oscillation. The periodicity of oscillation in the force differs from the structural periodicity of water. With a close proximity of the tip, the molecular dipoles align parallel to the surface, regardless of whether the tip and the surface are hydrophilic or hydrophobic.

On the convergence of multi-scale free energy simulations

In this work we employ simple model systems to evaluate the relative performance of two of the most important free energy methods: The Zwanzig equation (also known as "Free energy perturbation") and Bennett's acceptance ratio method (BAR). Although our examples should be transferable to other kinds of free energy simulations, we focus on applications of multi-scale free energy simulations. Such calculations are especially complex, since they connect two different levels of theory with very different requirements in terms of speed, accuracy, sampling and parallelizability. We try to reconcile all those different factors by developing some simple criteria to guide the early stages of the development of a free energy protocol. This is accomplished by quantifying how many λ intermediate steps and how many potential energy evaluations are necessary in order to reach a certain level of convergence.

A Multidimensional B-Spline Correction for Accurate Modeling Sugar Puckering in QM/MM Simulations

The computational efficiency of approximate quantum mechanical methods allows their use for the construction of multidimensional reaction free energy profiles. It has recently been demonstrated that quantum models based on the neglect of diatomic differential overlap (NNDO) approximation have difficulty modeling deoxyribose and ribose sugar ring puckers and thus limit their predictive value in the study of RNA and DNA systems. A method has been introduced in our previous work to improve the description of the sugar puckering conformational landscape that uses a multidimensional B-spline correction map (BMAP correction) for systems involving intrinsically coupled torsion angles. This method greatly improved the adiabatic potential energy surface profiles of DNA and RNA sugar rings relative to high-level ab initio methods even for highly problematic NDDO-based models. In the present work, a BMAP correction is developed, implemented, and tested in molecular dynamics simulations using the AM1/d-PhoT semiempirical Hamiltonian for biological phosphoryl transfer reactions. Results are presented for gas-phase adiabatic potential energy surfaces of RNA transesterification model reactions and condensed-phase QM/MM free energy surfaces for nonenzymatic and RNase A-catalyzed transesterification reactions. The results show that the BMAP correction is stable, efficient, and leads to improvement in both the potential energy and free energy profiles for the reactions studied, as compared with ab initio and experimental reference data. Exploration of the effect of the size of the quantum mechanical region indicates the best agreement with experimental reaction barriers occurs when the full CpA dinucleotide substrate is treated quantum mechanically with the sugar pucker correction.

Nucleotide Excision Repair Lesion-Recognition Protein Rad4 Captures a Pre-Flipped Partner Base in a Benzo[a]pyrene-Derived DNA Lesion: How Structure Impacts the Binding Pathway

The xeroderma pigmentosum C protein complex (XPC) recognizes a variety of environmentally induced DNA lesions and is the key in initiating their repair by the nucleotide excision repair (NER) pathway. When bound to a lesion, XPC flips two nucleotide pairs that include the lesion out of the DNA duplex, yielding a productively bound complex that can lead to successful lesion excision. Interestingly, the efficiencies of NER vary greatly among different lesions, influencing their toxicity and mutagenicity in cells. Though differences in XPC binding may influence NER efficiency, it is not understood whether XPC utilizes different mechanisms to achieve productive binding with different lesions. Here, we investigated the well-repaired 10R-(+)-cis-anti-benzo[a]pyrene-N²-dG (cis-B[a]P-dG) DNA adduct in a duplex containing normal partner C opposite the lesion. This adduct is derived from the environmental pro-carcinogen benzo[a]pyrene and is likely to be encountered by NER in the cell. We have extensively investigated its binding to the yeast XPC orthologue, Rad4, using umbrella sampling with restrained molecular dynamics simulations and free energy calculations. The NMR solution structure of this lesion in duplex DNA has shown that the dC complementary to the adducted dG is flipped out of the DNA duplex in the absence of XPC. However, it is not known whether the "pre-flipped" base would play a role in its recognition by XPC. Our results show that Rad4 first captures the displaced dC, which is followed by a tightly coupled lesion-extruding pathway for productive binding. This binding path differs significantly from the one deduced for the small cis-syn cyclobutane pyrimidine dimer lesion opposite mismatched thymines [ Mu , H. , ( 2015 ) Biochemistry , 54 ( 34 ), 5263 - 7 ]. The possibility of multiple paths that lead to productive binding to XPC is consistent with the versatile lesion recognition by XPC that is required for successful NER.

Charge-dependent many-body exchange and dispersion interactions in combined QM/MM simulations

Accurate modeling of the molecular environment is critical in condensed phase simulations of chemical reactions. Conventional quantum mechanical/molecular mechanical (QM/MM) simulations traditionally model non-electrostatic non-bonded interactions through an empirical Lennard-Jones (LJ) potential which, in violation of intuitive chemical principles, is bereft of any explicit coupling to an atom's local electronic structure. This oversight results in a model whereby short-ranged exchange-repulsion and long-ranged dispersion interactions are invariant to changes in the local atomic charge, leading to accuracy limitations for chemical reactions where significant atomic charge transfer can occur along the reaction coordinate. The present work presents a variational, charge-dependent exchange-repulsion and dispersion model, referred to as the charge-dependent exchange and dispersion (QXD) model, for hybrid QM/MM simulations. Analytic expressions for the energy and gradients are provided, as well as a description of the integration of the model into existing QM/MM frameworks, allowing QXD to replace traditional LJ interactions in simulations of reactive condensed phase systems. After initial validation against QM data, the method is demonstrated by capturing the solvation free energies of a series of small, chlorine-containing compounds that have varying charge on the chlorine atom. The model is further tested on the SN2 attack of a chloride anion on methylchloride. Results suggest that the QXD model, unlike the traditional LJ model, is able to simultaneously obtain accurate solvation free energies for a range of compounds while at the same time closely reproducing the experimental reaction free energy barrier. The QXD interaction model allows explicit coupling of atomic charge with many-body exchange and dispersion interactions that are related to atomic size and provides a more accurate and robust representation of non-electrostatic non-bonded QM/MM interactions.

Structured and Unstructured Binding of an Intrinsically Disordered Protein as Revealed by Atomistic Simulations

Intrinsically disordered proteins (IDPs) are a set of proteins that lack a definite secondary structure in solution. IDPs can acquire tertiary structure when bound to their partners; therefore, the recognition process must also involve protein folding. The nature of the transition state (TS), structured or unstructured, determines the binding mechanism. The characterization of the TS has become a major challenge for experimental techniques and molecular simulations approaches since diffusion, recognition, and binding is coupled to folding. In this work we present atomistic molecular dynamics (MD) simulations that sample the free energy surface of the coupled folding and binding of the transcription factor c-myb to the cotranscription factor CREB binding protein (CBP). This process has been recently studied and became a model to study IDPs. Despite the plethora of available information, we still do not know how c-myb binds to CBP. We performed a set of atomistic biased MD simulations running a total of 15.6 μs. Our results show that c-myb folds very fast upon binding to CBP with no unique pathway for binding. The process can proceed through both structured or unstructured TS's with similar probabilities. This finding reconciles previous seemingly different experimental results. We also performed Go-type coarse-grained MD of several structured and unstructured models that indicate that coupled folding and binding follows a native contact mechanism. To the best of our knowledge, this is the first atomistic MD simulation that samples the free energy surface of the coupled folding and binding processes of IDPs.

Multipolar Ewald methods, 2: applications using a quantum mechanical force field

A fully quantum mechanical force field (QMFF) based on a modified “divide-and-conquer” (mDC) framework is applied to a series of molecular simulation applications, using a generalized Particle Mesh Ewald method extended to multipolar charge densities. Simulation results are presented for three example applications: liquid water, p-nitrophenylphosphate reactivity in solution, and crystalline N,N-dimethylglycine. Simulations of liquid water using a parametrized mDC model are compared to TIP3P and TIP4P/Ew water models and experiment. The mDC model is shown to be superior for cluster binding energies and generally comparable for bulk properties. Examination of the dissociative pathway for dephosphorylation of p-nitrophenylphosphate shows that the mDC method evaluated with the DFTB3/3OB and DFTB3/OPhyd semiempirical models bracket the experimental barrier, whereas DFTB2 and AM1/d-PhoT QM/MM simulations exhibit deficiencies in the barriers, the latter for which is related, in part, to the anomalous underestimation of the p-nitrophenylate leaving group pKa. Simulations of crystalline N,N-dimethylglycine are performed and the overall structure and atomic fluctuations are compared with the experiment and the general AMBER force field (GAFF). The QMFF, which was not parametrized for this application, was shown to be in better agreement with crystallographic data than GAFF. Our simulations highlight some of the application areas that may benefit from using new QMFFs, and they demonstrate progress toward the development of accurate QMFFs using the recently developed mDC framework.

Impact of Ser17 Phosphorylation on the Conformational Dynamics of the Oncoprotein MDM2

MDM2 is an important oncoprotein that downregulates the activity of the tumor suppressor protein p53 via binding of its N-terminal domain to the p53 transactivation domain. The first 24 residues of the MDM2 N-terminal domain form an intrinsically disordered "lid" region that interconverts on a millisecond time scale between "open" and "closed" states in unliganded MDM2. While the former conformational state is expected to facilitate p53 binding, the latter competes in a pseudo-substrate manner with p53 for its binding site. Phosphorylation of serine 17 in the MDM2 lid region is thought to modulate the equilibrium between "open" and "closed" lid states, but contradictory findings on the favored lid conformational state upon phosphorylation have been reported. Here, the nature of the conformational states of MDM2 pSer17 and Ser17Asp variants was addressed by means of enhanced sampling molecular dynamics simulations. Detailed analyses of the computed lid conformational ensembles indicate that both lid variants stabilize a "closed" state, with respect to wild type. Nevertheless, the nature of the closed-state conformational ensembles differs significantly between the pSer17 and Ser17Asp variants. Thus, care should be applied in the interpretation of biochemical experiments that use phosphomimetic variants to model the effects of phosphorylation on the structure and dynamics of this disordered protein region.

Force Field for Mg(2+), Mn(2+), Zn(2+), and Cd(2+) Ions That Have Balanced Interactions with Nucleic Acids

Divalent metal ions are of fundamental importance to the function and folding of nucleic acids. Divalent metal ion-nucleic acid interactions are complex in nature and include both territorial and site specific binding. Commonly employed nonbonded divalent ion models, however, are often parametrized against bulk ion properties and are subsequently utilized in biomolecular simulations without considering any data related to interactions at specific nucleic acid sites. Previously, we assessed the ability of 17 different nonbonded Mg(2+) ion models to reproduce different properties of Mg(2+) in aqueous solution including radial distribution functions, solvation free energies, water exchange rates, and translational diffusion coefficients. In the present work, we depart from the recently developed 12-6-4 potential models for divalent metal ions developed by Li and Merz and tune the pairwise parameters for Mg(2+), Mn(2+), Zn(2+), and Cd(2+) binding dimethyl phosphate, adenosine, and guanosine in order to reproduce experimental site specific binding free energies derived from potentiometric pH titration data. We further apply these parameters to investigate a metal ion migration previously proposed to occur during the catalytic reaction of the hammerhead ribozyme. The new parameters are shown to be accurate and balanced for nucleic acid binding in comparison with available experimental data and provide an important tool for molecular dynamics and free energy simulations of nucleic acids where these ions may exhibit different binding modes.

Assessment of metal-assisted nucleophile activation in the hepatitis delta virus ribozyme from molecular simulation and 3D-RISM

The hepatitis delta virus ribozyme is an efficient catalyst of RNA 2'-O-transphosphorylation and has emerged as a key experimental system for identifying and characterizing fundamental features of RNA catalysis. Recent structural and biochemical data have led to a proposed mechanistic model whereby an active site Mg(2+) ion facilitates deprotonation of the O2' nucleophile, and a protonated cytosine residue (C75) acts as an acid to donate a proton to the O5' leaving group as noted in a previous study. This model assumes that the active site Mg(2+) ion forms an inner-sphere coordination with the O2' nucleophile and a nonbridging oxygen of the scissile phosphate. These contacts, however, are not fully resolved in the crystal structure, and biochemical data are not able to unambiguously exclude other mechanistic models. In order to explore the feasibility of this model, we exhaustively mapped the free energy surfaces with different active site ion occupancies via quantum mechanical/molecular mechanical (QM/MM) simulations. We further incorporate a three-dimensional reference interaction site model for the solvated ion atmosphere that allows these calculations to consider not only the rate associated with the chemical steps, but also the probability of observing the system in the presumed active state with the Mg(2+) ion bound. The QM/MM results predict that a pathway involving metal-assisted nucleophile activation is feasible based on the rate-controlling transition state barrier departing from the presumed metal-bound active state. However, QM/MM results for a similar pathway in the absence of Mg(2+) are not consistent with experimental data, suggesting that a structural model in which the crystallographically determined Mg(2+) is simply replaced with Na(+) is likely incorrect. It should be emphasized, however, that these results hinge upon the assumption of the validity of the presumed Mg(2+)-bound starting state, which has not yet been definitively verified experimentally, nor explored in depth computationally. Thus, further experimental and theoretical study is needed such that a consensus view of the catalytic mechanism emerges.

Recognition of Damaged DNA for Nucleotide Excision Repair: A Correlated Motion Mechanism with a Mismatched cis-syn Thymine Dimer Lesion

Mammalian global genomic nucleotide excision repair requires lesion recognition by XPC, whose detailed binding mechanism remains to be elucidated. Here we have delineated the dynamic molecular pathway and energetics of lesion-specific and productive binding by the Rad4/yeast XPC lesion recognition factor, as it forms the open complex [Min, J. H., and Pavletich, N. P. (2007) Nature 449, 570–575; Chen, X., et al. (2015) Nat. Commun. 6, 5849] that is required for excision. We investigated extensively a cis-syn cyclobutane pyrimidine dimer in mismatched duplex DNA, using high-level computational approaches. Our results delineate a preferred correlated motion mechanism, which provides for the first time an atomistic description of the sequence of events as Rad4 productively binds to the damaged DNA.

Free energy landscape and transition pathways from Watson-Crick to Hoogsteen base pairing in free duplex DNA

Houghton (HG) base pairing plays a central role in the DNA binding of proteins and small ligands. Probing detailed transition mechanism from Watson–Crick (WC) to HG base pair (bp) formation in duplex DNAs is of fundamental importance in terms of revealing intrinsic functions of double helical DNAs beyond their sequence determined functions. We investigated a free energy landscape of a free B-DNA with an adenosine–thymine (A–T) rich sequence to probe its conformational transition pathways from WC to HG base pairing. The free energy landscape was computed with a state-of-art two-dimensional umbrella molecular dynamics simulation at the all-atom level. The present simulation showed that in an isolated duplex DNA, the spontaneous transition from WC to HG bp takes place via multiple pathways. Notably, base flipping into the major and minor grooves was found to play an important role in forming these multiple transition pathways. This finding suggests that naked B-DNA under normal conditions has an inherent ability to form HG bps via spontaneous base opening events.

Elucidation of Ligand-Dependent Modulation of Disorder-Order Transitions in the Oncoprotein MDM2

Numerous biomolecular interactions involve unstructured protein regions, but how to exploit such interactions to enhance the affinity of a lead molecule in the context of rational drug design remains uncertain. Here clarification was sought for cases where interactions of different ligands with the same disordered protein region yield qualitatively different results. Specifically, conformational ensembles for the disordered lid region of the N-terminal domain of the oncoprotein MDM2 in the presence of different ligands were computed by means of a novel combination of accelerated molecular dynamics, umbrella sampling, and variational free energy profile methodologies. The resulting conformational ensembles for MDM2, free and bound to p53 TAD (17-29) peptide identify lid states compatible with previous NMR measurements. Remarkably, the MDM2 lid region is shown to adopt distinct conformational states in the presence of different small-molecule ligands. Detailed analyses of small-molecule bound ensembles reveal that the ca. 25-fold affinity improvement of the piperidinone family of inhibitors for MDM2 constructs that include the full lid correlates with interactions between ligand hydrophobic groups and the C-terminal lid region that is already partially ordered in apo MDM2. By contrast, Nutlin or benzodiazepinedione inhibitors, that bind with similar affinity to full lid and lid-truncated MDM2 constructs, interact additionally through their solubilizing groups with N-terminal lid residues that are more disordered in apo MDM2.

Life as we know it depends on interactions between proteins. There is substantial evidence that many interactions between proteins involve very flexible protein regions. These disordered regions may undergo disorder/order transitions upon forming an interaction with another protein. Many successful approaches to medicinal chemistry are based on mimicking the interactions of biological molecules with man-made small molecules. However how drug-like small-molecules may modulate protein disorder is currently poorly understood, largely because it is difficult to measure in details this type of interaction with experimental methods. Here we have used computer simulations to resolve with great details the process by which different small-molecules modulate the flexibility of a disordered region of the protein MDM2. This protein is overexpressed in many cancers and small-molecules that recognize MDM2 have been developed over the last decade as possible novel anti-cancer agents. We show that the flexible MDM2 “lid” region adopts different conformational states in the presence of different small-molecules. Our results suggest why some classes of small-molecules form favorable interactions with the lid region, whereas others do not. These findings may prove crucial to develop new and more effective MDM2 inhibitors, and more generally to help drug designers target disordered proteins regions with small-molecules.

Comparison of structural, thermodynamic, kinetic and mass transport properties of Mg(2+) ion models commonly used in biomolecular simulations

The prevalence of Mg(2+) ions in biology and their essential role in nucleic acid structure and function has motivated the development of various Mg(2+) ion models for use in molecular simulations. Currently, the most widely used models in biomolecular simulations represent a nonbonded metal ion as an ion-centered point charge surrounded by a nonelectrostatic pairwise potential that takes into account dispersion interactions and exchange effects that give rise to the ion's excluded volume. One strategy toward developing improved models for biomolecular simulations is to first identify a Mg(2+) model that is consistent with the simulation force fields that closely reproduces a range of properties in aqueous solution, and then, in a second step, balance the ion-water and ion-solute interactions by tuning parameters in a pairwise fashion where necessary. The present work addresses the first step in which we compare 17 different nonbonded single-site Mg(2+) ion models with respect to their ability to simultaneously reproduce structural, thermodynamic, kinetic and mass transport properties in aqueous solution. None of the models based on a 12-6 nonelectrostatic nonbonded potential was able to reproduce the experimental radial distribution function, solvation free energy, exchange barrier and diffusion constant. The models based on a 12-6-4 potential offered improvement, and one model in particular, in conjunction with the SPC/E water model, performed exceptionally well for all properties. The results reported here establish useful benchmark calculations for Mg(2+) ion models that provide insight into the origin of the behavior in aqueous solution, and may aid in the development of next-generation models that target specific binding sites in biomolecules.