Brownian dynamics simulations of the recognition of the scorpion toxin P05 with the small-conductance calcium-activated potassium channels

The recognition of the scorpion toxin P05 and the small-conductance, calcium-activated potassium (SK) channels, rsk1, rsk2, and rsk3, has been studied by means of the Brownian dynamics (BD) method. All of the 25 available structures of P05 in the RCSB Protein Data Bank determined by NMR were considered during the simulation, which indicated that the conformation of P05 affects both the recognition and binding between the two proteins significantly. Comparing the top four high-frequency structures of P05 binding to the SK channels, we found that the rsk2 channel, with high frequencies and lowest electrostatic interaction energies (E (int)(ES)), is the most favorable for P05 binding, while rsk3 is intermediate, and rsk1 is the least favorable. Among the 25 structures of P05, the 13th structure docks into the binding site of the rsk2 channel with the highest probability and most favorable electrostatic interactions. From the P05-rsk2 channel binding model, we identified the residues critical for the recognition of these two proteins through triplet contact analyses. P05 locates around the extracellular mouth of the SK channels and contacts the SK channels using its alpha-helix rather than beta-sheets. The critical triplet contacts for recognition between P05 and the rsk2 channel are Arg6 (P05)-Asp364 (SK), Arg7 (P05)-Asn368 (SK), and Arg13 (P05)-Asp341 (SK). The structure of the P05-rsk2 complex with the most favorable electrostatic interaction energy was further refined by molecular mechanics, showing that six hydrogen bonding interactions exist between P05 and the rsk2 channel: one hydrogen bond is formed between Arg6 (P05) and Asp364(D) (rsk2); Arg7 (P05) forms three hydrogen bonds with Asp341(B) (rsk2)) and Asp364(C) (rsk2); two hydrogen bonds are formed by Arg13 (P05) with Asp341(A) (rsk2) and Asp364(B) (rsk2). The simulation results are in good agreement with the previous molecular biological experiments and can explain the binding phenomena between P05 and SK channels at the level of molecular structure. The consistency between the results of the BD simulations and the experimental data indicated that our 3D model of the P05-rsk2 channel complex is reasonable and can be employed in further biological studies, such as rational design of the novel therapeutic agents blocking the small-conductance, calcium-activated and apamin-sensitive potassium channels, and for mutagenesis studies in both toxins and SK channels. In particular, both the BD simulations and the molecular mechanics refinements indicate that residue Asp364 of the rsk2 channel is critical for its recognition and binding functionality towards P05. This phenomenon has not been appreciated in the previous mutagenesis experiments, indicating that this might be a new clue for further functional study of SK channels.     

Membrane partitioning of peptide aggregates: coarse-grained molecular dynamics simulations

Abstract

Coarse-grained molecular dynamics (CGMD) simulation technique (MARTINI force field) is applied to monitor the aggregation of helical peptides representing the transmembrane sequence and its extension of bone marrow stromal cell antigen 2 (BST-2). One of the peptides is coupled with a protein transducing domain (PTD) of nine arginine residues (R9) at its N-terminal side as well as a peptide, pep11**, which has been shown to bind to human papilloma virus 16 (HPV16) E6 oncoprotein. A short hydrophobic stretch of the transmembrane domain (TMD) of BST-2 aggregates the fastest and inserts into a lipid membrane. An aggregate of R9-pep11** attaches to the membrane via simultaneous contact of many arginine residues. Monomers from the aggregates of the shortest of the hydrophobic TMDs dissolve into the opposing leaflet when the aggregate spans the bilayer. A ‘flipping’ of the individual monomeric peptides is not observed.

Communicated by Ramaswamy H. Sarma     

Brownian dynamics simulations of intramolecular energy transfer

A novel technique for modelling intramolecular energy transfer is presented. Brownian dynamics calculations are used to compute the trajectories of donor and acceptor species, and the instantaneous orientation factor is calculated during each temporal iteration. In this work, several model systems are considered. Trajectories were computed for energy transfer between a flexible donor and a rigidly fixed acceptor. We have considered configurations where the donor is, (1) tethered to a fixed point in space, but free to diffuse rotationally, and (2) constrained to wobble in a cone. The luminescence decay of the donor is ‘measured’, and a non-single-exponential decay is observed for configurations of efficient energy transfer. Luminescence anisotropy measurements of constrained and unconstrained donors reflect the contribution of both energy transfer and rotational diffusion to the shape of the anisotropy decay curve.     

Dissecting the energetics of an antibody-antigen interface by alanine shaving and molecular grafting.

