Dynamic bond constraints in protein Langevin dynamics

Bond constraint algorithms for molecular dynamics typically take, as the target constraint lengths, the values of the equilibrium bond lengths defined in the potential. In Langevin form, the equations of motion are temperature dependent, which gives the average value for the individual bond lengths a temperature dependence. In addition to this, locally constant force fields can shift the local equilibrium bond lengths. To restore the average bond lengths in constrained integration to their unconstrained values, we suggest changing the constraint length used by popular constraint methods such as RATTLE [H. C. Andersen, J. Comput. Phys. 52, 23 (1983)] at each step. This allows us to more accurately capture the equilibrium bond length changes (with respect to the potential) due to the local equilibration and temperature effects. In addition, the approximations to the unconstrained nonbonded energies are closer using the dynamic constraint method than a traditional fixed constraint algorithm. The mechanism for finding the new constrained lengths involves one extra calculation of the bonded components of the force, and therefore adds O(N) time to the constraint algorithm. Since most molecular dynamics calculations are dominated by the O(N2) nonbonded forces, this new method does not take significantly more time than a fixed constraint algorithm.     

Intramolecular vibrational coupling contribution to temperature dependence of vibrational mode frequencies

High-frequency vibrational modes in molecules in solution are sensitive to temperature and shift either to lower or higher frequencies with the temperature increase. These frequency shifts are often attributed to specific interactions of the molecule and to the solvent polarization effect. We found that a substantial and often dominant contribution to sensitivity of vibrational high-frequency modes to temperature originates from anharmonic interactions with other modes in the molecule. The temperature dependencies were measured for several modes in ortho-, meta-, and para-isomers of acetylbenzonitrile in solution and in a solid matrix and compared to the theoretical predictions originated from the intramolecular vibrational coupling (IVC) evaluated using anharmonic density functional theory calculations. It is found that the IVC contribution is essential for temperature dependencies of all high-frequency vibrational modes and is dominant for many modes. As such, the IVC contribution alone permits predicting the main trend in the temperature dependencies, especially for vibrational modes with smaller transition dipoles. In addition, an Onsager reaction field theory was used to describe the solvent contribution to the temperature dependencies.     

Low-frequency vibrational spectrum of water in the hydration layer of a protein: a molecular dynamics simulation study

An atomistic molecular dynamics simulation has been carried out to understand the low-frequency intermolecular vibrational spectrum of water present in the hydration layer of the protein villin headpiece subdomain or HP-36. An attempt is made to explore how the heterogeneous rigidity of the hydration layers of different segments (three alpha helices) of the protein, strength of the protein-water hydrogen bonds, and their differential relaxation behavior influence the distribution of the intermolecular vibrational density of states of water in the hydration layers. The calculations revealed that compared to bulk water these bands are nonuniformly blue-shifted for water near the helices, the extent of shifts being more pronounced for water molecules hydrogen bonded to the protein residues. It is further noticed that the larger blue shift observed for the water molecules hydrogen bonded to helix 2 residues correlates excellently with the slowest structural relaxation of these hydrogen bonds. These results can be verified by suitable experimental measurements.     

The structure of liquid alcohols and the temperature dependence of vibrational bandwidth

To obtain the information about the structure peculiarities of liquid alcohols the temperature dependence of Raman profiles widths for these objects was studied. The results of these investigations have shown that in the region 150–340 K the widths of Raman bands in liquid methanol and ethanol are constant. From the data on the constancy of the number of hydrogen bonds per one molecule (in the same temperature range) the conclusion on the dynamic stability of the cluster structure of alcohols was made. The broadening of vibrational bands is determined by inhomogeneous broadening and by the dephasing of intramolecular vibrations due to the hydrogen bond dissociation.     

The Langevin Hull: Constant pressure and temperature dynamics for non-periodic systems

We have developed a new isobaric-isothermal (NPT) algorithm which applies an external pressure to the facets comprising the convex hull surrounding the system. A Langevin thermostat is also applied to the facets to mimic contact with an external heat bath. This new method, the "Langevin Hull", can handle heterogeneous mixtures of materials with different compressibilities. These systems are problematic for traditional affine transform methods. The Langevin Hull does not suffer from the edge effects of boundary potential methods, and allows realistic treatment of both external pressure and thermal conductivity due to the presence of an implicit solvent. We apply this method to several different systems including bare metal nanoparticles, nanoparticles in an explicit solvent, as well as clusters of liquid water. The predicted mechanical properties of these systems are in good agreement with experimental data and previous simulation work.     

