Dynamics-stability relationships in apo- and holomyoglobin: a combined neutron scattering and molecular dynamics simulations study

The removal of the heme group from myoglobin (Mb) results in a destabilization of the protein structure. The dynamic basis of the destabilization was followed by comparative measurements on holo- (holo-Mb) and apomyoglobin (apo-Mb). Mean-squared displacements (MSD) and protein resilience on the picosecond-to-nanosecond timescale were measured by elastic incoherent neutron scattering. Differences in thermodynamic parameters, MSD, and resilience were observed for both proteins. The resilience of holo-Mb was significantly lower than that of apo-Mb, indicating entropic stabilization by a higher degree of conformational sampling in the heme-bound folded protein. Molecular dynamics simulations provided site-specific information. Averaged over the whole structure, the molecular dynamics simulations yielded similar MSD and resilience values for the two proteins. The mobility of residues around the heme group in holo-Mb showed a smaller MSD and higher resilience compared to the same residue group in apo-Mb. It is of interest that in holo-Mb, higher MSD values are observed for the residues outside the heme pocket, indicating an entropic contribution to protein stabilization by heme binding, which is in agreement with experimental results.     

Molecular dynamics simulations of human DNA methyltransferase 3B with selective inhibitor nanaomycin A

DNA methyltransferases (DNMTs) are involved in epigenetic regulation of the genome and are promising targets for therapeutic intervention in cancer and other diseases. Until now, very limited information is available concerning the molecular dynamics of DNMTs. The natural product nanaomycin A is the first selective inhibitor of DNMT3B that induce genomic demethylation. Herein we report long (>100 ns) molecular dynamics simulations for human DNMT3B bound to nanaomycin A with and without the presence of the cofactor S-adenosyl-l-methionine (SAM). We concluded that SAM favors the binding of nanaomycin A to DNMT3B. Key interactions of nanaomycin A with DNMT3B involve long lasting interactions with Arg731, Arg733, Arg832, and the catalytic Cys651. Results further support the previous hypothesis that nanaomycin A has key interactions with amino acid residues involved in the mechanism of methylation. This work represents one of the first molecular dynamics studies of DNMT3B. Results of this work shed light on the structure and binding recognition process of a key epigenetic enzyme with a small molecule inhibitor.     

Distal Val346Ile mutation in inducible NO synthase promotes substrate-dependent NO confinement

The function of inducible NO synthase (WT iNOS) depends on the release of NO from the ferric heme before the enzyme is reduced. Key parameters controlling ligand dynamics include the distal and proximal heme pocket amino acids, as well as the inner solvent molecules. In this work, we tested how a point mutation in the distal heme side of WT iNOS affected the geminate rebinding of NO by ultrafast kinetics and molecular dynamics simulations. The mutation sequestered much of the photodissociated NO close to the heme compared to WT iNOS, with a main picosecond phase accounting for 78% of the rebinding to the arginine-bound Val346Ile protein. Consequently, the probability of NO release from Val346Ile decreased as compared to that from WT iNOS, provided the substrate binding site is filled. These data are rationalized by a steric effect of the Ile methyl group inducing events mediated by the substrate, transmitted via the propionates to the NO and the protein. This model is consistent with the role of the H-bonding network involving the heme, the substrate, and the BH4 cofactor in controlling NO release, with a key role of the heme propionates [Gautier et al. (2006) Nitric Oxide 15, 312]. These data support the effect of Val346Ile mutation in decreasing NO release and slowing down NO synthesis compared to WT iNOS determined by single turnover catalysis [Wang et al. (2004) J. Biol. Chem. 279, 19018].     

Homology modeling of a heme protein, lignin peroxidase, from the crystal structure of cytochrome c peroxidase.

