Extracting the Energy Barrier Distribution of a Disordered System from the Instantaneous Normal Mode Density of States: Applications to Peptides and Proteins

Experimental Confirmation of the Universal Law for the Vibrational Density of States of Liquids

An analytical model describing the vibrational density of states (VDOS) of liquids has long been elusive, owing to the complexities of liquid dynamics. Nevertheless, Zaccone and Baggioli have recently developed such a model which was proposed to be the universal law for the vibrational density of states of liquids. Distinct from the Debye law, g(ω) ∝ ω², for solids, the universal law for liquids reveals a linear relationship, g(ω) ∝ ω, in the low-energy region. We have confirmed this universal law with experimental VDOS measured by inelastic neutron scattering on real liquid systems including water, liquid metal, and polymer liquids, and have applied this model to extract the effective relaxation rate for the short time dynamics for each liquid. The model has also been further evaluated in the prediction of the specific heat with comparison to existing experimental data as well as with values obtained by different approaches.

Finite-Size Effects and Optimal System Sizes in Simulations of Surfactant Micelle Self-Assembly

The spontaneous formation of micelles in aqueous solutions is governed by the amphipathic nature of surfactants and is practically interesting due to the regular use of micelles as membrane mimics, for the characterization of protein structure, and for drug design and delivery. We performed a systematic characterization of the finite-size effect observed in single-component dodecylphosphocholine (DPC) micelles with the coarse-grained MARTINI model. Of multiple coarse-grained solvent models investigated using large system sizes, the nonpolarizable solvent model was found to most accurately reproduce SANS spectra of 100 mM DPC in aqueous solution. We systematically investigated the finite-size effect at constant 100 mM concentration in 23 systems of sizes 40-150 DPC, confirming the finite-size effect to manifest as an oscillation in the mean micelle aggregation number about the thermodynamic aggregation number as the system size increases. The oscillations in aggregation number mostly diminish once the system supports the formation of three micelles. Similar oscillations were observed in the estimated critical micelle concentration with a mean value of 1.10 mM, which is in agreement with experiment to 0.1 mM. The accuracy of using a multiscale simulation approach to avoid finite-size effects in the micelle size distribution and SANS spectra using MARTINI and CHARMM36 was explored using multiple long time scale 500 DPC coarse-grained simulations, which were back-mapped to CHARMM36 all-atom systems. It was found that the MARTINI model generally occupies more volume than the all-atom model, leading to the formation of micelles that are of a reasonable radius of gyration but are smaller in aggregation number. The systematic characterization of the finite-size effect and exploration of multiscale modeling presented in this work provide guidance for the accurate modeling of micelles in simulations.

Persistent Protein Motions in a Rugged Energy Landscape Revealed by Normal Mode Ensemble Analysis

Proteins are allosteric machines that couple motions at distinct, often distant, sites to control biological function. Low-frequency structural vibrations are a mechanism of this long-distance connection and are often used computationally to predict correlations, but experimentally identifying the vibrations associated with specific motions has proved challenging. Spectroscopy is an ideal tool to explore these excitations, but measurements have been largely unable to identify important frequency bands. The result is at odds with some previous calculations and raises the question what methods could successfully characterize protein structural vibrations. Here we show the lack of spectral structure arises in part from the variations in protein structure as the protein samples the energy landscape. However, by averaging over the energy landscape as sampled using an aggregate 18.5 μs of all-atom molecular dynamics simulation of hen egg white lysozyme and normal-mode analyses, we find vibrations with large overlap with functional displacements are surprisingly concentrated in narrow frequency bands. These bands are not apparent in either the ensemble averaged vibrational density of states or isotropic absorption. However, in the case of the ensemble averaged anisotropic absorption, there is persistent spectral structure and overlap between this structure and the functional displacement frequency bands. We systematically lay out heuristics for calculating the spectra robustly, including the need for statistical sampling of the protein and inclusion of adequate water in the spectral calculation. The results show the congested spectrum of these complex molecules obscures important frequency bands associated with function and reveal a method to overcome this congestion by combining structurally sensitive spectroscopy with robust normal mode ensemble analysis.

Instantaneous normal modes and the protein glass transition

In the instantaneous normal mode method, normal mode analysis is performed at instantaneous configurations of a condensed-phase system, leading to modes with negative eigenvalues. These negative modes provide a means of characterizing local anharmonicities of the potential energy surface. Here, we apply instantaneous normal mode to analyze temperature-dependent diffusive dynamics in molecular dynamics simulations of a small protein (a scorpion toxin). Those characteristics of the negative modes are determined that correlate with the dynamical (or glass) transition behavior of the protein, as manifested as an increase in the gradient with T of the average atomic mean-square displacement at approximately 220 K. The number of negative eigenvalues shows no transition with temperature. Further, although filtering the negative modes to retain only those with eigenvectors corresponding to double-well potentials does reveal a transition in the hydration water, again, no transition in the protein is seen. However, additional filtering of the protein double-well modes, so as to retain only those that, on energy minimization, escape to different regions of configurational space, finally leads to clear protein dynamical transition behavior. Partial minimization of instantaneous configurations is also found to remove nondiffusive imaginary modes. In summary, examination of the form of negative instantaneous normal modes is shown to furnish a physical picture of local diffusive dynamics accompanying the protein glass transition.

Molecular dynamics study on the solvent dependent heme cooling following ligand photolysis in carbonmonoxy myoglobin

