TIP3P Water Model Research Articles

Interfacial Water Structure at Zwitterionic Membrane/Water Interface: The Importance of Interactions between Water and Lipid Carbonyl Groups.

In this work, atomistic molecular dynamics (MD) simulations of palmitoyl-oleoyl-phosphatidylcholine (POPC) bilayer were carried out to investigate the effect of water models on membrane dipole potential, which is primarily associated with the preferential orientation of molecular dipoles at the membrane–water interface. We discovered that the overestimation of the dipole potential by the TIPS3P water model can be effectively reduced by the TIP4P water model. On the one hand, the TIP4P water model decreases the negative contribution of lipid to the dipole potential through influencing the orientation of lipid headgroups. On the other hand, the TIP4P water model reduces the positive contribution of water to the dipole potential by increasing the preference of H-down orientation (the water dipole orients toward the bilayer center). Interestingly, the TIP4P water model affects the orientation of interfacial water molecules more obviously than that of lipid headgroups, leading to the decrease in the dipole potential. Furthermore, the MD results revealed that the water close to the positively charged choline (namely, N-associated water) prefers the H-down orientation while the water around the negatively charged phosphate (namely, P-associated water) favors the H-up orientation, in support of recent experimental and MD studies. However, interfacial water molecules are more strongly influenced by the phosphate groups than by the choline groups, resulting in the net H-up orientation (the water dipole orients toward the bilayer center) in the region of lipid headgroups. In addition, it is intriguing that the preference of H-up orientation decreases when water molecules penetrate more deeply into the lipid bilayer. This is attributed to the counteracting effect of lipid carbonyl groups, and the effect varies with the lipid chains (oleoyl and palmitoyl chains), suggesting the important role of lipid carbonyl groups.

Performance of CHARMM36m with modified water model in simulating intrinsically disordered proteins: a case study

Molecular dynamics simulations can be a powerful tool to complement experiments in the study of the structures and dynamics of intrinsically disordered proteins. Though the accuracy of the physics-based all-atom force fields has improved significantly in simulating structured proteins over the past twenty years, most of these force fields face a big challenge to simulate flexible proteins. Recently, CHARMM36m with modified TIP3P model was proposed as a possible solution to simulate intrinsically disordered proteins. Here, we tested the proposed solution using an extensively studied protein, namely NCBD, to explore the performance of CHARMM36m plus modified TIP3P water. Our results suggest that the modified TIP3P water model does enhance the sampling of conformational space compared to the standard TIP3P water model. However, the new CHARMM36m force field still leads to over-compact structures and over-stabilized helices.

The accuracy of force fields on the simulation of intrinsically disordered proteins: A benchmark test on the human p53 tumor suppressor

Intrinsically disordered proteins (IDPs) are a class of proteins without stable three-dimensional structures under physiological conditions. IDPs exhibit high dynamic nature and could be described by structural ensembles. As one of the most widely used tools, molecular dynamics (MD) simulation could provide the atomic descriptions of the structural ensemble of IDPs. However, the accuracy of the MD simulation largely depends on the accuracy of the force field. In this paper, we compared the structural ensembles of the activation domain 1 (AD1) in p53 tumor suppressor obtained from the widely used force fields, AMBER99SB-ILDN, CHARMM27, CHARMM36m with different water models. The results show that CHARMM36m generates more extended conformations than other force fields, while CHARMM27 prefers to sample the [Formula: see text]-helical structure. Moreover, the chemical shifts obtained by CHARMM36m are the closest to the experimental measurements. These results indicate that the CHARMM36m force field performs best in characterizing the structure properties of p53 AD1. Water models are also critical to describe the structural ensemble of IDPs. TIP4P water model can obtain more extended conformations and produce more local helical conformations than the TIP3P model in our simulation. In addition, we also compare the chemical shifts predicted by different chemical shift predicting programs with experimental measurements, the results show that SHIFTX2 obtains the best performance in the chemical shifts prediction.

Selectivity of 18-crown-6 ether to alkali ions by density functional theory and molecular dynamics simulation

