Molecular dynamics simulations of proteins in lipid bilayers

With recent advances in X-ray crystallography of membrane proteins promising many new high-resolution structures, molecular dynamics simulations will become increasingly valuable for understanding membrane protein function, as they can reveal the dynamic behavior concealed in the static structures. Dramatic increases in computational power, in synergy with more efficient computational methodologies, now allow us to carry out molecular dynamics simulations of any structurally known membrane protein in its native environment, covering timescales of up to 0.1 micros. At the frontiers of membrane protein simulations are ion channels, aquaporins, passive and active transporters, and bioenergetic proteins.     

Transmembrane helix-helix interactions: comparative simulations of the glycophorin a dimer

The glycophorin helix dimer is a paradigm for the exploration of helix-helix interactions in integral membrane proteins. Two NMR structures of the dimer are known, one in a detergent micelle and one in a lipid bilayer. Multiple (4 x 50 ns) molecular dynamics simulations starting from each of the two NMR structures, with each structure in either a dodecyl phosphocholine (DPC) micelle or a dimyristoyl phosphatidylcholine (DMPC) bilayer, have been used to explore the conformational dynamics of the helix dimer. Analysis of the helix-helix interaction, mediated by the GxxxG sequence motif, suggests convergence of the simulations to a common model. This is closer to the NMR structure determined in a bilayer than to micelle structure. The stable dimer interface in the final simulation model is characterized by (i) Gly/Gly packing and (ii) Thr/Thr interhelix H-bonds. These results demonstrate the ability of extended molecular dynamics simulations in a lipid bilayer environment to refine membrane protein structures or models derived from experimental data obtained in protein/detergent micelles.     

Interaction of anesthetics with open and closed conformations of a potassium channel studied via molecular dynamics and normal mode analysis

A variety of experiments suggest that membrane proteins are important targets of anesthetic molecules, and that ion channels interact differently with anesthetics in their open and closed conformations. The availability of an open and a closed structural model for the KirBac1.1 potassium channel has made it possible to perform a comparative analysis of the interactions of anesthetics with the same channel in its open and closed states. To this end, all-atom molecular dynamics simulations supplemented by normal mode analysis have been employed to probe the interactions of the inhalational anesthetic halothane with both an open and closed conformer of KirBac1.1 embedded in a lipid bilayer. Normal mode analysis on the closed and open channel, in the presence and absence of halothane, reveals that the anesthetic modulates the global as well as the local dynamics of both conformations differently. In the case of the open channel, the observed reduction of flexibility of residues in the inner helices suggests a functional modification action of anesthetics on ion channels. In this context, preferential quenching of the aromatic residue motion and modulation of global dynamics by halothane may be seen as steps toward potentiating or favoring open state conformations. These molecular dynamics simulations provide the first insights into possible specific interactions between anesthetic molecules and ion channels in different conformations.     

High-Performance Drug Discovery: Computational Screening by Combining Docking and Molecular Dynamics Simulations

Virtual compound screening using molecular docking is widely used in the discovery of new lead compounds for drug design. However, this method is not completely reliable and therefore unsatisfactory. In this study, we used massive molecular dynamics simulations of protein-ligand conformations obtained by molecular docking in order to improve the enrichment performance of molecular docking. Our screening approach employed the molecular mechanics/Poisson-Boltzmann and surface area method to estimate the binding free energies. For the top-ranking 1,000 compounds obtained by docking to a target protein, approximately 6,000 molecular dynamics simulations were performed using multiple docking poses in about a week. As a result, the enrichment performance of the top 100 compounds by our approach was improved by 1.6–4.0 times that of the enrichment performance of molecular dockings. This result indicates that the application of molecular dynamics simulations to virtual screening for lead discovery is both effective and practical. However, further optimization of the computational protocols is required for screening various target proteins.     

Ultrafast and downhill protein folding

