Experimentally verifying molecular dynamics simulations through fluorescence anisotropy measurements

Molecular dynamics simulations of Trp side-chain conformational flexibility in the gramicidin A channel

Gramicidin A (gA) is prototypical peptide antibiotic and a model ion channel former. Configured in the solid-state NMR beta(6.5)-helix channel conformation, gA was subjected to 1-ns molecular dynamics (MD) gas phase simulations using the all-atom charmm22 force field to ascertain the conformational stability of the Trp side chains as governed by backbone and neighboring side-chain contacts. Three microcanonical trajectories were computed using different initial atomic velocities for each of twenty different initial structures. For each set, one of the four Trp side chains in each monomer was initially positioned in one of the five non-native conformations (A. E. Dorigo et al., Biophysical Journal, 1999, Vol. 76, 1897-1908), the other Trps being positioned in the native state, o1. In three additional control simulations, all Trps were initiated in the native conformation. After equilibration, constraints were removed and subsequent conformational changes of the initially constrained Trp were measured. The chi(1) was more flexible than chi(2.1). The energetically optimal orientation, o1 (Dorigo et al., 1999), was the most stable in all four Trp positions (9, 11, 13, 15) and remained unchanged for the entire 1 ns simulation in 19 of 24 trials. Changes in chi(1) from each of the 5 suboptimal states occur readily. Two of the non-native conformations reverted readily to o1, whereas the other three converted to an intermediate state, i2. There were frequent interconversions between i2 and o1. We speculate that experimentally observed Trp stability is caused by interactions with the lipid-water interface, and that stabilization of one of the suboptimal conformations in gA, such as i2, by lipid headgroups could produce a secondary, metastable conformational state. This could explain recent experimental studies of differences in the channel conductance dispersity between gA and a Trp-to-Phe gA analog, gramicidin M (gM, J. C. Markham et al., Biochimica et Biophysica Acta, 2001, Vol. 1513, 185-192).

Trp42 rotamers report reduced flexibility when the inhibitor acetyl-pepstatin is bound to HIV-1 protease

The Q7K/L331/L631 HIV-1 protease mutant was expressed in Escherichia coli and the effect of binding a substrate-analog inhibitor, acetyl-pepstatin, was investigated by fluorescence spectroscopy and molecular dynamics. The dimeric enzyme has four intrinsic tryptophans, located at positions 6 and 42 in each monomer. Fluorescence spectra and acrylamide quenching experiments show two differently accessible Trp populations in the apoenzyme with k(q1) = 6.85 x 10(9) M(-1) s(-1) and k(q2) = 1.88 x 10(9) M(-1) s(-1), that merge into one in the complex with k(q) = 1.78 x 10(9) M(-1) s(-1). 500 ps trajectory analysis of Trp X1/X2 rotameric interconversions suggest a model to account for the observed Trp fluorescence. In the simulations, Trp6/Trp6B rotameric interconversions do not occur on this timescale for both HIV forms. In the apoenzyme simulations, however, both Trp42s and Trp42Bs are flipping between X1/X2 states; in the complexed form, no such interconverions occur. A detailed investigation of the local Trp environments sampled during the molecular dynamics simulation suggests that one of the apoenzyme Trp42B rotameric interconversions would allow indole-quencher contact, such as with nearby Tyr59. This could account for the short lifetime component. The model thus interprets the experimental data on the basis of the conformational fluctuations of Trp42s alone. It suggests that the rotameric interconversions of these Trps, located relatively far from the active site and at the very start of the flap region, becomes restrained when the apoenzyme binds the inhibitor. The model is thus consistent with associating components of the fluorescence decay in HIV-1 protease to ground state conformational heterogeneity.

X-ray structure and conformational dynamics of the HIV-1 protease in complex with the inhibitor SDZ283-910: agreement of time-resolved spectroscopy and molecular dynamics simulations