Alanine-scanning mutagenesis on human growth hormone (hGH) identified 5 primary determinants (Arg 8, Asn 12, Arg 16, Asp 112, and Asp 116) for binding to a monoclonal antibody (MAb 3) (Jin L, Fendly BM, Wells JA, 1992, J Mol Biol 226:851-865). To further analyze the energetic importance of residues surrounding these five, we mutated all neighboring residues to alanine in groups of 7-16 (a procedure we call alanine shaving). Even the most extremely mutated variant, with 16 alanine substitutions, caused less than a 10-fold reduction in binding affinity to MAb3. By comparison, mutating any 1 of the 5 primary determinants to alanine caused a 6- to > 500-fold reduction in affinity. Replacing any of the 4 charged residues (Arg 8, Arg 16, Asp 112, and Asp 116) with a homologous residue (i.e., Arg to Lys or Asp to Glu) caused nearly as large a reduction in affinity as the corresponding alanine replacement. It was possible to graft the 5 primary binding determinants onto a nonbinding homologue of hGH, human placental lactogen (hPL), which has 86% sequence identity to hGH. The grafted hPL mutant bound 10-fold less tightly than hGH to MAb3 but bound as well as hGH when 2 additional framework mutations were introduced. Attempts to recover binding affinity by grafting the MAb3 epitope onto more distantly related scaffolds having a similar 4-helix bundle motif, such as human prolactin (23% sequence identity) or granulocyte colony-stimulating factor, were unsuccessful.(ABSTRACT TRUNCATED AT 250 WORDS)     

Rapid binding of a cationic active site inhibitor to wild type and mutant mouse acetylcholinesterase: Brownian dynamics simulation including diffusion in the active site gorge

It is known that anionic surface residues play a role in the long-range electrostatic attraction between acetylcholinesterase and cationic ligands. In our current investigation, we show that anionic residues also play an important role in the behavior of the ligand within the active site gorge of acetylcholinesterase. Negatively charged residues near the gorge opening not only attract positively charged ligands from solution to the enzyme, but can also restrict the motion of the ligand once it is inside of the gorge. We use Brownian dynamics techniques to calculate the rate constant k_on for wild type and mutant acetylcholinesterase with a positively charged ligand. These calculations are performed by allowing the ligand to diffuse within the active site gorge. This is an extension of previously reported work in which a ligand was allowed to diffuse only to the enzyme surface. By setting the reaction criteria for the ligand closer to the active site, better agreement with experimental data is obtained. Although a number of residues influence the movement of the ligand within the gorge, Asp74 is shown to play a particularly important role in this function. Asp74 traps the ligand within the gorge, and in this way helps to ensure a reaction. © 1998 John Wiley & Sons, Inc. Biopoly 46: 465–474, 1998     

Deciphering the GPER/GPR30-agonist and antagonists interactions using molecular modeling studies,molecular dynamics,and docking simulations

The G-protein coupled estrogen receptor 1 GPER/GPR30 is a transmembrane seven-helix (7TM) receptor involved in the growth and proliferation of breast cancer. Due to the absence of a crystal structure of GPER/GPR30, in this work, molecular modeling studies have been carried out to build a three-dimensional structure, which was subsequently refined by molecular dynamics (MD) simulations (up to 120 ns). Furthermore, we explored GPER/GPR30’s molecular recognition properties by using reported agonist ligands (G1, estradiol (E2), tamoxifen, and fulvestrant) and the antagonist ligands (G15 and G36) in subsequent docking studies. Our results identified the E2 binding site on GPER/GPR30, as well as other receptor cavities for accepting large volume ligands, through GPER/GPR30 π–π, hydrophobic, and hydrogen bond interactions. Snapshots of the MD trajectory at 14 and 70 ns showed almost identical binding motifs for G1 and G15. It was also observed that C107 interacts with the acetyl oxygen of G1 (at 14 ns) and that at 70 ns the residue E275 interacts with the acetyl group and with the oxygen from the other agonist whereas the isopropyl group of G36 is oriented toward Met141, suggesting that both C107 and E275 could be involved in the protein activation. This contribution suggest that GPER1 has great structural changes which explain its great capacity to accept diverse ligands, and also, the same ligand could be recognized in different binding pose according to GPER structural conformations.     

Sequence dependence of amyloid fibril formation: insights from molecular dynamics simulations

The clarification of the physico-chemical determinants underlying amyloid deposition is critical for our understanding of misfolding diseases. With this purpose we have performed a systematic all-atom molecular dynamics (MD) study of a series of single point mutants of the de novo designed amyloidogenic peptide STVIIE. Sixteen different 50ns long simulations using explicit solvent have been carried out starting from four different conformations of a polymeric six-stranded beta-sheet. The simulations have provided evidence for the influence of a small number of site-specific hydrophobic interactions on the packing and stabilization of nascent aggregates, as well as the interplay between side-chain interactions and the net charge of the molecule on the strand arrangement of polymeric beta-sheets. This MD analysis has also shed light into the origin of the position dependence on mutation of beta-sheet polymerization that was found experimentally for this model system. Our results suggest that MD can be applied to detect critical positions for beta-sheet aggregation within a given amyloidogenic stretch. Studies similar to the one presented here can guide site-directed mutations or the design of drugs that specifically disrupt the key stabilizing interactions of beta-sheet aggregates.     

