Comparison of implicit solvent models for the simulation of protein-surface interactions

Empirical force field-based molecular simulations can provide valuable atomistic-level insights into protein-surface interactions in aqueous solution. While the implicit treatment of solvation effects is desired as a means of improving simulation efficiency, existing implicit solvent models were primarily developed for the simulation of peptide or protein behavior in solution alone, and thus may not be appropriate for protein interactions with synthetic material surfaces. The objective of this research was to calculate the change in free energy as a function of surface-separation distance for peptide-surface interactions using different empirical force field-based implicit solvation models (ACE, ASP, EEF1, and RDIE with the CHARMM 19 force field), and to compare these results with the same calculations conducted using density functional theory (DFT) combined with the self-consistent reaction field (SCRF) implicit solvation model. These comparisons show that distinctly different types of behavior are predicted with each implicit solvation method, with ACE providing the best overall agreement with DFT/SCRF calculations. These results also identify areas where ACE is in need of improvement for this application and provide a basis for subsequent parameter refinement.     

Microcalorimetric studies of the mechanism of interaction between designed peptides and hydrophobic adsorbents

To understand the mechanism of interaction between peptides and peptides with hydrophobic ligands, the oligomers (GWG, GWWG, GWWWG) were designed and synthesized to study adsorption behavior with octyl sepharose and CM-octyl sepharose. By batch equilibrium binding analysis and dilution heat of peptide solution measurement, the binding isotherm and adsorption enthalpy were obtained and the binding thermodynamics parameters were calculated and analyzed. In the isotherm analysis, we reveled that the affinity of GWG for both adsorbents is stronger than that of GWWG and GWWWG. The results demonstrate that the cation-pi interaction between the peptides and the buffer molecules is significant for solutions of peptides with tryptophan residues, and the solvation is competitive with the hydrophobic interaction between the peptides and the hydrophobic ligands. From the dilution heat measurements, we observed an endothermic dilution heat for GWG and exothermic for GWWG and GWWWG. All these results indicate that the increased tryptophan chain length can promote the solvation behavior of the peptides by the peptide-buffer interaction in this buffer system. Comparing the types of ligands reveals that the binding affinities of each peptide for the two adsorbents are similar. However, the mechanism of adsorption for peptides with hydrophobic ligands might be quite different with respect to the binding enthalpy between peptides and adsorbents. The adsorption of the peptides on octyl sepharose is an entropy-driven process for all the peptides. In contrast, the adsorption of CM-octyl sepharose with GWG and GWWG is an enthalpy-driven process, whereas that with GWWWG is entropy-driven. These findings indicate that the amount of tryptophan controls the characteristics of the peptides and the interaction mechanism in the binding procedure. This study of the adsorption mechanism of the designed peptide could provide fundamental information for peptide purification and amino acid residue behavior in peptide drug design.     

Molecular dynamics study of peptides in implicit water: ab initio folding of beta-hairpin, beta-sheet, and beta beta alpha-motif

In this communication, we have demonstrated that molecular dynamics simulations using a GB implicit solvation model with the all-atom based force field (CHARMM19) can describe the spontaneous folding of small peptides in aqueous solution. The native structures of peptides with various structural motifs (beta-hairpin, beta-sheet, and betabetaalpha-moiety) were successfully predicted within reasonable time scales by MD simulations at moderately elevated temperatures. It is expected that the present simulations provide further insight into mechanism/pathways of the peptide folding.     

Benchmarking implicit solvent folding simulations of the amyloid beta(10-35) fragment

