Convergence of folding free energy landscapes via application of enhanced sampling methods in a distributed computing environment

We have implemented the serial replica exchange method (SREM) and simulated tempering (ST) enhanced sampling algorithms in a global distributed computing environment. Here we examine the helix-coil transition of a 21 residue alpha-helical peptide in explicit solvent. For ST, we demonstrate the efficacy of a new method for determining initial weights allowing the system to perform a random walk in temperature space based on short trial simulations. These weights are updated throughout the production simulation by an adaptive weighting method. We give a detailed comparison of SREM, ST, as well as standard MD and find that SREM and ST give equivalent results in reasonable agreement with experimental data. In addition, we find that both enhanced sampling methods are much more efficient than standard MD simulations. The melting temperature of the Fs peptide with the AMBER99phi potential was calculated to be about 310 K, which is in reasonable agreement with the experimental value of 334 K. We also discuss other temperature dependent properties of the helix-coil transition. Although ST has certain advantages over SREM, both SREM and ST are shown to be powerful methods via distributed computing and will be applied extensively in future studies of complex bimolecular systems.     

Practical considerations for building GROMOS-compatible small-molecule topologies

Molecular dynamics simulations are being applied to increasingly complex systems, including those involving small endogenous compounds and drug molecules. In order to obtain meaningful and accurate data from these simulations, high-quality topologies for small molecules must be generated in a manner that is consistent with the derivation of the force field applied to the system. Often, force fields are designed for use with macromolecules such as proteins, making their transferability to other species challenging. Investigators are increasingly attracted to automated topology generation programs, although the quality of the resulting topologies remains unknown. Here we assess the applicability of the popular PRODRG server that generates small-molecule topologies for use with the GROMOS family of force fields. We find that PRODRG does not reproduce topologies for even the most well-characterized species in the force field due to inconsistent charges and charge groups. We assessed the effects of PRODRG-derived charges on several systems: pure liquids, amino acids at a hydrophobic-hydrophilic interface, and an enzyme-cofactor complex. We found that partial atomic charges generated by PRODRG are largely incompatible with GROMOS force fields, and the behavior of these systems deviates substantially from that of simulations using GROMOS parameters. We conclude by proposing several points as "best practices" for parametrization of small molecules under the GROMOS force fields.     

Converging free energies of binding in cucurbit[7]uril and octa-acid host–guest systems from SAMPL4 using expanded ensemble simulations

Molecular containers such as cucurbit[7]uril (CB7) and the octa-acid (OA) host are ideal simplified model test systems for optimizing and analyzing methods for computing free energies of binding intended for use with biologically relevant protein–ligand complexes. To this end, we have performed initially blind free energy calculations to determine the free energies of binding for ligands of both the CB7 and OA hosts. A subset of the selected guest molecules were those included in the SAMPL4 prediction challenge. Using expanded ensemble simulations in the dimension of coupling host–guest intermolecular interactions, we are able to show that our estimates in most cases can be demonstrated to fully converge and that the errors in our estimates are due almost entirely to the assigned force field parameters and the choice of environmental conditions used to model experiment. We confirm the convergence through the use of alternative simulation methodologies and thermodynamic pathways, analyzing sampled conformations, and directly observing changes of the free energy with respect to simulation time. Our results demonstrate the benefits of enhanced sampling of multiple local free energy minima made possible by the use of expanded ensemble molecular dynamics and may indicate the presence of significant problems with current transferable force fields for organic molecules when used for calculating binding affinities, especially in non-protein chemistries.     

Helix formation in a pentapeptide: experiment and force-field dependent dynamics

We used a combined approach of experiment and simulation to determine the helical population and folding pathway of a small helix forming blocked pentapeptide, Ac-(Ala)(5)-NH(2). Experimental structural characterization of this blocked peptide was carried out with far UV circular dichroism spectroscopy, FTIR, and NMR measurements. These measurements confirm the presence of the α-helical state in a buffer solution. Direct molecular dynamics and replica-exchange simulations of the pentapeptide were performed using several popular force fields with explicit solvent. The simulations yielded statistically reliable estimates of helix populations, melting curves, folding, and nucleation times. The distributions of conformer populations are used to measure folding cooperativity. Finally, a statistical analysis of the sample of helix-coil transition paths was performed. The details of the calculated helix populations, folding kinetics and pathways vary with the employed force field. Interestingly, the helix populations, folding, and unfolding times obtained from most of the studied force fields are in qualitative agreement with each other and with available experimental data, with the deviations corresponding to several kcal/mol in energy at 300 K. Most of the force fields also predict qualitatively similar transition paths, with unfolding initiated at the C-terminus. Accuracy of potential energy parameters, rather than conformational sampling may be the limiting factor in current molecular simulations.     

