Molecular dynamics studies of a DNA-binding protein: 2. An evaluation of implicit and explicit solvent models for the molecular dynamics simulation of the Escherichia coli trp repressor.

Although aqueous simulations with periodic boundary conditions more accurately describe protein dynamics than in vacuo simulations, these are computationally intensive for most proteins. Trp repressor dynamic simulations with a small water shell surrounding the starting model yield protein trajectories that are markedly improved over gas phase, yet computationally efficient. Explicit water in molecular dynamics simulations maintains surface exposure of protein hydrophilic atoms and burial of hydrophobic atoms by opposing the otherwise asymmetric protein-protein forces. This properly orients protein surface side chains, reduces protein fluctuations, and lowers the overall root mean square deviation from the crystal structure. For simulations with crystallographic waters only, a linear or sigmoidal distance-dependent dielectric yields a much better trajectory than does a constant dielectric model. As more water is added to the starting model, the differences between using distance-dependent and constant dielectric models becomes smaller, although the linear distance-dependent dielectric yields an average structure closer to the crystal structure than does a constant dielectric model. Multiplicative constants greater than one, for the linear distance-dependent dielectric simulations, produced trajectories that are progressively worse in describing trp repressor dynamics. Simulations of bovine pancreatic trypsin were used to ensure that the trp repressor results were not protein dependent and to explore the effect of the nonbonded cutoff on the distance-dependent and constant dielectric simulation models. The nonbonded cutoff markedly affected the constant but not distance-dependent dielectric bovine pancreatic trypsin inhibitor simulations. As with trp repressor, the distance-dependent dielectric model with a shell of water surrounding the protein produced a trajectory in better agreement with the crystal structure than a constant dielectric model, and the physical properties of the trajectory average structure, both with and without a nonbonded cutoff, were comparable.     

A molecular dynamics simulation of bacteriophage T4 lysozyme.

An analysis of a 400 ps molecular dynamics simulation of the 164 amino acid enzyme T4 lysozyme is presented. The simulation was carried out with all hydrogen atoms modeled explicitly, the inclusion of all 152 crystallographic waters and at a temperature of 300 K. Temporal analysis of the trajectory versus energy, hydrogen bond stability, r.m.s. deviation from the starting crystal structure and radius of gyration, demonstrates that the simulation was both stable and representative of the average experimental structure. Average structural properties were calculated from the enzyme trajectory and compared with the crystal structure. The mean value of the C alpha displacements of the average simulated structure from the X-ray structure was 1.1 +/- 0.1 A; differences of the backbone phi and psi angles between the average simulated structure and the crystal structure were also examined. Thermal-B factors were calculated from the simulation for heavy and backbone atoms and both were in good agreement with experimental values. Relationships between protein secondary structure elements and internal motions were studied by examining the positional fluctuations of individual helix, sheet and turn structures. The structural integrity in the secondary structure units was preserved throughout the simulation; however, the A helix did show some unusually high atomic fluctuations. The largest backbone atom r.m.s. fluctuations were found in non-secondary structure regions; similar results were observed for r.m.s. fluctuations of non-secondary structure phi and psi angles. In general, the calculated values of r.m.s. fluctuations were quite small for the secondary structure elements. In contrast, surface loops and turns exhibited much larger values, being able to sample larger regions of conformational space. The C alpha difference distance matrix and super-positioning analyses comparing the X-ray structure with the average dynamics structure suggest that a 'hinge-bending' motion occurs between the N- and C-terminal domains.     

Protein dynamics in solution and in a crystalline environment: a molecular dynamics study