Investigating hydration dependence of dynamics of confined water: monolayer, hydration water and Maxwell-Wagner processes

The dynamics of water confined in silica matrices MCM-41 C10 and C18, with pore diameter of 21 and 36 A, respectively, is examined by broadband dielectric spectroscopy (10(-2)-10(9) Hz) and differential scanning calorimetry for a wide temperature interval (110-340 K). The dynamics from capillary condensed hydration water and surface monolayer of water are separated in the analysis. Contrary to previous reports, the rotational dynamics are shown to be virtually independent on the hydration level and pore size. Moreover, a third process, also reported for other systems, and exhibiting a saddlelike temperature dependence is investigated. We argue that this process is due to a Maxwell-Wagner process and not to strongly bound surface water as previously suggested in the literature. The dynamics of this process is strongly dependent on the amount of hydration water in the pores. The anomalous temperature dependence can then easily be explained by a loss of hydration water at high temperatures in contradiction to previous explanations.     

Insight into the role of hydration on protein dynamics

The potential energy surface of a protein is rough. This intrinsic energetic roughness affects diffusion, and hence the kinetics. The dynamics of a system undergoing Brownian motion on this surface in an implicit continuum solvent simulation can be tuned via the frictional drag or collision frequency to be comparable to that of experiments or explicit solvent simulations. We show that the kinetic rate constant for a local rotational isomerization in stochastic simulations with continuum solvent and a collision frequency of 2 ps(-1) is about 10(4) times faster than that in explicit water and experiments. A further increase in the collision frequency to 60 ps(-1) slows down the dynamics, but does not fully compensate for the lack of explicit water. We also show that the addition of explicit water does not only slow down the dynamics by increasing the frictional drag, but also increases the local energetic roughness of the energy landscape by as much as 1.0 kcal/mol.     

Low-frequency vibrational modes and temperature dependence of NQR frequencies in 1,2,4,5-tetrabromobenzene

A temperature variation study of the nuclear quadrupole resonance (NQR) frequencies for the two resonance lines in 1,2,4,5-tetrabromobenzene was carried out in the range 77–300 K. The far-infrared spectrum of 1,2,4,5-tetrabromobenzene was recorded and the low-frequency vibrational modes noted. From the nine observed frequencies, the torsional modes contributing to the change in NQR frequencies with temperature are identified by correlating the calculated values of the NQR frequencies at different temperatures with the observed values.     

Molecular dynamics with the united-residue model of polypeptide chains. II. Langevin and Berendsen-bath dynamics and tests on model alpha-helical systems

The implementation of molecular dynamics (MD) with our physics-based protein united-residue (UNRES) force field, described in the accompanying paper, was extended to Langevin dynamics. The equations of motion are integrated by using a simplified stochastic velocity Verlet algorithm. To compare the results to those with all-atom simulations with implicit solvent in which no explicit stochastic and friction forces are present, we alternatively introduced the Berendsen thermostat. Test simulations on the Ala(10) polypeptide demonstrated that the average kinetic energy is stable with about a 5 fs time step. To determine the correspondence between the UNRES time step and the time step of all-atom molecular dynamics, all-atom simulations with the AMBER 99 force field and explicit solvent and also with implicit solvent taken into account within the framework of the generalized Born/surface area (GBSA) model were carried out on the unblocked Ala(10) polypeptide. We found that the UNRES time scale is 4 times longer than that of all-atom MD simulations because the degrees of freedom corresponding to the fastest motions in UNRES are averaged out. When the reduction of the computational cost for evaluation of the UNRES energy function is also taken into account, UNRES (with hydration included implicitly in the side chain-side chain interaction potential) offers about at least a 4000-fold speed up of computations relative to all-atom simulations with explicit solvent and at least a 65-fold speed up relative to all-atom simulations with implicit solvent. To carry out an initial full-blown test of the UNRES/MD approach, we ran Berendsen-bath and Langevin dynamics simulations of the 46-residue B-domain of staphylococcal protein A. We were able to determine the folding temperature at which all trajectories converged to nativelike structures with both approaches. For comparison, we carried out ab initio folding simulations of this protein at the AMBER 99/GBSA level. The average CPU time for folding protein A by UNRES molecular dynamics was 30 min with a single Alpha processor, compared to about 152 h for all-atom simulations with implicit solvent. It can be concluded that the UNRES/MD approach will enable us to carry out microsecond and, possibly, millisecond simulations of protein folding and, consequently, of the folding process of proteins in real time.     