Molecular dynamics (MD) simulation was employed to explore the structural and thermodynamic characteristics of 18-crown-6 ether selectivity to alkali metal ions in aqueous solution. The radial distribution functions (RDFs), potentials of mean force (PMFs) and binding free energies between Li+, Na+, K+, Rb+, Cs+ and 18-crown-6 ether with SPC/E and TIP4P water models in aqueous solutions at 278.15 K, 298.15 K and 318.15 K were obtained. It shows that K+ is selected over other alkali metal ions by 18-crown-6 ether in aqueous solution at all water models and temperatures we considered, and increasing temperature is beneficial for the binding of 18-crown-6 ether to alkali metal ions. The calculated selective sequence (Li+ ≪ Na+ < K+ > Rb+ > Cs+) in aqueous solution is close to the experimental results. The gas-phase binding energies of 18-crown-6 to alkali metal ions and the minimum potential energies with corresponding minimum positions of alkali ions complexes were got by density functional theory (DFT). The results reveal that the selectivity of 18-crown-6 to alkali metal ions in gas-phase follows the sequence, Li+ > Na+ > K+ > Rb+ > Cs+, and follows their ionic potentials. The effect of solvent on the selectivity of 18-crown-6 ether to alkali ions is confirmed by comparing the selective sequences in gas-phase and aqueous solution. The noncovalent interactions sites and types between alkali ions and 18-crown-6 ether were determined by electrostatic potential (ESP) and the average reduced density gradient (aRDG), which shows that the interactions between alkali ions and 18-crown-6 are electrostatic dominated.

Binding Insight of Anti-HIV Phytocompounds with Prime Targets of HIV: A Molecular Dynamics Simulation Analysis.

Despite intense efforts, AIDS is difficult to tackle by current anti-retroviral therapy (ART) due to its side effects; therefore, there is an urgent need to discover potential, multitarget and low-cost anti-HIV compounds. We have shown that few phytocompounds can potentially inhibit the prime targets of HIV namely GP120 envelope protein, reverse transcriptase, protease, integrase and ribonulcease. In this study, top ranked prioritized compounds were subjected to Molecular Dynamics (MD) simulation in order to study the conformational dynamics and integrity of crucial interaction in the receptor sites. The system was built for selected protein-ligand complex using TIP3P water model and OPLS_2005 force field. Trajectories were recorded up to 20 ns simulation time in Desmond module of Schrödinger software. As a result of a comprehensive analysis of molecular properties and dynamics of the complexes, it has been concluded that Chebulic acid, Curcumin and Mulberroside C could be developed as envelope glycoprotein GP120 inhibitor, reverse transcriptase inhibitor and protease inhibitor respectively. However, the fluctuation of Chebulic acid with respect to integrase and ribonuclease protein was higher during the simulation. These findings can aid in the designing of the structural properties for more effective anti-HIV compounds against the given targets.

Choice of Force Field for Proteins Containing Structured and Intrinsically Disordered Regions

Biomolecular force fields optimized for globular proteins fail to properly reproduce properties of intrinsically disordered proteins. In particular, parameters of the water model need to be modified to improve applicability of the force fields to both ordered and disordered proteins. Here, we compared performance of force fields recommended for intrinsically disordered proteins in molecular dynamics simulations of three proteins differing in the content of ordered and disordered regions (two proteins consisting of a well-structured domain and of a disordered region with and without a transient helical motif and one disordered protein containing a region of increased helical propensity). The obtained molecular dynamics trajectories were used to predict measurable parameters, including radii of gyration of the proteins and chemical shifts, residual dipolar couplings, paramagnetic relaxation enhancement, and NMR relaxation data of their individual residues. The predicted quantities were compared with experimental data obtained within this study or published previously. The results showed that the NMR relaxation parameters, rarely used for benchmarking, are particularly sensitive to the choice of force-field parameters, especially those defining the water model. Interestingly, the TIP3P water model, leading to an artificial structural collapse, also resulted in unrealistic relaxation properties. The TIP4P-D water model, combined with three biomolecular force-field parameters for the protein part, significantly improved reliability of the simulations. Additional analysis revealed only one particular force field capable of retaining the transient helical motif observed in NMR experiments. The benchmarking protocol used in our study, being more sensitive to imperfections than the commonly used tests, is well suited to evaluate the performance of newly developed force fields.