Ultrafast folding proteins have served an important role in benchmarking molecular dynamics simulations and testing protein folding theories. These proteins are simple enough and fold fast enough that realistic simulations are possible, which facilitates the direct comparison of absolute folding rates and folding mechanisms with those observed experimentally. Such comparisons have achieved remarkable success, but have also revealed the shortcomings that remain in experiment, theory and simulation alike. Some ultrafast folding proteins may fold without encountering an activation barrier (downhill folding), allowing the exploration of the molecular timescale of folding and the roughness of the energy landscape. The biological significance of ultrafast folding remains uncertain, but its practical significance is crucial to progress in understanding how proteins fold.     

Collective motions in proteins: a covariance analysis of atomic fluctuations in molecular dynamics and normal mode simulations.

A method is described for identifying collective motions in proteins from molecular dynamics trajectories or normal mode simulations. The method makes use of the covariances of atomic positional fluctuations. It is illustrated by an analysis of the bovine pancreatic trypsin inhibitor. Comparison of the covariance and cross-correlation matrices shows that the relative motions have many similar features in the different simulations. Many regions of the protein, especially regions of secondary structure, move in a correlated manner. Anharmonic effects, which are included in the molecular dynamics simulations but not in the normal analysis, are of some importance in determining the larger scale collective motions, but not the more local fluctuations. Comparisons of molecular dynamics simulations in the present and absence of solvent indicate that the environment is of significance for the long-range motions.     

High pressure simulations of biomolecules

Pressure is a thermodynamic variable which is particularly suitable for exploration of the properties of biological macromolecules. For proteins, in particular, denaturation induced by pressure is different from that induced by temperature or denaturants. The response of proteins to pressure changes can provide information on properties of their native and non-native states. This review focuses on molecular dynamics studies of the effect of pressure on detailed atomic models of proteins. It also reports on other theoretical approaches, such as Monte Carlo simulations, which have been used to study simplified models. Another purpose of this review is to try to point out potential future studies that may be both interesting and feasible, with constantly increasing computing power.     

The MUMO (minimal under-restraining minimal over-restraining) method for the determination of native state ensembles of proteins

While reliable procedures for determining the conformations of proteins are available, methods for generating ensembles of structures that also reflect their flexibility are much less well established. Here we present a systematic assessment of the ability of ensemble-averaged molecular dynamics simulations with ensemble-averaged NMR restraints to simultaneously reproduce the average structure of proteins and their associated dynamics. We discuss the effects that under-restraining (overfitting) and over-restraining (underfitting) have on the structures generated in ensemble-averaged molecular simulations. We then introduce the MUMO (minimal under-restraining minimal over-restraining) method, a procedure in which different observables are averaged over a different number of molecules. As both over-restraining and under-restraining are significantly reduced in the MUMO method, it is possible to generate ensembles of conformations that accurately characterize both the structure and the dynamics of native states of proteins. The application of the MUMO method to the protein ubiquitin yields a high-resolution structural ensemble with an RDC Q-factor of 0.19.     

Molecular dynamics simulations of membrane proteins

Membrane proteins control the traffic across cell membranes and thereby play an essential role in cell function from transport of various solutes to immune response via molecular recognition. Because it is very difficult to determine the structures of membrane proteins experimentally, computational methods have been increasingly used to study their structure and function. Here we focus on two classes of membrane proteins—ion channels and transporters—which are responsible for the generation of action potentials in nerves, muscles, and other excitable cells. We describe how computational methods have been used to construct models for these proteins and to study the transport mechanism. The main computational tool is the molecular dynamics (MD) simulation, which can be used for everything from refinement of protein structures to free energy calculations of transport processes. We illustrate with specific examples from gramicidin and potassium channels and aspartate transporters how the function of these membrane proteins can be investigated using MD simulations.     

Predicting large-scale conformational changes in proteins using energy-weighted normal modes