A Robust and Cost-Efficient Scheme for Accurate Conformational Energies of Organic Molecules

Several standard semiempirical methods as well as the MMFF94 force field approximation have been tested in reproducing 8 DLPNO-CCSD(T)/cc-pVTZ level conformational energies and spatial structures for 37 organic molecules representing pharmaceuticals, drugs, catalysts, synthetic precursors, industry-related chemicals (37conf8 database). All contemporary semiempirical methods surpass their standard counterparts resulting in more reliable conformational energies and spatial structures, even though at significantly higher computational costs. However, even these methods show unexpected failures in reproducing energy differences between several conformers of the crown ether 1,4,7,10,13,16-hexaoxacyclooctadecane (18-crown-6). Inexpensive force field MMFF94 approximation groups with contemporary semiempirical methods in reproducing the correct order of conformational energies and spatial structures, although the performance in predicting absolute conformational energies compares to standard semiempirical methods. Based on these findings, we suggest a two-step strategy for reliable yet feasible conformational search and sampling in realistic-size flexible organic molecules: i) geometry optimization/preselection of relevant conformers using the MMFF94 force field; ii) single-point energy evaluations using a contemporary semiempirical method. We expect that developed database 37conf8 is going to be useful for development of semiempirical methods.     

Globular structures of a helix-coil copolymer: self-consistent treatment

A self-consistent-field theory was developed in the grand canonical ensemble formulation to study transitions in a helix-coil multiblock globule. Helical and coil parts are treated as stiff rods and self-avoiding walks of variable lengths correspondingly. The resulting field theory takes, in addition to the conventional Zimm-Bragg, [J. Chem. Phys. 31, 526 (1959)] parameters, also three-dimensional interaction terms into account. The appropriate differential equations which determine the self-consistent fields were solved numerically with finite element method. Three different phase states are found: open chain, amorphous globule, and nematic liquid-crystalline (LC) globule. The LC-globule formation is driven by the interplay between the hydrophobic helical segment attraction and the anisotropic globule surface energy of an entropic nature. The full phase diagram of the helix-coil copolymer was calculated and thoroughly discussed. The suggested theory shows a clear interplay between secondary and tertiary structures in globular homopolypeptides.     

Are Protein Force Fields Getting Better? A Systematic Benchmark on 524 Diverse NMR Measurements

Recent hardware and software advances have enabled simulation studies of protein systems on biophysically-relevant timescales, often revealing the need for improved force fields. Although early force field development was limited by the lack of direct comparisons between simulation and experiment, recent work from several labs has demonstrated direct calculation of NMR observables from protein simulations. Here we quantitatively evaluate recent molecular dynamics force fields against a suite of 524 chemical shift and J coupling ((3)JH(N)H(α), (3)JH(N)C(β), (3)JH(α)C', (3)JH(N)C', and (3)JH(α)N) measurements on dipeptides, tripeptides, tetra-alanine, and ubiquitin. Of the force fields examined (ff96, ff99, ff03, ff03*, ff03w, ff99sb*, ff99sb-ildn, ff99sb-ildn-phi, ff99sb-ildn-nmr, CHARMM27, OPLS-AA), two force fields (ff99sb-ildn-phi, ff99sb-ildn-nmr) combining recent side chain and backbone torsion modifications achieve high accuracy in our benchmark. For the two optimal force fields, the calculation error is comparable to the uncertainty in the experimental comparison. This observation suggests that extracting additional force field improvements from NMR data may require increased accuracy in J coupling and chemical shift prediction. To further investigate the limitations of current force fields, we also consider conformational populations of dipeptides, which were recently estimated using vibrational spectroscopy.     

Comparison of charge models for fixed-charge force fields: small-molecule hydration free energies in explicit solvent