SAMPL6 Octanol-water partition coefficients from alchemical free energy calculations with MBIS atomic charges.

In molecular modeling the description of the interactions between molecules forms the basis for a correct prediction of macroscopic observables. Here, we derive atomic charges from the implicitly polarized electron density of 11 molecules in the SAMPL6 challenge using the Hirshfeld-I and Minimal Basis Set Iterative Stockholder (MBIS) partitioning method. These atomic charges combined with other parameters in the GAFF force field and different water/octanol models were then used in alchemical free energy calculations to obtain hydration and solvation free energies, which after correction for the polarization cost, result in the blind prediction of the partition coefficient. From the tested partitioning methods and water models the S-MBIS atomic charges with the TIP3P water model presented the smallest deviation from the experiment. Conformational dependence of the free energies and the energetic cost associated with the polarization of the electron density are discussed.

Molecular dynamics study of effect of wettability on rapid boiling over copper nanoparticle in water pool

Many special phenomena occur during rapid boiling of a liquid when it is largely superheated near the thermodynamic critical temperature which is much higher than the saturation temperature. However, although there are many experimental work reported the observation of rapid ejection near the surface region of the materials, they are limited in either length or time scales due to the limitation of classical macroscopic theory. In this work, molecular dynamics simulation of heat transfer from nanoparticle to a surrounding liquid pool is carried out to study the effect of solid-liquid interfacial wettability on the boiling phenomena of water around the nanoparticle. The interactions among copper atoms are described by the embedded atom method (EAM) potentials, and the TIP4P water model is used to describe water atoms. The results show that the interfacial thermal conductance is influenced by the interfacial wettability. Increasing the wettability between the particle and the fluid reduces the interfacial resistance. The nucleation of bubbles and formation of vapor occur more quickly for increased solid-liquid interfacial wettability. The results might change due to the uncertainties of the potential functions, therefore further investigation is worthwhile to study rapid boiling of water on nanoparticles using different water models.

Hydrophobicity of an isobutane dimer in water, methanol and acetonitrile as solvents — A classical molecular dynamics study

The potential of mean forces (PMFs) has been determined for an isobutane dimer in various solvent environments such as water, methanol and acetonitrile at a temperature of 298 K and pressure of 1 bar using GROMACS software. All the molecular dynamics (MD) simulations are carried out using a TIP3P water model under a CHARMM36 forcefield. Following Umbrella Sampling technique, PMFs are calculated and analyzed using Weighted Histogram Analysis Method (WHAM) and coordination number of first solvation shell is extracted for all solvents using radial distribution function. The shape of PMFs contains contact minima, solvent-separated minima and desolvation maxima. The values of contact minima are not affected much by solvent environment and found to be at 0.5377, 0.5480 and 0.5495 nm for water, methanol and acetonitrile respectively. The corresponding energy depths are found −0.9134, −0.7080 and −0.5295 kcalmol[Formula: see text]. The variation observed at solvent-separated minima is noticeable and found at 0.9012, 0.9721 and 0.9151 nm for water, methanol and acetonitrile, respectively. The coordination number of the first solvation shell by taking an isobutane molecule as a reference from their center of mass is found to be 28.1, 16.9 and 14.8 for water, methanol and acetonitrile, respectively. There is a soft hydrophobic interaction between isobutane dimer and solvents like methanol and acetonitrile relative to water, might be due to the presence of competitive methyl group of methanol and acetonitrile in the solvent medium.

Comprehensive evaluation of the MM-GBSA method on bromodomain-inhibitor sets.

