Computing ensembles of transitions from stable states: Dynamic importance sampling

There is an increasing dataset of solved biomolecular structures in more than one conformation and increasing evidence that large-scale conformational change is critical for biomolecular function. In this article, we present our implementation of a dynamic importance sampling (DIMS) algorithm that is directed toward improving our understanding of important intermediate states between experimentally defined starting and ending points. This complements traditional molecular dynamics methods where most of the sampling time is spent in the stable free energy wells defined by these initial and final points. As such, the algorithm creates a candidate set of transitions that provide insights for the much slower and probably most important, functionally relevant degrees of freedom. The method is implemented in the program CHARMM and is tested on six systems of growing size and complexity. These systems, the folding of Protein A and of Protein G, the conformational changes in the calcium sensor S100A6, the glucose-galactose-binding protein, maltodextrin, and lactoferrin, are also compared against other approaches that have been suggested in the literature. The results suggest good sampling on a diverse set of intermediates for all six systems with an ability to control the bias and thus to sample distributions of trajectories for the analysis of intermediate states.     

Quantum mechanical calculation of the effects of stiff and rigid constraints in the conformational equilibrium of the alanine dipeptide

If constraints are imposed on a macromolecule, two inequivalent classical models may be used: the stiff and the rigid one. This work studies the effects of such constraints on the conformational equilibrium distribution (CED) of the model dipeptide HCO-L-Ala-NH(2)without any simplifying assumption. We use ab initio quantum mechanics calculations including electron correlation at the MP2 level to describe the system, and we measure the conformational dependence of all the correcting terms to the naive CED based in the potential energy surface that appear when the constraints are considered. These terms are related to mass-metric tensors determinants and also occur in the Fixman's compensating potential. We show that some of the corrections are non-negligible if one is interested in the whole Ramachandran space. On the other hand, if only the energetically lower region, containing the principal secondary structure elements, is assumed to be relevant, then, all correcting terms may be neglected up to peptides of considerable length. This is the first time, as far as we know, that the analysis of the conformational dependence of these correcting terms is performed in a relevant biomolecule with a realistic potential energy function.     

Analysis of the adaptive multilevel splitting method on the isomerization of alanine dipeptide

We apply the adaptive multilevel splitting (AMS) method to the C _eq → C _ax transition of alanine dipeptide in vacuum. Some properties of the algorithm are numerically illustrated, such as the unbiasedness of the probability estimator and the robustness of the method with respect to the reaction coordinate. We also calculate the transition time obtained via the probability estimator, using an appropriate ensemble of initial conditions. Finally, we show how the AMS method can be used to compute an approximation of the committor function. © 2019 Wiley Periodicals, Inc.     

Conformational sampling with Poisson–Boltzmann forces and a stochastic dynamics/Monte Carlo method: Application to alanine dipeptide

We apply a combination of stochastic dynamics and Monte Carlo methods (MC/SD) to alanine dipeptide, with solvation forces derived from a Poisson–Boltzmann model supplemented with apolar terms. Our purpose is to study the effects of the model parameters, such as the friction constant and the size of the electrostatic finite difference grid, on the rate of conformational sampling and on the accuracy of the resulting free energy map. For dialanine, a converged Ramachandran map is produced in significantly less time than what is required by stochastic dynamics or Monte Carlo alone. MC/SD is also shown to be faster, per timestep, than explicit methods. © 1997 John Wiley & Sons, Inc. J Comput Chem 18 : 1750–1759, 1997     

Effects of external electromagnetic fields on the conformational sampling of a short alanine peptide

