Electric Fields at the Active Site of an Enzyme: Direct Comparison of Experiment with Theory

See allHide authors and affiliations

Science  14 Jul 2006:
Vol. 313, Issue 5784, pp. 200-204
DOI: 10.1126/science.1127159


The electric fields produced in folded proteins influence nearly every aspect of protein function. We present a vibrational spectroscopy technique that measures changes in electric field at a specific site of a protein as shifts in frequency (Stark shifts) of a calibrated nitrile vibration. A nitrile-containing inhibitor is used to deliver a unique probe vibration to the active site of human aldose reductase, and the response of the nitrile stretch frequency is measured for a series of mutations in the enzyme active site. These shifts yield quantitative information on electric fields that can be directly compared with electrostatics calculations. We show that extensive molecular dynamics simulations and ensemble averaging are required to reproduce the observed changes in field.

The organization of charged and polar groups in the folded state of proteins produces large electric fields that influence nearly every aspect of protein function. Electrostatics calculations suggest that these fields vary markedly from site to site in magnitude and direction (13). The consequence of such variations can be appreciated by considering a dipolar transition state that separates a unit charge over a distance of 1 Å. If this dipole were parallel to a field of 10 MV/cm, not atypical of fields that are estimated to be present in proteins, the energy of the transition state would be lowered by 9.6 kJ/mol; thus, the magnitude and direction of the field could have a substantial effect on the rate of reaction. Simulations on a large number of enzymes support this hypothesis and have shown that preorganized electric fields in active sites contribute substantially to transition-state stabilization (4).

Colorful maps of electrostatic potentials are routinely included in papers describing or analyzing protein structures and are often used to speculate on the electrostatic contributions to a variety of protein functions. However, there are relatively few experiments that can be used to test these calculations directly. Common benchmarks related to changes in free energies, such as pKa shifts (where Ka is the acid dissociation constant), redox potential shifts, or binding constants, depend on a convolution of factors including electrostatic interactions, making direct comparisons to calculated potentials difficult. Similarly, electrostatic interactions can contribute substantially to the effect of mutations on the rates of enzyme-catalyzed reactions, but it is difficult to isolate their contribution.

In contrast, spectroscopic observables that relate directly to electric fields provide a straight-forward connection to calculated potentials (5, 6). The vibrational Stark effect, which describes the effect of an electric field on a molecular vibration, provides a particularly straightforward approach. The vibrational Stark tuning rate gives the sensitivity of a probe vibration to an electric field and can be calibrated by measuring the vibrational Stark spectrum in a known external electric field (7, 8). Once calibrated, a probe vibration acts as a local reporter of its electrostatic environment; its frequency shifts in response to changes in nearby amino acids or protonation states, offering a direct measurement of the electric field change at specific sites. For many vibrations, the effect of an electric field on the stretching frequency is dominated by the projection of Math, the change in dipole moment between the ground and excited states of the probe transition, along the electric field vector (812) Math(1) where Δν̄obs is the observed frequency shift in cm–1, h is Planck's constant, c is the speed of light, and Math is the change in field at the probe location between two states of the protein—for example, two conformational states, two protonation states, or two mutants.

Previous vibrational Stark effect measurements in proteins have focused almost exclusively on the special case of diatomic ligands bound to heme cofactors, most notably carbon monoxide–bound myoglobin (7, 13). Bound CO is an ideal vibrational probe because its frequency is isolated, the transition is intense, and one can obtain high-quality light-minus-dark difference spectra through CO photolysis. We sought to identify vibrational probes that could be applied more generally. The nitrile Embedded Image stretching mode has nearly ideal properties for a vibrational probe, because its frequency is distinct from those of any protein absorption, it is quite intense, Math is relatively large, and Math lies along the nitrile bond in mononitriles, oriented from the nitrogen toward the carbon (810, 12). The well-defined relationship between molecular structure and the direction of Math provides a straightforward connection between observed frequency shifts, Δν̄obs, and changes in field, Math, along the nitrile bond. Aromatic nitriles can be incorporated into proteins through chemical and biosynthetic strategies, and individual nitrile vibrations can be observed in the background of a protein (14, 15).