MM-PB/GBSA methods represent a higher-level scoring theory than docking. This study reports an extensive testing of different MM-GBSA scoring schemes on two bromodomain (BRD) datasets. The first set is composed of 24 BRPF1 complexes, and the second one is a nonredundant set constructed from the PDBbind and composed of 28 diverse BRD complexes. A variety of MM-GBSA schemes were analyzed to evaluate the performance of four protocols with different numbers of minimization and MD steps, 10 different force fields and three different water models. Results showed that neither additional MD steps nor unfixing the receptor atoms improved scoring or ranking power. On the contrary, our results underscore the advantage of fixing receptor atoms or limiting the number of MD steps not only for a reduction in the computational costs but also for boosting the prediction accuracy. Among Amber force fields tested, ff14SB and its derivatives rather than ff94 or polarized force fields provided the most accurate scoring and ranking results. The TIP3P water model yielded the highest scoring and ranking power compared to the others. Posing power was further evaluated for the BRPF1 set. A slightly better posing power for the protocol which uses both minimization and MD steps with a fixed receptor than the one which uses only minimization with a fully flexible receptor-ligand system was observed. Overall, this study provides insights into the usage of the MM-GBSA methods for screening of BRD inhibitors, substantiating the benefits of shorter protocols and latest force fields and maintaining the crystal waters for accuracy.

Effects of ion–water Lennard-Jones potentials on the hydration dynamics around a monovalent atomic ion in molecular dynamics simulations

ABSTRACTThe molecular dynamics (MD) method is a promising technique to dissect the atomistic details of water dynamics around a solute. However, the quantitative predictions of experimentally measured kinetic properties, e.g. translational diffusion coefficient (D) and rotational relaxation time (τ), are not straightforward. Current water models have failed to reproduce these properties quantitatively; therefore, the fine-tuning of water models is required. In this study, we examined the effects of ion–water Lennard-Jones (LJ) potentials on the water dynamics around a monovalent atomic ion. For the TIP5P water model, we introduced new LJ potentials for the ion–hydrogen and ion–pseudoatom interactions, which were zero in the original TIP5P model. The hydration properties, i.e. D, τ, and the radius of the first hydration shell (rMO), were examined for various parameter settings. As a result, the new parameters certainly improved the reproducibility of the hydration properties in correspondence with experimental values. However, it is still difficult to reproduce faster rotational relaxation of hydration water than that of bulk water. In addition, we found that the three hydration properties (D, τ, and rMO) were artificially correlated in the MD simulations.

From a Liquid to a Crystal without Going through a First-Order Phase Transition: Determining the Free Energy of Melting with Glassy Intermediates.

The excess free energy of a liquid relative to an Einstein crystal reference state is calculated without going through a first-order phase transition. This is accomplished by going through an arrested glassy state to avoid a direct liquid to gas or liquid to crystal transition. The method is demonstrated by calculating the free energy difference between liquid water and ice Ih using the TIP4P and WAIL water models. TIP4P ice Ih melts at 232 ± 1 K, in close agreement with other estimates in the literature. WAIL ice melts at 272 ± 1 K, in good agreement with that of real water, which serves as a good validation of the quality of the WAIL model. The glassy intermediate method is easy to implement and amicable to parallel executions. We expect this method to have broad applications for calculating the liquid excess free energies for other materials.

Dynamical aspects of supercooled TIP3P-water in the grooves of DNA.

We investigate by molecular dynamics simulations the mobility of the water located at the DNA minor and major grooves. We employ the TIP3P water model, and our system is analyzed for a range of temperatures 190-300 K. For high temperatures, the water at the grooves shows an Arrhenius behavior similar to that observed in the bulk water. At lower temperatures, a departure from the bulk behavior is observed. This slowing down in the dynamics is compared with the dynamics of the hydrogen of the DNA at the grooves and with the autocorrelation functions of the water hydrogen bonds. Our results indicate that the hydrogen bonds of the water at the minor grooves are highly correlated, which suggests that this is the mechanism for the slow dynamics at this high confinement.

QM/MM Calculations for the Cl- + CH3Cl SN2 Reaction in Water Using CM5 Charges and Density Functional Theory.