In molecular simulations with fixed-charge force fields, the choice of partial atomic charges influences numerous computed physical properties, including binding free energies. Many molecular mechanics force fields specify how nonbonded parameters should be determined, but various choices are often available for how these charges are to be determined for arbitrary small molecules. Here, we compute hydration free energies for a set of 44 small, neutral molecules in two different explicit water models (TIP3P and TIP4P-Ew) to examine the influence of charge model on agreement with experiment. Using the AMBER GAFF force field for nonbonded parameters, we test several different methods for obtaining partial atomic charges, including two fast methods exploiting semiempirical quantum calculations and methods deriving charges from the electrostatic potentials computed with several different levels of ab initio quantum calculations with and without a continuum reaction field treatment of solvent. We find that the best charge sets give a root-mean-square error from experiment of roughly 1 kcal/mol. Surprisingly, agreement with experimental hydration free energies does not increase substantially with increasing level of quantum theory, even when the quantum calculations are performed with a reaction field treatment to better model the aqueous phase. We also find that the semiempirical AM1-BCC method for computing charges works almost as well as any of the more computationally expensive ab initio methods and that the root-mean-square error reported here is similar to that for implicit solvent models reported in the literature. Further, we find that the discrepancy with experimental hydration free energies grows substantially with the polarity of the compound, as does its variation across theory levels.     

Prediction of adsorption of small molecules in porous materials based on ab initio force field method

Computational prediction of adsorption of small molecules in porous materials has great impact on the basic and applied research in chemical engineering and material sciences. In this work,we report an approach based on grand canonical ensemble Monte Carlo(GCMC) simulations and ab initio force fields. We calculated the adsorption curves of ammonia in ZSM-5 zeolite and hydrogen in MOF-5(a metal-organic-framework material). The predictions agree well with experimental data. Because the predictions are based on the first principle force fields,this approach can be used for the adsorption prediction of new molecules or materials without experimental data as guidance.     

Foldamer simulations: novel computational methods and applications to poly-phenylacetylene oligomers

We apply several methods to probe the ensemble kinetic and structural properties of a model system of poly-phenylacetylene (pPA) oligomer folding trajectories. The kinetic methods employed included a brute force accounting of conformations, a Markovian state matrix method, and a nonlinear least squares fit to a minimalist kinetic model used to extract the folding time. Each method gave similar measures for the folding time of the 12-mer chain, calculated to be on the order of 7 ns for the complete folding of the chain from an extended conformation. Utilizing both a linear and a nonlinear scaling relationship between the viscosity and the folding time to correct for a low simulation viscosity, we obtain an upper and a lower bound for the approximate folding time within the range 70 ns相似文献   

Using massively parallel simulation and Markovian models to study protein folding: examining the dynamics of the villin headpiece

We report on the use of large-scale distributed computing simulation and novel analysis techniques for examining the dynamics of a small protein. Matters addressed include folding rate, very long time scale kinetics, ensemble properties, and interaction with water. The target system for the study, the villin headpiece, has been of great interest to experimentalists and theorists both. Sampling totaled nearly 500 mus-the most extensive published to date for a system of villin's size in explicit solvent with all atom detail-and was in the form of tens of thousands of independent molecular dynamics trajectories, each several tens of nanoseconds in length. We report on kinetics sensitivity analyses that, using a set of short simulations, probed the role of water in villin's folding and sensitivity to the simulation's electrostatics treatment. By constructing Markovian state models (MSMs) from the collected data, we were able to propagate dynamics to times far beyond those directly simulated and to rapidly compute mean first passage times, long time kinetics (tens of microseconds), and evolution of ensemble property distributions over long times, otherwise currently impossible. We also tested our MSM by using it to predict the structure of villin de novo.     

Toward ab initio refinement of protein X-ray crystal structures: interpreting and correlating structural fluctuations

The refinement of protein crystal structures currently involves the use of empirical restraints and force fields that are known to work well in many situations but nevertheless yield structural models with some features that are inconsistent with detailed chemical analysis and therefore warrant further improvement. Ab initio electronic structure computational methods have now advanced to the point at which they can deliver reliable results for macromolecules in realistic times using linear-scaling algorithms. The replacement of empirical force fields with ab initio methods in a final refinement stage could allow new structural features to be identified in complex structures, reduce errors and remove computational bias from structural models. In contrast to empirical approaches, ab initio refinements can only be performed on models that obey basic qualitative chemical rules, imposing constraints on the parameter space of existing refinements, and this in turn inhibits the inclusion of unlikely structural features. Here, we focus on methods for determining an appropriate ensemble of initial structural models for an ab initio X-ray refinement, modeling as an example the high-resolution single-crystal X-ray diffraction data reported for the structure of lysozyme (PDB entry “2VB1”). The AMBER force field is used in a Monte Carlo calculation to determine an ensemble of 8 structures that together embody all of the partial atomic occupancies noted in the original refinement, correlating these variations into a set of feasible chemical structures while simultaneously retaining consistency with the X-ray diffraction data. Subsequent analysis of these results strongly suggests that the occupancies in the empirically refined model are inconsistent with protein energetic considerations, thus depicting the 2VB1 structure as a deep-lying minimum in its optimized parameter space that actually embodies chemically unreasonable features. Indeed, density-functional theory calculations for one specific nitrate ion with an occupancy of 62% indicate that water replaces this ion 38% of the time, a result confirmed by subsequent crystallographic analysis. It is foreseeable that any subsequent ab initio refinement of the whole structure would need to locate a globally improved structure involving significant changes to 2VB1 which correct these identified local structural inconsistencies.     