The effect of a solvent and a crystalline environment on the dynamics of proteins is investigated by the method of computer simulation. Three 25-ps molecular dynamics simulations at 300 K of the bovine pancreatic trypsin inhibitor (BPTI), consisting of 454 heavy atoms, are compared: one of BPTI in vacuo, one of BPTI in a box with 2647 spherical nonpolar solvent atoms, and one of BPTI surrounded by fixed crystal image atoms. Both average and time-dependent molecular properties are examined to determine the effect of the environment on the behavior of the protein. The dynamics of BPTI in solution or in the crystal environment are found to be very similar to that found in the vacuum calculation. The primary difference in the average properties is that the equilibrium structure in the presence of solvent or the crystal field is significantly closer to the X-ray structure than is the vacuum result; concomitantly, the more realistic environment leads to a number density closer to experiment. The presence of solvent has a negligible effect on the overall magnitude of the positional or dihedral angle fluctuations in the interior of the protein; however, there are changes in the decay times of the fluctuations of interior atoms. For surface residues, both the magnitude and the time course of the motions are significantly altered by the solvent. There tends to be an increase in the displacements of long side chains and the flexible parts of the main chain that protrude into the solvent. Further, these motions tend to have a more diffusive character with longer relaxation times than in vacuo. The crystal environment has a specific effect on a number of side chains which are held in relatively fixed positions through hydrogen-bond and electric interactions with the neighboring protein atoms. Most of the effects of the solution environment seem to be sufficiently nonspecific that it may be possible to model them by applying a mean field and stochastic dynamic methods.     

A 175-psec molecular dynamics simulation of camphor-bound cytochrome P-450cam.

The structure and internal motions of cytochrome P-450cam, a monooxygenase heme enzyme with 414 amino acid residues, with camphor bound at the active site have been evaluated on the basis of a 175-psec molecular dynamics simulation carried out at 300 K. All hydrogen atoms were explicitly modeled, and 204 crystallographic waters were included in the simulation. Based on an analysis of the time course of the trajectory versus potential energy, root mean square deviation, radius of gyration, and hydrogen bonding, the simulation was judged to be stable and representative of the average experimental structure. The averaged structural properties of the enzyme were evaluated from the final 135 psec of the simulation. The average atomic displacement from the X-ray structure was 1.39 A for all heavy atoms and 1.17 A for just C-alpha atoms. The average root-mean-square (rms) fluctuations of all heavy atoms and backbone atoms were 0.42 and 0.37 A, respectively. The computed rms fluctuations were in reasonable agreement with the experimentally determined temperature factors. All 13 segments of alpha-helix and 5 segments of beta-sheet were well preserved with the exception of the N-terminal half of helix F which alternated between an alpha-helix and a 310-helix. In addition there were in general only small variations in the relative orientation of adjacent alpha-helices. The rms fluctuations of the backbone dihedral angles in the secondary structure elements were almost uniformly smaller, with the fluctuation in alpha-helices and beta-sheets, 31 and 10% less, respectively, than those in nonsecondary structure regions. The reported crystal structure contains kinks in both helices C and I. In the simulation, both of these regions showed high mobility and large deviations from their starting positions. Since the kink in the I helix is at the oxygen binding site, these motions may have mechanistic implications.     

Temperature-dependent molecular dynamics and restrained X-ray refinement simulations of a Z-DNA hexamer