The prototypical SN2 reaction of chloride ion with methyl chloride has been reinvestigated in aqueous solution using QM/MM methodology featuring MO6-2X/6-31+G(d) calculations with the TIP4P water model, and partial charges were computed with the CM5 method. Though the DFT method yields excellent gas-phase energetics for the reaction, the QM/MM approach is found to yield overestimation of the activation barrier by ca. 12 kcal/mol. The discrepancy is traced to underestimate of the magnitude of the partial charges on the chlorine atoms in the transition structure. When CM1 or CM3 charges based on semiempirical wave functions are used instead, the agreement with experiment is much improved. The findings emphasize the sensitivity of the results of QM/MM calculations to the choice of QM method, the MM force field, and implementation of the QM/MM interface.

Interactions of Water and Alkanes: Modifying Additive Force Fields to Account for Polarization Effects.

Atomistic biomolecular simulations predominantly utilize additive force fields (FF), where the electrostatic potential is modeled by fixed point charges. Among other consequences, the lack of polarizability in these models undermines the balance of hydrophilic/hydrophobic nonbonded interactions. Simulations of water/alkane systems using the TIP3P water model and CHARMM36 parameters reveal a 1 kcal/mol overestimate of the experimental transfer free energy of water to hexadecane; more recent optimized water models (SPC/E, TIP4P/2005, TIP4P-Ew, TIP3P-FB, TIP4P-FB, OPC, TIP4P-D) overestimate this transfer free energy by approximately 2 kcal/mol. In contrast, the polarizable SWM4-NDP and SWM6 water models reproduce experimental values to within statistical error. As an alternative to explicitly modeling polarizability, this paper develops an efficient automated workflow to optimize pair-specific Lennard-Jones parameters within an additive FF. Water/hexadecane is used as a prototype and the free energy of water transfer to hexadecane as a target. The optimized model yields quantitative agreement with the experimental transfer free energy and improves the water/hexadecane interfacial tension by 6%. Simulations of five different lipid bilayers show a strong increase of water permeabilities compared to the unmodified CHARMM36 lipid FF which consistently improves match with experiment: the order-of-magnitude underestimate for monounsaturated bilayers is rectified and the factor of 2.8-4 underestimate for saturated bilayers is turned into a factor of 1.5-3 overestimate. While agreement with experiment is decreased for the diffusion constant of water in hexadecane, alkane transfer free energies, and the bilayers' area per lipid, the method provides a permeant-specific route to achieve a wide range of heterogeneous observables via rapidly optimized pairwise parameters.

Parallel Implementation of a Sequential Markov Chain in Monte Carlo Simulations of Physical Systems with Pairwise Interactions.

In molecular simulations performed by Markov Chain Monte Carlo (typically employing the Metropolis criterion), each state of a system is obtained by a small random modification of the previous state. Therefore, the process consists of an immense number of small, quick to calculate steps, which are inherently sequential and hence considered to be very hard to parallelise. Here, we present a novel protocol for efficient calculation of multiple sequential steps in parallel. To this end, we first precompute in parallel energy components of all states achievable in a sequence of steps. Then we select a single path through all achievable states, which is identical with the path obtained with the sequential algorithm. As an example, we carried out simulations of the TIP5P water model with the new protocol and compared results with those obtained using the standard Metropolis Monte Carlo scheme. The implementation on the Titan X (Pascal) graphic processor (GPU) architectures allows for a 30-fold speedup in comparison with a simulation on a single core of a multicore CPU. The protocol is general and not limited to the GPU; it can also be used on multicore CPU when the longest possible length of the single simulation is required.

Calculation of the water-octanol partition coefficient of cholesterol for SPC, TIP3P, and TIP4P water.