Non‐equilibrium molecular dynamics simulations of a solvated 21‐residue polyalanine (A21) peptide, featuring a high propensity for helix formation, have been performed at 300 K and 1 bar in the presence of external electromagnetic (e/m) fields in the microwave region (2.45 GHz) and an r.m.s. electric field intensity range of 0.01–0.05 V/Å. To investigate how the field presence affects transitions between the conformational states of a protein, we report 16 independent 40 ns‐trajectories of A21 starting from both extended and fully folded states. We observe folding‐behavior of the peptide consistent with prior simulation and experimental studies. The peptide displays a natural tendency to form stable elements of secondary structure which are stabilized by tertiary interactions with proximate regions of the peptide. Consistent with our earlier work, the presence of external e/m fields disrupts this behavior, involving a mechanism of localized dipolar alignment which serves to enhance intra‐protein perturbations in hydrogen bonds (English, et al., J. Chem. Phys. 2010 , 133, 091105), leading to more frequent transitions between shorter‐lifetime states. © 2012 Wiley Periodicals, Inc.     

Stochastic dynamics simulation of alanine dipeptide: Including solvation interaction determined by boundary element method

A merger of the Poisson–Boltzmann equation and stochastic dynamics simulation is examined using illustrative calculations of alanine dipeptide. The boundary element method (BEM) is used to calculate the hydration forces acting on the solute molecule based on the surroundings. Computational efficiency is achieved by the use of a simple hydration model and a coarse boundary element. Nonetheless, the conformational distribution obtained from this new method is reasonable compared with other theoretical and computational results. Detailed analysis has been accomplished in terms of the hydration interactions and solvation energies. The results indicate that the new simulation method provides an obvious improvement over the conventional stochastic dynamics simulation technique. The further improvement of the hydration model and future application to large molecules are also discussed. © 1997 John Wiley & Sons, Inc. J Comput Chem 18 : 1440–1449, 1997     

A structure-based design approach for the identification of novel inhibitors: application to an alanine racemase

