Valence selectivity of the gramicidin channel: a molecular dynamics free energy perturbation study.

The valence selectivity of the gramicidin channel is examined using computer simulations based on atomic models. The channel interior is modeled using a gramicidin-like periodic poly (L,D)-alanine beta-helix. Free energy perturbation calculations are performed to obtain the relative affinity of K+ and Cl- for the channel. It is observed that the interior of the gramicidin channel provides an energetically favorable interaction site for a cation but not for an anion. Relative to solvation in bulk water, the carbonyl CO oxygens can provide a favorable interaction to stabilize K+, whereas the amide NH hydrogens are much less effective in stabilizing Cl-. The results of the calculations demonstrate that, as a consequence of the structural asymmetry of the backbone charge distribution, a K+ cation can partition spontaneously from bulk water to the interior of the gramicidin channel, whereas a Cl- anion cannot.     

Photoproducts of bacteriorhodopsin mutants: a molecular dynamics study.

Molecular dynamics simulations of wild-type bacteriorhodopsin (bR) and of its D85N, D85T, D212N, and Y57F mutants have been carried out to investigate possible differences in the photoproducts of these proteins. For each mutant, a series of 50 molecular dynamics simulations of the photoisomerization and subsequent relaxation process were completed. The photoproducts can be classified into four distinct classes: 1) 13-cis retinal, with the retinal N-H+ bond oriented toward Asp-96; 2) 13-cis retinal, with the N-H+ oriented toward Asp-85 and hydrogen-bonded to a water molecule; 3) 13,14-di-cis retinal; 4) all-trans retinal. Simulations of wild-type bR and of its Y57F mutant resulted mainly in class 1 and class 2 products; simulations of D85N, D85T, and D212N mutants resulted almost entirely in class 1 products. The results support the suggestion that only class 2 products initiate a functional pump cycle. The formation of class 1 products for the D85N, D85T, and D212N mutants can explain the reversal of proton pumping under illumination by blue and yellow light.     

Effect of methylation on the stability and solvation free energy of amylose and cellulose fragments: a molecular dynamics study

Molecular dynamics (MD) simulations were used to study the stability and solvation of amylose and cellulose fragments. The recently developed gromos carbohydrate force field was further tested by simulating maltose, cellobiose, and maltoheptaose. The MD simulations reproduced fairly well the favorable conformations of disaccharides defined by the torsional angles related with the glycosidic bond and the radius gyration of maltoheptaose. The effects of methylation at different hydroxyl groups on the stability of amylose and cellulose fragments were investigated. The methylations of O-2 and O-3 reduce the stability of a single helix more than methylation at O-6, while the latter reduces the stability of a double helix more. Solvation free-energy differences between the unsubstituted amylose and cellulose fragments and the methylated species were studied using the single-step perturbation method. It was found that methylation at O-2 has the biggest effect, in agreement with experiment.     

Minor groove deformability of DNA: a molecular dynamics free energy simulation study

The conformational deformability of nucleic acids can influence their function and recognition by proteins. A class of DNA binding proteins including the TATA box binding protein binds to the DNA minor groove, resulting in an opening of the minor groove and DNA bending toward the major groove. Explicit solvent molecular dynamics simulations in combination with the umbrella sampling approach have been performed to investigate the molecular mechanism of DNA minor groove deformations and the indirect energetic contribution to protein binding. As a reaction coordinate, the distance between backbone segments on opposite strands was used. The resulting deformed structures showed close agreement with experimental DNA structures in complex with minor groove-binding proteins. The calculated free energy of minor groove deformation was approximately 4-6 kcal mol(-1) in the case of a central TATATA sequence. A smaller equilibrium minor groove width and more restricted minor groove mobility was found for the central AAATTT and also a significantly ( approximately 2 times) larger free energy change for opening the minor groove. The helical parameter analysis of trajectories indicates that an easier partial unstacking of a central TA versus AT basepair step is a likely reason for the larger groove flexibility of the central TATATA case.     

The stability of transmembrane helices: A molecular dynamics study on the isolated helices of bacteriorhodopsin