Long range behavior of high-rank topological multipole moments

The construction of a high-rank multipolar force field (for peptides) is a complex task, leading to several intermediate questions in need of a clear answer. Here we focus on the convergence of the (electrostatic) multipolar expansion at medium and long range. Using molecular electron densities, quantum chemical topology (QCT) defines the atoms as finite volumes, each endowed with multipole moments. The terms in the multipole expansion are grouped according to powers of the internuclear distance, R(-L). Given two atom types at a given distance, we determine which rank (L) is necessary for the electrostatic energy to converge to the exact interaction energy within a certain error. With this information, the rank of the expansion for each interaction can be adapted to the required accuracy and the available computing power.     

Paramfit: Automated optimization of force field parameters for molecular dynamics simulations

The generation of bond, angle, and torsion parameters for classical molecular dynamics force fields typically requires fitting parameters such that classical properties such as energies and gradients match precalculated quantum data for structures that scan the value of interest. We present a program, Paramfit, distributed as part of the AmberTools software package that automates and extends this fitting process, allowing for simplified parameter generation for applications ranging from single molecules to entire force fields. Paramfit implements a novel combination of a genetic and simplex algorithm to find the optimal set of parameters that replicate either quantum energy or force data. The program allows for the derivation of multiple parameters simultaneously using significantly fewer quantum calculations than previous methods, and can also fit parameters across multiple molecules with applications to force field development. Paramfit has been applied successfully to systems with a sparse number of structures, and has already proven crucial in the development of the Assisted Model Building with Energy Refinement Lipid14 force field. © 2014 Wiley Periodicals, Inc.     

Computer simulations of ligand–protein binding with ensembles of protein conformations: A Monte Carlo study of HIV‐1 protease binding energy landscapes

We present the results of molecular docking simulations with HIV‐1 protease for the sb203386 and skf107457 inhibitors by Monte Carlo simulated annealing. A simplified piecewise linear energy function, the standard AMBER force field, and the AMBER force field with solvation and a soft‐core smoothing component are employed in simulations with a single‐protein conformation to determine the relationship between docking simulations with a simple energy function and more realistic force fields. The temperature‐dependent binding free energy profiles of the inhibitors interacting with a single protein conformation provide a detailed picture of relative thermodynamic stability and a distribution of ligand binding modes in agreement with experimental crystallographic data. Using the simplified piecewise linear energy function, we also performed Monte Carlo docking simulations with an ensemble of protein conformations employing preferential biased sampling of low‐energy protein conformations, and the results are analyzed in connection with the free energy profiles. ©1999 John Wiley & Sons, Inc. Int J Quant Chem 72: 73–84, 1999     

An approach to developing a force field for molecular simulation of martensitic phase transitions between phases with subtle differences in energy and structure

d,l-Norleucine is one of only a few molecules whose crystals exhibit a martensitic or displacive-type phase transformation where the emerging phase shows a topotaxial relationship with the parent phase. The molecular mechanism for such phase transformations, particularly in molecular crystals, is not well understood. Crystalline phases that exhibit displacive phase transitions tend to be very similar in structure and energy. Consequently, the development of a force field for such phases is challenging as the phase behavior is determined by subtle differences in their lattice energies and entropies. We report an approach for developing a force field for such phases with an application to d,l-norleucine. The proposed procedure includes calculation of the phase diagram of the crystalline phases as a function of temperature to identify the best force field. d,l-Norleucine also presents an additional problem since in the solid state it exists as a zwitterion that is unstable in vacuo and therefore cannot be characterized using high-level ab initio calculations in the gas phase. However, a stable zwitterion could be obtained using Onsager's reaction-field continuum model for a solvent (SCRF) using both Hartree-Fock and density functional theory. A number of force fields and the various sets of partial charges obtained from the SCRF calculations were screened for their ability to reproduce the crystal structures of the two known phases, alpha and beta, of d,l-norleucine. Selected parameter sets were then employed in free energy minimizations to identify the best set on the basis of a correct prediction of the alpha-beta phase transition. The Williams' nonbonded parameters combined with partial charges from SCRF-Polarized Continuum Model calculation were found to reproduce the structures of the phases accurately and also maintained their stability in extended molecular dynamics simulations in the Parrinello-Rahman constant stress ensemble. Moreover, we were also able to successfully simulate the phase transformation of the beta- to the alpha-phase. The identified force field should enable detailed studies of the phase transformations exhibited by crystals of d,l-norleucine and hence enhance our understanding of martensitic-type transformations in molecular crystals.     