A pathogenetic feature of Alzhemier disease is the aggregation of monomeric beta-amyloid proteins (Abeta) to form oligomers. Usually these oligomers of long peptides aggregate on time scales of microseconds or longer, making computational studies using atomistic molecular dynamics models prohibitively expensive and making it essential to develop computational models that are cheaper and at the same time faithful to physical features of the process. We benchmark the ability of our implicit solvent model to describe equilibrium and dynamic properties of monomeric Abeta(10-35) using all-atom Langevin dynamics (LD) simulations, since Alphabeta(10-35) is the only fragment whose monomeric properties have been measured. The accuracy of the implicit solvent model is tested by comparing its predictions with experiment and with those from a new explicit water MD simulation, (performed using CHARMM and the TIP3P water model) which is approximately 200 times slower than the implicit water simulations. The dependence on force field is investigated by running multiple trajectories for Alphabeta(10-35) using the CHARMM, OPLS-aal, and GS-AMBER94 force fields, whereas the convergence to equilibrium is tested for each force field by beginning separate trajectories from the native NMR structure, a completely stretched structure, and from unfolded initial structures. The NMR order parameter, S2, is computed for each trajectory and is compared with experimental data to assess the best choice for treating aggregates of Alphabeta. The computed order parameters vary significantly with force field. Explicit and implicit solvent simulations using the CHARMM force fields display excellent agreement with each other and once again support the accuracy of the implicit solvent model. Alphabeta(10-35) exhibits great flexibility, consistent with experiment data for the monomer in solution, while maintaining a general strand-loop-strand motif with a solvent-exposed hydrophobic patch that is believed to be important for aggregation. Finally, equilibration of the peptide structure requires an implicit solvent LD simulation as long as 30 ns.     

Secondary structure bias in generalized Born solvent models: comparison of conformational ensembles and free energy of solvent polarization from explicit and implicit solvation

The effects of the use of three generalized Born (GB) implicit solvent models on the thermodynamics of a simple polyalanine peptide are studied via comparing several hundred nanoseconds of well-converged replica exchange molecular dynamics (REMD) simulations using explicit TIP3P solvent to REMD simulations with the GB solvent models. It is found that when compared to REMD simulations using TIP3P the GB REMD simulations contain significant differences in secondary structure populations, most notably an overabundance of alpha-helical secondary structure. This discrepancy is explored via comparison of the differences in the electrostatic component of the free energy of solvation (DeltaDeltaG(pol)) between TIP3P (via thermodynamic Integration calculations), the GB models, and an implicit solvent model based on the Poisson equation (PE). The electrostatic components of the solvation free energies are calculated using each solvent model for four representative conformations of Ala10. Since the PE model is found to have the best performance with respect to reproducing TIP3P DeltaDeltaG(pol) values, effective Born radii from the GB models are compared to effective Born radii calculated with PE (so-called perfect radii), and significant and numerous deviations in GB radii from perfect radii are found in all GB models. The effect of these deviations on the solvation free energy is discussed, and it is shown that even when perfect radii are used the agreement of GB with TIP3P DeltaDeltaG(pol) values does not improve. This suggests a limit to the optimization of the effective Born radius calculation and that future efforts to improve the accuracy of GB models must extend beyond such optimizations.     

Calculation of the water-cyclohexane transfer free energies of neutral amino acid side-chain analogs using the OPLS all-atom force field

We calculated the free energy of solvation of the neutral analogs of 18 amino acid side-chains (not including glycine and proline) using the OPLS all-atom force field in TIP4P water, SPC water, and cyclohexane by molecular dynamics simulation and thermodynamic integration. The average unsigned errors in the free energies of solvation in TIP4P, SPC, and cyclohexane are 4.4, 4.9, and 2.1 kJ/mol respectively. Most of the calculated hydration free energies are not favorable enough compared to experiment. The largest errors are found for tryptophan, histidine, glutamic acid, and glutamine. The average unsigned errors in the free energy of transfer from TIP4P to cyclohexane and from SPC to cyclohexane are 4.0 and 4.1 kJ/mol, respectively. The largest errors, of more than 7.5 kJ/mol, are found for histidine, glutamine, and glutamatic acid.     

New-generation amber united-atom force field