Here, we used vibrational Stark shifts to measure changes in electric fields caused by mutations in the active site of the enzyme human aldose reductase (hALR2) (Fig. 1). Aldose reductase catalyzes the first step in the polyol pathway, an alternative route of glucose metabolism that has been implicated in a number of diabetic complications (16). hALR2 has been the subject of numerous inhibition and structural studies, culminating recently in a collection of inhibitor-bound structures solved at ultrahigh resolution (17, 18). This wealth of structural data provides an excellent starting point for the design of mutations anticipated to alter the field in the active site and electrostatics calculations at various levels of theory (19). In addition, the nitrile containing inhibitor IDD743 (indexed by the Institute for Diabetes Discovery) (Fig. 1C) had been described and shown to inhibit hALR2 with a median inhibitory concentration (IC50) of 7 nM (18), providing a convenient method of delivering a nitrile probe vibration to the active site in a well-defined orientation. Ultrahigh-resolution crystal structures have been reported for hALR2 with two closely related inhibitors, IDD393 and IDD594 (Fig. 1C), bound in the active site (17, 20). The orientations of IDD393 and IDD594 in these structures are very similar, and the IDD594 structure is of sufficient resolution (0.66 Å) to identify 77% of active-site hydrogens, clearly defining the protonation state of several titratable residues in the active site.

Fig. 1.

Model of the hALR2/NADP+/IDD743 structure. (A) Location of IDD743 in the tertiary structure. (B) Location of side chains mutated in this study relative to the nitrile of IDD743 (green). (C) Chemical structure of IDD743 and closely related inhibitors, IDD393 and IDD594, for which high-resolution structures are available (17, 20).

The vibrational Stark spectrum of the IDD743 nitrile stretch is similar to those obtained for other aromatic nitriles (8, 10), and a fit to the vibrational Stark spectrum yields a magnitude for Math of 0.041 D (fig. S1). Defining the nitrile axis as pointing from carbon toward nitrogen, Math, leading to Eq. 2: Math(2) where Δν̄obs is the cm–1 change in nitrile stretching frequency observed in a mutant of hALR2 relative to the wild-type (WT) protein, ν̄mutantν̄WT, and ΔF is the MV/cm change in the electric field along the nitrile axis between the mutant and wild type (F∥,mutantF∥,WT). According to this calibration, a mutation that causes a change in field of +1.0 MV/cm along the IDD743 nitrile would cause a shift in the nitrile stretching frequency of +0.69 cm–1.

To demonstrate the ability of a nitrile vibration to report on electric field changes caused by mutation, we designed a series of mutants predicted to produce large changes in the field along the nitrile of IDD743. Mutations were initially screened by performing continuum electrostatics calculations on the wild-type protein and a large number of potential mutants [input structures were modeled from the available IDD393 and IDD594 structures (21)]. For the V47→D47 (V47D) mutation (22), continuum calculations performed with the internal dielectric set to 4 predicted large changes in both electric potential at the nitrile and electric field projection along the nitrile (Fig. 2), giving rise to a predicted Δν̄ of –6.8 cm–1. From the mutations screened, eight were chosen for actual mutagenesis and infrared experiments, two of which made direct contacts with the inhibitor but were predicted to minimally perturb the field at the nitrile (Y48F and H110A), three of which were predicted to cause a measurable positive change in field at the nitrile (W20Y, V47N, and F121E), and three of which were predicted to cause a measurable negative change in field at the nitrile (V47D, Q49R, and K77M) (table S3). On the basis of the IDD393 and IDD594 structures, each mutated residue is expected to lie within 11 Å of the IDD743 nitrile, which is positioned at the active-site cleft of the enzyme. The selected mutations included buried (Y48F, K77M, and H110A), partially exposed (W20Y, V47D, and V47N), and fully solvated (Q49R and F121E) sites, and each either added or removed substantial charge or polarity.

Fig. 2.

(A) Result of continuum electrostatics calculation, with the interior dielectric set to 4, on the initial model of wild-type hALR2/NADP+/IDD743. Electrostatic potentials are plotted on the surface of the enzyme from –10 kT/e (red) to +10 kT/e (blue), where k is the Boltzmann constant, T is temperature, and e is the charge of the electron. (B) Potentials from (A) plotted on a two-dimensional slice running through the nitrile of IDD743 (green). (C) Expanded view of slice in (B) in the vicinity of the nitrile probe. (D) Result of continuum electrostatics calculation, with the interior dielectric set to 4, on the initial model of V47D hALR2/NADP+/IDD743. Electrostatic potentials are plotted on the same slice as in (C). The calculated change in field along the nitrile bond for the V47D mutation is –3.86 kT/eÅ relative to the wild-type protein (1 kT/eÅ = 2.57 MV/cm at 298 K).

