Molecular dynamics study of Desulfovibrio africanus cytochrome c3 in oxidized and reduced forms

A 5-ns molecular dynamics study of a tetraheme cytochrome in fully oxidized and reduced forms was performed using the CHARMM molecular modeling program, with explicit water molecules, Langevin dynamics thermalization, Particle Mesh Ewald long-range electrostatics, and quantum mechanical determination of heme partial charges. The simulations used, as starting points, crystallographic structures of the oxidized and reduced forms of the acidic cytochrome c(3) from Desulfovibrio africanus obtained at pH 5.6. In this paper we also report structures for the two forms obtained at pH 8. In contrast to previous cytochrome c(3) dynamics simulations, our model is stable. The simulation structures agree reasonably well with the crystallographic ones, but our models show higher flexibility and the water molecules are more labile. We have compared in detail the differences between the simulated and experimental structures of the two redox states and observe that the hydration structure is highly dependent on the redox state. We have also analyzed the interaction energy terms between the hemes, the protein residues, and water. The direct electrostatic interaction between hemes is weak and nearly insensitive to the redox state, but the remaining terms are large and contribute in a complex way to the overall potential energy differences that we see between the redox states.     

Molecular dynamics simulations of rubredoxin from Clostridium pasteurianum: Changes in structure and electrostatic potential duringredox reactions

Molecular dynamics simulations of Clostridium pasteurianum rubredoxin in the oxidized and reduced forms have been performed. Good agreement between both forms and crystal data has been obtained (rms deviation of backbone atoms of 1.06 and 1.42 Å, respectively), which was due in part to the use of explicit solvent and counterions. The reduced form exhibits an unexpected structural change: the redox site becomes much more solvent-accessible, so that water enters a channel between the surface and the site, but with little actual structural rearrangement (the rms deviation of backbone atoms between the oxidized and reduced is 0.77 Å). The increase in solvent accessibility is also seen, although to a much lesser extent, between the oxidized and reduced crystal structures of Pyrococcus furiosus rubredoxin, but no high resolution crystal or nuclear magnetic resonance solution data exist for reduced C. pasteurianum rubredoxin. The electrostatic potential at the iron site and fluctuations in the potential, which contribute to both the redox and electron transfer properties, have also been evaluated for both the oxidized and the reduced simulations. These results show that the backbone plays a significant role (62–70 kcall/mol/e) and the polar sidechains contribute relatively little (0–4 kcal/mol/e) to the absolute electrostatic potential at the iron of rubredoxin for both forms. However, both groups contribute significantly to the change in redox state by becoming more polarized and more densely packed around the redox site upon reduction. Furthermore, these results show that the solvent becomes much more polarized in the reduced form than in the oxidized form, even excluding the penetrating water. Finally, the simulation indicates that the contribution of the charged side chains to the electrostatic potential is largely canceled by that of the counterions. © 1995 Wiley-Liss, Inc.     

Reorganization and conformational changes in the reduction of tetraheme cytochromes

Molecular dynamics simulation (MD) constitutes an alternative to time-consuming experiments for studying conformational changes. We apply MD on a redox system where experimental information exists for the fully oxidized and fully reduced states: tetraheme cytochrome c3. Instead of doing one simulation for each state, we apply 10 4-ns replicas for both states, which provides robust statistics to characterize the redox changes. Besides these long simulations, we perform 120 short ones (50 ps), where an equilibrated oxidized state is perturbed to a reduced state. This allows the application of a nonequilibrium method, the subtraction technique, which makes it possible to characterize the different timescales of conformational changes. Reduction induces conformational changes in the N-terminus and on the loops spanning residues 36-42 and 88-93, which correlate very well with experiments, demonstrating the applicability of this methodology. We also analyze the effect of reduction on hydrogen bonds, solvent accessible surface and bound water, the changes being found to involve the hemes and propionate groups. Redox-induced protonation is also investigated, by protonating the propionates D from hemes I and IV. Although this change in the former does not have major conformational consequences, it induces in the latter conformational changes beyond the ones obtained with reduction.     