Nanocontraction flows of short-chain polyethylene via molecular dynamics simulations

Nanocontraction flows of liquid short-chain polyethylene ([CH₂]₅₀) that were uniformly extruded by a constant-speed piston into a surrounding vacuum from a reservoir through an abrupt contraction nozzle were performed by employing molecular dynamics simulations. The extrudate exhibits a similar die swell phenomenon around the exit of the nozzle. In addition, numerous molecular chains are strongly adsorbed on the external surface of the nozzle. At high extrusion speeds, the velocity and temperature profiles in the nozzle show convex and concave parabolic curves, respectively, whereas the profiles are relatively flat at lower speeds. Near the internal boundary of the nozzle, the wall slip is inspected. Significantly, during the flow, the molecular chains undergo structural deformation, including compressed, stretched and shrunk motions. Comparisons with related experimental observations show that the simulated probability distributions of the bending and dihedral angles, and variations of the squared radius of gyration and orientations, are in reasonable agreement.     

Symmetric K+ and Mg2+ ion-binding sites in the 5S rRNA loop E inferred from molecular dynamics simulations

Potassium binding to the 5 S rRNA loop E motif has been studied by molecular dynamics at high (1.0 M) and low (0.2 M) concentration of added KCl in the presence and absence of Mg²⁺. A clear pattern of seven deep groove K⁺ binding sites or regions, in all cases connected with guanine N7/O6 atoms belonging to GpG, GpA, and GpU steps, was identified, indicating that the LE deep groove is significantly more ionophilic than the equivalent groove of regular RNA duplexes. Among all, two symmetry-related sites (with respect to the central G·A pair) were found to accommodate K⁺ ions with particularly long residence times. In a preceding molecular dynamics study by Auffinger et al. in the year 2003, these two sites were described as constituting important Mg²⁺ binding locations. Altogether, the data suggest that these symmetric sites correspond to the loop E main ion binding regions. Indeed, they are located in the deep groove of an important ribosomal protein binding motif associated with a fragile pattern of non-Watson-Crick pairs that has certainly to be stabilized by specific Mg²⁺ ions in order to be efficiently recognized by the protein. Besides, the other sites accommodate monovalent ions in a more diffuse way pointing out their lesser significance for the structure and function of this motif. Ion binding to the shallow groove and backbone atoms was generally found to be of minor importance since, at the low concentration, no well defined binding site could be characterized while high K⁺ concentration promoted mostly unspecific potassium binding to the RNA backbone. In addition, several K⁺ binding sites were located in positions equivalent to water molecules from the first hydration shell of divalent ions in simulations performed with magnesium, indicating that ion binding regions are able to accommodate both mono- and divalent ionic species. Overall, the simulations provide a more precise but, at the same time, a more intricate view of the relations of this motif with its ionic surrounding.     

Mixing of nanofluids: molecular dynamics simulations and modelling

In this paper we study a mixing scheme, which has recently been proposed for microfluids, on the nanoscale. We do this by performing a series of nonequilibrium molecular dynamics simulations. On the nanoscale the chaotic mixing regime is captured. We discover a new phenomenon where the two mixing fluids exchange positions after leaving the mixing intersection. The results from the molecular dynamics simulations also reveal complex spatio-temporal stream velocity profiles generated by the mixing device. We find that these profiles can be modelled through an approximate analytical solution to the Navier–Stokes equation.     

Molecular dynamics simulations of hydrophobic collapse of ubiquitin.

Nine nonnative conformations of ubiquitin, generated during two different thermal denaturation trajectories, were simulated under nearly native conditions (62 degrees C). The simulations included all protein and solvent atoms explicitly, and simulation times ranged from 1-2.4 ns. The starting structures had alpha-carbon root-mean-square deviations (RMSDs) from the crystal structure of 4-12 A and radii of gyration as high as 1.3 times that of the native state. In all but one case, the protein collapsed when the temperature was lowered and sampled conformations as compact as those reached in a control simulation beginning from the crystal structure. In contrast, the protein did not collapse when simulated in a 60% methanol:water mixture. The behavior of the protein depended on the starting structure: during simulation of the most native-like starting structures (<5 A RMSD to the crystal structure) the RMSD decreased, the number of native hydrogen bonds increased, and the secondary and tertiary structure increased. Intermediate starting structures (5-10 A RMSD) collapsed to the radius of gyration of the control simulation, hydrophobic residues were preferentially buried, and the protein acquired some native contacts. However, the protein did not refold. The least native starting structures (10-12 A RMSD) did not collapse as completely as the more native-like structures; instead, they experienced large fluctuations in radius of gyration and went through cycles of expansion and collapse, with improved burial of hydrophobic residues in successive collapsed states.     

