Molecular mechanics is powerful for its speed in atomistic simulations, but an accurate force field is required. The Amber ff99SB force field improved protein secondary structure balance and dynamics from earlier force fields like ff99, but weaknesses in side chain rotamer and backbone secondary structure preferences have been identified. Here, we performed a complete refit of all amino acid side chain dihedral parameters, which had been carried over from ff94. The training set of conformations included multidimensional dihedral scans designed to improve transferability of the parameters. Improvement in all amino acids was obtained as compared to ff99SB. Parameters were also generated for alternate protonation states of ionizable side chains. Average errors in relative energies of pairs of conformations were under 1.0 kcal/mol as compared to QM, reduced 35% from ff99SB. We also took the opportunity to make empirical adjustments to the protein backbone dihedral parameters as compared to ff99SB. Multiple small adjustments of φ and ψ parameters were tested against NMR scalar coupling data and secondary structure content for short peptides. The best results were obtained from a physically motivated adjustment to the φ rotational profile that compensates for lack of ff99SB QM training data in the β-ppII transition region. Together, these backbone and side chain modifications (hereafter called ff14SB) not only better reproduced their benchmarks, but improved secondary structure content in small peptides, and reproduction of NMR χ1 scalar coupling measurements for proteins in solution. We also discuss the Amber ff12SB parameter set, a preliminary version of ff14SB that includes most of its improvements.

Molecular dynamics (MD) simulations have become increasingly popular in studying the motions and functions of biomolecules. The accuracy of the simulation, however, is highly determined by the molecular mechanics (MM) force field (FF), a set of functions with adjustable parameters to compute the potential energies from atomic positions. However, the overall quality of the FF, such as our previously published ff99SB and ff14SB, can be limited by assumptions that were made years ago. In the updated model presented here (ff19SB), we have significantly improved the backbone profiles for all 20 amino acids. We fit coupled φ/ψ parameters using 2D φ/ψ conformational scans for multiple amino acids, using as reference data the entire 2D quantum mechanics (QM) energy surface. We address the polarization inconsistency during dihedral parameter fitting by using both QM and MM in aqueous solution. Finally, we examine possible dependency of the backbone fitting on side chain rotamer. To extensively validate ff19SB parameters, and to compare to results using other Amber models, we have performed a total of ∼5 ms MD simulations in explicit solvent. Our results show that after amino-acid-specific training against QM data with solvent polarization, ff19SB not only reproduces the differences in amino-acid-specific Protein Data Bank (PDB) Ramachandran maps better but also shows significantly improved capability to differentiate amino-acid-dependent properties such as helical propensities. We also conclude that an inherent underestimation of helicity is present in ff14SB, which is (inexactly) compensated for by an increase in helical content driven by the TIP3P bias toward overly compact structures. In summary, ff19SB, when combined with a more accurate water model such as OPC, should have better predictive power for modeling sequence-specific behavior, protein mutations, and also rational protein design. Of the explicit water models tested here, we recommend use of OPC with ff19SB.

<p>Molecular dynamics (MD) simulations have become increasingly popular in studying the motions and functions of biomolecules. The accuracy of the simulation, however, is highly determined by the molecular mechanics (MM) force field (FF), a set of functions with adjustable parameters to compute the potential energies from atomic positions. However, the overall quality of the FF, such as our previously published ff99SB and ff14SB, can be limited by assumptions that were made years ago. In the updated model presented here (ff19SB), we have significantly improved the backbone profiles for all 20 amino acids. We fit coupled ϕ/ψ parameters using 2D ϕ/ψ conformational scans for multiple amino acids, using as reference data the entire 2D quantum mechanics (QM) energy surface. We address the polarization inconsistency during dihedral parameter fitting by using both QM and MM in solution. Finally, we examine possible dependency of the backbone fitting on side chain rotamer. To extensively validate ff19SB parameters, we have performed a total of ~5 milliseconds MD simulations in explicit solvent. Our results show that after amino-acid specific training against QM data with solvent polarization, ff19SB not only reproduces the differences in amino acid specific Protein Data Bank (PDB) Ramachandran maps better, but also shows significantly improved capability to differentiate amino acid dependent properties such as helical propensities. We also conclude that an inherent underestimation of helicity is present in ff14SB, which is (inexactly) compensated by an increase in helical content driven by the TIP3P bias toward overly compact structures. In summary, ff19SB, when combined with a more accurate water model such as OPC, should have better predictive power for modeling sequence-specific behavior, protein mutations, and also rational protein design. </p>