Mean field stochastic boundary molecular dynamics simulation of a phospholipid in a membrane

Computer simulations of phospholipid membranes have been carried out by using a combined approach of molecular and stochastic dynamics and a mean field based on the Marcelja model. First, the single-chain mean field simulations of Pastor et al. [(1988) J. Chem. Phys. 89, 1112-1127] were extended to a complete dipalmitoylphosphatidylcholine molecule; a 102-ns Langevin dynamics simulation is presented and compared with experiment. Subsequently, a hexagonally packed seven-lipid array was simulated with Langevin dynamics and a mean field at the boundary and with molecular dynamics (and no mean field) in the center. This hybrid method, mean field stochastic boundary molecular dynamics, reduces bias introduced by the mean field and eliminates the need for periodic boundary conditions. As a result, simulations extending to tens of nanoseconds may be carried out by using a relatively small number of molecules to model the membrane environment. Preliminary results of a 20-ns simulation are reported here. A wide range of motions, including overall reorientation with a nanosecond decay time, is observed in both simulations, and good agreement with NMR, IR, and neutron diffraction data is found.     

Effects of environment on the structure of Pyrococcus furiosus rubredoxin: a molecular dynamics study

Molecular dynamics simulations based on a 0.95-A resolution crystal structure of Pyrococcus furiosus have been performed to elucidate the effects of the environment on the structure of rubredoxin, and proteins in general. Three 1-ns simulations are reported here: two crystalline state simulations at 123 and 300 K, and a solution state simulation at 300 K. These simulations show that temperature has a greater impact on the protein structure than the close molecular contacts of the crystal matrix in rubredoxin, although both have an effect on its dynamic properties. These results indicate that differences between NMR solution structures and X-ray crystal structures will be relatively minor if they are done at similar temperatures. In addition, the crystal simulations appears to mimic previous crystallographic experiments on the effects of cryo-temperature on temperature factors, and might provide a useful tool in the structural analysis of protein structures solved at cryo-temperatures.     

Structure and dynamics of calcium-activated calmodulin in solution.

Two 4-ns molecular dynamics simulations of calcium loaded calmodulin in solution have been performed, using both standard nonbonded cutoffs and Ewald summation to treat electrostatic interactions. Our simulation results are generally consistent with solution experimental studies of calmodulin structure and dynamics, including NMR, cross-linking, fluorescence and x-ray scattering. The most interesting result of the molecular dynamics simulations is the detection of large-scale structural fluctuations of calmodulin in solution. The globular N- and C-terminal domains tend to move approximately like rigid bodies, with fluctuations of interdomain distances within a 7 A range and of interdomain angles by up to 60 deg. Essential dynamics analysis indicates that the three dominant types of motion involve bending of the central helix in two perpendicular planes and a twist in which the domains rotate in opposite directions around the central helix. In the more realistic Ewald trajectory the protein backbone remains mostly within a 2-3 A root-mean-square distance from the crystal structure, the secondary structure within the domains is conserved and middle part of the central helix becomes disordered. The central helix itself exhibits limited fluctuations, with its bend angle exploring the 0-50 degrees range and the end-to-end distance falling in 39-43 A. The results of the two simulations were similar in many respects. However, the cutoff trajectory exhibited a larger deviation from the crystal, loss of several helical hydrogen bonds in the N-terminal domain and lack of structural disorder in the central helix.     

Leucine 41 is a gate for water entry in the reduction of Clostridium pasteurianum rubredoxin