Bacteriorhodopsin (bR) continues to be a proven testing ground for the study of integral membrane proteins (IMPs). It is important to study the stability of the individual helices of bR, as they are postulated to exist as independently stable transmembrane helices (TMHs) and also for their utility as templates for modeling other IMPs with the postulated seven-helix bundle topology. Toward this purpose, the seven helices of bR have been studied by molecular dynamics simulation in this study. The suitability of using the backbone-dependent rotamer library of side-chain conformations arrived at from the data base of globular protein structures in the case TMHs has been tested by another set of 7 helix simulations with the side-chain orientations taken from this library. The influence of the residue's net charge on the helix stability was examined by simulating the helices III, IV, and VI (from both of the above sets of helices) with zero net charge on the side chains. The results of these 20 simulations demonstrate in general the stability of the isolated helices of bR in conformity with the two-stage hypothesis of IMP folding. However, the helices I, II, V, and VII are more stable than the other three helices. The helical nature of certain regions of III, IV, and VI are influenced by factors such as the net charge and orientation of several residues. It is seen that the residues Arg, Lys, Asp, and Glu (charged residues), and Ser, Thr, Gly, and Pro, play a crucial role in the stability of the helices of bR. The backbone-dependent rotamer library for the side chains is found to be suitable for the study of TMHs in IMP. © 1996 John Wiley & Sons, Inc.     

Modeling of the structure of bacteriorhodopsin. A molecular dynamics study.

Secondary structure predictions for membrane proteins are relatively reliable and permit the construction of model structures that may serve as initial conformations for molecular dynamics simulations. This might provide a scheme to predict the three-dimensional structures of membrane proteins. The feasibility of such an approach is tested for bacteriorhodopsin. We were not able to fully predict the kidney-shaped structure of bacteriorhodopsin. However, features compatible with this structure developed in a simulation starting from a circular arrangement of the seven predicted helices. When instead we started from the kidney shape, assigning the seven predicted helices in different ways to those on the structure, we could distinguish between the different assignments on the basis of energy and tilt of the helices. In this way we could select the correct assignment from a few others. For the correct assignment, the helices spontaneously adopted a tilt that agrees remarkably well with the experimental model structure derived by others. The root-mean-square deviation between our best molecular dynamics structure and the experimental model structure is 3.8 A, caused mainly by deviations in the internal degrees of freedom of the helices.     

Ab initio molecular dynamics study of proton transfer in a polyglycine analog of the ion channel gramicidin A.

Proton transfer in biological systems is thought to often proceed through hydrogen-bonded chains of water molecules. The ion channel, gramicidin A (gA), houses within its helical structure just such a chain. Using the density functional theory based ab initio molecular dynamics Car-Parrinello method, the structure and dynamics of proton diffusion through a polyglycine analog of the gA ion channel has been investigated. In the channel, a proton, which is initially present as hydronium (H3O+), rapidly forms a strong hydrogen bond with a nearest neighbor water, yielding a transient H5O2+ complex. As in bulk water, strong hydrogen bonding of this complex to a second neighbor solvation shell is required for proton transfer to occur. Within gA, this second neighbor shell included not only a channel water molecule but also a carbonyl of the channel backbone. The present calculations suggest a transport mechanism in which a priori carbonyl solvation is a requirement for proton transfer.     

Proton transfer via the external proton channel of bacteriorhodopsin

The process of proton transfer across the membrane via the external proton channel in bacteriorhodopsin is considered. A possible amino acid composition of the channel is suggested and the step-by-step mechanism of proton transfer is proposed which agrees with the experimental data. The rate of proton transfer between fixed centers at several chains of the channel was estimated for which the spectroscopic data are available.     

Proton transfer along the external proton channel of bacteriorhodopsin

The process of proton transfer along a proton channel is considered using bacteriorhodopsin as a model system, for which a large body of experimental data is available. The possible amino acid composition of the external proton half-channel of bacteriorhodopsin and the stepwise scheme of proton transfer consistent with experimental data are proposed. The rate of proton transfer between fixed centers is assessed for certain regions of this channel for which spectroscopic data are available.     

Application of molecular dynamics and free energy perturbation methods to metalloporphyrin-ligand systems II: CO and dioxygen binding to myoglobin.

The protein contribution to the relative binding affinity of the ligands CO and O2 toward myoglobin (Mb) has been simulated using free energy perturbation calculations. The tautomers of the His E7 residue are different for the oxymyoglobin (MbO2) and carboxymyoglobin (MbCO) systems. This was modeled by performing two-step calculations that mutate the ligand and mutate the His E7 tautomers in separate steps. Differences in hydrogen bonding to the O2 and CO ligands were incorporated into the model. The O2 complex was calculated to be 2-3 kcal/mol more stable than the corresponding CO complex when compared to the same difference in an isolated heme control. This value agrees well with the experimental value of 2.0 kcal/mol. In qualitative agreement with experiments, the Fe-C-O bond is found to be bent (theta = 159.8 degrees) with a small tilt (theta = 6.2 degrees). The contributions made by each of the 29 residues--within the 9.0-A radius of the iron atom--to the free energy difference are separated into van der Waals and electrostatic contributions; the latter contributions are dominant. Aside from the proximal histidine and the heme group, the residues having the largest difference in free energy in mutating MbO2-->MbCO are His E7, Phe CD1, Phe CD4, Val E11, and Thr E10.     

