A coupling of homology modeling with multiple molecular dynamics simulation for identifying representative conformation of GPCR structures: a case study on human bombesin receptor subtype-3

In this study, a computational pipeline was therefore devised to overcome homology modeling (HM) bottlenecks. The coupling of HM with molecular dynamics (MD) simulation is useful in that it tackles the sampling deficiency of dynamics simulations by providing good-quality initial guesses for the native structure. Indeed, HM also relaxes the severe requirement of force fields to explore the huge conformational space of protein structures. In this study, the interaction between the human bombesin receptor subtype-3 and MK-5046 was investigated integrating HM, molecular docking, and MD simulations. To improve conformational sampling in typical MD simulations of GPCRs, as in other biomolecules, multiple trajectories with different initial conditions can be employed rather than a single long trajectory. Multiple MD simulations of human bombesin receptor subtype-3 with different initial atomic velocities are applied to sample conformations in the vicinity of the structure generated by HM. The backbone atom conformational space distribution of replicates is analyzed employing principal components analysis. As a result, the averages of structural and dynamic properties over the twenty-one trajectories differ significantly from those obtained from individual trajectories.     

A molecular dynamics study of the C-terminal fragment of the L7/L12 ribosomal protein

The crystallographic dimer of the C-terminal fragment (CTF) of the L7/L12 ribosomal protein has been subjected to molecular dynamics (MD) simulations. A 90 picosecond (ps) trajectory for the protein dimer, 19 water molecules and two counter ions has been calculated at constant temperature. Effects of intermolecular interactions on the structure and dynamics have been studied. The exact crystallographic symmetry is lost and the atomic fluctuations differ from one monomer to the other. The average MD structure is more stable than the X-ray one, as judged by accessible surface area and energy calculations. Crystal (non-dimeric) interactions have been simulated in another 40 ps trajectory by using harmonic restraints to represent intermolecular hydrogen bonds. The conformational changes with respect ot the X-ray structure are then virtually suppressed.The unrestrained dimer trajectory has been scanned for cooperative motions involving secondary structure elements. The intrinsic collective motions of the monomer are transmitted via intermolecular contacts to the dimer structure.The existence of a stable dimeric form of CTF, resembling the crystallographic one, has been documented. At the cost of fairly small energy expenditure the dimer has considerable conformational flexibility. This flexibility may endow the dimer with some functional potential as an energy transducer.     

Conformational sampling by NMR solution structures calculated with the program DIANA evaluated by comparison with long-time molecular dynamics calculations in explicit water

The NMR solution structure of bovine pancreatic trypsin inhibitor (BPTI) obtained by distance geometry calculations with the program DIANA is compared with groups of conformers generated by molecular dynamics (MD) simulations in explicit water at ambient temperature and pressure. The MD simulations started from a single conformer and were free or restrained either by the experimental NOE distance restraints or by time-averaged restraints; the groups of conformers were collected either in 10 ps intervals during 200 ps periods of simulation, or in 50 ps intervals during a 1 ns period of simulation. Overall, these comparisons show that the standard protein structure determination protocol with the program DIANA provides a picture of the protein structure that is in agreement with MD simulations using “realistic” potential functions over a nanosecond timescale. For well-constrained molecular regions there is a trend in the free MD simulation of duration 1 ns that the sampling of the conformation space is slightly increased relative to the DIANA calculations. In contrast, for surface-exposed side-chains that are less extensively constrained by the NMR data, the DIANA conformers tend to sample larger regions of conformational space than conformers selected from any of the MD trajectories. Additional insights into the behavior of surface side-chains come from comparison of the MD runs of 200 ps or 1 ns duration. In this time range the sampling of conformation space by the protein surface depends strongly on the length of the simulation, which indicates that significant side-chain transitions occur on the nanosecond timescale and that much longer simulations will be needed to obtain statistically significant data on side-chain dynamics.     

Trajectory analysis of NMR structure calculations