Molecular dynamics simulations of the Z-DNA hexamer 5BrdC-dG-5BrdC-dG-5BrdC-dG were performed at several temperatures between 100 K and 300 K. Above 250 K, a strong sequence-dependent flexibility in the nucleic acid is observed, with the guanine sugar and the phosphate of GpC sequences much more mobile than the cytosine sugar and phosphate of CpG sequences. At 300 K, the hexamer is in dynamic equilibrium between several Z forms, including the crystallographically determined ZI and ZII forms. The local base-pair geometry, however, is not very variable, except for the roll of the base-pairs. The hexamer molecular dynamics trajectories have been used to test the restrained parameter crystallographic refinement model for nucleic acids. X-ray diffraction intensities corresponding to observed diffraction data were computed. The average structures obtained from the simulations were then refined against the calculated intensities, using a restrained least-squares program developed for nucleic acids in order to analyse the effects of the refinement model on the derived quantities. In general, the temperature dependence of the atomic fluctuations determined directly from the refined Debye-Waller factors is in reasonably good agreement with the results obtained by calculating the atomic fluctuations directly from the Z-DNA molecular dynamics trajectories. The agreement is best for refinement of temperature factors without restraints. At the highest temperature studied (300 K), the effect of the refinement on the most mobile atoms (phosphates) is to significantly reduce the mean-square atomic fluctuations estimated from the refined Debye-Waller factors below the actual values (less than (delta r)2 greater than congruent to 0.5 A2). Analysis of the temperature-dependence of the mean-square atomic fluctuations provides information concerning the conformational potential within which the atoms move. The calculated temperature-dependence and anharmonicity of the Z-DNA helix are compared with the results observed for proteins. The average structures from the simulations were refined against the experimental X-ray intensities. It is found that low-temperature molecular dynamics simulations provide a useful tool for optimizing the refinement of X-ray structures.     

Molecular dynamics simulation of carboxy-myoglobin embedded in a trehalose-water matrix

We report on a molecular dynamics (MD) simulation of carboxy-myoglobin (MbCO) embedded in a water-trehalose system. The mean square fluctuations of protein atoms, calculated at different temperatures in the 100-300 K range, are compared with those from a previous MD simulation on an H2O-solvated MbCO and with experimental data from M?ssbauer spectroscopy and incoherent elastic neutron scattering on trehalose-coated MbCO. The results show that, for almost all the atomic classes, the amplitude of the nonharmonic motions stemming from the interconversion among the protein's conformational substates is reduced with respect to the H2O-solvated system, and their onset is shifted toward higher temperature. Moreover, our simulation shows that, at 300 K, the heme performs confined diffusive motions as a whole, leaving the underlying harmonic vibrations unaltered.     

Locally accessible conformations of proteins: multiple molecular dynamics simulations of crambin.

Multiple molecular dynamics (MD) simulations of crambin with different initial atomic velocities are used to sample conformations in the vicinity of the native structure. Individual trajectories of length up to 5 ns sample only a fraction of the conformational distribution generated by ten independent 120 ps trajectories at 300 K. The backbone atom conformational space distribution is analyzed using principal components analysis (PCA). Four different major conformational regions are found. In general, a trajectory samples only one region and few transitions between the regions are observed. Consequently, the averages of structural and dynamic properties over the ten trajectories differ significantly from those obtained from individual trajectories. The nature of the conformational sampling has important consequences for the utilization of MD simulations for a wide range of problems, such as comparisons with X-ray or NMR data. The overall average structure is significantly closer to the X-ray structure than any of the individual trajectory average structures. The high frequency (less than 10 ps) atomic fluctuations from the ten trajectories tend to be similar, but the lower frequency (100 ps) motions are different. To improve conformational sampling in molecular dynamics simulations of proteins, as in nucleic acids, multiple trajectories with different initial conditions should be used rather than a single long trajectory.     

Combination of Neutron Scattering and Molecular Dynamics to Determine Internal Motions in Biomolecules

Abstract

The forms and frequencies of atomic dynamics on the pico- and nanosecond timescales are accessible experimentally using incoherent neutron scattering. Molecular dynamics simulations cover the same space and time domains and neutron scattering intensities can be calculated from the simulations for direct comparison with experiment. To illustrate the complementarity of neutron scattering and molecular dynamics we examine measured and simulation-derived elastic incoherent scattering profiles from myoglobin and from the crystalline alanine dipeptide. Elastic incoherent scattering gives information on the geometry of the volume accessible to the atoms in the samples. The simulation-derived dipeptide elastic scattering profiles are in reasonable accord with experiment, deviations being due to the sampling limitations in the simulations and experimental detector normalisation procedures. The simulated dynamics is decomposed, revealing characteristic profiles due to rotational diffusional and translational vibrational motions of the methyl groups. In myoglobin, for which the timescale of the simulation matches more closely that accessible to the experiment, good agreement is seen for the elastic incoherent structure factor. This indicates that the space sampled by the hydrogen atoms in the protein on the timescale <100 ps is well represented by the simulation. Part of the helix atom fluctuations can be described in terms of rigid helix motions.     