The time scale and mechanism of vibrational energy relaxation of the heme moiety in myoglobin was studied using molecular dynamics simulation. Five different solvent models, including normal water, heavy water, normal glycerol, deuterated glycerol and a nonpolar solvent, and two forms of the heme, one native and one lacking acidic side chains, were studied. Structural alteration of the protein was observed in native myoglobin glycerol solution and native myoglobin water solution. The single-exponential decay of the excess kinetic energy of the heme following ligand photolysis was observed in all systems studied. The relaxation rate depends on the solvent used. However, this dependence cannot be explained using bulk transport properties of the solvent including macroscopic thermal diffusion. The rate and mechanism of heme cooling depends upon the detailed microscopic interaction between the heme and solvent. Three intermolecular energy transfer mechanisms were considered: (i) energy transfer mediated by hydrogen bonds, (ii) direct vibration-vibration energy transfer via resonant interaction, and (iii) energy transfer via vibration-translation or vibration-rotation interaction, or in other words, thermal collision. The hydrogen bond interaction and vibration-vibration interaction between the heme and solvent molecules dominates the energy transfer in native myoglobin aqueous solution and native myoglobin glycerol solutions. For modified myoglobin, the vibration-vibration interaction is also effective in glycerol solution, different from aqueous solution. Thermal collisions form the dominant energy transfer pathway for modified myoglobin in water solution, and for both native myoglobin and modified myoglobin in a nonpolar environment. For native myoglobin in a nonpolar solvent solution, hydrogen bonds between heme isopropionate side chains and nearby protein residues, absent in the modified myoglobin nonpolar solvent solution, are key interactions influencing the relaxation pathways.

Vibrational population relaxation of carbon monoxide in the heme pocket of photolyzed carbonmonoxy myoglobin: comparison of time-resolved mid-IR absorbance experiments and molecular dynamics simulations

The vibrational energy relaxation of carbon monoxide in the heme pocket of sperm whale myoglobin was studied by using molecular dynamics simulation and normal mode analysis methods. Molecular dynamics trajectories of solvated myoglobin were run at 300 K for both the delta- and epsilon-tautomers of the distal His-64. Vibrational population relaxation times of 335 +/- 115 ps for the delta-tautomer and 640 +/- 185 ps for the epsilon-tautomer were estimated by using the Landau-Teller model. Normal mode analysis was used to identify those protein residues that act as the primary "doorway" modes in the vibrational relaxation of the oscillator. Although the CO relaxation rates in both the epsilon- and delta-tautomers are similar in magnitude, the simulations predict that the vibrational relaxation of the CO is faster in the delta-tautomer with the distal His playing an important role in the energy relaxation mechanism. Time-resolved mid-IR absorbance measurements were performed on photolyzed carbonmonoxy hemoglobin (Hb(13)CO). From these measurements, a T(1) time of 600 +/- 150 ps was determined. The simulation and experimental estimates are compared and discussed.

Energy landscape of a native protein: jumping-among-minima model

We have investigated energy landscape of human lysozyme in its native state by using principal component analysis and a model, jumping-among-minima (JAM) model. These analyses are applied to 1 nsec molecular dynamics trajectory of the protein in water. An assumption embodied in the JAM model allows us to divide protein motions into intra-substate and inter-substate motions. By examining intra-substate motions, it is shown that energy surfaces of individual conformational substates are nearly harmonic and mutually similar. As a result of principal component analysis and JAM model analysis, protein motions are shown to consist of three types of collective modes, multiply hierarchical modes, singly hierarchical modes, and harmonic modes. Multiply hierarchical modes, the number of which accounts only for 0.5% of all modes, dominate contributions to total mean-square atomic fluctuation. Inter-substate motions are observed only in a small-dimensional subspace spanned by the axes of multiplyhierarchical and singly hierarchical modes. Inter-substate motions have two notable time components: faster component seen within 200 psec and slower component. The former involves transitions among the conformational substates of the low-level hierarchy, whereas the latter involves transitions of the higher level substates observed along the first four multiply hierarchical modes. We also discuss dependence of the subspace, which contains conformational substates, on time duration of simulation.

Ligand-induced conformational changes in ras p21: a normal mode and energy minimization analysis

A normal mode and energy minimization of ras p21 is used to determine the flexibility of the protein and the origin of the conformational differences between GTP and GDP-bound forms. To preserve the integrity of the structures, a hydration shell of water molecules was included as part of the system. Certain low-frequency modes were found to have high involvement coefficients with the conformational transition between the GTP and GDP-bound structures; the involvement coefficients of some of the modes increase when the gamma-phosphate group is removed. Two unstable modes that appear in the GTP-bound structure upon deletion of the gamma-phosphate group were determined and shown to have dominant contributions in the regions of switch I and switch II; there was also a significant displacement of loop 1. The initial motion in these regions is predicted by the modes to be approximately perpendicular to the direction of the transition from the GTP-bound state to the GDP-bound state. The overall conformational change in the switch I and II regions involves rearrangements of the protein backbone within these regions, rather than rigid body motion. Differences in the low-frequency modes of the GTP and GDP-bound forms appear to play a role in ligand binding. A coupling between the helix alpha3 position and the deletion of the gamma-phosphate group may be involved in the interaction with GAP. The oncogenic mutation G12D leads to a global increase in the rigidity of the protein. Thus, the mutant is likely to have a higher barrier for the conformational change to the inactive form; this would slow the transition and could be related to its oncogenic properties.