Summary NMR as well as X-ray crystallography are used to determine the three-dimensional structures of macromolecules at atomic resolution. Structure calculation generates coordinates that are compatible with NMR data from randomly generated initial structures. We analyzed the trajectory taken by structures during NMR structure calculation in conformational space, assuming that the distance between two structures in conformational space is the root-mean-square deviation between the two structures. The coordinates of a structure in conformational space were obtained by applying the metric multidimensional scaling method. As an example, we used a 22-residue peptide, -Conotoxin GIIIA, and a simulated annealing protocol of XPLOR. We found that the three-dimensional solution of the multidimensional scaling analysis is sufficient to describe the overall configuration of the trajectories in conformational space. By comparing the trajectories of the entire calculation with those of the converged calculation, random sampling of conformational space is readily discernible. Trajectory analysis can also be used for optimization of protocols of NMR structure calculation, by examining individual trajectories.Abbreviations MD molecular dynamics - MDS multidimensional scaling - rmsd root-mean-square deviation - armsd angular rmsd - R multiple correlation coefficient - YASAP yet another simulated annealing protocol - PCA principal component analysis     

Brownian and essential dynamics studies of the HIV-1 integrase catalytic domain

The three-dimensional structure of the active site region of the enzyme HIV-1 integrase is not unambiguously known. This region includes a flexible peptide loop that cannot be well resolved in crystallographic determinations. Here we present two different computational approaches with different levels of resolution and on different time-scales to understand this flexibility and to analyze the dynamics of this part of the protein. We have used molecular dynamics simulations with an atomic model to simulate the region in a realistic and reliable way for 1 ns. It is found that parts of the loop wind up after 300 ps to extend an existing helix. This indicates that the helix is longer than in the earlier crystal structures that were used as basis for this study. Very recent crystal data confirms this finding, underlining the predictive value of accurate MD simulations. Essential dynamics analysis of the MD trajectory yields an anharmonic motion of this loop. We have supplemented the MD data with a much lower resolution Brownian dynamics simulation of 600 ns length. It provides ideas about the slow-motion dynamics of the loop. It is found that the loop explores a conformational space much larger than in the MD trajectory, leading to a "gating"-like motion with respect to the active site.     

Brownian and Essential Dynamics Studies of the HIV-1 Integrase Catalytic Domain

Abstract

The three-dimensional structure of the active site region of the enzyme HIV-1 integrase is not unambiguously known. This region includes a flexible peptide loop that cannot be well resolved in crystallographic determinations. Here we present two different computional approaches with different levels of resolution and on different time-scales to understand this flexibility and to analyze the dynamics of this part of the protein. We have used molecular dynamics simulations with an atomic model to simulate the region in a realistic and reliable way for 1 ns. It is found that parts of the loop wind up after 300 ps to extend an existing helix. This indicates that the helix is longer than in the earlier crystal structures that were used as basis for this study. Very recent crystal data confirms this finding, underlining the predictive value of accurate MD simulations. Essential dynamics analysis of the MD trajectory yields an anharmonic motion of this loop. We have supplemented the MD data with a much lower resolution Brownian dynamics simulation of 600 ns length. It provides ideas about the slow-motion dynamics of the loop. It is found that the loop explores a conformational space much larger than in the MD trajectory, leading to a “gating”-like motion with respect to the active site.     

Molecular dynamics of myoglobin at 298 degrees K. Results from a 300-ps computer simulation.

We have carried out a very long (300 ps) molecular dynamics simulation of the protein myoglobin. This trajectory is approximately three times longer than the longest previous molecular dynamics simulation of a protein, and ten times longer than protein simulations of comparable size (1,423 atoms in our model). Here we report results from this long simulation concerning the average structure, the mean square fluctuations of atoms about the average structure, and the nuclear magnetic resonance order parameters for various groups in myoglobin. The results demonstrate that the average coordinates change very slowly during the simulation. The relative atomic mobilities are well described by the simulation. For both the mean square atomic fluctuations and the order parameters, however, there are significant quantitative differences when values calculated using shorter portions of the trajectory are compared with results obtained for the entire 300-ps simulation. The implications of this result for obtaining converged properties from protein molecular dynamics simulations for comparison with experiment are discussed.     

Conformational dynamics of cytochrome c: correlation to hydrogen exchange.