Based on the X-ray structure of the human immunodeficiency virus type-1 (HIV-1) protease in complex with the statine-derived inhibitor SDZ283-910, a 542 ps molecular dynamics trajectory was computed. For comparison with the 805 ps trajectory obtained for the uncomplexed enzyme, the theoretical fluorescence anisotropy decay of the unliganded protease and the inhibitor complex was calculated from the trajectories of the Trp6A/Trp6B and Trp42A/Trp42B transition dipole moments. This enabled us to directly compare the simulated data with the experimental picosecond time-resolved fluorescence data. Fitting both experimental and simulated data to the Kohlrausch-Williams-Watts (KWW) function exp(-t/tauk)beta revealed a very good agreement for the uncomplexed protease as well as for the SDZ283-910 complex. Binding of the inhibitor induced a faster decay of both the experimental and the computed protease fluorescence anisotropy decay. By this integrative approach, the atomic detail of inhibitor-induced changes in the conformational dynamics of the HIV-1 protease was experimentally verified and will be used for further inhibitor optimisation.

Simulated dipeptide recognition by vancomycin

The antimicrobial activity of vancomycin and related glycopeptide antibiotics is due to stereospecific recognition of polypeptide components in bacterial cell walls. To better understand how these antibiotics recognize polypeptide determinants, we have developed dynamic models of the complexes formed by the vancomycin aglycon and two different dipeptide ligands, Ac-D-ala-D-ala and Ac-D-ala-gly. Molecular dynamics simulations of the two complexes, initially conditioned with distance constraints derived from two-dimensional nuclear magnetic resonance (NMR) studies, are conformationally stable and propagate in a manner consistent with the NMR-derived constraints after the constraints are removed. Free energy calculations accurately predict the relative binding affinity of these two complexes and help validate the simulation models for detailed structural analysis. Although the two ligands adopt similar conformations when bound to the antibiotic, there are clear differences in the configuration of intermolecular hydrogen bonds, the overall shape of the antibiotic, and other structural features of the two complexes. This analysis illustrates how complex structural and dynamic factors interrelate and contribute to differences in binding affinity.

Time-resolved fluorescence spectroscopy

Time-resolved fluorescence spectroscopy is used to monitor molecular interactions and motions that occur in the picosecond-nanosecond time range, and is especially useful in the analysis of biomolecular structure and dynamics. Recent advances in the application of time-resolved fluorescence spectroscopy to biological systems have led to a better understanding of the origin of nonexponential fluorescence decay in proteins, the use of tryptophan analogs as unique spectroscopic probes of protein-protein interactions, the detailed characterization of protein-folding processes and intermediates, and the development of new approaches to the study of DNA-protein interactions.

Molecular dynamics simulations of peptides from BPTI: a closer look at amide-aromatic interactions

Molecular dynamics (MD) simulations of short peptides in water were performed to establish whether it is possible to reproduce experimental data from chemical shift measurements by nuclear magnetic resonance spectroscopy. Three different tetrapeptides were studied. The first, YTGP (Tyr-Thr-Gly-Pro), shows an electrostatic interaction between the aromatic ring of Tyr and the backbone amide hydrogen atom of Gly. The second, YTAP (Tyr-Thr-Ala-Pro), cannot make such an interaction because of steric hindrance of the Ala side chain and hence does not show a well-defined conformation. The third, FTGP (Phe-Thr-Gly-Pro), is shown to alternate between two different conformations. It is demonstrated that small differences in chemical shift, corresponding to these slightly different conformations, can be quantitatively modeled in MD simulations when using the proper force-field parameters and water model Explicit inclusion of hydrogen atoms o the aromatic rings is essential for a proper description of electrostatic interactions, but the choice of the water model is equally important. We found that a combination of the SPC/E water model and a revised GROMOS87 force field gives close agreement with experiment, while the same and other force fields in combination with SPC or TIP3P water did not reproduce the NMR data at all. Simulations of a longer peptide from bovine pancreatic trypsin inhibitor, containing the YTGP sequence, did show the interaction between the aromatic ring and the amide hydrogen, but not as pronounced as the simulations of shorter periods.

Molecular dynamics simulation of the rare amino acid LL-dityrosine and a dityrosine-containing peptide: comparison with time-resolved fluorescence