Simulation analysis of the temperature dependence of lignin structure and dynamics

Lignins are hydrophobic, branched polymers that regulate water conduction and provide protection against chemical and biological degradation in plant cell walls. Lignins also form a residual barrier to effective hydrolysis of plant biomass pretreated at elevated temperatures in cellulosic ethanol production. Here, the temperature-dependent structure and dynamics of individual softwood lignin polymers in aqueous solution are examined using extensive (17 μs) molecular dynamics simulations. With decreasing temperature the lignins are found to transition from mobile, extended to glassy, compact states. The polymers are composed of blobs, inside which the radius of gyration of a polymer segment is a power-law function of the number of monomers comprising it. In the low temperature states the blobs are interpermeable, the polymer does not conform to Zimm/Stockmayer theory, and branching does not lead to reduction of the polymer size, the radius of gyration being instead determined by shape anisotropy. At high temperatures the blobs become spatially separated leading to a fractal crumpled globule form. The low-temperature collapse is thermodynamically driven by the increase of the translational entropy and density fluctuations of water molecules removed from the hydration shell, thus distinguishing lignin collapse from enthalpically driven coil-globule polymer transitions and providing a thermodynamic role of hydration water density fluctuations in driving hydrophobic polymer collapse. Although hydrophobic, lignin is wetted, leading to locally enhanced chain dynamics of solvent-exposed monomers. The detailed characterization obtained here provides insight at atomic detail into processes relevant to biomass pretreatment for cellulosic ethanol production and general polymer coil-globule transition phenomena.     

Solvent viscosity dependence of the protein folding dynamics

Solvent viscosity has been frequently adopted as an adjustable parameter in various computational studies (e.g., protein folding simulations) with implicit solvent models. A common approach is to use low viscosities to expedite simulations. While using viscosities lower than that of aqueous is unphysical, such treatment is based on observations that the viscosity affects the kinetics (rates) in a well-defined manner as described by Kramers' theory. Here, we investigate the effect of viscosity on the detailed dynamics (mechanism) of protein folding. On the basis of a simple mathematical model, we first show that viscosity may indeed affect the dynamics in a complex way. By applying the model to the folding of a small protein, we demonstrate that the detailed dynamics is affected rather pronouncedly especially at unphysically low viscosities, cautioning against using such viscosities. In this regard, our model may also serve as a diagnostic tool for validating low-viscosity simulations. It is also suggested that the viscosity dependence can be further exploited to gain information about the protein folding mechanism.     

On the temperature dependence of electronically non-adiabatic vibrational energy transfer in molecule-surface collisions

Here we extend a recently introduced state-to-state kinetic model describing single- and multi-quantum vibrational excitation of molecular beams of NO scattering from a Au(111) metal surface. We derive an analytical expression for the rate of electronically non-adiabatic vibrational energy transfer, which is then employed in the analysis of the temperature dependence of the kinetics of direct overtone and two-step sequential energy transfer mechanisms. We show that the Arrhenius surface temperature dependence for vibrational excitation probability reported in many previous studies emerges as a low temperature limit of a more general solution that describes the approach to thermal equilibrium in the limit of infinite interaction time and that the pre-exponential term of the Arrhenius expression can be used not only to distinguish between the direct overtone and sequential mechanisms, but also to deduce their relative contributions. We also apply the analytical expression for the vibrational energy transfer rates introduced in this work to the full kinetic model and obtain an excellent fit to experimental data, the results of which show how to extract numerical values of the molecule-surface coupling strength and its fundamental properties.     

A generalized Langevin dynamics approach to model solvent dynamics effects on proteins via a solvent-accessible surface. The carboxypeptidase A inhibitor protein as a model

