A three-point method for evaluations of AMBER force field parameters: an application to copper-based artificial nucleases

We present the theoretical evaluation of new AMBER force field parameters for 12 copper-based nucleases with bis(2-pyridylmethyl) amine, 2,2′-dipyridylamine, imidazole, N,N-bis(2-benzimidazolylmethyl) amine and their derivative ligands based on first-principles electronic structure calculations at the B3LYP level of theory. A three-point approach was developed to accurately and efficiently evaluate the force field parameters for the copper-based nucleases with the ligands. The protocol of RESP atomic charges has been used to calculate the atomic charge distributions of the studied copper-based nucleases. The evaluated force field parameters and RESP atomic charges have been successfully applied in the testing molecular mechanics calculations and molecular dynamics simulations on the nucleases and the nuclease–DNA complexes, respectively. It has been demonstrated that the developed force field parameters and atomic charges can consistently reproduce molecular geometries and conformations in the available X-ray crystal structures and can reasonably predict the interaction properties of the nucleases with DNA. The developed force field parameters in this work provide an extension of the AMBER force field for its application to computational modeling and simulations of the copper-based artificial nucleases associated with DNA.     

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.     

Lennard-Jones parameters for the combined QM/MM method using the B3LYP/6-31G*/AMBER potential

A combined DFT quantum mechanical and AMBER molecular mechanical potential (QM/MM) is presented for use in molecular modeling and molecular simulations of large biological systems. In our approach we evaluate Lennard-Jones parameters describing the interaction between the quantum mechanical (QM) part of a system, which is described at the B3LYP/6-31+G* level of theory, and the molecular mechanical (MM) part of the system, described by the AMBER force field. The Lennard-Jones parameters for this potential are obtained by calculating hydrogen bond energies and hydrogen bond geometries for a large set of bimolecular systems, in which one hydrogen bond monomer is described quantum mechanically and the other is treated molecular mechanically. We have investigated more than 100 different bimolecular systems, finding very good agreement between hydrogen bond energies and geometries obtained from the combined QM/MM calculations and results obtained at the QM level of theory, especially with respect to geometry. Therefore, based on the Lennard-Jones parameters obtained in our study, we anticipate that the B3LYP/6-31+G*/AMBER potential will be a precise tool to explore intermolecular interactions inside a protein environment.     

Free energy perturbation and molecular dynamics calculations of copper binding to azurin

Free energy perturbation/molecular dynamics simulations have been carried out on copper/azurin systems calculating the binding affinities of copper (II) ion to azurin either in the native or in the unfolded state. In order to test the validity of the strategy adopted for the calculations and to establish what force field is suitable for these kinds of calculations, three different force fields, AMBER, CVFF, and CFF, have been alternatively used for the calculations and the results have been compared with experimental data obtained by spectroscopic titrations of copper (II)/azurin solutions and denaturation experiments. Our findings have pointed out that only CFF gives satisfactory results, thus providing a reliable tool for copper binding simulations in copper protein.     

Elucidating the Origin of Conformational Energy Differences in Substituted 1,3-Dioxanes: A Combined Theoretical and Experimental Study

13C NMR spectroscopy, ab initio quantum mechanics, and molecular mechanics have been used to investigate the trans-4-(trifluoromethyl)-2,2,6-trimethyl-1,3-dioxane chair/twist-boat equilibrium. The molecular mechanics calculations were based upon the MM3 and AMBER force fields. A 6-31G basis set was used for the ab initio calculations, and MP2 correlation corrections were applied. Both the ab initio and AMBER molecular mechanics calculations are consistent with the (13)C NMR chemical shift differences for the trans-4-(trifluoromethyl)-2,2,6-trimethyl-1,3-dioxane conformers. The predicted chair to twist-boat equilibrium suggested by the MM3 calculations is not consistent with the experimental data. These results support the suggestion by Howard et al. (Howard, A. E.; Cieplak, P.; Kollman, P. A. J. Comput.Chem. 1995, 16, 243-261) on the critical role of electrostatic interactions in determining the chair/twist-boat equilibrium.     

Calculations of the free energy of interaction of the c-Fos-c-Jun coiled coil: effects of the solvation model and the inclusion of polarization effects