Wild-type and mutant proteins were expressed and purified as described in the literature (23). Proper folds were confirmed by kinetic assay or by circular dichroism spectroscopy, and cofactor binding was confirmed by kinetic assay or by the spectral shift in the 340-nm absorption band of reduced nicotinamide adenine dinucleotide phosphate (NADPH) upon binding (21). Samples for infrared spectroscopy were prepared as 1:1:1 complexes of hALR2:NADP+:IDD743, the state of the enzyme present in the IDD393 and IDD594 crystal structures [for details of mutant characterization, sample preparation, and infrared absorption experiments, see (21)]. Although the nitrile absorption band widths (Fig. 3A) were larger than the observed frequency shift for the mutations studied, frequency shifts measurable with our experimental resolution (0.5 cm–1) were observed for several mutants (table S2). The V47D mutation shifts the nitrile stretching frequency by –1.6 cm–1, corresponding to a change in field along the nitrile bond of –2.3 MV/cm relative to the wild-type protein. This change in field is much smaller than that predicted from continuum electrostatics calculations on the initial model for V47D (–9.9 MV/cm). In general, the Stark shifts observed in the mutant spectra correlated poorly with predictions based on substitution of the calculated fields into Eq. 2 (Fig. 4).

Fig. 3.

(A) Infrared absorption spectra of IDD743 bound to wild-type (black) and V47D (red) hALR2 in the nitrile stretch region. Absorption spectra were scaled to an absorbance of 1.0 at the absorption maximum (typical absorption maxima were between 0.001 and 0.002). Nitrile stretching frequencies for the mutants studied were as follows (in cm–1): wild type (2241.6), W20Y (2240.9), V47D (2240.0), V47N (2238.5), Y48F (2240.2), Q49R (2241.2), K77M (2239.4), H110A (2242.1), and F121E (2240.7) (B) Average trajectory of electric field along the nitrile bond of IDD7434 bound to wild-type (black) and V47D (red) hALR2. Each point is the average of continuum electrostatics calculations, with the interior dielectric set to 2, performed on multiple MD-generated structures. The predicted change in field relaxes from greater than –5 to –1.75 kT/eÅ in the first 200 ps of the average trajectory. (C) Normalized distribution of electric field along the nitrile of IDD743 for wild-type (black) and V47D (red). Electric field values were obtained from continuum electrostatics calculations, with the interior dielectric set to 2, performed on structures obtained from MD trajectories at time points greater than 200 ps. Dotted lines in (B) and (C) represent the mean values of the distributions provided in (C).

Fig. 4.

Correlation between observed changes in frequency and calculated changes in electric field along the nitrile bond of IDD743, ΔF, for all mutants relative to wild-type: H110A (◯), Q49R (▢), W20Y (△), F121E (♢), Y48F (⚫), V47D (◼), K77M (▲), and V47N (♦). (A) ΔF calculated using initial structure models, with the interior dielectric set to 4. (B) ΔF calculated using energy-minimized structure models, with the interior dielectric set to 4. (C) ΔF calculated from the mean of electric field distributions. Distributions were generated by performing continuum electrostatics calculations on structures obtained from MD trajectories at time points greater than 200 ps. Each individual continuum electrostatics calculation was performed with the internal dielectric set to 4. (D) ΔF calculated from the mean of field distributions as in (C), but with the internal dielectric of each continuum electrostatics calculation set to 2. The dashed lines in (A) to (D) represent the linear correlation expected from Eq. 2.

Because the discrepancy between observed and expected frequency shifts was greatest for the two mutations closest to the nitrile, V47D and V47N, we explored the effect of side chain conformation on the calculated fields for these mutations. Conformational flexibility is a common problem in using continuum electrostatics calculations to predict changes in potential or field caused by mutations or protonation events (24). Uncertainty in side chain conformations is often modeled by choosing a slightly higher value for the internal dielectric, but such dielectric scaling is not expected to remove local errors (25, 26). This problem was most evident for the V47N mutation. Here, reversing the orientation of the side chain amide led to a sign reversal of the predicted change in field, a particularly notable example of the importance of side chain orientation.

To evaluate the consequences of uncertainty and appropriate averaging in side chain positions, we used a large set of independent molecular dynamics (MD) simulations to relax the initial models toward an ensemble of structures more representative of equilibrium conformations. Techniques that combine MD simulations with electrostatics calculations have become increasingly popular and range in complexity from combined quantum mechanical and molecular mechanical methods to approximations of the continuum model such as the generalized Born method (2). The choice of method depends on the quality of the structural data available, the size of the system investigated, and the time scale of the motions of interest. Both structural reorganization and the consideration of a diverse thermal ensemble are crucial. Because structural data are not yet available for the mutations studied (and even with structures, the data are often not of sufficient resolution, especially for side chain orientations) we chose a simulation method that allows the long time-scale dynamics that might be required to observe side chain reorganization to the equilibrium ensemble. Dynamics simulations also allow an estimate of the distribution in electric field produced by thermal fluctuations about an equilibrium structure. This distribution of calculated electric fields should better represent the time-averaged environment reported in the infrared spectra.