We present a numerical study of the relative solubility of cholesterol in octanol and water. Our calculations allow us to compare the accuracy of the computed values of the excess chemical potential of cholesterol for several widely used water models (SPC, TIP3P, and TIP4P). We compute the excess solvation free energies by means of a cavity-based method [L. Li et al., J. Chem. Phys. 146(21), 214110 (2017)] which allows for the calculation of the excess chemical potential of a large molecule in a dense solvent phase. For the calculation of the relative solubility ("partition coefficient," log10 ) of cholesterol between octanol and water, we use the OPLS/AA force field in combination with the SPC, TIP3P, and TIP4P water models. For all water models studied, our results reproduce the experimental observation that cholesterol is less soluble in water than in octanol. While the experimental value for the partition coefficient is log10 = 3.7, SPC, TIP3P, and TIP4P give us a value of log10 = 4.5, 4.6, and 2.9, respectively. Therefore, although the results for the studied water models in combination with the OPLS/AA force field are acceptable, further work to improve the accuracy of current force fields is needed.

Revealing the interaction mechanism of peptide with specific platinum facets by stochastic tunneling-basin hopping method

This study uses the stochastic tunneling- basin hopping (STUN-BH) method to determine the most stable configurations for the peptides S7 (Ac-Ser-Ser-Phe-Pro-Gln-Pro-Asn-CONH2) and T7 (Ac-Thr-Leu-Thr-Thr-Leu-Thr-Asn-CONH2) on both Pt(1 1 1) and Pt(1 0 0) surfaces. The Amber99sb-ILDN force field was used for S7 and T7, and the TIP3P water model was used for the water environment. The CHARMM-METAL force field was adopted to model the S7 and T7 interactions with Pt(1 1 1) and Pt(1 0 0). The simulation results show that the lock-and-key (LAK) matches of both six-membered and five-membered rings of S7 on Pt(1 1 1) cause a more stable adsorption of S7 on Pt(1 1 1) than on Pt(1 0 0). For T7, the STUN-BH simulation results for multiple T7 molecules on Pt(1 1 1) and Pt(1 0 0) show the configurations of T7 molecules on Pt(1 0 0) are relatively more extended than those on Pt(1 1 1). Moreover, the order parameter of T7 molecules on Pt(1 0 0) is higher than that of T7 molecules on Pt(1 1 1), indicating that it is relatively easier for the T7 molecules to form a complete T7 layer on Pt(1 0 0) than on Pt(1 1 1).

The temperature-dependent structure, hydrogen bonding and other related dynamic properties of the standard TIP3P and CHARMM-modified TIP3P water models

It was found in the literature that the use of CHARMM-modified TIP3P (mTIP3P) water model in molecular dynamics (MD) simulations resulted in a less realistic peptide folding, as compared to the standard TIP3P (sTIP3P) water model. It is hypothesised that the hydrogen bonding (HB) network in both sTIP3P and mTIP3P could play a role in the discrepancy of the MD simulation results obtained. The lack of direct comparison between the structure, HB and other dynamic properties of sTIP3P and mTIP3P in the existing literature is likely due to the perception that these two water models should share similar structural and dynamic properties, despite a slight difference in their force-field potentials. In the results obtained from our MD simulations that were run from 10 °C to 100 °C at 1 atm, distinct differences are observed in their equilibrated densities, O⋯O radial pair distribution function, HB and other dynamic properties, such as dynamic viscosity and self-diffusion coefficient. The mTIP3P water model shows a more extensive HB network than sTIP3P by having a higher number of HB per molecule, higher dynamic viscosity and higher latent heat of vaporization.

Heat transfer around copper nanoparticle with high superheat in water pool: A molecular dynamics simulation

Rapid boiling of a liquid occurs when it is highly superheated near the thermodynamic critical temperature that is much higher than the saturation temperature. Many experimental investigation have been conducted to observe the rapid ejection phenomenon near the heated surface. However, these study were limited by either length or time scale because of the limitation of classical macroscopic theory. In this work, a molecular dynamics simulation of heat transfer from a nanoparticle to its surrounding liquid pool was carried out to study the effect of nanoparticle diameter on the heat transfer of water around it. The interactions among copper atoms were described by the Embedded Atom Method (EAM) potentials, and the TIP4P water model was used to model the water atoms. The results show that the temperature of water increases more quickly for nanoparticle with larger diameter, the heat flux increases apparently when the nanoparticle diameter increases, and the interfacial thermal resistance is relatively larger for bigger nanoparticle due to the vibrational mismatch present at the interface.