A 3-dimensional model of lignin peroxidase (LiP) was constructed based on its sequence homology with other peroxidases, particularly cytochrome c peroxidase, the only protein with a known crystal structure in the peroxidase family. The construction of initial conformations of insertions and deletions was assisted by secondary structure predictions, amphipathic helix predictions, and consideration of the specific protein environment. A succession of molecular dynamics simulations of these regions with surrounding residues as constraints were carried out to relax the bond lengths and angles. Full protein molecular dynamics simulations with explicit consideration of bound waters were performed to relax the geometry and to identify dynamically flexible regions of the successive models for further refinement. Among the important functionally relevant structural features predicted are: (i) four disulfide bonds are predicted to be formed between Cys3 and Cys15, Cys14 and Cys285, Cys34 and Cys120 and Cys249 and Cys317; (ii) a glycosylation site, Asn257, was located on the surface; (iii) Glu40 was predicted to form a salt bridge with Arg43 on the distal side of the heme and was considered as a possible origin for the pH dependence of compound I formation; and (iv) two candidate substrate binding sites with a cluster of surface aromatic residues and flexible backbones were found in the refined model, consistent with the nature of known substrates of LiP. Based on these predicted structural features of the model, further theoretical and experimental studies are proposed to continue to elucidate the structure and function of LiP.     

Tethered-Hopping Model for Protein-DNA Binding and Unbinding Based on Sox2-Oct1-Hoxb1 Ternary Complex Simulations

The sliding and hopping models encapsulate the essential protein-DNA binding process for binary complex formation and dissociation. However, the effects of a cofactor protein on the protein-DNA binding process that leads to the formation of a ternary complex remain largely unknown. Here we investigate the effect of the cofactor Sox2 on the binding and unbinding of Oct1 with the Hoxb1 control element. We simulate the association of Oct1 with Sox2-Hoxb1 using molecular dynamics simulations, and the dissociation of Oct1 from Sox2-Hoxb1 using steered molecular dynamics simulations, in analogy to a hopping event of Oct1. We compare the kinetic and thermodynamic properties of three model complexes (the wild-type and two mutants) in which the Oct1-DNA base-specific interactions or the Sox2-Oct1 protein-protein interactions are largely abolished. We find that Oct1-DNA base-specific interactions contribute significantly to the total interaction energy of the ternary complex, and that nonspecific Oct1-DNA interactions are sufficient for driving the formation of the protein-DNA interface. The Sox2-Oct1 protein-protein binding interface is largely hydrophobic, with remarkable shape complementarity. This interface promotes the formation of the ternary complex and slows the dissociation of Oct1 from its DNA-binding site. We propose a simple two-step reaction model of protein-DNA binding, called the tethered-hopping model, that explains the importance of the cofactor Sox2 and may apply to similar ternary protein-DNA complexes.     

β‐strand interactions at the domain interface critical for the stability of human lens γD‐crystallin

Human age‐onset cataracts are believed to be caused by the aggregation of partially unfolded or covalently damaged lens crystallin proteins; however, the exact molecular mechanism remains largely unknown. We have used microseconds of molecular dynamics simulations with explicit solvent to investigate the unfolding process of human lens γD‐crystallin protein and its isolated domains. A partially unfolded folding intermediate of γD‐crystallin is detected in simulations with its C‐terminal domain (C‐td) folded and N‐terminal domain (N‐td) unstructured, in excellent agreement with biochemical experiments. Our simulations strongly indicate that the stability and the folding mechanism of the N‐td are regulated by the interdomain interactions, consistent with experimental observations. A hydrophobic folding core was identified within the C‐td that is comprised of a and b strands from the Greek key motif 4, the one near the domain interface. Detailed analyses reveal a surprising non‐native surface salt‐bridge between Glu135 and Arg142 located at the end of the ab folded hairpin turn playing a critical role in stabilizing the folding core. On the other hand, an in silico single E135A substitution that disrupts this non‐native Glu135‐Arg142 salt‐bridge causes significant destabilization to the folding core of the isolated C‐td, which, in turn, induces unfolding of the N‐td interface. These findings indicate that certain highly conserved charged residues, that is, Glu135 and Arg142, of γD‐crystallin are crucial for stabilizing its hydrophobic domain interface in native conformation, and disruption of charges on the γD‐crystallin surface might lead to unfolding and subsequent aggregation.     