The leucine zipper region of activator protein-1 (AP-1) comprises the c-Jun and c-Fos proteins and constitutes a well-known coiled coil protein-protein interaction motif. We have used molecular dynamics (MD) simulations in conjunction with the molecular mechanics/Poisson-Boltzmann generalized-Born surface area [MM/PB(GB)SA] methods to predict the free energy of interaction of these proteins. In particular, the influence of the choice of solvation model, protein force field, and water potential on the stability and dynamic properties of the c-Fos-c-Jun complex were investigated. Use of the AMBER polarizable force field ff02 in combination with the polarizable POL3 water potential was found to result in increased stability of the c-Fos-c-Jun complex. MM/PB(GB)SA calculations revealed that MD simulations using the POL3 water potential give the lowest predicted free energies of interaction compared to other nonpolarizable water potentials. In addition, the calculated absolute free energy of binding was predicted to be closest to the experimental value using the MM/GBSA method with independent MD simulation trajectories using the POL3 water potential and the polarizable ff02 force field, while all other binding affinities were overestimated.     

New force field for molecular simulation of guanidinium-based ionic liquids

An all-atom force field was proposed for a new class of room temperature ionic liquids (RTILs), N,N,N',N'-tetramethylguanidinium (TMG) RTILs. The model is based on the AMBER force field with modifications on several parameters. The refinements include (1) fitting the vibration frequencies for obtaining force coefficients of bonds and angles against the data obtained by ab initio calculations and/or by experiments and (2) fitting the torsion energy profiles of dihedral angles for obtaining torsion parameters against the data obtained by ab initio calculations. To validate the force field, molecular dynamics (MD) simulations at different temperatures were performed for five kinds of RTILs, where TMG acts as a cation and formate, lactate, perchlorate, trifluoroacetate, and trifluoromethylsulfonate act as anions. The predicted densities were in good agreement with the experimental data. Radial distribution functions (RDFs) and spatial distribution functions (SDFs) were investigated to depict the microscopic structures of the RTILs.     

Force field development for cofactors in the photosystem II

We present a set of force field (FF) parameters compatible with the AMBER03 FF to describe five cofactors in photosystem II (PSII) of oxygenic photosynthetic organisms: plastoquinone‐9 (three redox forms), chlorophyll‐a, pheophytin‐a, heme‐b, and β‐carotene. The development of a reliable FF for these cofactors is an essential step for performing molecular dynamics simulations of PSII. Such simulations are important for the calculation of absorption spectrum and the further investigation of the electron and energy transfer processes. We have derived parameters for partial charges, bonds, angles, and dihedral‐angles from solid theoretical models using systematic quantum mechanics (QM) calculations. We have shown that the developed FF parameters are in good agreement with both ab initio QM and experimental structural data in small molecule crystals as well as protein complexes. © 2012 Wiley Periodicals, Inc.     

Mechanism of proteolysis in matrix metalloproteinase‐2 revealed by QM/MM modeling

The mechanism of enzymatic peptide hydrolysis in matrix metalloproteinase‐2 (MMP‐2) was studied at atomic resolution through quantum mechanics/molecular mechanics (QM/MM) simulations. An all‐atom three‐dimensional molecular model was constructed on the basis of a crystal structure from the Protein Data Bank (ID: 1QIB), and the oligopeptide Ace‐Gln‐Gly～Ile‐Ala‐Gly‐Nme was considered as the substrate. Two QM/MM software packages and several computational protocols were employed to calculate QM/MM energy profiles for a four‐step mechanism involving an initial nucleophilic attack followed by hydrogen bond rearrangement, proton transfer, and C? N bond cleavage. These QM/MM calculations consistently yield rather low overall barriers for the chemical steps, in the range of 5–10 kcal/mol, for diverse QM treatments (PBE0, B3LYP, and BB1K density functionals as well as local coupled cluster treatments) and two MM force fields (CHARMM and AMBER). It, thus, seems likely that product release is the rate‐limiting step in MMP‐2 catalysis. This is supported by an exploration of various release channels through QM/MM reaction path calculations and steered molecular dynamics simulations. © 2015 Wiley Periodicals, Inc.     

An evaluation of force-field treatments of aromatic interactions

Experimental measurements of edge-to-face aromatic interactions have been used to test a series of molecular mechanics force fields. The experimental data were determined for a range of differently substituted aromatic rings using chemical double mutant cycles on hydrogen-bonded zipper complexes. These complexes were truncated for the purposes of the molecular mechanics calculations so that problems of conformational searching and the optimisation of large structures could be avoided. Double-mutant cycles were then carried out in silico using these truncated systems. Comparison of the experimental aromatic interaction energies and the X-ray crystal structures of these truncated complexes with the calculated data show that conventional molecular mechanics force fields (MM2, MM3, AMBER and OPLS) do not perform well. However, the XED force field which explicitly represents electron anisotropy as an expansion of point charges around each atom reproduces the trends in interaction energy and the three-dimensional structures exceedingly well. Collapsing the XED charges onto atom centres or the use of semi-empirical atom-centred charges within the XED force field gives poor results. Thus the success of XED is not related to the methods used to assign the atomic charge distribution but can be directly attributed to the use of off-atom centre charges.     