The fluorescence of the rare amino acid LL-dityrosine, which is found in insoluble biological materials with structural features, was recently shown to decay non-exponentially (Kungl et al. (1992) J. Fluorescence 2, 63-74). Here we investigated the time-resolved fluorescence of a dityrosine-containing peptide (DCP) to study the influence of side chains on the fluorescence decay of the chromophore. The fluorescence decay of DCP was best fitted by three exponential terms including a sub-nanosecond rise term, the values of which are quite similar to the parameters obtained for the decay of free dityrosine. They were found to depend on the pH of the aqueous solution but not on the temperature. Analysis by an exponential series method revealed broad fluorescence lifetime distributions for DCP. Compared to the corresponding analysis of dityrosine transients, similar lifetime centers were found whereas the widths of the distributions were found broader for DCP. Molecular dyamics (MD) simulations of dityrosine at 300 K show that chi 1 and chi 2 side chain conformers (rotamers) of both tyrosine subunits interconvert on a picosecond timescale. The rates of interconversion were shown to depend critically upon the MD technique applied: in vacuo simulations yielded lower interconversion rates compared to stochastic dynamics (SD) and full MD (water explicitly included). However, MD simulations of the dityrosine-containing peptide revealed no interconversions of the chi 1 and chi 2 side chain rotamers of both tyrosine subunits within a 400 ps trajectory. Interconversions could be induced by raising the temperature of the system (DCP plus solvent) to 340 K. Side chain rotamers of dityrosine are not stable on a fluorescence time scale but are stable when a dityrosine-containing peptide is regarded. Nevertheless both molecules yield similar fluorescence decay patterns. We therefore conclude that the rotamer model proposed for the fluorescence decay of tyrosine and tryptophan cannot be applied to the fluorescence decay of dityrosine and peptides containing this chromophore. This should be of future interest when dityrosine is used as an intrinsic sensor to study complex dityrosine-containing macromolecules by fluorescence spectroscopy.

Helix-capping interaction in lambda Cro protein: a free energy simulation analysis

The stability mutant Tyr-26-->Asp was studied in the Cro protein from bacteriophage lambda using free energy molecular dynamics simulations. The mutant was calculated to be more stable than the wild type by 3.0 +/- 1.7 kcal/mol/monomer, in reasonable agreement with experiment (1.4 kcal/mol/monomer). Moreover, the aspartic acid in the mutant was found to form a capping interaction with the amino terminus of the third alpha-helix of Cro. The simulations were analyzed to understand better the source of the stability of this helix-capping interaction and to examine the results in light of previous explanations of stabilizing helix caps--namely, a model of local unsatisfied hydrogen bonds at the helix termini and the helix macrodipole model. Analysis of the simulations shows that the stabilizing effect of this charged helical cap is due both to favorable hydrogen bonds with backbone NH groups at the helix terminus and to favorable electrostatic interactions (but not hydrogen bonds) with their carbonyls (effectively the next row of local dipoles in the helix). However, electrostatic interactions are weak or negligible with backbone dipolar groups in the helix further away from the terminus. Moreover, the importance of other local electrostatic interactions with polar side chains near the helix terminus, which are neglected in most treatments of this effect, are shown to be important. Thus, the results support a model that is intermediate between the two previous explanations: both unsatisfied hydrogen bonds at the helix terminus and other, local preoriented dipolar groups stabilize the helix cap. These findings suggest that similar interactions with preoriented dipolar groups may be important for cooperativity in other charge-dipole interactions and may be employed to advantage for molecular design.

13C NMR and fluorescence analysis of tryptophan dynamics in wild-type and two single-Trp variants of Escherichia coli thioredoxin