Biological electron transfer is an efficient process even though the distances between the redox moieties are often quite large. It is therefore of great interest to gain an understanding of the physical basis of the rates and driving forces of these reactions. The structural relaxation of the protein that occurs upon change in redox state gives rise to the reorganizational energy, which is important in the rates and the driving forces of the proteins involved. To determine the structural relaxation in a redox protein, we have developed methods to hold a redox protein in its final oxidation state during crystallization while maintaining the same pH and salt conditions of the crystallization of the protein in its initial oxidation state. Based on 1.5 A resolution crystal structures and molecular dynamics simulations of oxidized and reduced rubredoxins (Rd) from Clostridium pasteurianum (Cp), the structural rearrangements upon reduction suggest specific mechanisms by which electron transfer reactions of rubredoxin should be facilitated. First, expansion of the [Fe-S] cluster and concomitant contraction of the NH...S hydrogen bonds lead to greater electrostatic stabilization of the extra negative charge. Second, a gating mechanism caused by the conformational change of Leucine 41, a nonpolar side chain, allows transient penetration of water molecules, which greatly increases the polarity of the redox site environment and also provides a source of protons. Our method of producing crystals of Cp Rd from a reducing solution leads to a distribution of water molecules not observed in the crystal structure of the reduced Rd from Pyrococcus furiosus. How general this correlation is among redox proteins must be determined in future work. The combination of our high-resolution crystal structures and molecular dynamics simulations provides a molecular picture of the structural rearrangement that occurs upon reduction in Cp rubredoxin.     

Influence of protein flexibility on the redox potential of rubredoxin: Energy minimization studies

A theoretical investigation of the protein contribution to the redox potential of the iron–sulfur protein rubredoxin is presented. Structures of the oxidized and reduced forms of the protein were obtained by energy minimizing the oxidized crystal structure of Clostridium pasteurianum rubredoxin with appropriate charges and parameters. By including 102 crystal waters, structures close to the original crystal structure were obtained (rms difference of 1.16 Å), even with extensive minimization, thus allowing accurate calculations of comparative energies. Our calculations indicate an energy change of about –60 kcal/mol (2.58 eV) in the protein alone upon reduction. This energy change was due to both the change in charge of the redox site and the subsequent relaxation of the protein. An energy minimization procedure for the relaxation gives rms differences between the oxidized and reduced states of about 0.2 Å. The changes were small and occurred in both the backbone and sidechain mainly near the Fe–S center but contributed about – 16 kcal/mol (0.69 eV) to the total protein contribution. Although the neglect of certain effects such as electronic polarization may make the relaxation energies calculated an upper limit, the results indicate that protein relaxation contributes substantially to the redox potential. © 1993 Wiley-Liss, Inc.     

Differences between apo and three holo forms of the intestinal fatty acid binding protein seen by molecular dynamics computer calculations

It is commonly believed that binding affinity can be estimated by consideration of local changes of ligand and protein. This paper discusses a set of molecular dynamics simulations of intestinal fatty acid binding protein addressing the protein's response to presence or absence of different ligands. A 5-ns simulation was performed of the protein without a ligand, and three simulations (one 5-ns and two 2-ns) were performed with different fatty acids bound. The results indicate that, although the basic protein structure is unchanged by the presence of the ligand, other properties are significantly affected by ligand binding. For example, zero-time covariance patterns between protein, bound waters, and ligand vary between the different simulations. Moreover, the interaction energies between ligand and specific residues indicate that different ligands are stabilized in different ways. In sum, the results suggest that binding thermodynamics within this system will need to be calculated not from a subset of nearby protein:ligand interactions, but will depend on a knowledge of the motions coupling together water, protein, and ligand.     

Filter flexibility and distortion in a bacterial inward rectifier K+ channel: simulation studies of KirBac1.1