To obtain statistically meaningful values of the time-averaged fields, we simulated multiple MD trajectories for up to 500 ps. Each initial model was energy minimized, and 50 independent implicit-solvent MD simulations were initiated. Trajectories were calculated with the distributed computing platform Folding@Home (27) with structures saved every 0.5 ps (21). After the completion of these simulations, continuum electrostatics calculations were performed on each saved structure, providing the electric potentials and fields as a function of time along each trajectory. The potentials and fields for all structures were calculated using an interior dielectric value of 4 and repeated using an internal dielectric of 2 to test the dependence of the fields on this parameter. For each minimized mutant structure, we performed more than 5000 continuum calculations at time steps along multiple trajectories and then calculated the time-dependent ensemble average projections of the electric field along the nitrile bond (Fig. 3B).

The importance of relaxing initial models is evident from the predicted change in field as a function of time along these trajectories. Ensemble average field trajectories of this type were generated for each mutant (fig. S3) and the predicted field was stable within 200 ps for all mutants except V47N. For the V47N mutation, a very slow relaxation toward more negative field values was observed through 500 ps. Inspection of individual trajectories identified the flip of the Asn47 side chain amide as the origin of this slow relaxation, with the conformation aligning the side chain amide dipole with the nitrile bond vector accumulating in time. Comparative analysis of five additional V47N ensembles, each starting from a different initial conformation, suggested that the low-field Asn rotamer is favored at equilibrium. MD simulations initiated from models with the Asn47 side chain amide in the favored conformation led to stable values of the field by 100 ps (fig. S3). The orientation of side chain amides is a major source of uncertainty in electrostatics calculations because the vast majority of structures are of insufficient resolution to define their orientation. The Asn47 trajectories suggest that the barrier to rotation for side chain amides can be high, and that the orientation of side chain amides can have a profound effect on local fields. Thus, the observed negative frequency shift for the V47N mutation provides a strong prediction for the relative orientation of the Asn47 amide when IDD743 is bound.

To obtain relaxed values of the field along the IDD743 nitrile, we collected more than 2000 structures from trajectory time points longer than 200 ps for the wild-type protein and for each mutant; we then performed continuum calculations, first with an interior dielectric of 4 and then with an internal dielectric of 2, for each structure, generating field distributions past 200 ps (21) (Fig. 3C). The mean values of these distributions were then used to calculate the change in field between mutants and the wild-type protein (28). Using bootstrap analysis, we estimated that the statistical error for changes in field calculated from these distributions is at most 0.05 MV/cm, although we expect additional systematic error in the calculated fields to arise from MD and continuum electrostatics parameters that have not been optimized to reproduce electric field values.

The substantial improvement in the correlation relative to that obtained from the initial models demonstrates the importance of relaxing side chain conformations (a process that requires hundreds of picoseconds or more). Perhaps the most important aspect of the calculated distributions of field is that the width of these distributions is greater than the difference in means for any two mutants (fig. S4). These widths demonstrate that the sensitivity of the field to small structural changes can be much greater than the average change in field caused by mutation. For several of the mutants studied, the fields predicted from initial models and energy-minimized structures derived from those models fell in the far wings of the later calculated distribution, resulting in the poor frequency-to-field correlation (Fig. 4, A and B). Attempts to accurately calculate the time-averaged field from a single structure rely on obtaining a structure that represents the mean of these distributions, a challenge even when good structural data are available.

Continuum models are one of the most popular methods of calculating electrostatic properties of macromolecules. For many applications, the property of interest is a state of the system for which complete structural information is not available—for example, the effect of a mutation, side chain protonation, or the modification of functional groups in a bound ligand or inhibitor. Our results reinforce the importance of MD in calculating average electric fields or potentials. The observed frequency shifts also demonstrate the utility of nitrile vibrations as probes of electrostatic fields in proteins. Because obtaining high-quality nitrile absorption spectra in protein solutions is straightforward (21), this functional group should be widely applicable as a probe of electrostatics. For hALR2, the approach could be readily extended to measure changes in field at solvent-excluded sites in the inhibitor binding pocket, where several mutations are predicted to produce much larger frequency shifts, including mutations at residues believed to be important for inhibitor specificity. Nitriles are a relatively common substituent in inhibitors and marketed drugs, and this method provides a strategy to probe the active sites of a large number of enzymes. Techniques to incorporate nitriles by means of nonnatural amino acids or chemical modification of amino acids have also been demonstrated and should allow similar measurements in a variety of protein systems. Direct comparison between experimental and computational field distributions should lead to improved electrostatics calculation methodology.

Supporting Online Material

Materials and Methods

Figs. S1 to S4

Tables S1 to S3


Databases S1 to S3

References and Notes

View Abstract

Navigate This Article