Molecular dynamics simulations of the bacterial periplasmic heme binding proteins ShuT and PhuT

ShuT and PhuT are two periplasmic heme binding proteins that shuttle heme between the outer and inner membranes of the Gram-negative bacteria. Periplasmic binding proteins (PBPs) generally exhibit considerable conformational changes during the ligand binding process, whereas ShuT and PhuT belong to a class of PBPs that do not show such behavior based on their apo and holo crystal structures. By employing a series of molecular dynamic simulations on the ShuT and the PhuT, the dynamics and functions of the two PBPs were investigated. Through monitoring the distance changes between the two conserved glutamates of ShuT and PhuT, it was found the two PBPs were more flexible than previously assumed, exhibiting obvious opening-closing motions which were more remarkable in the apo runs of ShuT. Based on the results of the domain motion analysis, large scale conformational transitions were found in all apo runs of ShuT and PhuT, hinting that the domain motions of the two PBPs may be intrinsic. On the basis of the results of the principle component analysis, distinct opening-closing and twisting motion tendencies were observed not only in the apo, but also in the holo simulations of the two PBPs. The Gaussian network model was applied in order to analyze the hinge bending regions. The most important bending regions of ShuT and PhuT are located around the midpoints of their respective connecting helixes. Finally, the flexibilities and the details of the simulations of ShuT and PhuT were discussed. Characterized by the remarkably large flexibilities, the loop constituted by Ala 169, Gly170 and Gly171 of ShuT and the beta-turn constituted by Ala176, Gly177 and Gly178 of PhuT may be important for the functions of the two PBPs. Furthermore, the Asn254 of ShuT and the Arg228 of PhuT may be indispensable for the binding or unbinding of heme, since it is involved in the important hydrogen bonding to the propionate side-chains of heme.     

Recognition mechanism of Wilms’ tumour suppressor protein and DNA triplets: insights from molecular dynamics simulation and free energy analysis

The Wilms’ tumour suppressor protein (WT1) plays a multifaceted role in human cancer processes. Mutations on its DNA recognition domain could lead to Denys–Drash syndrome, and alternate splicing results in insertion of the tripeptide Lys–Thr–Ser (KTS) between the third and fourth zinc fingers (ZFs), leading to changes in the DNA-binding function. However, detailed recognition mechanisms of the WT1–DNA complex have not been explored. To clarify the mutational effects upon WT1 towards DNA binding at the atomic level, molecular dynamics simulations and the molecular mechanics/Poisson Boltzmann surface area (MM/PBSA) method were employed. The simulation results indicate that mutations in ZF domains (E427Q and Q369H) may weaken the binding affinity, and the statistical analyses of the hydrogen bonds and hydrophobic interactions show that eight residues (Lys351, Arg366, Arg375, Arg376, Lys399, Arg403, Arg424 and Arg430) have a significant influence on recognition and binding to DNA. Insertion of the tripeptide KTS could form an immobilized hydrogen-bonding network with Arg403, affecting the flexibility and angle of the linker between ZF3 and ZF4, thus influencing the recognition between the protein and the DNA triplet at its 5′ terminus. These results represent the first step towards a thorough characterization of the WT1 recognition mechanisms, providing a better understanding of the structure–function relationship of WT1 and its mutants.     

Conformational Coupling between Receptor and Kinase Binding Sites through a Conserved Salt Bridge in a Signaling Complex Scaffold Protein