The bacterial channel KirBac1.1 provides a structural homolog of mammalian inward rectifier potassium (Kir) channels. The conformational dynamics of the selectivity filter of Kir channels are of some interest in the context of possible permeation and gating mechanisms for this channel. Molecular dynamics simulations of KirBac have been performed on a 10-ns timescale, i.e., comparable to that of ion permeation. The results of five simulations (total simulation time 50 ns) based on three different initial ion configurations and two different model membranes are reported. These simulation data provide evidence for limited (<0.1 nm) filter flexibility during the concerted motion of ions and water molecules within the filter, such local changes in conformation occurring on an approximately 1-ns timescale. In the absence of K(+) ions, the KirBac selectivity filter undergoes more substantial distortions. These resemble those seen in comparable simulations of other channels (e.g., KcsA and KcsA-based homology models) and are likely to lead to functional closure of the channel. This suggests filter distortions may provide a mechanism of K-channel gating in addition to changes in the hydrophobic gate formed at the intracellular crossing point of the M2 helices. The simulation data also provide evidence for interactions of the "slide" (pre-M1) helix of KirBac with phospholipid headgroups.     

A tree-based algorithm for determining the effects of solvation on the structure of salivary gland tripeptide NH3+-D-PHE-D-GLU-GLY-COO-

A D-enantiomeric analog of the submandibular gland rat-1 tripeptide FEG (Seq: NH(3)(+)-Phe-Glu-Gly-COO(-)) called feG (Seq: NH(3)(+)-D-Phe-D-Glu-Gly-COO(-)) was examined by molecular dynamics simulations in water. Previous in vacuo simulations suggested a conformation consisting predominantly of interactions between the Phe side chain and glutamyl-carboxyl group and a carboxyl/amino termini interaction. The solvated peptide was simulated using two approaches which were compared-a single 400-ns simulation and a "simulation tree." The "tree" approach utilized 45 10-ns simulations with different conformations used as initial structures for given trajectories. We demonstrate that multiple short duration simulations are able to describe the same conformational space as that described by longer simulations. Furthermore, previously described in vacuo interactions were confirmed with amendments: the previously described head-to-tail arrangement of the amino and carboxyl termini, was not observed; the interaction between the glutamyl carboxyl and Phe side chain describes only one of a continuum of conformations present wherein the aromatic residue remains in close proximity to the glutamyl carbonyl group, and also interacts with either of the two available carboxyl groups. Finally, utilizing only two separate 10-ns trajectories, we were able to better describe the conformational space than a single 60-ns trajectory, realizing a threefold decrease in the computational complexity of the problem.     

HIV-1 integrase catalytic core: molecular dynamics and simulated fluorescence decays.

Two molecular dynamics simulations have been carried out on the HIV-1 integrase catalytic core starting from fully determined crystal structures. During the first one, performed in the absence of divalent cation (6-ns long), the catalytic core took on two main conformations. The conformational transition occurs at approximately 3.4 ns. In contrast, during the second one, in the presence of Mg(2+) (4-ns long), there were no such changes. The molecular dynamics simulations were used to compute the fluorescence intensity decays emitted by the four tryptophan residues considered as the only chromophores. The decay was computed by following, frame by frame, the amount of chromophores that remained excited at a certain time after light absorption. The simulation took into account the quenching through electron transfer to the peptide bond and the fluorescence resonance energy transfer between the chromophores. The fit to the experimental intensity decays obtained at 5 degrees C and at 30 degrees C is very good. The fluorescence anisotropy decays were also simulated. Interestingly, the fit to the experimental anisotropy decay was excellent at 5 degrees C and rather poor at 30 degrees C. Various hypotheses such as dimerization and abnormal increase of uncorrelated internal motions are discussed.     

Molecular dynamics study of the dominant-negative E219K polymorphism in human prion protein

Human prion diseases are associated with misfolding or aggregation of the Human Prion Protein (HuPrP). Missense mutations in the HuPrP gene, contribute to conversion of HuPrP^C to HuPrP^Sc and amyloid formation. Based on our previous comprehensive study, three missense mutations, from two different functional groups, i.e. disease-related mutations, and protective mutations, were selected and extensive molecular dynamics simulations were performed on these three mutants to compare their dynamics and conformations with those of the wildtype HuPrP. In addition to simulations of monomeric forms of mutants, in order to study the dominant-negative effect of protective mutation (E219K), 30-ns simulations were performed on E219K-wildtype and wildtype-wildtype dimeric forms. Our results indicate that, although after 30-ns simulations the global three-dimensional structure of models remain fairly intact, the disease-related mutations (V210I and Q212P) introduce local structural changes, i.e. close contact changes and secondary structure changes, in addition to global flexibility changes. Furthermore, our results support the loss of hydrophobic interaction due to the mutations in hydrophobic core that has been reported by previous NMR and computational studies. On the other hand, this protective mutation (E219K) results in helix elongation, and significant increases of overall flexibility of E219K mutant during 30-ns simulation. In conclusion, the simulations of dimeric forms suggest that the dominant-negative effect of this protective mutation (E219K) is due to the incompatible structures and dynamics of allelic variants during conversion process.     