Evidence of reduced flexibility in disulfide bridge-depleted azurin: a molecular dynamics simulation study.

Two molecular dynamics simulations have been performed for 2 ns, at room temperature, on fully hydrated wild type and Cys3Ala/Cys26Ala double-mutant azurin, to investigate the role of the unique disulfide bridge on the structure and dynamics of the protein. The results show that the removal of the [bond]SS[bond] bond does not affect the structural features of the protein, whereas alterations of the dynamical properties are observed. The root mean square fluctuations of the atomic positions are, on average, considerably reduced in the azurin mutant with respect to the wild type form. The number of intramolecular hydrogen bonds between protein backbone atoms that are lost during the simulation, with respect to the starting configuration, are reduced in the absence of the disulfide bond. The analysis of the dynamical cross-correlation map, characterising the protein co-ordinated internal motions, demonstrates in the mutated azurin a significant decrease in anti-correlated displacements between protein residues, with the only exception occurring in the region of the mutation sites. The overall findings show a relevant reduction in flexibility as a consequence of the disulfide bridge depletion in azurin, suggesting that the [bond]SS[bond] bond is a structural element which significantly contributes to the dynamic properties of the native protein.     

Translational hydration water dynamics drives the protein glass transition

Experimental and computer simulation studies have revealed the presence of a glass-like transition in the internal dynamics of hydrated proteins at approximately 200 K involving an increase of the amplitude of anharmonic dynamics. This increase in flexibility has been correlated with the onset of protein activity. Here, we determine the driving force behind the protein transition by performing molecular dynamics simulations of myoglobin surrounded by a shell of water. A dual heat bath method is used with which, in any given simulation, the protein and solvent are held at different temperatures, and sets of simulations are performed varying the temperature of the components. The results show that the protein transition is driven by a dynamical transition in the hydration water that induces increased fluctuations primarily in side chains in the external regions of the protein. The water transition involves activation of translational diffusion and occurs even in simulations where the protein atoms are held fixed.     

Visualization of energetics and conformations from molecular computer simulations

The quantity of data generated from molecular dynamics simulations and energy minimizations of macromolecules is overwhelming. It is an arduous task to extract the relevant and interesting information from the numerous coordinate sets produced. To help solve this problem, the authors have developed a method to aid the visualization of the relevant information from the simulations. This approach combines animation of the results on a high performance graphics device, such as the PS300, with colour-coded atoms based on changes in energy or conformation. The method will be illustrated using as examples: the molecular mechanics minimization of a nonapeptide, the molecular dynamics simulation of the protein myoglobin, including the analysis of the motion of helices during a 300ps trajectory, and changes in sugar puckering that occur during the molecular dynamics simulation of a DNA oligomer. The method is also applicable for analysing energy components and conformational properties of a fixed conformation.     

Structure and internal mobility of proteins: a molecular dynamics study of hen egg white lysozyme

The structure and internal motions of the protein hen egg white lysozyme are studied by analysis of simulation and experimental data. A molecular dynamics simulation and an energy minimization of the protein in vacuum have been made and the results compared with high-resolution structures and temperature factors of hen egg white lysozyme in two different crystal forms and of the homologous protein human lysozyme. The structures obtained from molecular dynamics and energy minimization have root-mean-square deviations for backbone atoms of 2.3 Å and 1.1–1.3 Å, respectively, relative to the crystal structures; the different crystal structures have root-mean-square deviations of 0.73–0.81 Å for the backbone atoms. In comparing the backbone dihedral angles, the difference between the dynamics and the crystal structure on which it is based is the same as that between any two crystal structures. The internal fluctuations of atomic positions calculated from the molecular dynamics trajectory agree well with the temperature factors from the three structures. Simulation and crystal results both show that there are large motions for residues involved in exposed turns of the backbone chain, relatively smaller motions for residues involved in the middle of helices or β-sheet structures, and relatively small motions of residues near disulfide bridges. Also, both the simulation and crystal data show that side-chain atoms have larger fluctuations than main-chain atoms. Moreover, the regions that have large deviations among the x-ray crystal structures, which indicates flexibility, are found to have large fluctuations in the simulation.     