We report a new structure-based strategy for the identification of novel inhibitors. This approach has been applied to Bacillus stearothermophilus alanine racemase (AlaR), an enzyme implicated in the biosynthesis of the bacterial cell wall. The enzyme catalyzes the racemization of l- and d-alanine using pyridoxal 5-phosphate (PLP) as a cofactor. The restriction of AlaR to bacteria and some fungi and the absolute requirement for d-alanine in peptidoglycan biosynthesis make alanine racemase a suitable target for drug design. Unfortunately, known inhibitors of alanine racemase are not specific and inhibit the activity of other PLP-dependent enzymes, leading to neurological and other side effects.This article describes the development of a receptor-based pharmacophore model for AlaR, taking into account receptor flexibility (i.e. a `dynamic' pharmacophore model). In order to accomplish this, molecular dynamics (MD) simulations were performed on the full AlaR dimer from Bacillus stearothermophilus (PDB entry, 1sft) with a d-alanine molecule in one active site and the non-covalent inhibitor, propionate, in the second active site of this homodimer. The basic strategy followed in this study was to utilize conformations of the protein obtained during MD simulations to generate a dynamic pharmacophore model using the property mapping capability of the LigBuilder program. Compounds from the Available Chemicals Directory that fit the pharmacophore model were identified and have been submitted for experimental testing.The approach described here can be used as a valuable tool for the design of novel inhibitors of other biomolecular targets.     

Critical assessment of side-chain conformational space sampling procedures designed for quantifying the effect of side-chain environment

We introduce a family of procedures designed to sample side-chain conformational space at particular locations in protein structures. These procedures (CRSP) use intensive cycles of random assignment of side-chain conformations followed by minimization to determine all the conformations that a group of side-chains can adopt simultaneously. First, we consider a procedure evolving in the dihedral space (dCRSP). Our results suggest that it can accurately map low-energy conformations adopted by clusters of side-chains of a protein. dCRSP is relatively insensitive to various important parameters, and it is sufficiently accurate to capture efficiently the constraint induced by the environment on the conformations a particular side-chain can adopt. Our results show that dCRSP, compared with molecular dynamics (MD), can overcome the problem of the limited set of conformations reached in a reasonable amount of simulations. Next, we introduce procedures (vCRSP) in which valence angles are relaxed, and we assess how efficiently they quantify the conformational entropy of side-chains in the protein native state. For simple peptides, entropies obtained with vCRSP are fully compatible with those obtained with a Monte Carlo procedure. For side-chains in a protein environment, however, vCRSP appears of limited use. Finally, we consider a two-step procedure that combines dCRSP and vCRSP. Our tests suggest that it is able to overcome the limitations of vCRSP. We also note that dCRSP provides a reasonable initial approximation. This family of procedures offers promise in quantifying the contribution of conformational entropy to the energetics of protein structures.     

Sparsity‐weighted outlier FLOODing (OFLOOD) method: Efficient rare event sampling method using sparsity of distribution

As an extension of the Outlier FLOODing (OFLOOD) method [Harada et al., J. Comput. Chem. 2015, 36, 763], the sparsity of the outliers defined by a hierarchical clustering algorithm, FlexDice, was considered to achieve an efficient conformational search as sparsity‐weighted “OFLOOD.” In OFLOOD, FlexDice detects areas of sparse distribution as outliers. The outliers are regarded as candidates that have high potential to promote conformational transitions and are employed as initial structures for conformational resampling by restarting molecular dynamics simulations. When detecting outliers, FlexDice defines a rank in the hierarchy for each outlier, which relates to sparsity in the distribution. In this study, we define a lower rank (first ranked), a medium rank (second ranked), and the highest rank (third ranked) outliers, respectively. For instance, the first‐ranked outliers are located in a given conformational space away from the clusters (highly sparse distribution), whereas those with the third‐ranked outliers are nearby the clusters (a moderately sparse distribution). To achieve the conformational search efficiently, resampling from the outliers with a given rank is performed. As demonstrations, this method was applied to several model systems: Alanine dipeptide, Met‐enkephalin, Trp‐cage, T4 lysozyme, and glutamine binding protein. In each demonstration, the present method successfully reproduced transitions among metastable states. In particular, the first‐ranked OFLOOD highly accelerated the exploration of conformational space by expanding the edges. In contrast, the third‐ranked OFLOOD reproduced local transitions among neighboring metastable states intensively. For quantitatively evaluations of sampled snapshots, free energy calculations were performed with a combination of umbrella samplings, providing rigorous landscapes of the biomolecules. © 2015 Wiley Periodicals, Inc.     

Analyzing the robustness of the MM/PBSA free energy calculation method: application to DNA conformational transitions

The ability to predict and characterize free energy differences associated with conformational equilibria or the binding of biomolecules is vital to understanding the molecular basis of many important biological functions. As biological studies focus on larger molecular complexes and properties of the genome, proteome, and interactome, the development and characterization of efficient methods for calculating free energy becomes increasingly essential. The aim of this study is to examine the robustness of the end-point free energy method termed the molecular mechanics Poisson-Boltzmann solvent accessible surface area (MM/PBSA) method. Specifically, applications of MM/PBSA to the conformational equilibria of nucleic acid (NA) systems are explored. This is achieved by comparing A to B form DNA conformational free energy differences calculated using MM/PBSA with corresponding free energy differences determined with a more rigorous and time-consuming umbrella sampling algorithm. In addition, the robustness of NA MM/PBSA calculations is also evaluated in terms of the sensitivity towards the choice of force field and the choice of solvent model used during conformational sampling. MM/PBSA calculations of the free energy difference between A-form and B-form DNA are shown to be in very close agreement with the PMF result determined using an umbrella sampling approach. Further, it is found that the MM/PBSA conformational free energy differences were also in agreement using either the CHARMM or AMBER force field. The influence of ionic strength on conformational stability was particularly insensitive to the choice of force field. Finally, it is also shown that the use of a generalized Born implicit solvent during conformational sampling results in free energy estimates that deviate slightly from those obtained using explicitly solvated MD simulations in these NA systems.     

Enhanced conformational sampling method for proteins based on the TaBoo SeArch algorithm: Application to the folding of a mini‐protein,chignolin

The conformational samplings are indispensible for obtaining reliable canonical ensembles, which provide statistical averages of physical quantities such as free energies. However, the samplings of vast conformational space of biomacromolecules by conventional molecular dynamics (MD) simulations might be insufficient, due to their inadequate accessible time‐scales for investigating biological functions. Therefore, the development of methodologies for enhancing the conformational sampling of biomacromolecules still remains as a challenging issue in computational biology. To tackle this problem, we newly propose an efficient conformational search method, which is referred as TaBoo SeArch (TBSA) algorithm. In TBSA, an inverse energy histogram is used to select seeds for the conformational resampling so that states with high frequencies are inhibited, while states with low frequencies are efficiently sampled to explore the unvisited conformational space. As a demonstration, TBSA was applied to the folding of a mini‐protein, chignolin, and automatically sampled the native structure (C_α root mean square deviation < 1.0 Å) with nanosecond order computational costs started from a completely extended structure, although a long‐time 1‐µs normal MD simulation failed to sample the native structure. Furthermore, a multiscale free energy landscape method based on the conformational sampling of TBSA were quantitatively evaluated through free energy calculations with both implicit and explicit solvent models, which enable us to find several metastable states on the folding landscape. © 2015 Wiley Periodicals, Inc.     

Finding the transition state without initial guess: the growing string method for Newton trajectory to isomerization and enantiomerization reaction of alanine dipeptide and poly(15)alanine

We report a new, high-dimensional application of a method for finding a transition state (TS) between a reactant and a product on the potential energy surface: the search of a growing string along a reaction path defined by any Newton trajectory in combination with the Berny method (Quapp, J Chem Phys (2005), 122, 174106; we have provided this algorithm on a web page). Two given minima are connected by a one-dimensional, but usually curvilinear reaction coordinate. It leads to the TS region. The application of the method to alanine dipeptide finds the TS of the isomerisation C(7 ax) --> C(5), some TSs of the enantiomerisation of C(7 ax) from L-form to quasi-D-form, and it finds the TS region of a transition of a partly unfolded, bent structure which turns back into a mainly alpha-helix in the Ac(Ala)(15)NHMe polyalanine (all at the quantum mechanical level B3LYP/6-31G: the growing string calculation is interfaced with the Gaussian03 package). The formation or dissolvation of some alpha- or 3(10)-hydrogen bonds of the helix are discussed along the TS pathway, as well as the case of an enantiomer at the central residue of the helix.     

Calculation of rotational partition functions by an efficient Monte Carlo importance sampling technique

The evaluation of the classical rotational partition function represented by a configuration integral over all external and internal rotational degrees of freedom of nonrigid chain polyatomic molecules is described. The method of Pitzer and Gwinn is used to correct the classical partition function for quantum mechanical effects at low temperatures. The internal rotor hindrance and all coupling arising from the external and internal rotational degrees of freedom are explicitly taken into account. Importance sampling Monte Carlo based on the adaptive VEGAS algorithm to perform multidimensional integration is implemented within the TINKER program package. A multidimensional potential energy hypersurface is calculated with the MM3(2000) molecular mechanics force field. Numerical tests are performed on a number of small n-alkanes (from ethane to octane), for which the absolute entropies calculated at three different temperatures are compared both with the experimental values and with the previous theoretical results. The application of a more efficient importance sampling technique developed here results in a substantial reduction of statistical errors in the evaluation of the configuration integral for a given number of Monte Carlo steps. Error estimates for the calculated entropies are given, and possible sources of systematic errors, and their importance for a reliable prediction of the absolute entropy, are discussed.     

A method for conformational sampling of loops in proteins based on adiabatic decoupling and temperature or force scaling

A method for conformational Boltzmann sampling of loops in proteins in aqueous solution is presented that is based on adiabatic decoupling molecular dynamics (MD) simulation with temperature or force scaling. To illustrate the enhanced sampling, the loop from residues 33 to 43 in the bovine protein ribonuclease A is adiabatically decoupled from the rest of the protein and the solvent with a mass scaling factor s(m) =1000 and the sampling is enhanced with a scaling of the temperature using s(T) =2 or of the force using s(V) =0.667. Over 5 ns of simulation the secondary structure of the protein remains unaltered while a combined dihedral-angle conformational cluster analysis shows an increase of conformations outside the first most populated cluster of loop conformations for adiabatic decoupling MD with temperature scaling using s(T) =2 or force scaling using s(V) =0.667 compared to the standard MD simulation. The atom-positional root-mean-square fluctuations of the C(α) atoms of the loop show an increase in the movement of the loop as well, indicating that adiabatic decoupling MD with upscaling of the temperature or downscaling of the force is a promising method for conformational Boltzmann sampling.     

Transition pathways in complex systems: Application of the finite-temperature string method to the alanine dipeptide

The finite-temperature string method proposed by E, et al. [W. E, W. Ren, and E. Vanden-Eijnden, Phys. Rev. B 66, 052301 (2002)] is a very effective way of identifying transition mechanisms and transition rates between metastable states in systems with complex energy landscapes. In this paper, we discuss the theoretical background and algorithmic details of the finite-temperature string method, as well as the application to the study of isomerization reaction of the alanine dipeptide, both in vacuum and in explicit solvent. We demonstrate that the method allows us to identify directly the isocommittor surfaces, which are approximated by hyperplanes, in the region of configuration space where the most probable transition trajectories are concentrated. These results are verified subsequently by computing directly the committor distribution on the hyperplanes that define the transition state region.     

Synthesis of a new quaternary phosphonium salt: NMR study of the conformational structure and dynamics

A novel phosphonium salt based on pyridoxine was synthesized. Conformational analysis of the compound in solution was performed using dynamic NMR experiments and calculations. The obtained results revealed some differences in the conformational transitions and the energy parameters of the conformational exchange of the studied compound in comparison to previously reported data for other phosphorus‐containing pyridoxine derivatives. It was shown that increasing the substituent at the C‐11 carbon leads to greater differences in the populations of stable states and the corresponding equilibrium energies. Copyright © 2015 John Wiley & Sons, Ltd.     

Nonadiabatic effects in the two-dimensional infrared spectra of peptides: application to alanine dipeptide

A method of simulating two-dimensional infrared spectra accounting for nonadiabatic effects is presented. The method is applied to the amide I modes of a dipeptide. The information necessary to construct the time-dependent Hamiltonian for the system is extracted from molecular dynamics simulations using a recently published ab initio-based model. It is shown that the linear absorption spectrum agrees with experiment only if the nonadiabatic effects are accounted for. The two-dimensional infrared spectrum is predicted for a range of mixing times. It is shown that population transfer between the amide I site vibrations affects the anisotropy at longer mixing times. It is also demonstrated that the population transfer can, to a good approximation, be extracted from the simulated spectra using a procedure that should also be applicable to experimental spectra.     

The use of a generalized born model for the analysis of protein conformational transitions: a comparative study with explicit solvent simulations for chemotaxis Y protein (CheY)

To investigate whether implicit solvent models are appropriate for mechanistic studies of conformational transition in proteins, a recently developed generalized Born model (GBSW) was applied to a small signaling protein, chemotaxis protein Y (CheY), with different combinations of the phosphorylation state and conformation of the system; the results were compared to explicit solvent simulations using a stochastic boundary condition. The subtle but distinct conformational transitions involved in CheY activation makes the system ideally suited for comparing implicit and explicit solvent models because these conformational transitions are potentially accessible in both types of simulations. The structural and dynamical properties analyzed include not only those localized to the active site region but also throughout the protein, such as sidechain methyl group order parameters, backbone hydrogen bonding lifetime and occupancy as well as principal components of the trajectories. Overall, many properties were well reproduced by the GBSW simulations when compared with the explicit solvent calculations, although a number of observations consistently point to the suggestion that the current parameterization of the GBSW model tends to overestimate hydrogen-bonding interactions involving both charged groups and (charge-neutral) backbone atoms. This deficiency led to overstabilization of certain secondary structural motifs and more importantly, qualitatively different behaviors for the active site groups (Thr 87, Ala 88, the beta4-alpha4 loop) in response to phosphorylation, when compared with explicit solvent simulations. The current study highlights the value of carrying out both explicit and implicit solvent simulations for complementary mechanistic insights in the analysis of conformational transition in biomolecules.