Key structural features of the actin filament Arp2/3 complex branch junction revealed by molecular simulation

We investigated the structure, properties and dynamics of the actin filament branch junction formed by actin-related protein (Arp) 2/3 complex using all-atom molecular dynamics (MD) simulations based on a model fit to a reconstruction from electron tomograms. Simulations of the entire structure consisting of 31 protein subunits together with solvent molecules containing ～3 million atoms were performed for an aggregate time of 175 ns. One 75-ns simulation of the original reconstruction was compared to two 50-ns simulations of alternate structures, showing that the hypothesized branch junction structure is very stable. Our simulations revealed that the interface between Arp2/3 complex and the mother actin filament features a large number of salt bridges and hydrophobic contacts, many of which are dynamic and formed/broken on the timescale of the simulation. The simulations suggest that the DNase binding loops in Arp3, and possibly Arp2, form stabilizing contacts with the mother filament. Unbiased comparison of models sampled from the MD simulation trajectory with the primary experimental electron tomography data identified regions were snapshots from the simulation provide atomic details of the model structures and also pinpoints regions where the initial modeling based on the electron tomogram reconstruction may be suboptimal.     

Probing the flexibility of the DsbA oxidoreductase from Vibrio cholerae--a 15N - 1H heteronuclear NMR relaxation analysis of oxidized and reduced forms of DsbA

We have determined the structure of the reduced form of the DsbA oxidoreductase from Vibrio cholerae. The reduced structure shows a high level of similarity to the crystal structure of the oxidized form and is typical of this class of enzyme containing a thioredoxin domain with an inserted alpha-helical domain. Proteolytic and thermal stability measurements show that the reduced form of DsbA is considerably more stable than the oxidized form. NMR relaxation data have been collected and analyzed using a model-free approach to probe the dynamics of the reduced and oxidized states of DsbA. Akaike's information criteria have been applied both in the selection of the model-free models and the diffusion tensors that describe the global motions of each redox form. Analysis of the dynamics reveals that the oxidized protein shows increased disorder on the pico- to nanosecond and micro- to millisecond timescale. Many significant changes in dynamics are located either close to the active site or at the insertion points between the domains. In addition, analysis of the diffusion data shows there is a clear difference in the degree of interdomain movement between oxidized and reduced DsbA with the oxidized form being the more rigid. Principal components analysis has been employed to indicate possible concerted movements in the DsbA structure, which suggests that the modeled interdomain motions affect the catalytic cleft of the enzyme. Taken together, these data provide compelling evidence of a role for dynamics in the catalytic cycle of DsbA.     

A probabilistic approach for estimating water permeability in pressure-driven membranes

A probabilistic approach is proposed to estimate water permeability in a cellulose triacetate (CTA) membrane. Water transport across the membrane is simulated in reverse osmosis mode by means of non-equilibrium molecular dynamics (MD) simulations. Different membrane configurations obtained by an annealing MD simulation are considered and simulation results are analyzed by using a hierarchical Bayesian model to obtain the permeability of the different membranes. The estimated membrane permeability is used to predict full-scale water flux by means of a process-level Monte Carlo simulation. Based on the results, the parameters of the model are observed to converge within 5-ns total simulation time. The results also indicate that the use of unique structural configurations in MD simulations is essential to capture realistic membrane properties at the molecular scale. Furthermore, the predicted full-scale water flux based on the estimated permeability is within the same order of magnitude of bench-scale experimental measurement of 1.72×10^?5 m/s.     