We report the development of a method to improve the sampling of protein conformational space in molecular simulations. It is shown that a principal component analysis of energy-weighted normal modes in Cartesian coordinates can be used to extract vectors suitable for describing the dynamics of protein substructures. The method can operate with either atomistic or user-defined coarse-grained models of protein structure. An implicit reverse coarse-graining allows the dynamics of all-atoms to be recovered when a coarse-grained model is used. For an external test set of four proteins, it is shown that the new method is more successful than normal mode analysis in describing the large-scale conformational changes observed on ligand binding. The method has potential applications in protein-ligand and protein-protein docking and in biasing molecular dynamics simulations.     

Projection of Monte Carlo and molecular dynamics trajectories onto the normal mode axes: human lysozyme

A method is presented to describe the internal motions of proteins obtained from molecular dynamics or Monte Carlo simulations as motions of normal mode variables. This method calculates normal mode variables by projecting trajectories of these simulations onto the axes of normal modes and expresses the trajectories as a linear combination of normal mode variables. This method is applied to the result of the molecular dynamics and the Monte Carlo simulations of human lysozyme. The motion of the lowest frequency mode extracted from the simulations represents the hinge bending motion very faithfully. Analysis of the obtained motions of the normal mode variables provides an explanation of the anharmonic aspects of protein dynamics as due first to the anharmonicity of the actual potential energy surface near a minimum and second to trans-minimum conformational changes.     

Markov state modeling of membrane transport proteins

The flux of ions and molecules in and out of the cell is vital for maintaining the basis of various biological processes. The permeation of substrates across the cellular membrane is mediated through the function of specialized integral membrane proteins commonly known as membrane transporters. These proteins undergo a series of structural rearrangements that allow a primary substrate binding site to be accessed from either side of the membrane at a given time. Structural insights provided by experimentally resolved structures of membrane transporters have aided in the biophysical characterization of these important molecular drug targets. However, characterizing the transitions between conformational states remains challenging to achieve both experimentally and computationally. Though molecular dynamics simulations are a powerful approach to provide atomistic resolution of protein dynamics, a recurring challenge is its ability to efficiently obtain relevant timescales of large conformational transitions as exhibited in transporters. One approach to overcome this difficulty is to adaptively guide the simulation to favor exploration of the conformational landscape, otherwise known as adaptive sampling. Furthermore, such sampling is greatly benefited by the statistical analysis of Markov state models. Historically, the use of Markov state models has been effective in quantifying slow dynamics or long timescale behaviors such as protein folding. Here, we review recent implementations of adaptive sampling and Markov state models to not only address current limitations of molecular dynamics simulations, but to also highlight how Markov state modeling can be applied to investigate the structure–function mechanisms of large, complex membrane transporters.     

Comparison of molecular dynamics simulation methods for amyloid β(1-42) monomers containing D-aspartic acid residues for predicting retention times in chromatography

Molecular dynamics simulations of amyloid β(1-42) containing D-aspartic acid residues were performed using several continuous solvent models to investigate the usefulness of simulation methods for D-amino acid-containing proteins and peptides. Normal molecular dynamics simulations and replica exchange molecular dynamics simulations, which are one of the generalized-ensemble algorithms, were performed. Because the β-structure contents of amyloid β(1-42) peptides obtained by replica exchange molecular dynamics simulations with Onufriev-Bashford-Case generalized Born implicit solvent were qualitatively consistent with experimental data, replica exchange molecular dynamics rather than other methods appeared to be more reasonable for calculations of amyloid β(1-42) containing D-aspartic acid residues. Computational results revealed that peptides with stereoinversion of Asp23 tend to form β-sheet structures by themselves, in contrast to the wild-type peptides that form β-sheet structures only after aggregation. These results are expected to be useful for computational investigations of proteins and peptides such as prediction of retention time of peptides and proteins containing D-aspartic acid residues.     

Pressure-induced water transport in membrane channels studied by molecular dynamics