The rotational motion of tryptophan side chains in oxidized and reduced wild-type (WT) Escherichia coli thioredoxin and in two single-tryptophan variants of E. coli thioredoxin was studied in solution in the temperature range 20-50 degrees C from 13C-NMR relaxation rate measurements at 75.4 and 125.7 MHz and at 20 degrees C from steady-state and time-resolved trp fluorescence anisotropy measurements. Tryptophan enriched with 13C at the delta 1 and epsilon 3 sites of the indole ring was incorporated into WT thioredoxin and into two single-trp mutants, W31F and W28F, in which trp-28 or trp-31 of WT thioredoxin was replaced, respectively, with phenylalanine. The NMR relaxation data were interpreted using the Lipari and Szabo "model-free" approach (G. Lipari and A. Szabo. 1982. J. Amer. Chem. Soc. 104:4546-4559) with trp steady-state anisotropy data included for the variants at 20 degrees C. Values for the correlation time for the overall rotational motion (tau m) from NMR of oxidized and reduced WT thioredoxin at 35 degrees C agree well with those given by Stone et al. (Stone, M. J., K. Chandrasekhar, A. Holmgren, P. E. Wright, and H. J. Dyson. 1993. Biochemistry. 32:426-435) from 15N NMR relaxation rates, and the dependence of tau m on viscosity and temperature was in accord with the Stokes-Einstein relationship. Order parameters (S2) near 1 were obtained for the trp side chains in the WT proteins even at 50 degrees C. A slight increase in the amplitude of motion (decrease in S2) of trp-31, which is near the protein surface, but not of trp-28, which is partially buried in the protein matrix, was observed in reduced relative to oxidized WT thioredoxin. For trp-28 in W31F, order parameters near 1 (S2 > or = 0.8) at 20 degrees C were found, whereas trp-31 in W28F yielded the smallest order parameters (S2 approximately 0.6) of any of the cases. Analysis of time-resolved anisotropy decays in W28F and W31F yielded S2 values in good agreement with NMR, but gave tau m values about 60% smaller. Generally, values of tau e, the effective correlation time for the internal motion, were < or = 60 ps from NMR, whereas somewhat longer times were obtained from fluorescence. The ability of NMR and fluorescence techniques to detect subnanosecond motions in proteins reliably is examined.

Structure and dynamics of the active site gorge of acetylcholinesterase: synergistic use of molecular dynamics simulation and X-ray crystallography

The active site of acetylcholinesterase (AChE) from Torpedo californica is located 20 A from the enzyme surface at the bottom of a narrow gorge. To understand the role of this gorge in the function of AChE, we have studied simulations of its molecular dynamics. When simulations were conducted with pure water filling the gorge, residues in the vicinity of the active site deviated quickly and markedly from the crystal structure. Further study of the original crystallographic data suggests that a bis-quaternary decamethonium (DECA) ion, acquired during enzyme purification, residues in the gorge. There is additional electron density within the gorge that may represent small bound cations. When DECA and 2 cations are placed within the gorge, the simulation and the crystal structure are dramatically reconciled. The small cations, more so than DECA, appear to stabilize part of the gorge wall through electrostatic interactions. This part of the gorge wall is relatively thin and may regulate substrate, product, and water movement through the active site.

Quaternary ligand binding to aromatic residues in the active-site gorge of acetylcholinesterase

Binding sites of Torpedo acetylcholinesterase (EC 3.1.1.7) for quaternary ligands were investigated by x-ray crystallography and photoaffinity labeling. Crystal structures of complexes with ligands were determined at 2.8-A resolution. In a complex with edrophonium, and quaternary nitrogen of the ligand interacts with the indole of Trp-84, and its m-hydroxyl displays bifurcated hydrogen bonding to two members of the catalytic triad, Ser-200 and His-440. In a complex with tacrine, the acridine is stacked against the indole of Trp-84. The bisquaternary ligand decamethonium is oriented along the narrow gorge leading to the active site; one quaternary group is apposed to the indole of Trp-84 and the other to that of Trp-279, near the top of the gorge. The only major conformational difference between the three complexes is in the orientation of the phenyl ring of Phe-330. In the decamethonium complex it lies parallel to the surface of the gorge; in the other two complexes it is positioned to make contact with the bound ligand. This close interaction was confirmed by photoaffinity labelling by the photosensitive probe 3H-labeled p-(N,N-dimethylamino)benzenediazonium fluoroborate, which labeled, predominantly, Phe-330 within the active site. Labeling of Trp-279 was also observed. One mole of label is incorporated per mole of AcChoEase inactivated, indicating that labeling of Trp-279 and that of Phe-330 are mutually exclusive. The structural and chemical data, together, show the important role of aromatic groups as binding sites for quaternary ligands, and they provide complementary evidence assigning Trp-84 and Phe-330 to the "anionic" subsite of the active site and Trp-279 to the "peripheral" anionic site.

A molecular dynamics study of solvent behavior around a protein