A generalized Langevin dynamics (GLD) scheme is derived for (bio)macromolecules having internal structure, arbitrary shapes and a size larger than solvent molecules (i.e. proteins). The concept of solvent-accessible surface area (SASA) is used to incorporate solvent effects via external forces thereby avoiding its explicit molecular representation. A simulation algorithm is implemented in the GROMOS molecular dynamics (MD) program including random forces and memory effects, while solvation effects enter via derivatives of the surface area. The potato carboxypeptidase inhibitor (PCI), a small protein, is used to numerically test the approach. This molecule has N- and C-terminal tails whose structure and fluctuations are solvent dependent. A 1-ns MD trajectory was analyzed in depth. X-ray and NMR structures are used in conjunction with MD simulations with and without explicit solvent to gauge the quality of the results. All the analyses showed that the GLD simulation approached the results obtained for the MD simulation with explicit simple-point-charge-model water molecules. The SASAs of the polar atoms show a natural exposure towards the solvent direction. A FLS solvent simulation was completed in order to sense memory effects. The approach and results presented here could be of great value for developing alternatives to the use of explicit solvent molecules in the MD simulation of proteins, expanding its use and the time-scale explored. Received: 2 February 2000 / Revised: 12 March 2000 / Accepted: 26 May 2000 / Published online: 2 November 2000     

Molecular dynamics study of the hydration of the hydroxyl radical at body temperature

Classical molecular dynamics (MD) simulation of ˙OH in liquid water at 37 °C has been performed using flexible models of the solute and solvent molecules. We derived the Morse function describing the bond stretching of the radical and the potential for ˙OH-H(2)O interactions, including short-range interactions of hydrogen atoms. Scans of the potential energy surface of the ˙OH-H(2)O complex have been performed using the DFT method with the B3LYP functional and the 6-311G(d,p) basis set. The DFT-derived partial charges, ±0.375e, and the equilibrium bond-length, 0.975 ?, of ˙OH resulted in the dipole moment of 1.76 D. The radical-water radial distribution functions revealed that ˙OH is not built into the solvent structure but it rather occupies distortions or cavities in the hydrogen-bonded network. The solvent structure at 37 °C has been found to be the same as that of pure water. The hydration cage of the radical comprises 13-14 water molecules. The estimated hydration enthalpy -42 ± 5 kJ mol(-1) is comparable with the experimental value -39 ± 6 kJ mol(-1) for 25 °C. Inspection of hydrogen bonds showed the importance of short-range interaction of hydrogen atoms and indicated that neglect of the angular condition greatly overestimates the number of the H-acceptor radical-water bonds. The mean number ?n = 0.85 of radical-water H-bonds has been calculated using geometric definition of H-bond and ?n = 0.62 has been obtained when the energetic condition, E(da)≤-8 kJ mol(-1), was additionally considered. The continuous lifetimes of 0.033 ps and 0.024 ps have been estimated for the radical H-donor and the H-acceptor bonds, respectively. Within statistical uncertainty the radical self-diffusion coefficient, (2.9 ± 0.6) × 10(-9) m(2) s(-1), is the same as (3.1 ± 0.5) × 10(-9) m(2) s(-1) calculated for water in solution and in pure solvent. To the best of our knowledge, this is the first study of the ˙OH(aq) properties at a biologically relevant body temperature.     

How protein surfaces induce anomalous dynamics of hydration water

Water around biomolecules slows down with respect to pure water, and both rotation and translation exhibit anomalous time dependence in the hydration shell. The origin of such behavior remains elusive. We use molecular dynamics simulations of water dynamics around several designed protein models to establish the connection between the appearance of the anomalous dynamics and water-protein interactions. For the first time we quantify the separate effect of protein topological and energetic disorder on the hydration water dynamics. When a static protein structure is simulated, we show that both types of disorder contribute to slow down water diffusion, and that allowing for protein motion, increasing the spatial dimensionality of the interface, reduces the anomalous character of hydration water. The rotation of water is, instead, altered by the energetic disorder only; indeed, when electrostatic interactions between the protein and water are switched off, water reorients even faster than in the bulk. The dynamics of water is also related to the collective structure--à voir the hydrogen bond (H-bond) network--formed by the solvent enclosing the protein surface. We show that, as expected for a full hydrated protein, when the protein surface offers pinning sites (charged or polar sites), the superficial water-water H-bond network percolates throughout the whole surface, hindering the water diffusion, whereas it does not when the protein surface lacks electrostatic interactions with water and the water diffusion is enhanced.     