A method is proposed to measure the water permeability of membrane channels by means of molecular dynamics simulations. By applying a constant force to the bulk water molecules and a counter force on the complementary system, a hydrostatic pressure difference across the membrane can be established, producing a net directional water flow. The hydraulic or osmotic permeability can then be determined by the ratio of the water flux and the pressure difference. The method is applied and tested on an aquaglyceroporin channel through a series of simulations totaling 5 ns in duration.     

Quantification and visualization of molecular surface flexibility

Two new methods for the quantification and visualization of the flexibility of molecular surfaces are presented. Both methods rely on results of molecular dynamics (MD) simulations. Whereas method I is based on a simple but fast grid-counting algorithm, method II uses a mapping function that allows for a sharp and clear visualization of atomic RMS fluctuations on a molecular surface. To demonstrate the scope of the methods, MD simulations of two proteins, PTI and ubiquitin, were performed. The flexibility data are mapped onto the molecular surfaces of the proteins and visualized using texture mapping technology available on modern workstations.     

Molecular dynamics simulations of a K+ channel blocker: Tc1 toxin from Tityus cambridgei

Toxins that block voltage-gated potassium (Kv) channels provide a possible template for improved homology models of the Kv pore. In assessing the interactions of Kv channels and their toxins it is important to determine the dynamic flexibility of the toxins. Multiple 10 ns duration molecular dynamics simulations combined with essential dynamics analysis have been used to explore the flexibility of four different Kv channel-blocking toxins. Three toxins (Tc1, AgTx and ChTx) share a common fold. They also share a common pattern of conformational dynamics, as revealed by essential dynamics analysis of the simulation results. This suggests that some aspects of dynamic behaviour are conserved across a single protein fold class. In each of these three toxins, the residue exhibiting minimum flexibility corresponds to a conserved lysine residue that is suggested to interact with the filter domain of the channel. Thus, comparative simulations reveal functionally important conservation of molecular dynamics as well as protein fold across a family of related toxins.     

Effective energy function for proteins in solution

A Gaussian solvent-exclusion model for the solvation free energy is developed. It is based on theoretical considerations and parametrized with experimental data. When combined with the CHARMM 19 polar hydrogen energy function, it provides an effective energy function (EEF1) for proteins in solution. The solvation model assumes that the solvation free energy of a protein molecule is a sum of group contributions, which are determined from values for small model compounds. For charged groups, the self-energy contribution is accounted for primarily by the exclusion model. Ionic side-chains are neutralized, and a distance-dependent dielectric constant is used to approximate the charge-charge interactions in solution. The resulting EEF1 is subjected to a number of tests. Molecular dynamics simulations at room temperature of several proteins in their native conformation are performed, and stable trajectories are obtained. The deviations from the experimental structures are similar to those observed in explicit water simulations. The calculated enthalpy of unfolding of a polyalanine helix is found to be in good agreement with experimental data. Results reported elsewhere show that EEF1 clearly distinguishes correctly from incorrectly folded proteins, both in static energy evaluations and in molecular dynamics simulations and that unfolding pathways obtained by high-temperature molecular dynamics simulations agree with those obtained by explicit water simulations. Thus, this energy function appears to provide a realistic first approximation to the effective energy hypersurface of proteins.     

Conformational landscapes of membrane proteins delineated by enhanced sampling molecular dynamics simulations

The expansion of computational power, better parameterization of force fields, and the development of novel algorithms to enhance the sampling of the free energy landscapes of proteins have allowed molecular dynamics (MD) simulations to become an indispensable tool to understand the function of biomolecules. The temporal and spatial resolution of MD simulations allows for the study of a vast number of processes of interest. Here, we review the computational efforts to uncover the conformational free energy landscapes of a subset of membrane proteins: ion channels, transporters and G-protein coupled receptors. We focus on the various enhanced sampling techniques used to study these questions, how the conclusions come together to build a coherent picture, and the relationship between simulation outcomes and experimental observables.