Prediction of adsorption of small molecules in porous materials based on <Emphasis Type="Italic">ab initio</Emphasis> force field method

Computational prediction of adsorption of small molecules in porous materials has great impact on the basic and applied research in chemical engineering and material sciences. In this work, we report an approach based on grand canonical ensemble Monte Carlo (GCMC) simulations and ab initio force fields. We calculated the adsorption curves of ammonia in ZSM-5 zeolite and hydrogen in MOF-5 (a metal-organic-framework material). The predictions agree well with experimental data. Because the predictions are based on the first principle force fields, this approach can be used for the adsorption prediction of new molecules or materials without experimental data as guidance. Supported by the National “973” Program (Grant No. 2007CB209701, 2003CB615804)     

A review of Cloud computing technologies for comprehensive microRNA analyses

Cloud computing revolutionized many fields that require ample computational power. Cloud platforms may also provide huge support for microRNA analysis mainly through disclosing scalable resources of different types. In Clouds, these resources are available as services, which simplifies their allocation and releasing. This feature is especially useful during the analysis of large volumes of data, like the one produced by next generation sequencing experiments, which require not only extended storage space but also a distributed computing environment. In this paper, we show which of the Cloud properties and service models can be especially beneficial for microRNA analysis. We also explain the most useful services of the Cloud (including storage space, computational power, web application hosting, machine learning models, and Big Data frameworks) that can be used for microRNA analysis. At the same time, we review several solutions for microRNA and show that the utilization of the Cloud in this field is still weak, but can increase in the future when the awareness of their applicability grows.     

Electrohydrodynamics of a drop under nonaxisymmetric electric fields

We consider the electrohydrodynamics of a spherical drop in a nonaxisymmetric electric field, which can be approximated by the sum of a uniform field and a linear straining field. We obtain the analytic solution of the three-dimensional flow fields inside and outside a drop for the Stokes flow regime by using Lamb's general solution and the leaky dielectric model. With the analytic solution, the dielectrophoretic migration velocity of a drop is obtained as a function of the type and the frequency of the imposed electric field. The direction of drop motion is found to be parallel to the dielectrophoretic force. The analytic solution is also used to investigate the characteristics of the interfacial flow under various nonaxisymmetric electric fields. While investigating the interfacial flow, we find a surface vortex structure under certain nonaxisymmetric electric fields, which is found to be related to the chaotic mixing inside the drop. Finally, we consider the chaotic features of three-dimensional flows inside the drop under static nonaxisymmetric electric fields.     

Empirical force fields for complex hydrated calcio-silicate layered materials

The use of empirical force fields is now a standard approach in predicting the properties of hydrated oxides which are omnipresent in both natural and engineering applications. Transferability of force fields to analogous hydrated oxides without rigorous investigations may result in misleading property predictions. Herein, we focus on two common empirical force fields, the simple point charge ClayFF potential and the core-shell potential to study tobermorite minerals, the most prominent family of Calcium-Silicate-Hydrates that are complex hydrated oxides. We benchmark the predictive capabilities of these force fields against first principles results. While the structural information seem to be in close agreement with DFT results, we find that for higher order properties such as elastic constants, the core-shell potential quantitatively improves upon the simple point charge model, and shows a larger degree of transferability to complex materials. In return, to remedy the deficiencies of the simple point charge potential for hydrated calcio-silicates, we suggest using both structural data and elasticity data for potential calibration, a new force field potential, CSH-FF. This re-parameterized version of ClayFF is then applied to simulating an atomistic model of cement (Pellenq et al., PNAS, 2009). We demonstrate that this force field improves the predictive capabilities of ClayFF, being considerably less computational intensive than the core-shell model.