We have developed a new-generation Amber united-atom force field for simulations involving highly demanding conformational sampling such as protein folding and protein-protein binding. In the new united-atom force field, all hydrogens on aliphatic carbons in all amino acids are united with carbons except those on Calpha. Our choice of explicit representation of all protein backbone atoms aims at minimizing perturbation to protein backbone conformational distributions and to simplify development of backbone torsion terms. Tests with dipeptides and solvated proteins show that our goal is achieved quite successfully. The new united-atom force field uses the same new RESP charging scheme based on B3LYP/cc-pVTZ//HF/6-31g** quantum mechanical calculations in the PCM continuum solvent as that in the Duan et al. force field. van der Waals parameters are empirically refitted starting from published values with respect to experimental solvation free energies of amino acid side-chain analogues. The suitability of mixing new point charges and van der Waals parameters with existing Amber covalent terms is tested on alanine dipeptide and is found to be reasonable. Parameters for all new torsion terms are refitted based on the new point charges and the van der Waals parameters. Molecular dynamics simulations of three small globular proteins in the explicit TIP3P solvent are performed to test the overall stability and accuracy of the new united-atom force field. Good agreements between the united-atom force field and the Duan et al. all-atom force field for both backbone and side-chain conformations are observed. In addition, the per-step efficiency of the new united-atom force field is demonstrated for simulations in the implicit generalized Born solvent. A speedup around two is observed over the Duan et al. all-atom force field for the three tested small proteins. Finally, the efficiency gain of the new united-atom force field in conformational sampling is further demonstrated with a well-known toy protein folding system, an 18 residue polyalanine in distance-dependent dielectric. The new united-atom force field is at least a factor of 200 more efficient than the Duan et al. all-atom force field for ab initio folding of the tested peptide.     

Surveying implicit solvent models for estimating small molecule absolute hydration free energies

Implicit solvent models are powerful tools in accounting for the aqueous environment at a fraction of the computational expense of explicit solvent representations. Here, we compare the ability of common implicit solvent models (TC, OBC, OBC2, GBMV, GBMV2, GBSW, GBSW/MS, GBSW/MS2 and FACTS) to reproduce experimental absolute hydration free energies for a series of 499 small neutral molecules that are modeled using AMBER/GAFF parameters and AM1-BCC charges. Given optimized surface tension coefficients for scaling the surface area term in the nonpolar contribution, most implicit solvent models demonstrate reasonable agreement with extensive explicit solvent simulations (average difference 1.0-1.7 kcal/mol and R(2)=0.81-0.91) and with experimental hydration free energies (average unsigned errors=1.1-1.4 kcal/mol and R(2)=0.66-0.81). Chemical classes of compounds are identified that need further optimization of their ligand force field parameters and others that require improvement in the physical parameters of the implicit solvent models themselves. More sophisticated nonpolar models are also likely necessary to more effectively represent the underlying physics of solvation and take the quality of hydration free energies estimated from implicit solvent models to the next level.     

Role of solvent in determining conformational preferences of alanine dipeptide in water