Partial hessian fitting for determining force constant parameters in molecular mechanics

We present a new protocol for deriving force constant parameters that are used in molecular mechanics (MM) force fields to describe the bond‐stretching, angle‐bending, and dihedral terms. A 3 × 3 partial matrix is chosen from the MM Hessian matrix in Cartesian coordinates according to a simple rule and made as close as possible to the corresponding partial Hessian matrix computed using quantum mechanics (QM). This partial Hessian fitting (PHF) is done analytically and thus rapidly in a least‐squares sense, yielding force constant parameters as the output. We herein apply this approach to derive force constant parameters for the AMBER‐type energy expression. Test calculations on several different molecules show good performance of the PHF parameter sets in terms of how well they can reproduce QM‐calculated frequencies. When soft bonds are involved in the target molecule as in the case of secondary building units of metal‐organic frameworks, the MM‐optimized geometry sometimes deviates significantly from the QM‐optimized one. We show that this problem is rectified effectively by use of a simple procedure called Katachi that modifies the equilibrium bond distances and angles in bond‐stretching and angle‐bending terms. © 2016 Wiley Periodicals, Inc.     

Classical and quantum mechanical/molecular mechanical molecular dynamics simulations of alanine dipeptide in water: comparisons with IR and vibrational circular dichroism spectra

We have implemented the combined quantum mechanical (QM)/molecular mechanical (MM) molecular dynamics (MD) simulations of alanine dipeptide in water along with the polarizable and nonpolarizable classical MD simulations with different models of water. For the QM/MM MD simulation, the alanine dipeptide is treated with the AM1 or PM3 approximations and the fluctuating solute dipole moment is calculated by the Mulliken population analysis. For the classical MD simulations, the solute is treated with the polarizable or nonpolarizable AMBER and polarizable CHARMM force fields and water is treated with the TIP3P, TIP4P, or TIP5P model. It is found that the relative populations of right-handed alpha-helix and extended beta and P(II) conformations in the simulation trajectory strongly depend on the simulation method. For the QM/MM MD simulations, the PM3/MM shows that the P(II) conformation is dominant, whereas the AM1/MM predicts that the dominant conformation is alpha(R). Polarizable CHARMM force field gives almost exclusively P(II) conformation and other force fields predict that both alpha-helical and extended (beta and P(II)) conformations are populated with varying extents. Solvation environment around the dipeptide is investigated by examining the radial distribution functions and numbers and lifetimes of hydrogen bonds. Comparing the simulated IR and vibrational circular dichroism spectra with experimental results, we concluded that the dipeptide adopts the P(II) conformation and PM3/MM, AMBER03 with TIP4P water, and AMBER polarizable force fields are acceptable for structure determination of the dipeptide considered in this paper.     

Development of a molecular mechanics (MM2) force field for α-chlorosilanes

Molecular mechanics (MM2) parameters for silanes which have a Si-C-Cl fragment have been developed based on available experimental data and ab initio molecular orbital (MO) calculations. Molecular properties, mainly rotational barriers and geometries, of α-chlorosilanes have been studied using our new MM2 parameter set. Changes in the Si-C bond lengths and several bond angles of α-chlorosilanes due to the additional attachment of polar atom(s) have been investigated utilizing ab initio calculations. An electronegativity correction to both bond lengths and angles helps MM2 to reproduce results from ab initio calculations. The new force field has been applied to the conformational analysis of l-(chloromethyl)-1,2-dimethylsilacyclopentane, a model used in our studies of rearrangements of α-halosilanes.     

Development of polyphosphate parameters for use with the AMBER force field