A delicate balance of different types of intramolecular interactions makes the folded states of proteins marginally more stable than the unfolded states. Experiments use thermal, chemical, or mechanical stress to perturb the folding equilibrium for examining protein stability and the protein folding process. Elucidation of the mechanism by which chemical denaturants unfold proteins is crucial; this study explores the nature of urea−aromatic interactions relevant in urea-assisted protein denaturation. Free energy profiles corresponding to the unfolding of Trp-cage miniprotein in the presence and absence of urea at three different temperatures demonstrate the distortion of the hydrophobic core to be a crucial step. Exposure of the Trp6 residue to the solvent is found to be favored in the presence of urea. Previous experiments showed that urea has a high affinity for aromatic groups of proteins. We show here that this is due to the remarkable ability of urea to form stacking and NH−π interactions with aromatic groups of proteins. Urea−nucleobase stacking interactions have been shown to be crucial in urea-assisted RNA unfolding. Examination of these interactions using microsecond-long unrestrained simulations shows that urea−aromatic stacking interactions are stabilizing and long lasting. Further MD simulations, thermodynamic integration, and quantum mechanical calculations on aromatic model systems reveal that such interactions are possible for all the aromatic amino acid side-chains. Finally, we validate the ubiquitous nature of urea−aromatic stacking interactions by analyzing experimental structures of urea transporters and proteins crystallized in the presence of urea or urea derivatives.

Temperature replica exchange molecular dynamics (REMD) is a widely used enhanced sampling method for accelerating biomolecular simulations. During the past 2 decades, several variants of REMD have been developed to further improve the rate of conformational sampling of REMD. One such variant, reservoir REMD (RREMD), was shown to improve the rate of conformational sampling by around 5−20×. Despite the significant increase in the sampling speed, RREMD methods have not been widely used because of the difficulties in building the reservoir and also because of the code not being available on the graphics processing units (GPUs). In this work, we ported the Amber RREMD code onto GPUs making it 20× faster than the central processing unit code. Then, we explored protocols for building Boltzmann-weighted reservoirs as well as non-Boltzmann reservoirs and tested how each choice affects the accuracy of the resulting RREMD simulations. We show that, using the recommended protocols outlined here, RREMD simulations can accurately reproduce Boltzmann-weighted ensembles obtained by much more expensive conventional temperature-based REMD simulations, with at least 15× faster convergence rates even for larger proteins (>50 amino acids) compared to conventional REMD.

The recently developed AlphaFold2 (AF2) algorithm predicts proteins’ 3D structures from amino acid sequences. The open AlphaFold protein structure database covers the complete human proteome. Using an industry-leading molecular docking method (Glide), we investigated the virtual screening performance of 37 common drug targets, each with an AF2 structure and known holo and apo structures from the DUD-E data set. In a subset of 27 targets where the AF2 structures are suitable for refinement, the AF2 structures show comparable early enrichment of known active compounds (avg. EF 1%: 13.0) to apo structures (avg. EF 1%: 11.4) while falling behind early enrichment of the holo structures (avg. EF 1%: 24.2). With an induced-fit protocol (IFD-MD), we can refine the AF2 structures using an aligned known binding ligand as the template to improve the performance in structure-based virtual screening (avg. EF 1%: 18.9). Glide-generated docking poses of known binding ligands can also be used as templates for IFD-MD, achieving similar improvements (avg. EF 1% 18.0). Thus, with proper preparation and refinement, AF2 structures show considerable promise for in silico hit identification.

The recently developed AlphaFold2 (AF2) algorithm predicts proteins’ 3D structures from amino acid sequences. The open AlphaFold Protein Structure Database covers the complete human proteome. It shows great potential to provide structural information to enable and enhance existing and new drug discovery projects. Using an industry-leading molecular docking method (Glide), we benchmarked the virtual screening performance of 28 common drug targets each with an AF2 structure and known holo and apo structures from the DUD-E dataset. The AF2 structures show comparable early enrichment of known active compounds (avg. EF 1%: 13.16) to apo structures (avg. EF 1%: 11.56), while falling behind early enrichment of the holo structures (avg. EF 1%: 24.81). We also demonstrated that with the IFD-MD induced-fit docking approach, we can refine the AF2 structures using a known binding ligand to improve the performance in structure-based virtual screening (avg. EF 1%: 19.25). Thus, with proper preparation and refinement, AF2 structures show considerable promise for in silico hit identification.

Urea has long been used to investigate protein folding and, more recently, RNA folding. Studies have proposed that urea denatures RNA by participating in stacking interactions and hydrogen bonds with nucleic acid bases. In this study, the ability of urea to form unconventional stacking interactions with RNA bases is investigated using ab initio calculations (RI-MP2 and CCSD(T) methods with the aug-cc-pVDZ basis set). A total of 29 stable nucleobase-urea stacked complexes are identified in which the intermolecular interaction energies (up to −14 kcal/mol) are dominated by dispersion effects. Natural bond orbital (NBO) and atoms in molecules (AIM) calculations further confirm strong interactions between urea and nucleobases. Calculations on model systems with multiple urea and water molecules interacting with a guanine base lead to a hypothesis that urea molecules along with water are able to form cage-like structures capable of trapping nucleic acid bases in extrahelical states by forming both hydrogen bonded and dispersion interactions, thereby contributing to the unfolding of RNA in the presence of urea in aqueous solution.