Evidence from a variety of spectroscopic probes indicates that (phi, psi) values corresponding to the left-handed polyproline II helix (P(II)) are preferred for short alanine-based peptides in water. On the basis of results from theoretical studies, it is believed that the observed preference is dictated by favorable peptide-solvent interactions, which are realized through formation of optimal hydrogen-bonding water bridges between peptide donor and acceptor groups. In the present study, we address this issue explicitly by analyzing the hydration structure and thermodynamics of 16 low-energy conformers of the alanine dipeptide (N-acetylalanine-N'-methylamide) in liquid water. Monte Carlo simulations in the canonical ensemble were performed under ambient conditions with all-atom OPLS parameters for the alanine dipeptide and the TIP5P model for water. We find that the number of hydrogen-bonded water molecules connecting the peptide group donor and acceptor atoms has no effect on the solvation thermodynamics. Instead, the latter are determined by the work done to fully hydrate the peptide. This work is minimal for conformations that are characterized by a minimal overlap of the primary hydration shells around the peptide donor and acceptor atoms. As a result, peptide-solvent interactions favor "compact" conformations that do not include P(II)-like geometries. Our main conclusion is that the experimentally observed preference for P(II) does not arise due to favorable direct interactions between the peptide and water molecules. Instead, the latter act to unmask underlying conformational preferences that are a consequence of minimizing intrapeptide steric conflicts.     

Molecular dynamics simulations of peptide-surface interactions

Proteins, which are bioactive molecules, adsorb on implants placed in the body through complex and poorly understood mechanisms and directly influence biocompatibility. Molecular dynamics modeling using empirical force fields provides one of the most direct methods of theoretically analyzing the behavior of complex molecular systems and is well-suited for the simulation of protein adsorption behavior. To accurately simulate protein adsorption behavior, a force field must correctly represent the thermodynamic driving forces that govern peptide residue-surface interactions. However, since existing force fields were developed without specific consideration of protein-surface interactions, they may not accurately represent this type of molecular behavior. To address this concern, we developed a host-guest peptide adsorption model in the form of a G(4)-X-G(4) peptide (G is glycine, X is a variable residue) to enable determination of the contributions to adsorption free energy of different X residues when adsorbed to functionalized Au-alkanethiol self-assembled monolayers (SAMs). We have previously reported experimental results using surface plasmon resonance (SPR) spectroscopy to measure the free energy of peptide adsorption for this peptide model with X = G and K (lysine) on OH and COOH functionalized SAMs. The objectives of the present research were the development and assessment of methods to calculate adsorption free energy using molecular dynamics simulations with the GROMACS force field for these same peptide adsorption systems, with an oligoethylene oxide (OEG) functionalized SAM surface also being considered. By comparing simulation results to the experimental results, the accuracy of the selected force field to represent the behavior of these molecular systems can be evaluated. From our simulations, the G(4)-G-G(4) and G(4)-K-G(4) peptides showed minimal to no adsorption to the OH SAM surfaces and the G(4)-K-G(4) showed strong adsorption to the COOH SAM surface, which is in agreement with our SPR experiments. Contrary to our experimental results, however, the simulations predicted a relatively strong adsorption of G(4)-G-G(4) peptide to the COOH SAM surface. In addition, both peptides were unexpectedly predicted to adsorb to the OEG surface. These findings demonstrate the need for GROMACS force field parameters to be rebalanced for the simulation of peptide adsorption behavior on SAM surfaces. The developed methods provide a direct means of assessing, modifying, and validating force field performance for the simulation of peptide and protein adsorption to surfaces, without which little confidence can be placed in the simulation results that are generated with these types of systems.     

NMR solvent shifts of adenine in aqueous solution from hybrid QM/MM molecular dynamics simulations

We present first principles calculations of the NMR solvent shift of adenine in aqueous solution. The calculations are based on snapshots sampled from a molecular dynamics simulation, which were obtained via a hybrid quantum-mechanical/mechanical modeling approach, using an all-atom force field (TIP3P). We find that the solvation via the strongly fluctuating hydrogen bond network of water leads to nontrivial changes in the NMR spectra of the solutes regarding the ordering of the resonance lines. Although there are still sizable deviations from experiment, the overall agreement is satisfactory for the 1H and 15N NMR shifts. Our work is another step toward a realistic first-principles prediction of NMR chemical shifts in complex chemical environments.     

Adsorption of antimicrobial indolicidin-derived peptides on hydrophobic surfaces

The hydrophobic interaction between antimicrobial peptides and membrane hydrophobic cores is usually related to their cytotoxicity. In this study, the adsorption mechanism of five plasma membrane-associated peptides, indolicidin (IL) and its four derivatives, with hydrophobic ligands was investigated to understand the relationship between peptide hydrophobicity and bioactivity. The hydrophobic adsorption mechanisms of IL and its derivatives were interpreted thermodynamically and kinetically by reversed-phase chromatography (RPC) analysis and surface plasmon resonance (SPR) measurement, respectively. IL and its derivatives possess a similar random coil structure in both aqueous and organic solvents. Thermodynamic analysis showed that the binding enthalpy of peptides with higher electropositivity was lower than those with lower electropositivity and exhibited unfavorable binding entropy. Higher electropositivity peptides adsorbed to the hydrophobic surface arising from the less bound solvent on the peptide surface. A comparison with the kinetic analysis showed that IL and its derivatives adopt a two-state binding model (i.e., adsorption onto and self-association on the hydrophobic acyl chain) to associate with the hydrophobic surface, and the binding affinity of peptide self-association correlates well with peptide hemolysis. Consequently, this study provided a novel concept for understanding the action of plasma membrane-associated peptides.     

Attractive surface force in the presence of dissolved gas: a molecular approach

Despite widespread evidence of the influence of dissolved air on hydrophobic interaction, the mechanisms of observed effects are still unknown. Although some experiments indicate that adsorbed gases can modify the structure of water next to hydrophobic surfaces, gas effects on measured forces have been observed only at large surface separations. Gas-specific depletion of water at a hydrophobic surface has been detected but was not reproduced in subsequent measurements. We use computer simulations to study short-ranged hydrophobic attraction in the absence and presence of dissolved gas and monitor gas adsorption at molecular resolution inaccessible in experiments. Although we observe a significant accumulation of dissolved gases at hydrophobic surfaces, even in supersaturated gas solutions surface concentrations remain too low to induce any significant change in the local structure of water and short-range surface forces. We present direct calculations of the hydrophobic force between model hydrocarbon plates at separations between 1.5 and 4 nm. Although stronger, the calculated solvation force has a similar decay rate as deduced from recent surface force apparatus measurements at a somewhat lower contact angle. Within the statistical uncertainty, short-range attraction is not affected by the presence of dissolved nitrogen, even in supersaturated solution with a gas fugacity as high as 30 atm. Comparisons of the adsorption behavior of N2, O2, CO2, and Ar reveal similar features in contrast to the peculiar suppression of water depletion reported for an Ar solution in a neutron reflectivity experiment. Our calculations reveal a notable difference between pathways to the capillary evaporation of pure water and gas-phase nucleation in confined supersaturated gas solutions.     

Assessing the quality of absolute hydration free energies among CHARMM‐compatible ligand parameterization schemes

Multipurpose atom‐typer for CHARMM (MATCH), an atom‐typing toolset for molecular mechanics force fields, was recently developed in our laboratory. Here, we assess the ability of MATCH‐generated parameters and partial atomic charges to reproduce experimental absolute hydration free energies for a series of 457 small neutral molecules in GBMV2, Generalized Born with a smooth SWitching (GBSW), and fast analytical continuum treatment of solvation (FACTS) implicit solvent models. The quality of hydration free energies associated with small molecule parameters obtained from ParamChem, SwissParam, and Antechamber are compared. Given optimized surface tension coefficients for scaling the surface area term in the nonpolar contribution, these automated parameterization schemes with GBMV2 and GBSW demonstrate reasonable agreement with experimental hydration free energies (average unsigned errors of 0.9–1.5 kcal/mol and R² of 0.63–0.87). GBMV2 and GBSW consistently provide slightly more accurate estimates than FACTS, whereas Antechamber parameters yield marginally more accurate estimates than the current generation of MATCH, ParamChem, and SwissParam parameterization strategies. Modeling with MATCH libraries that are derived from different CHARMM topology and parameter files highlights the importance of having sufficient coverage of chemical space within the underlying databases of these automated schemes and the benefit of targeting specific functional groups for parameterization efforts to maximize both the breadth and the depth of the parameterized space. © 2013 Wiley Periodicals, Inc.     

Critical importance of length-scale dependence in implicit modeling of hydrophobic interactions

The existence of length-scale dependence of hydrophobic solvation has important implications in the equilibrium of disordered, partially folded, and folded protein conformations. Neglecting this dependence, such as in popular solute surface-area based implicit solvent models with fixed surface tension coefficients, severely limits the ability to accurately model protein conformational equilibrium. We illustrate such fundamental limitations by examining the potentials of mean force of forming dimeric and trimeric nonpolar clusters and propose a new empirical model that effectively captures the context dependence of the local effective surface tension. Further optimization of the new model with other components of the implicit solvent force fields provides promise to significantly improve one's ability to simulate protein folding and conformational transitions. The existence of length-scale dependence of hydrophobic solvation has important implications in the equilibrium of disordered, partially folded, and folded protein conformations. Neglecting this dependence, such as in popular solute surface-area based implicit solvent models with fixed surface tension coefficients, severely limits the ability to accurately model protein conformational equilibrium. We illustrate such fundamental limitations by examining the potentials of mean force of forming dimeric and trimeric nonpolar clusters and propose a new empirical model that effectively captures the context dependence of the local effective surface tension. Further optimization of the new model with other components of the implicit solvent force fields provides promise to significantly improve one's ability to simulate protein folding and conformational transitions.