Comparison of different schemes to treat long-range electrostatic interactions in molecular dynamics simulations of a protein crystal

Eight molecular dynamics simulations of a ubiquitin crystal unit cell were performed to investigate the effect of different schemes to treat the long-range electrostatic interactions as well as the need to include counter ions. A crystal system was chosen as the test system, because the higher charge density compared with a protein in solution makes it more sensitive to the way of treating the electrostatic interactions. Three different schemes of treating the long-range interactions were compared: straight cutoff, reaction-field approximation, and a lattice-sum method (P3M). For each of these schemes, two simulations were performed, one with and one without the counter ions. Two additional simulations with a reaction-field force and different initial placements of the counter ions were performed to examine the effect of the initial positions of the ions. The inclusion of long-range electrostatic interactions using either a reaction-field or a lattice-sum method proved to be necessary for the simulation of crystals. These two schemes did not differ much in their ability to reproduce the crystallographic structure. The inclusion of counter ions, on the other hand, seems not necessary for obtaining a stable simulation. The initial positions of the ions have a visible but small effect on the simulation.     

Position-resolved free energy of solvation for amino acids in lipid membranes from molecular dynamics simulations

Studies of insertion and interactions of amino acids in lipid membranes are pivotal to our understanding of membrane protein structure and function. Calculating the insertion cost as a function of transmembrane helix sequence is thus an important step towards improved membrane protein prediction and eventually drug design. Here, we present position-dependent free energies of solvation for all amino acid analogs along the membrane normal. The profiles cover the entire region from bulk water to hydrophobic core, and were produced from all-atom molecular dynamics simulations. Experimental differences corresponding to mutations and costs for entire segments match experimental data well, and in addition the profiles provide the spatial resolution currently not available from experiments. Polar side-chains largely maintain their hydration and assume quite ordered conformations, which indicates the solvation cost is mainly entropic. The cost of solvating charged side-chains is not only significantly lower than for implicit solvation models, but also close to experiments, meaning these could well maintain their protonation states inside the membrane. The single notable exception to the experimental agreement is proline, which is quite expensive to introduce in vivo despite its hydrophobicity--a difference possibly explained by kinks making it harder to insert helices in the translocon.     

Brownian dynamics simulations of ion atmospheres around polyalanine and B-DNA: effects of biomolecular dielectric

We have extended an earlier Brownian dynamics simulation algorithm for simulating the structural dynamics of ions around biomolecules to accommodate dielectric inhomogeneity. The electrostatic environment of a biomolecule immersed in water was obtained by numerically solving the Poisson equation with the biomolecule treated as a low dielectric region and the solvent treated as a high dielectric region. Instead of using the mean-field type approximations of ion interactions as in the Poisson-Boltzmann model, the ions were treated explicitly by allowing them to evolve dynamically under the electrostatic field of the biomolecule. This model thus accounts for ion-ion correlations and the finite-size effects of the ions. For a 13-residue alpha-helical polyalanine and a 12-base-pair bp B-form DNA, we found that the choice of the dielectric constant of the biomolecule has much larger effects on the mean ionic structure around the biomolecule than on the fluctuational and dynamical properties of the ions surrounding the biomolecule.     

Investigations of the thermostability of rubredoxin models using molecular dynamics simulations.

The affects of differences in amino acid sequence on the temperature stability of the three-dimensional structure of the small beta-sheet protein, rubredoxin (Rd), was revealed when a set of homology models was subjected to molecular dynamics simulations at relatively high temperatures. Models of Rd from the hyperthermophile, Pyrococcus furiosus (Pf), an organism that grows optimally at 100 degrees C, were compared to three mesophilic Rds of known X-ray crystal structure. Simulations covering the limits of known Rd thermostabilities were carried out at temperatures of 300 K, 343 K, 373 K, and 413 K. They suggest that Rd stability is correlated with structural dynamics. Because the dynamic behavior of three Pf Rd models was consistently different from the dynamic behavior of the three mesophilic Rd structures, detailed analysis of the temperature-dependent dynamic behavior was carried out. The major differences between the models of the protein from the hyperthermophile and the others were: (1) an obvious temperature-dependent transition in the mesophilic structures not seen with the Pf Rd models, (2) consistent AMBER energy for the Pf Rd due to differences in nonbonded interaction terms, (3) less variation in the average conformations for the Pf Rd models with temperature, and (4) the presence of more extensive secondary structure for the Pf Rd models. These unsolvated dynamics simulations support a simple, general hypothesis to explain the hyperthermostability of Pf Rd. Its structure simplifies the conformational space to give a single minimum accessible over an extreme range of temperatures, whereas the mesophilic proteins sample a more complex conformational space with two or more minima over the same temperature range.     

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.