Structural stability of wild type and mutated alpha-keratin fragments: molecular dynamics and free energy calculations

The objective of this study is to investigate the influence of point mutations on the structural stability of coiled coil fragments of the human hair intermediate filament by molecular dynamics simulations and free energy calculations. Mutations in the helix termination motif of human hair keratin gene hHb6 seem to be connected to the hereditary hair dystrophy Monilethrix. The most common mutations reported are Glu413Lys and Glu413Asp, located at the C-terminal end of the coiled coil 2B rod domain of the IF. According to our simulations, significant conformational changes of the side chains at the mutation and neighboring sites occur due to the Glu413Lys mutation. Furthermore, the differences in electrostatic interactions cause a large change in free energy during transformation of Glu413 to Lys calculated by the thermodynamic integration approach. It is speculated that the structural rearrangement necessary to adapt the interactions in the mutated coiled coil leads to changes in the IF assembly or its stability. The second mutation, Glu413Asp, only leads to a small value of the calculated free energy difference that is within the error limits of the simulations. Thus, it has to be concluded that this mutation does not affect the coiled coil stability.     

Determination of the number of water molecules in the proton pathway of bacteriorhodopsin using neutron diffraction data

It has been shown that water molecules participate in the proton pathway of bacteriorhodopsin. Large efforts have been made to determine with various biophysical methods the number of water molecules involved. Neutron diffraction H2O/D2O exchange experiments have been often used to reveal the position of water even with low-resolution diffraction data. With this technique, care must be taken with the limitations of the difference Fourier method which are commonly applied to analyze the data. In this paper we compare the results of the difference Fourier method applied to measured diffraction data (not presented here) and models with those from alternative methods introduced here: (1) a computer model calculation procedure to determine a label's scattering length density based on a comparison of intensity differences derived from models and intensity differences from our measurements; (2) a method based on the Parseval formula. Both alternative methods have been evaluated and tested using results of neutron diffraction experiments on purple membranes (Hauss et al. 1994). Our findings indicate that the difference Fourier method applied to low-resolution diffraction data can successfully determine the position of localized water molecules but underestimates their integrated scattering length density in the presence of labels in other positions. Furthermore, we present the results of neutron diffraction experiments on purple membranes performed to determine the number of water molecules in the projected area of the Schiff base at 86%, 75% and 57% relative humidity (r.h.). We found 19 +/- 2 exchangeable protons at 75% r.h., which means at least 8-9 water molecules are indispensable for normal pump function.     

Conformation and fluctuations of free stefin B: a molecular dynamics study.

Molecular dynamics study was performed on the cysteine proteinase inhibitor stefin B. Structure of inhibitor from the complex with papain was used as a starting point. Amino terminal "trunk" of the inhibitor which lies extended along the cleft of the enzyme in the complex, folded onto the body of inhibitor during MD simulation, thereby reducing the total and particularly hydrophobic surface exposed to the solvent. This effect counterbalances hydrophobic contribution of the "trunk" and explains why its deletion in stefin B and related inhibitors doesn't reduce the dissociation constant. The rest of stefin B conformation is conserved together with main chain hydrogen bonds. Fluctuations of C alpha atoms resembles crystallographic B factors with exception of residues in contact with enzyme.     

Calculation of the relative binding free energy of 2'GMP and 2'AMP to ribonuclease T1 using molecular dynamics/free energy perturbation approaches