The solvent structure and behavior around a protein were examined by analyzing a trajectory of molecular dynamics simulation of the trp-holorepressor in a periodic box of water. The calculated self diffusion coefficient indicated that the solvent within 10 A of the protein had lower mobility. Examination of the solvent diffusion around different atoms of different kinds of residues showed no general tendency. This fact suggested that the solvent mobility is not influenced significantly by the kind of the atom or residue they solvated. Distribution analysis around the protein revealed two peaks of water oxygen: a sharp one at 2.8 A around polar and charged atoms and a broad one at approximately 3.4 A around apolar atoms. The former was stabilized by water-protein hydrogen bonds, and the latter was stabilized by water-water hydrogen bonds, suggesting the existence of a hydrophobic shell. An analysis of protein atom-water radial distribution functions confirmed these shell structures around polar or charged atoms and apolar ones.

Improved molecular dynamics simulations for the determination of peptide structures

In this article a few methods or modifications proven to be useful in the conformational examination of peptides and related molecules by molecular dynamics are illustrated. The first is the explicit use of organic solvents in the simulations. For many cases such solvents are appropriate since the nmr measurements (or other experimental observations) were carried out in the same solvent. Here, the use of dimethylsulfoxide and chloroform in molecular dynamics is described, with some advantages of the use of these solvents high-lighted. A constant allowing for the scaling of the nonbonded interactions of the force field, an idea previously employed in distance geometry and simulated annealing, has been implemented. The usefulness of this method is that when the nonbonded term is turned to zero, atoms can pass through each other, while the connectivity of the molecule is maintained. It will be shown that such simulations, if a sufficient driving force is present (i.e., nuclear Overhauser effects restraints), can produce the correct stereoconfiguration (i.e., chiral center) as well as configurational isomer (i.e., cis/trans isomers). Lastly, a penalty term for coupling constants directly related to the Karplus curve has been implemented into the potential energy force field. The advantages of this method over the commonly used dihedral angle restraining are discussed. In particular, it is shown that with more than one coupling constant about a dihedral angle a great reduction of the allowed conformational space is obtained.

Effect of antibody binding on protein motions studied by hydrogen-exchange labeling and two-dimensional NMR

We have used hydrogen-exchange labeling detected by 2D NMR to study antibody-protein interactions for two monoclonal antibodies raised against horse cytochrome c. The data show that these antibodies bind mainly to the large 37-59 omega-loop of the cytochrome c molecule. In addition, the results provide some suggestive evidence concerning units of local structural flexibility in cytochrome c.

Binding of an antiviral agent to a sensitive and a resistant human rhinovirus. Computer simulation studies with sampling of amino acid side-chain conformation. I. Mapping the rotamers of residue 188 of viral protein 1

The mutation of valine 188 to leucine in the viral protein 1 of human rhinovirus 14 renders the virus resistant to certain antiviral compounds. Thermodynamic-cycle perturbation theory provides a means of calculating the difference in the binding free energies of an antiviral compound to the wild-type virus and to the mutant virus. In calculating the relevant free-energy differences in molecular dynamics simulations, it is important to sample the multiple rotational isomers of residue 188 correctly. In general, these rotamers will not be fully sampled during a single molecular dynamics simulation. However, the contributions of all the rotamers to the free-energy differences associated with mutation of residue 188 may be considered explicitly once they have been identified and their relative free energies determined. Therefore, we describe here the mapping of the rotamers of residue 188 by steric-bump search and energy minimization techniques, and by the computation of potentials of mean force (p.m.f.s.) using umbrella sampling. The usefulness, validity and efficiency of these methods of examining rotameric states is discussed. Adiabatic mapping by energy minimization was found to be unreliable for this residue due to the small magnitude of its interactions with the surrounding protein atoms. Ambiguities in the adiabatic maps were resolved by computing p.m.f.s. The p.m.f. for valine 188 in the unliganded wild-type virus shows a minimum corresponding to the crystallographically observed conformation of valine 188. The p.m.f.s. for valine 188 in the liganded virus and for leucine 188 in the unliganded mutant virus suggest that the experimentally observed conformations may be interpreted as averages of a number of conformations corresponding to those at the minima in the p.m.f.s. The calculations suggest also that the conformation of leucine 188 may change when the ligand binds. The use of the calculated p.m.f.s. to compute the difference in the free energy of binding of an antiviral compound to the wild-type and mutant rhinoviruses is described in the accompanying article.