Justification for the continuous model of the structure of liquid water from an analysis of the temperature dependence of the vibrational spectra

Conclusion The method proposed by Zhukovskii for calculating the temperature transformation of the contours [Eq. (6)], which, as we have shown here, effectively describes the behavior of the infrared absorption and Raman spectra of water over a wide range of temperatures, is based entirely on the Boltzmann law of the distribution of the energies of the hydrogen bonds within the limits of the continuous contour P(E). The success of this procedure, in our view, proves convincingly the existence of this distribution P(E). This provides direct experimental justification for the basic postulates of the continuous model of the structure of liquid water. The view that there is a continuous statistical distribution of the frequencies of the OH (OD) oscillators perturbed by H-bonding makes it possible to develop a systematic theory of the form of the valence bands both in pure water [10] and in its complexes with organic bases [9], in good agreement with the experiment. Thus the continuous model of the structure of water provides a constructive method for calculating the form of the broad vibrational bands of water and their temperature dependence.Institute of Chemical Kinetics and Combustion, Siberian Branch, Academy of Sciences of the USSR. Translated from Zhurnal Strukturnoi Khimii, Vol. 21, No. 3, pp. 95–105, May–June, 1980.     

Composition and temperature dependence of cesium-borate glasses by molecular dynamics

The structural aspects of xCs2O-(1-x)B2O3 glasses have been investigated by molecular dynamics as functions of Cs2O content (x=0.2, 0.3, and 0.4) and temperature (T=300 and 1250 K). The tetrahedral (B?4-) and triangular (B?3,B?2O-, and B?O2 (2-)) short-range order borate units were found to be the structure-building entities of the simulated glasses [?=bridging oxygen (BO) and O-=nonbridging oxygen (NBO) atom]. The increase of Cs2O content results in the progressive increase of the NBO-containing triangle population at the expense of the BO4- tetrahedral units. The same effect is caused by temperature increase at a fixed Cs2O content, and this was associated with the "fragile" characteristics of alkali borate glasses. A comparison of simulated Cs and Li borates showed very similar structures at x=0.2, but dissimilar ones when the alkali content exceeds this composition. In particular, for x>0.2 Cs borates exhibit a preference for NBO formation relative to Li borates. Differences in the microstructure of sites hosting Cs ions were found, and this permits their classification into bridging (b type) and nonbridging type (nb type) of sites. b-type sites consist exclusively of BO atoms, while both BO and NBO atoms participate in nb-type sites. These differences in Cs-site local bonding characteristics were found to be reflected on the Cs-O(site) vibration frequencies. Also, the computed Cs-O vibrational responses for simulated Cs borates were found to compare well with experimental far-infrared spectra.     

Temperature dependence of anisotropic protein backbone dynamics

The measurement of (15)N NMR spin relaxation, which reports the (15)N-(1)H vector reorientational dynamics, is a widely used experimental method to assess the motion of the protein backbone. Here, we investigate whether the (15)N-(1)H vector motions are representative of the overall backbone motions, by analyzing the temperature dependence of the (15)N-(1)H and (13)CO-(13)C(alpha) reorientational dynamics for the small proteins binase and ubiquitin. The latter dynamics were measured using NMR cross-correlated relaxation experiments. The data show that, on average, the (15)N-(1)H order parameters decrease only by 2.5% between 5 and 30 degrees C. In contrast, the (13)CO-(13)C(alpha) order parameters decrease by 10% over the same temperature trajectory. This strongly indicates that there are polypeptide-backbone motions activated at room temperature that are not sensed by the (15)N-(1)H vector. Our findings are at variance with the common crank-shaft model for protein backbone dynamics, which predicts the opposite behavior. This study suggests that investigation of the (15)N relaxation alone would lead to underestimation of the dynamics of the protein backbone and the entropy contained therein.     

Design of quasisymplectic propagators for Langevin dynamics

A vector field splitting approach is discussed for the systematic derivation of numerical propagators for deterministic dynamics. Based on the formalism, a class of numerical integrators for Langevin dynamics are presented for single and multiple time step algorithms.