Picosecond timescale rigid-helix and side-chain motions in deoxymyoglobin

The contribution of rigidbody motions to the atomic trajectories in a 100 ps molecular dynamics simulation of deoxymyoglobin is examined. Two typesof rigid-body motions are considered: one in which the helices are rigid units and one in which the side-chains are rigid units. Using a quaternionbased algorithm, fits of the rigid reference structures are made to each time frame of the simulation to derive trajectories of the rigid-body motions. The fitted trajectories are analysed in terms of atomic position fluctuations, mean-square displacements as a function of time, velocity autocorrelation functions and densities of states. The results are compared with the corresponding quantities calculated from the full trajectory. The relative contribution of the rigid helix motions to the helix atom dynamics depends on which quantity is examined and on which subset of atoms is chosen: rigid-helix motions contribute 86% of the rms helix backbone atomic position fluctuations, but 30% of the helix,: atom (backbone and side-chain) mean square displacements and only 1.1% of total kinetic energy. Only very low-frequency motions contribute to the rigid-helix dynamics; the rigid-body analysis allows characteristic rigid-helix vibrations to be identified and described. Treating the side-chains as rigid bodies is foundto be an excellent approximation to both their diffusive and vibrationalmean-square displacements: 96% of side-chain atom mean-square displacements originate from rigid side-Chain motions. However, the errors in theside-chain atomic positional fits are not always small. An analysis is madeof factors contributing to the positional error for different types of side-chain. © Wiley-Liss, Inc.     

A simple way to compute protein dynamics without a mechanical model

We found that in proteins the average atomic fluctuation is linearly related to the square of the atomic distance from the center of mass of the protein. Using this simple relation, we can accurately compute the temperature factors of proteins of a wide range of sizes and folds, and the correlation of the fluctuations in proteins. This simple relation provides a direct link between protein dynamics and the static protein's geometrical shape and offers a simple way to compute protein dynamics without either long time trajectory integration or any matrix operations.     

Long-term molecular dynamics simulation of copper plastocyanin in water

A long molecular dynamics simulation (1.1 ns) of fully hydrated plastocyanin has been performed and analysed to relate protein dynamics to structural elements and functional properties. The solvated structure is described in detail by the analysis of H-bond network. During all the simulation, the crystal H-bond network is maintained in the beta-sheet regions, while several H-bonds are broken or formed on the external surface of the protein. To evaluate whether such changes could be due to conformational rearrangements or to solvent competition, we have examined the average number of H-bonds between protein atoms and water molecules, and the root mean square deviations from crystal structure as a function of protein residues. Protein mobility and flexibility have been examined by positional and dihedral angle rms fluctuations. Finally, cross-correlation maps have revealed the existence of correlated motions among residues connected by hydrogen bonds.     

Convergence of molecular dynamics simulations of membrane proteins

The central question in evaluating almost any result from a molecular dynamics simulation is whether the calculation has converged. Unfortunately, assessing the ergodicity of a single trajectory is very difficult to do. In this work, we assess the sampling of molecular dynamics simulations of the membrane protein rhodopsin by comparing the results from 26 independent trajectories, each run for 100 ns. By examining principal components and cluster populations, we show that rhodopsin's fluctuations are not well described by 100 ns of dynamics, and that the sampling is not fully converged even for individual loops. The results serve as a reminder of the caution required when interpreting molecular dynamics simulations of macromolecules.     