We study the dynamical fluctuations of horse heart cytochrome c by molecular dynamics (MD) simulations in aqueous solution, at four temperatures: 300 K, 360 K, 430 K, and 550 K. Each simulation covers a production time of at least 1.5 nanoseconds (ns). The conformational dynamics of the system is analyzed in terms of collective motions that involve the whole protein, and local motions that involve the formation and breaking of intramolecular hydrogen bonds. The character of the MD trajectories can be described within the framework of rugged energy landscape dynamics. The MD trajectories sample multiple conformational minima, with basins in protein conformational space being sampled for a few hundred picoseconds. The trajectories of the system in configurational space can be described in terms of diffusion of a particle in real space with a waiting time distribution due to partial trapping in shallow minima. As a consequence of the hierarchical nature of the dynamics, the mean square displacement autocorrelation function, <|x(t) - x(0)|2>, exhibits a power law dependence on time, with an exponent of around 0.5 for times shorter than 100 ps, and an exponent of 1.75 for longer times. This power law behavior indicates that the system exhibits suppressed diffusion (sub-diffusion) in sampling of configurational space at time scales shorter than 100 ps, and enhanced (super-diffusion) at longer time scales. The multi-basin feature of the trajectories is present at all temperatures simulated. Structural changes associated with inter-basin displacements correspond to collective motions of the Omega loops and coiled regions and relative motions of the alpha-helices as rigid bodies. Similar motions may be involved in experimentally observed amide hydrogen exchange. However, some groups showing large correlated motions do not expose the amino hydrogens to the solvent. We show that large fluctuations are not necessarily correlated to hydrogen exchange. For example, regions of the proteins forming alpha helices and turns show significant fluctuations, but as rigid bodies, and the hydrogen bonds involved in the formation of these structures do not break in proportion to these fluctuations. Proteins 1999;36:175-191. Published 1999 Wiley-Liss, Inc.     

Molecular dynamics of sickle and normal hemoglobins

Molecular dynamics (MD) simulations have been carried out for 62.5 ps on crystal structures of deoxy sickle cell hemoglobin (HbS) and normal deoxy hemoglobin (HbA) using the CHARMM MD algorithm, with a time step of 0.001 ps. In the trajectory analysis of the 12.5–62.5 (50 ps) simulation, oscillations of the radius of gyration and solvent-accessible surface area were calculated. HbS exhibited a general contraction during the simulation, while HbA exhibited a nearly constant size. The average deviations of simulated structures from the starting structures were found to be 1.8 Å for HbA and 2.3 Å for HbS. The average rms amplitudes of atomic motions (atomic flexibility) were about 0.7 Å for HbA and about 1.0 Å for HbS. The amplitudes of backbone motion correlate well with temperature factors derived from x-ray crystallography. A comparison of flexibility between the α- and β-chains in both HbA and HbS indicates that the β-chains generally exhibited greater flexibility than the α-chains, and that the HbS β-chains exhibit greater flexibility in the N-terminal and D- and F-helix regions than do those of HbA. The average amplitude of backbone torsional oscillations was about 9°, a value comparable with that of other simulations, with enhanced torsional oscillation occurring primarily at the ends of helices or in loop regions between helices. Comparison of atomic flexibility and torsional oscillation results suggests that the increased β-chain flexibility results from relatively concerted motions of secondary structure elements. The increased flexibility may play an important role in HbS polymerization. Time course analysis of conformational energy of association, hydrogen bonding and hydrophobic bonding (as calculated from solvent accessibility) shows that all three of these factors contribute to the stability of subunit association for both hemoglobins. © 1993 John Wiley & Sons, Inc.     

Mass-weighted molecular dynamics simulation and conformational analysis of polypeptide.

Atomic motions in protein molecules have been studied by molecular dynamics (MD) simulations; dynamics simulation methods have also been employed in conformational studies of polypeptide molecules. It was found that when atomic masses are weighted, the molecular dynamics method can significantly increase the sampling of dihedral conformation space in such studies, compared to a conventional MD simulation of the same total simulation time length. Herein the theoretical study of molecular conformation sampling by the molecular dynamics-based simulation method in which atomic masses are weighted is reported in detail; moreover, a numerical scheme for analyzing the extensive conformational sampling in the simulation of a tetrapeptide amide molecule is presented. From numerical analyses of the mass-weighted molecular dynamics trajectories of backbone dihedral angles, low-resolution structures covering the entire backbone dihedral conformation space of the molecule were determined, and the distribution of rotationally stable conformations in this space were analyzed quantitatively. The theoretical analyses based on the computer simulation and numerical analytical methods suggest that distinctive regimes in the conformational space of the peptide molecule can be identified.     