Bacterial chemotaxis is one of the best studied signal transduction pathways. CheW is a scaffold protein that mediates the association of the chemoreceptors and the CheA kinase in a ternary signaling complex. The effects of replacing conserved Arg62 of CheW with other residues suggested that the scaffold protein plays a more complex role than simply binding its partner proteins. Although R62A CheW had essentially the same affinity for chemoreceptors and CheA, cells expressing the mutant protein are impaired in chemotaxis. Using a combination of molecular dynamics simulations (MD), NMR spectroscopy, and circular dichroism (CD), we addressed the role of Arg62. Here we show that Arg62 forms a salt bridge with another highly conserved residue, Glu38. Although this interaction is unimportant for overall protein stability, it is essential to maintain the correct alignment of the chemoreceptor and kinase binding sites of CheW. Computational and experimental data suggest that the role of the salt bridge in maintaining the alignment of the two partner binding sites is fundamental to the function of the signaling complex but not to its assembly. We conclude that a key feature of CheW is to maintain the specific geometry between the two interaction sites required for its function as a scaffold.     

A relationship between heme binding and protein stability in cytochrome b5

Conformational changes and long-range effects are often observed in proteins when they associate with their ligands. In many cases, these structural perturbations are essential to function, and they are the result of complex networks of interactions. Here we used cytochrome b(5), a protein that undergoes extensive structural rearrangement upon heme binding, to seek a relationship between affinity for the cofactor and extent of refolding induced by its binding. Three variants of the water-soluble domain of the rat microsomal protein were chosen to affect the stability of the apoprotein or the holoprotein. Sequence alterations were introduced in the heme binding loop (type I mutations, D60R and (55)TENFED --> (55)TEPFEED, or PE), which is largely unstructured in the apoprotein state, and in the folded core of the apoprotein (type II mutation, P81A). Thermal and chemical denaturation experiments and heme transfer experiments were performed on these proteins. Type I mutations left the thermodynamic stability of the apoprotein unchanged. The first mutation (D60R) stabilized the holoprotein in a probable manifestation of enhanced helical propensity or improved electrostatic interactions. The second mutation (PE) decreased heme affinity and holoprotein stability in concert. For this protein, heme transfer experiments could be used to estimate the rate constant of heme loss from each of the heme orientational isomers. In contrast, the type II mutation resulted in a marked destabilization of the apoprotein but an intermediate effect on the holoprotein stability and heme affinity. These data supported that heme affinity could be modulated by the apoprotein stability and by specific residues remote from the heme binding site.     

Molecular dynamics simulations of early steps in RNA‐mediated conversion of prions

The rate‐limiting step in prion diseases is the initial transition of a prion protein from its native form into a mis‐folded state in which the protein not only forms cell‐toxic aggregates but also becomes infectious. Recent experiments implicate polyadenosine RNA as a possible agent for generating the initial seed. In order to understand the mechanism of RNA‐mediated mis‐folding and aggregation of prions, we dock polyadenosine RNA to mouse and human prion models. Changes in stability and secondary structure of the prions upon binding to polyadenosine RNA are evaluated by comparing molecular dynamics simulations of these complexes with that of the unbound prions.     

Energy Landscapes of Protein Aggregation and Conformation Switching in Intrinsically Disordered Proteins

The protein folding problem was apparently solved recently by the advent of a deep learning method for protein structure prediction called AlphaFold. However, this program is not able to make predictions about the protein folding pathways. Moreover, it only treats about half of the human proteome, as the remaining proteins are intrinsically disordered or contain disordered regions. By definition these proteins differ from natively folded proteins and do not adopt a properly folded structure in solution. However these intrinsically disordered proteins (IDPs) also systematically differ in amino acid composition and uniquely often become folded upon binding to an interaction partner. These factors preclude solving IDP structures by current machine-learning methods like AlphaFold, which also cannot solve the protein aggregation problem, since this meta-folding process can give rise to different aggregate sizes and structures. An alternative computational method is provided by molecular dynamics simulations that already successfully explored the energy landscapes of IDP conformational switching and protein aggregation in multiple cases. These energy landscapes are very different from those of ‘simple’ protein folding, where one energy funnel leads to a unique protein structure. Instead, the energy landscapes of IDP conformational switching and protein aggregation feature a number of minima for different competing low-energy structures. In this review, I discuss the characteristics of these multifunneled energy landscapes in detail, illustrated by molecular dynamics simulations that elucidated the underlying conformational transitions and aggregation processes.     