We present a calculation of the relative changes in binding free energy between the complex of ribonuclease T1 (RNase Tr) with its inhibitor 2'-guanosine monophosphate (2'GMP) and that of RNase T1-2'-adenosine monophosphate (2'AMP) by means of a thermodynamic perturbation method implemented with molecular dynamics. Using the available crystal structure of the RNase T1-2'GMP complex, the structure of the RNase T1-2'AMP complex was obtained as a final structure of the perturbation calculation. The calculated difference in the free energy of binding (delta delta Gbind) was 2.76 kcal/mol. This compares well with the experimental value of 3.07 kcal/mol. The encouraging agreement in delta delta Gbind suggests that the interactions of inhibitors with the enzyme are reasonably represented. Energy component analyses of the two complexes reveal that the active site of RNase T1 electrostatically stabilizes the binding of 2'GMP more than that of 2'AMP by 44 kcal/mol, while the van der Waals' interactions are similar in the two complexes. The analyses suggest that the mutation from Glu46 to Gln may lead to a preference of RNase T1 for adenine in contrast to the guanine preference of the wild-type enzyme. Although the molecular dynamics equilibration moves the atoms of the RNase T1-2'GMP system about 0.9 A from their X-ray positions and the mutation of the G to A in the active site increases the deviation from the X-ray structure, the mutation of the A back to G reduces the deviation. This and the agreement found for delta delta Gbind suggest that the molecular dynamics/free energy perturbation method will be useful for both energetic and structural analysis of protein-ligand interactions.     

Role of internal water molecules in bacteriorhodopsin

Internal water molecules are considered to play a crucial role in the functional processes of proton pump proteins. They may participate in hydrogen-bonding networks inside proteins that constitute proton pathways. In addition, they could participate in the switch reaction by mediating an essential proton transfer at the active site. Nevertheless, little has been known about the structure and function of internal water molecules in such proteins. Recent progress in infrared spectroscopy and X-ray crystallography provided new information on water molecules inside bacteriorhodopsin, the light-driven proton pump. The accumulated knowledge on bacteriorhodopsin in the last decade of the 20th century will lead to a realistic picture of internal water molecules at work in the 21st century. In this review, I describe how the role of water molecules has been studied in bacteriorhodopsin, and what should be known about the role of water molecules in the future.     

Azide reduces the hydrophobic barrier of the bacteriorhodopsin proton channel

The sensitivity of a nitroxide spin label to the polarity of its environment has been used to estimate the hydrophobic barrier of the proton channel of the transmembrane proton pump bacteriorhodopsin. By means of site-specific mutagenesis, single cysteine residues were introduced at 10 positions located at the protein surface, in the protein interior, and along the proton pathway. After reaction with a methanethiosulfonate spin label, the principle values of the hyperfine tensor A and the g-tensor were determined from electron paramagnetic resonance spectra measured at 170 K. The shape of the hydrophobic barrier of the proton channel is characterized in terms of a polarity index, DeltaA, determined from the variation of the hyperfine coupling constant Azz. The maximum of the hydrophobic barrier is found to be close to the retinal chromophore in the proton uptake pathway. The effect of the asymmetric distribution of charged and polar residues in the proton release and uptake pathways is clearly reflected in the behavior of the hydrophobic barrier. The presence of azide reduces the barrier height of both the cytoplasmic and extracellular channels. This finding supports the view of azide and other weakly acidic anions as catalysts for the formation of hydrogen-bonded networks in proton pathways of proteins.     

Thermodynamic coupling between activation and inactivation gating in potassium channels revealed by free energy molecular dynamics simulations

The amount of ionic current flowing through K(+) channels is determined by the interplay between two separate time-dependent processes: activation and inactivation gating. Activation is concerned with the stimulus-dependent opening of the main intracellular gate, whereas inactivation is a spontaneous conformational transition of the selectivity filter toward a nonconductive state occurring on a variety of timescales. A recent analysis of multiple x-ray structures of open and partially open KcsA channels revealed the mechanism by which movements of the inner activation gate, formed by the inner helices from the four subunits of the pore domain, bias the conformational changes at the selectivity filter toward a nonconductive inactivated state. This analysis highlighted the important role of Phe103, a residue located along the inner helix, near the hinge position associated with the opening of the intracellular gate. In the present study, we use free energy perturbation molecular dynamics simulations (FEP/MD) to quantitatively elucidate the thermodynamic basis for the coupling between the intracellular gate and the selectivity filter. The results of the FEP/MD calculations are in good agreement with experiments, and further analysis of the repulsive, van der Waals dispersive, and electrostatic free energy contributions reveals that the energetic basis underlying the absence of inactivation in the F103A mutation in KcsA is the absence of the unfavorable steric interaction occurring with the large Ile100 side chain in a neighboring subunit when the intracellular gate is open and the selectivity filter is in a conductive conformation. Macroscopic current analysis shows that the I100A mutant indeed relieves inactivation in KcsA, but to a lesser extent than the F103A mutant.