Molecular dynamics of the anticodon domain of yeast tRNA(Phe): codon-anticodon interaction

We have studied the effect of codon-anticodon interaction on the structure and dynamics of transfer RNAs using molecular dynamics simulations over a nanosecond time scale. From our molecular dynamical investigations of the solvated anticodon domain of yeast tRNA(Phe) in the presence and absence of the codon trinucleotides UUC and UUU, we find that, although at a gross level the structures are quite similar for the free and the bound domains, there are small but distinct differences in certain parts of the molecule, notably near the Y37 base. Comparison of the dynamics in terms of interatomic or inter-residual distance fluctuation for the free and the bound domains showed regions of enhanced rigidity in the loop region in the presence of codons. Because fluorescence experiments suggested the existence of multiple conformers of the anticodon domain, which interconvert on a much larger time scale than our simulations, we probed the conformational space using five independent trajectories of 500 ps duration. A generalized ergodic measure analysis of the trajectories revealed that at least for this time scale, all the trajectories populated separate parts of the conformational space, indicating a need for even longer simulations or enhanced sampling of the conformational space to give an unequivocal answer to this question.     

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.     

Prediction of titration properties of structures of a protein derived from molecular dynamics trajectories.

This paper explores the dependence of the molecular dynamics (MD) trajectory of a protein molecule on the titration state assigned to the molecule. Four 100-ps MD trajectories of bovine pancreatic trypsin inhibitor (BPTI) were generated, starting from two different structures, each of which was held in two different charge states. The two starting structures were the X-ray crystal structure and one of the solution structures determined by NMR, and the charge states differed only in the ionization state of N terminus. Although it is evident that the MD simulations were too short to sample fully the equilibrium distribution of structures in each case, standard Poisson-Boltzmann titration state analysis of the resulting configurations shows general agreement between the overall titration behavior of the protein and the charge state assumed during MD simulation: at pH 7, the total net charge of the protein resulting from the titration analysis is consistently lower for the protein with the N terminus assumed to be neutral than for the protein with the N terminus assumed to be charged. For most of the ionizable residues, the differences in the calculated pKaS among the four trajectories are statistically negligible and remain in good agreement with the data obtained by crystal structure titration and by experiment. The exceptions include the N terminus, which responds directly to the change of its imposed charge; the C terminus, which in the NMR structure interacts strongly with the former; and a few other residues (Arg 1, Glu 7, Tyr 35, and Arg 42) whose pKaS reflect the initial structure and the limited trajectory lengths. This study illustrates the importance of the careful assignment of protonation states at the start of MD simulations and points to the need for simulation methods that allow for the variation of the protonation state in the calculation of equilibrium properties.     

A statistical approach to the interpretation of molecular dynamics simulations of calmodulin equilibrium dynamics

A sample of 35 independent molecular dynamics (MD) simulations of calmodulin (CaM) equilibrium dynamics was prepared from different but equally plausible initial conditions (20 simulations of the wild-type protein and 15 simulations of the D129N mutant). CaM's radius of gyration and backbone mean-square fluctuations were analyzed for the effect of the D129N mutation, and simulations were compared with experiments. Statistical tests were employed for quantitative comparisons at the desired error level. The computational model predicted statistically significant compaction of CaM relative to the crystal structure, consistent with the results of small-angle X-ray scattering (SAXS) experiments. This effect was not observed in several previously reported studies of (Ca2+)(4)-CaM, which relied on a single MD run. In contrast to radius of gyration, backbone mean-square fluctuations showed a distinctly non-normal and positively skewed distribution for nearly all residues. Furthermore, the D129N mutation affected the backbone dynamics in a complex manner and reduced the mobility of Glu123, Met124, Ile125, Arg126, and Glu127 located in the adjacent alpha-helix G. The implications of these observations for the comparisons of MD simulations with experiments are discussed. The proposed approach may be useful in studies of protein equilibrium dynamics where MD simulations fall short of properly sampling the conformational space, and when the comparison with experiments is affected by the reproducibility of the computational model.     