Binding interactions between the core central domain of 16S rRNA and the ribosomal protein S15 determined by molecular dynamics simulations

The goal of the current study is to utilize molecular dynamic (MD) simulations to investigate the dynamic behavior of 16S rRNA in the presence and absence of S15 and to identify the binding interactions between these two molecules. The simulations show that: (i) 16S rRNA remains in a highly folded structure when it is bound to S15; (ii) in the absence of S15, 16S rRNA significantly alters its conformation and transiently forms conformations that are similar to the bound structure that make it available for binding with S15; (iii) the unbound rRNA spends the majority of its time in extended conformations. The formation of the extended conformations is a result of the molecule reaching a lower electrostatic energy and the formation of the highly folded, crystal-like conformation is a result of achieving a lower solvation energy. In addition, our MD simulations show that 16S rRNA and S15 bind across the major groove of helix 22 (H22) via electrostatic interactions. The negatively charged phosphate groups of G658, U740, G741 and G742 bind to the positively charged S15 residues Lys7, Arg34 and Arg37. The current study provides a dynamic view of the binding of 16S rRNA with S15.     

Dynamics of heme complexed with human serum albumin: a theoretical approach

Human serum albumin (HSA) is the most abundant protein in the blood serum. It binds several ligands and has an especially strong affinity for heme, hence becoming a natural candidate for oxygen transport. In order to analyze the interaction of HSA-heme, molecular dynamics simulations of HSA with bound heme were performed. Based on the results of X-ray diffraction, the binding site of the heme, localized in subdomain IB, was considered. We analyzed the fluctuations and their correlations along trajectories to detect collective motions. The role of H bonds and salt bridges in the stabilization of heme in its pocket was also investigated. Complementarily, the localization of water molecules in the hydrophobic pocket and the interaction with heme were discussed.     

Protein modeling of human prorenin using the molecular dynamics method

To study the activation-inactivation mechanism of the renin zymogen, prorenin, a tertiary structural model of human prorenin was constructed using computer graphics and molecular dynamics calculations, based on the pepsinogen structure. This prorenin model shows that the folded prosegment polypeptide can fit into the substrate binding cleft of the renin moiety. The three positively charged residues, Arg 10, Arg 15, and Arg 20, in the prosegment make salt bridges with Asp 225, Glu 331, and Asp 60, respectively, in renin. Arg 43, which is in the processing site, forms salt bridges with the catalytic residues of Asp 81 and Asp 269. These ionic interactions between the prosegment and the renin may contribute to keeping the prorenin structure as an inactive form.     

Relaxation kinetics and the glassiness of native proteins: coupling of timescales

We provide evidence that the onset of functional dynamics of folded proteins with elevated temperatures is associated with the effective sampling of its energy landscape under physiological conditions. The analysis is based on data describing the relaxation phenomena governing the backbone dynamics of bovine pancreatic trypsin inhibitor derived from molecular dynamics simulations, previously reported by us. By representing the backbone dynamics of the folded protein by three distinct regimes, it is possible to decompose its seemingly complex dynamics, described by a stretch exponential decay of the backbone motions. Of these three regimes, one is associated with the slow timescales due to the activity along the envelope of the energy surface defining the folded protein. Another, with fast timescales, is due to the activity along the pockets decorating the folded-state envelope. The intermediate regime emerges at temperatures where jumps between the pockets become possible. It is at the temperature window where motions corresponding to all three timescales become operative that the protein becomes active.     

Normal mode analysis of the horseradish peroxidase collective motions: correlation with spectroscopically observed heme distortions

Horseradish peroxidase C is a class III peroxidase whose structure is stabilized by the presence of two endogenous calcium atoms. Calcium removal has been shown to decrease the enzymatic activity of the enzyme and significantly affect the spectroscopically detectable properties of the heme, such as the spin state of the iron, heme normal modes, and distortions from planarity. In this work, we report on normal mode analysis (NMA) performed on models subjected to 2 ns of molecular dynamics simulations to describe the effect of calcium removal on protein collective motions and to investigate the correlation between active site (heme) and protein matrix fluctuations. We show that in the native peroxidase model, heme fluctuations are correlated to matrix fluctuations while they are not in the calcium-depleted model.     