Accurate force fields are essential for reproducing the conformational and dynamic behavior of condensed-phase systems. The popular AMBER force field has parameters for monophosphates, but they do not extend well to polyphorylated molecules such as ADP and ATP. This work presents parameters for the partial charges, atom types, bond angles, and torsions in simple polyphosphorylated compounds. The parameters are based on molecular orbital calculations of methyldiphosphate and methyltriphosphate at the RHF/6-31+G* level. The new parameters were fit to the entire potential energy surface (not just minima) with an RMSD of 0.62 kcal/mol. This is exceptional agreement and a significant improvement over the current parameters that produce a potential surface with an RMSD of 7.8 kcal/mol to that of the ab initio calculations. Testing has shown that the parameters are transferable and capable of reproducing the gas-phase conformations of inorganic diphosphate and triphosphate. Also, the parameters are an improvement over existing parameters in the condensed phase as shown by minimizations of ATP bound in several proteins. These parameters are intended for use with the existing AMBER 94/99 force field, and they will permit users to apply AMBER to a wider variety of important enzymatic systems.     

On the effect of a variation of the force field, spatial boundary condition and size of the QM region in QM/MM MD simulations

During the past years, the use of combined quantum-classical, QM/MM, methods for the study of complex biomolecular processes, such as enzymatic reactions and photocycles, has increased considerably. The quality of the results obtained from QM/MM calculations is largely dependent on five aspects to be considered when setting up a molecular model: the QM Hamiltonian, the MM Hamiltonian or force field, the boundary and coupling between the QM and MM regions, the size of the QM region and the boundary condition for the MM region. In this study, we systematically investigate the influence of a variation of the molecular mechanics force field and the size of the QM region in QM/MM MD simulations on properties of the photoactive part of the blue light photoreceptor protein AppA. For comparison, we additionally performed classical MD simulations and studied the effect of a variation of the type of spatial boundary condition. The classical boundary conditions and the force field used in a QM/MM MD simulation are shown to have non-neglegible effects upon the structural and energetic properties of the protein which makes it advisable to minimize computational artifacts in QM/MM MD simulations by application of periodic boundary conditions and a thermodynamically calibrated force field. A comparison of the structural and energetic properties of MD simulations starting from two alternative, different X-ray structures for the blue light utilizing flavin protein in its dark state indicates a slight preference of the two force fields used for the so-called Anderson structure over the Jung structure.     

Simulation of mesogenic diruthenium tetracarboxylates: Development of a force field for coordination polymers of the MMX type

A classical molecular mechanics force field, able to simulate coordination polymers (CP) based on ruthenium carboxylates (Ru₂(O₂CR_eq)₄L_ax) (eq = equatorial group containing aliphatic chains, L_ax= axial ligand), has been developed. New parameters extracted from experimental data and quantum calculations on short aliphatic chains model systems were included in the generalized AMBER force field. The proposed parametrization was evaluated using model systems with known structure, containing either short or long aliphatic chains; experimental results were reproduced satisfactorily. This modified force field, although in a preliminary stage, could then be applied to long chain liquid crystalline compounds. The resulting atomistic simulations allowed assessing the relative influence of the factors determining the CP conformation, determinant for the physical properties of these materials. © 2013 Wiley Periodicals, Inc.     

Force‐field parameters of the Ψ and Φ around glycosidic bonds to oxygen and sulfur atoms

The Ψ and Φ torsion angles around glycosidic bonds in a glycoside chain are the most important determinants of the conformation of a glycoside chain. We determined force‐field parameters for Ψ and Φ torsion angles around a glycosidic bond bridged by a sulfur atom, as well as a bond bridged by an oxygen atom as a preparation for the next study, i.e., molecular dynamics free energy calculations for protein‐sugar and protein‐inhibitor complexes. First, we extracted the Ψ or Φ torsion energy component from a quantum mechanics (QM) total energy by subtracting all the molecular mechanics (MM) force‐field components except for the Ψ or Φ torsion angle. The Ψ and Φ energy components extracted (hereafter called “the remaining energy components”) were calculated for simple sugar models and plotted as functions of the Ψ and Φ angles. The remaining energy component curves of Ψ and Φ were well represented by the torsion force‐field functions consisting of four and three cosine functions, respectively. To confirm the reliability of the force‐field parameters and to confirm its compatibility with other force‐fields, we calculated adiabatic potential curves as functions of Ψ and Φ for the model glycosides by adopting the Ψ and Φ force‐field parameters obtained and by energetically optimizing other degrees of freedom. The MM potential energy curves obtained for Ψ and Φ well represented the QM adiabatic curves and also these curves' differences with regard to the glycosidic oxygen and sulfur atoms. Our Ψ and Φ force‐fields of glycosidic oxygen gave MM potential energy curves that more closely represented the respective QM curves than did those of the recently developed GLYCAM force‐field. © 2009 Wiley Periodicals, Inc., J Comput Chem, 2009