Dramatic differences in the motions of the mouth of open and closed cytochrome P450BM-3 by molecular dynamics simulations

Molecular dynamics trajectories were calculated separately for each of the two molecules in the asymmetric unit of the crystal structure of the hemoprotein domain of cytochrome P450BM-3. Each simulation was 200 ps in length and included a 10 Å layer of explicit solvent. The simulated time-average structure of each P450BM-3 molecule is closer to its crystal structure than the two molecular dynamics time-averaged structures are to each other. In the crystal structure, molecule 2 has a more accessible substrate binding pocket than molecule 1, and this difference is maintained throughout the simulations presented here. In particular, the substrate docking regions of molecule 1 and molecule 2 diverge in the solution state simulations. The mouth of the substrate binding pocket is significantly more mobile in the simulation of molecule 2 than in the simulation of molecule 1. For molecule 1, the width of the mouth is only slightly larger than its X-ray value of 8.7 Å and undergoes fluctuations of about 1 Å. However, in molecule 2, the mouth of the substrate binding pocket is dramatically more open in the time-average molecular dynamics structure (14.7 Å) than in the X-ray structure (10.9 Å). Furthermore, this region of the protein undergoes large amplitude motions during the trajectory that are not seen in the trajectory of molecule 1, repeatedly opening and closing up to 7 Å. Presumably, the binding of different substrates will induce the mouth region to adopt different conformations from within the wide range of structures that are accessible. © 1995 Wiley-Liss, Inc.     

Study of the Villin headpiece folding dynamics by combining coarse-grained Monte Carlo evolution and all-atom molecular dynamics

The folding mechanism of the Villin headpiece (HP36) is studied by means of a novel approach which entails an initial coarse-grained Monte Carlo (MC) scheme followed by all-atom molecular dynamics (MD) simulations in explicit solvent. The MC evolution occurs in a simplified free-energy landscape and allows an efficient selection of marginally-compact structures which are taken as viable initial conformations for the MD. The coarse-grained MC structural representation is connected to the one with atomic resolution through a "fine-graining" reconstruction algorithm. This two-stage strategy is used to select and follow the dynamics of seven different unrelated conformations of HP36. In a notable case the MD trajectory rapidly evolves towards the folded state, yielding a typical root-mean-square deviation (RMSD) of the core region of only 2.4 A from the closest NMR model (the typical RMSD over the whole structure being 4.0 A). The analysis of the various MC-MD trajectories provides valuable insight into the details of the folding and mis-folding mechanisms and particularly about the delicate influence of local and nonlocal interactions in steering the folding process.     

Weight average approaches for predicting dynamical properties of biomolecules

Recent advances in atomistic molecular dynamics (MD) simulations of biomolecules allow us to explore their conformational spaces widely, observing large-scale conformational fluctuations or transitions between distinct structures. To reproduce or refine experimental data using MD simulations, structure ensembles, which are characterized by multiple structures and their statistical weights on the rugged free-energy landscapes, are often used. Here, we summarize weight average approaches for various experimental measurements. Weight average approaches are now applied to hybrid quantum mechanics/molecular mechanics MD simulations to predict fast vibrational motions in a protein with a high accuracy for better understanding of molecular functions from atomic structures.     

Structure and fluctuations of bacteriophage T4 glutaredoxin modelled by molecular dynamics

A 120 ps molecular dynamics (MD) trajectory for bacteriophage T4 glutaredoxin was calculated including non-inertial solvent effects. The potential energy attains an equilibrated regime after the first 20 ps. The r.m.s. difference of all non-hydrogen atoms between X-ray and average MD structures for the regular secondary structure is 0.99A which shows that the MD simulation reproduces the essentials of the structure with high accuracy. Loop displacements are detected, shown by the larger full structure all non-hydrogen atom r.m.s. difference of 1.2A. The fluctuation pattern derived from MD agrees fairly well with that derived from X-ray isotropic temperature factors. The active site is a stable structural region in this MD modellization. Structural changes are put in context with the protein's function.