Identifying the protein folding nucleus using molecular dynamics

Molecular dynamics simulations of folding in an off-lattice protein model reveal a nucleation scenario, in which a few well-defined contacts are formed with high probability in the transition state ensemble of conformations. Their appearance determines folding cooperativity and drives the model protein into its folded conformation. Amino acid residues participating in those contacts may serve as "accelerator pedals" used by molecular evolution to control protein folding rate.     

Molecular dynamics simulations and functional characterization of the interactions of the PAR2 ectodomain with factor VIIa

Signaling of the tissue factor‐FVIIa complex regulates angiogenesis, tumor growth, and inflammation. TF‐FVIIa triggers cell signaling events by cleavage of protease activated receptor (PAR2) at the Arg36‐Ser37 scissile bond. The recognition of PAR2 by the FVIIa protease domain is poorly understood. We perform molecular modeling and dynamics simulations to derive the PAR2‐FVIIa interactions. Docking of the PAR2 Arg36‐Ser37 scissile bond to the S1 site and subsequent molecular dynamics leads to interactions of the PAR2 ectodomain with P and P′ sites of the FVIIa catalytic cleft as well as to electrostatic interactions between a stably folded region of PAR2 and a cluster of basic residues remote from the catalytic cleft of FVIIa. To address the functional significance of this interaction for PAR2 cleavage, we employed two antibodies with epitopes previously mapped to this cluster of basic residues. Although these antibodies do not block the catalytic cleft, both antibodies completely abrogated PAR2 activation by TF‐FVIIa. Our simulations indicate a conformation of the PAR2 ectodomain that limits the cleavage site to no more than 33 Å from its membrane proximal residue. Since the active site of FVIIa in the TF‐FVIIa complex is ～75 Å above the membrane, cleavage of the folded conformation of PAR2 would require tilting of the TF‐FVIIa complex toward the membrane, indicating that additional cellular factors may be required to properly align the scissile bond of PAR2 with TF‐FVIIa. Proteins 2009. © 2009 Wiley‐Liss, Inc.     

Molecular dynamics simulations of human LRH-1: the impact of ligand binding in a constitutively active nuclear receptor

The liver receptor homologue 1 (LRH-1 (NR5A2)) belongs to the orphan nuclear receptor family, indicating that initially no ligand was known. Although recent studies have shown that ligand binding can be obtained, the biological relevance remains elusive. Here, we modify the observed X-ray ligand into a biologically more significant phospholipid (phosphatidylserine, PS) present in human, to study, by molecular dynamics (MD) simulations, the impact of the ligand on the receptor and the interaction with different cofactor peptides. Furthermore, we characterize the interactions between receptor and the cofactor peptides of DAX-1 (NR0B1), Prox1 and SHP LXXLL box 1 and 2 (NR0B2) in terms of specificity. Our MD simulation results show different interaction patterns for the SHP box2 compared to DAX-1, PROX1 and SHP box1. SHP box2 shows specific interactions at its more C-terminal end while the other investigated peptides show specific interactions at several positions but particularly at the +2 site. The peptide +2 side chain interacts with a charged amino acid of the receptor, in hLRH-1 Asp372. Together with the charge clamp residues Arg361 and Glu534, Asp372 forms a triangle shaped charge clamp responsible for peptide orientation and increased affinity. The binding of the PS ligand causes no overall structural changes of the receptor but affects the interactions with cofactor peptides. The cofactor peptides from SHP decrease its interaction with the receptor upon ligand binding while DAX-1 and PROX1 are unchanged or increase. The diverse ligand binding response of the cofactor provides an opportunity for drug design with the possibility to create agonist ligands to modify cofactor interaction.