A molecular dynamics study of the C-terminal fragment of the L7/L12 ribosomal protein. Secondary structure motion in a 150 picosecond trajectory

A 150 picosecond molecular dynamics computer simulation of the C-terminal fragment of the L7/L12 ribosomal protein from Escherichia coli is reported. The molecular dynamics results are compared with the available high-resolution X-ray data in terms of atomic positions, distances and positional fluctuations. Good agreement is found between the molecular dynamics results and the X-ray data. The form and parameters of the interaction potential energy function and the procedures for deriving it are discussed. Some current misunderstandings concerning the ways of evaluating the efficiency of molecular dynamics algorithms and of application of bond-length constraints in protein simulations are cleared up. The 150 picosecond trajectory has been scanned in a search for correlated motions within and between secondary structure elements. The beta-strands have diffusional stretching modes, and uncorrelated transversal displacements. The dynamic analysis of alpha-helices shows a variety of features. The atomic fluctuations differ between the helix ends; this effect reflects long time-scale motions. Two alpha-helices, alpha A and alpha C, show diffusive longitudinal stretching modes. The third helix, alpha B, has a correlated asymmetric longitudinal stretching; the N-terminal part dominates this behaviour. Furthermore, alpha B presents a librational motion with respect to the other parts of the molecule with a frequency of approximately 5 cm-1. This motion is coupled to helix stretching. Interestingly, the regions of highly conserved residues contain the most mobile parts of the molecule.     

Residence times of water molecules in the hydration sites of myoglobin

Hydration sites are high-density regions in the three-dimensional time-averaged solvent structure in molecular dynamics simulations and diffraction experiments. In a simulation of sperm whale myoglobin, we found 294 such high-density regions. Their positions appear to agree reasonably well with the distributions of waters of hydration found in 38 x-ray and 1 neutron high-resolution structures of this protein. The hydration sites are characterized by an average occupancy and a combination of residence time parameters designed to approximate a distribution of residence times. It appears that although the occupancy and residence times of the majority of sites are rather bulk-like, the residence time distribution is shifted toward the longer components, relative to bulk. The sites with particularly long residence times are located only in the cavities and clefts of the protein. This indicates that other factors, such as hydrogen bonds and hydrophobicity of underlying protein residues, play a lesser role in determining the residence times of the longest-lived sites.     

Effect of anisotropy and anharmonicity on protein crystallographic refinement. An evaluation by molecular dynamics

Molecular dynamics simulations are employed to determine the errors introduced by anharmonicity and anisotropy in the structure and temperature factors obtained for proteins by refinement of X-ray diffraction data. Simulations (25 ps and 300 ps) of metmyoglobin are used to generate time-averaged diffraction data at 1.5 A resolution. The crystallographic restrained-parameter least-squares refinement program PROLSQ is used to refine models against these simulated data. The resulting atomic positions and isotropic temperature factors are compared with the average structure and fluctuations calculated directly from the simulations. It is found that significant errors in the atomic positions and fluctuations are introduced by the refinement, and that the errors increase with the magnitude of the atomic fluctuations. Of particular interest is the fact that the refinement generally underestimates the atomic motions. Moreover, while the actual fluctuations go up to a mean-square value of about 5 A2, the X-ray results never go above approximately 2 A2. This systematic deviation in the motional parameters appears to be due to the use of a single-site isotropic model for the atomic fluctuations. Many atoms have multiple peaks in their probability distribution functions. For some atoms, the multiple peaks are seen in difference electron density maps and it is possible to include these in the refinement as disordered residues. However, for most atoms the refinement fits only one peak and neglects the rest, leading to the observed errors in position and temperature factor. The use of strict stereochemical restraints is inconsistent with the average dynamical structure; nevertheless, refinement with tight restraints results in structures that are comparable to those obtained with loose restraints and better than those obtained with no restraints. The results support the use of tight stereochemical restraints, but indicate that restraints on the variation of temperature factors are too restrictive.