Analysis and evaluation of channel models: simulations of alamethicin

Alamethicin is an antimicrobial peptide that forms stable channels with well-defined conductance levels. We have used extended molecular dynamics simulations of alamethicin bundles consisting of 4, 5, 6, 7, and 8 helices in a palmitoyl-oleolyl-phosphatidylcholine bilayer to evaluate and analyze channel models and to link the models to the experimentally measured conductance levels. Our results suggest that four helices do not form a stable water-filled channel and might not even form a stable intermediate. The lowest measurable conductance level is likely to correspond to the pentamer. At higher aggregation numbers the bundles become less symmetrical. Water properties inside the different-sized bundles are similar. The hexamer is the most stable model with a stability comparable with simulations based on crystal structures. The simulation was extended from 4 to 20 ns or several times the mean passage time of an ion. Essential dynamics analyses were used to test the hypothesis that correlated motions of the helical bundles account for high-frequency noise observed in open channel measurements. In a 20-ns simulation of a hexameric alamethicin bundle, the main motions are those of individual helices, not of the bundle as a whole. A detailed comparison of simulations using different methods to treat long-range electrostatic interactions (a twin range cutoff, Particle Mesh Ewald, and a twin range cutoff combined with a reaction field correction) shows that water orientation inside the alamethicin channels is sensitive to the algorithms used. In all cases, water ordering due to the protein structure is strong, although the exact profile changes somewhat. Adding an extra 4-nm layer of water only changes the water ordering slightly in the case of particle mesh Ewald, suggesting that periodicity artifacts for this system are not serious.     

Monitoring of root growth and redox conditions in paddy soil rhizotrons by redox electrodes and image analysis

The root-zone of wetland rice was monitored in a paddy soil throughout a vegetation period with the aid of a rhizotron experiment. For this purpose (a) digital images of the root-zone were taken daily, and (b) the redox potential was measured in situ every day. The images were processed by image analysis in order to display areas of oxidation and reduction in the soil. Therefore, thresholds were set to simplify the localization and quantification of discrete areas which were colourized due to the redox potential. Both, images and measured redox potentials, provide the basis for the visualization of the root and redox dynamics in the root-zone. The anaerobic root-zone of flooded paddy soils is significantly influenced by the aerenchymal transport of oxygen to rice roots. The release of oxygen into the rhizosphere, which causes different patterns of oxidized and reduced areas in the course of the vegetation period, also affects microbial communities such as methane producing archaea or methane oxidizing bacteria. The visualization of redox dynamics may, therefore, be useful to localize potential hotspots for the microorganisms in the root-zone of paddy soils. The reduced and oxidized conditions changed spatiotemporally. Oxidized areas were mostly found in the surrounding of active roots and in a distinct layer next to the soil surface. Reduced areas shifted from beneath the oxidized surface layer into sparsely-rooted soil. The ratio of the analyzed oxidized and reduced areas was oscillating with increasing intensity throughout the monitored vegetation period.     

Redox state of peroxy and ferryl intermediates in cytochrome c oxidase catalysis.

The redox states of the "peroxy" (P) and "ferryl" (F) intermediates formed during reoxidation of reduced bovine cytochrome c oxidase have been probed by reduction with both ferrocytochrome c and acetylpyridine NADH under anaerobic conditions using optical spectroscopy. The reduction of the P and F forms revealed that both are in very similar redox states. One-electron reduction of either the P or F form yields an optical spectrum primarily due to oxidized enzyme implying that the heme iron of cytochrome a3 is in the ferryl state in both forms. The F and P forms were found to be 1 and less than 1.3 oxidizing equiv, respectively, above the oxidized enzyme. The slightly higher oxidation state in the P form is interpreted as being due to an optically undetectable redox center presumably located in the binuclear cavity.