Comparison of Coupled Motions in Escherichia c oli and Bacillus s

Comparison of Coupled Motions in Escherichia coli and Bacillus subtilis Dihydrofolate ... A common set of residues that play a significant role in the...
0 downloads 0 Views 1MB Size
10130

J. Phys. Chem. B 2006, 110, 10130-10138

Comparison of Coupled Motions in Escherichia coli and Bacillus subtilis Dihydrofolate Reductase James B. Watney and Sharon Hammes-Schiffer* Department of Chemistry, 104 Chemistry Building, PennsylVania State UniVersity, UniVersity Park, PennsylVania 16802 ReceiVed: January 27, 2006

Hybrid quantum/classical molecular dynamics simulations are used to compare the role of protein motion in the hydride transfer reaction catalyzed by Escherichia coli and Bacillus subtilis dihydrofolate reductase (DHFR). These two enzymes have 44% sequence identity, and the experimentally determined structures and hydride transfer rates are similar. The simulations indicate that the tertiary structures of both enzymes evolve in a similar manner during the hydride transfer reaction. In both enzymes, the donor-acceptor distance decreases to ∼2.7 Å at the transition state configurations to enable hydride transfer. Zero point energy and hydrogen tunneling effects are found to be significant for both enzymes. Covariance and rank correlation analyses of motions throughout the protein and ligands illustrate that E. coli and B. subtilis DHFR exhibit both similarities and differences in the equilibrium fluctuations and the conformational changes along the collective reaction coordinate for hydride transfer. A common set of residues that play a significant role in the network of coupled motions leading to configurations conducive to hydride transfer for both E. coli and B. subtilis DHFR was identified. These results suggest a balance between conservation and flexibility in the thermal motions and conformational changes during hydride transfer. Homologous protein structures, in conjunction with conformational sampling, enable enzymes with different sequences to catalyze the same hydride transfer reaction with similar efficiency.

I. Introduction The enzyme dihydrofolate reductase (DHFR) is responsible for the catalytic reduction of 7,8-dihydrofolate (H2F) to 5,6,7,8tetrahydrofolate (H4F) in prokaryotes and eukaryotes. The reaction is carried out in the presence of the coenzyme, nicotinamide adenine dinucleotide phosphate (NADPH). The reaction proceeds via protonation of the N5 nitrogen of the H2F pterin ring and hydride transfer from the nicotinamide ring of NADPH to the C6 of the pterin ring. Tetrahydrofolate is a necessary precursor for the biosynthesis of purines, pyrimadines, and amino acids. As a result, DHFR has been of interest for pharmacological purposes. In addition to its practical applications, DHFR has also been used as a representative system to gain insight into fundamental aspects of enzymes. The role of protein motion in DHFR catalysis has been the subject of interest in a wide range of experimental and theoretical studies.1-23 X-ray crystal structures of Escherichia coli DHFR representing various stages in the reaction pathway show that regions of the protein undergo distinct conformational changes throughout the catalytic process.6 In particular, the Met20 loop (residues 9-24) is known to adopt three conformations, denoted open, closed, and occluded, depending on the stage of catalysis. NMR relaxation experiments have shown that binding of the cofactor and substrate alters the motion of the protein in regions both near and far from the active site.3,8,19 Furthermore, kinetic measurements of single, double, and triple site-specific mutants have identified nonadditive effects that indicate coupling between distal regions in DHFR.4,5,9,23 Measurements of the primary and secondary kinetic isotope effects

and their temperature dependences have also been used as a probe for coupled motion.12,18 Theory and simulation have also contributed to the understanding of protein motion in DHFR catalysis. Classical molecular dynamics simulations of E. coli DHFR have identified regions exhibiting highly correlated and anticorrelated motions in the reactant complex that are not observed in the product complex.7 These regions of the protein were found to coincide with the dynamic regions identified in NMR relaxation experiments.8,19 Mixed quantum/classical molecular dynamics simulations, in conjunction with a genomic analysis across 36 species of DHFR, were used to identify and characterize a network of coupled motions in E. coli DHFR.10,11 These equilibrium motions represent thermally averaged conformational changes of the protein and ligands along the collective reaction coordinate for hydride transfer. Simulations of hydride transfer in E. coli DHFR with other computational methods also exhibited some of these conformational changes.14,16 A detailed analysis of hydrogen bonding, electrostatic potentials, and correlated motions along the collective reaction coordinate for hydride transfer indicated significant changes in these properties throughout the protein.20 In addition, mixed quantum/classical molecular dynamics simulations of mutant E. coli DHFR enzymes illustrated that mutations far from the active site can impact the rate of hydride transfer by introducing nonlocal structural perturbations and altering the conformational sampling of configurations conducive to the chemical reaction.15,23 Other simulations of mutant E. coli DHFR enzymes have also identified nonlocal structural effects21 and alterations in the conformational sampling.17

10.1021/jp0605956 CCC: $33.50 © 2006 American Chemical Society Published on Web 04/29/2006

Coupled Motions in Dihydrofolate Reductase

J. Phys. Chem. B, Vol. 110, No. 20, 2006 10131 to study the slower thermally averaged conformational changes along the collective reaction coordinate. Section II describes the methodology for the hybrid quantum/classical molecular dynamics simulations and the analyses, section III presents the results, and section IV summarizes the conclusions derived from these simulations and analyses. II. Methods

Figure 1. X-ray crystal structures of E. coli (red) and B subtilis (blue) DHFR. The structures were obtained from ref 6 (PDB code 1rx2) and ref 24, respectively, and were aligned with PyMOL.34 The Met20, βFβG, and βG-βH loops, as well as the H4F substrate and NADPH cofactor, are labeled.

These previous theoretical studies have concentrated on DHFR from E. coli. Recently, however, the crystal structure of DHFR from the Bacillus subtilis has been solved.24 B. subtilis DHFR was included in an earlier genomic analysis of 36 species of DHFR, which identified several regions of high sequence conservation.10 E. coli and B. subtilis DHFR have a sequence identity of ∼44%. Figure 1 shows the superposition of the crystal structures of E. coli and B. subtilis DHFR with a CR root-mean-square deviation (rmsd) of 1.53 Å. The tertiary structures are very similar, and the substrate and cofactor are aligned nearly identically. Kinetic measurements of B. subtilis DHFR provide forward and backward pH-independent rates for hydride transfer of 720 and 8 s-1, respectively.24 The corresponding rates for forward and backward hydride transfer in E. coli DHFR are 950 and 0.55 s-1, respectively.25 In light of these similarities in structure and rates, along with the significant differences in sequence, E. coli and B. subtilis DHFR are ideal for a comparative study of coupled motions in enzyme catalysis. To our knowledge, this work represents the first simulation of B. subtilis DHFR in the literature. In this paper, we use hybrid quantum/classical molecular dynamics simulations to study the role of protein motion in B. subtilis DHFR. Our goal is to determine the similarities and differences in the network of coupled motions for E. coli and B. subtilis DHFR and to gain insight into the transferability of these types of motions between different species. We use a covariance correlation analysis to study the fast equilibrium motions about average structures and a rank correlation analysis

A. Hybrid Quantum/Classical Molecular Dynamics. As in our previous simulations of E. coli DHFR,10,11 the N5 position of the substrate pterin ring is assumed to be protonated prior to hydride transfer. Hence, the simulations described in this work focus on the transfer of the pro-R hydride from the NADPH coenzyme to the C6 position of the H3F+ substrate, as depicted in Figure 2. To study this hydride transfer reaction, we used a hybrid quantum/classical molecular dynamics method26,27 that allows for the inclusion of both electronic and nuclear quantum effects, as well as motions from the entire solvated enzyme. Since the details of this method are described elsewhere, we provide only a brief summary in this subsection. To treat the breakage and formation of bonds, the electronic potential energy surface is described using a two state empirical valence bond (EVB) method.28 In VB state 1, the transferring hydride is bonded to the donor, and in VB state 2, it is bonded to the acceptor. A slightly modified version of the GROMOS 43A1 force field29,30 is used for the diagonal elements of the 2 × 2 EVB Hamiltonian matrix. The modifications to the force field include explicit representation of all hydrogens bound to the donor as well as the use of Morse potentials to describe the donor-hydrogen and acceptor-hydrogen bonds. The parameters for the Morse potentials were acquired from previous work.11 The coupling V12 between the two EVB states and a constant energy shift ∆12 between the two states are represented as constant, adjustable parameters. These constants are chosen so that the forward and backward rates calculated from the quantum free energy profiles are consistent with the experimentally determined forward and backward hydride transfer rates. In the case of B. subtilis, the experimentally determined rate constants for forward and backward hydride transfer are 720 and 8 s-1, respectively. Using these rates, the values of V12 and ∆12 were determined to be 113.50 and 264.50 kJ/mol, respectively. Diagonalization of the 2 × 2 EVB Hamiltonian yields the ground-state electronic potential energy surface. The hybrid quantum/classical molecular dynamics method also includes the nuclear quantum effects of the transferring hydride nucleus. The classical subsystem, which includes all nuclei except the transferring hydrogen, is propagated according to Newton’s equations of motion. At each sampled configuration of classical nuclei, the ground-state three-dimensional vibrational wave function for the hydrogen is determined by solving the time-independent nuclear Schro¨dinger equation for the hydrogen moving on the ground-state electronic potential energy surface provided by the EVB method. This approach allows the inclusion of important nuclear quantum effects such as tunneling and zero point energy in the free energy profiles.

Figure 2. Hydride transfer reaction from the NADPH cofactor to the H3F+ substrate.

10132 J. Phys. Chem. B, Vol. 110, No. 20, 2006

Watney and Hammes-Schiffer

The free energy profiles are calculated along a collective reaction coordinate. When the transferring hydrogen nucleus is treated classically, the reaction coordinate is defined as the energy gap between VB states 1 and 2. When the transferring hydrogen nucleus is treated quantum mechanically, the reaction coordinate is defined as the expectation value of the energy gap between VB states 1 and 2 with respect to the ground-state hydrogen vibrational wave function. In the context of transition state theory, the choice of this collective reaction coordinate is justified by previous simulations in which the transmission coefficient was shown to be near unity for E. coli DHFR.11 To enhance sampling of the reaction coordinate, we used a mapping potential defined as28

Vmap(r,R;λ) ) (1 - λ)V11(r,R) + λV22(r,R)

(1)

where r denotes the coordinate of the transferring hydrogen and R denotes the coordinates of all other nuclei. As the mapping parameter, λ, is varied from 0 to 1, the system is driven from reactants to products. For each window specified by λ, molecular dynamics trajectories governed by the mapping potential are propagated. A portion of the total free energy profile is generated for each window using a standard binning procedure. The segments of the free energy profile are combined using the weighted histogram analysis method (WHAM).31 The nuclear quantum effects of the transferring hydrogen are included using a free energy perturbation formula.26 The initial coordinates for the protein, substrate, and cofactor were obtained from the crystal structure.24 The entire system was solvated with 5368 SPCE water molecules in a truncated octahedral box with a distance of 71.9 Å between opposing faces. The equilibration procedure was described in ref 26. All trajectories were propagated in a canonical ensemble at 298 K using a Berendsen thermostat and a time step of 1 fs. A cutoff radius of 14 Å in conjunction with a reaction field was used for all nonbonded interactions. Trajectories for a total of 10 values of the mapping parameter were propagated (λ ) 0.000, 0.050, 0.125, 0.250, 0.375, 0.500, 0.625, 0.750, 0.875, 0.950). Data were collected from each trajectory for 2.1 ns after 100 ps of additional equilibration, yielding a total of 21 ns of data. B. Analysis Tools. Previously, we performed both a covariance analysis and a rank correlation analysis to study motions in E. coli DHFR. These tools provide information about different types of motions in the system. The covariance analysis identifies correlations between fluctuations of pairs of atoms about an average structure. These fluctuations are local to the average structure and typically occur on the femtosecond to picosecond time scale. The rank correlation analysis identifies correlations between thermally averaged interatomic distances along the reaction coordinate and a target model. These motions represent conformational changes that may occur on the millisecond time scale of the hydride transfer reaction. Since these methods have been described in detail elsewhere, here we provide only a brief summary. The covariance analysis of the data is based on the procedure described in ref 20. The correlation of the pairwise atomic fluctuations is represented by an N × N covariance matrix, where N is the number of atomic sites under consideration. The matrix elements of the covariance matrix are given by32,33

Cjk ) 〈(rj - 〈rj〉)•(rk - 〈rk〉)〉

(2)

where the brackets 〈...〉 denote an ensemble average over configurations. Each matrix element is a scalar quantity indicating the average correlation between pairwise atomic fluctuations.

Figure 3. Classical (blue) and adiabatic quantum (red) free energy profiles as functions of a collective reaction coordinate for the hydride transfer reaction catalyzed by B. subtilis DHFR. The minimum reactant free energies of the classical and adiabatic quantum profiles are set to zero, and the difference in free energy barriers is indicated. The dotted lines denote experimental values determined within the context of transition state theory from the forward and backward hydride transfer rates.

Matrix elements that are positive indicate pairwise fluctuations that are in phase (i.e. moving in the same direction), whereas matrix elements that are negative indicate pairwise fluctuations that are out of phase (i.e. moving in opposite directions). Note that a consequence of the dot product in the definition of Cjk is that completely uncorrelated pairwise fluctuations are indistinguishable from orthogonal pairwise fluctuations. In this paper, we provide figures of the covariance cross-correlation matrices with elements given by

Sjk )

Cjk (CjjCkk)1/2

(3)

In our covariance analysis of E. coli and B. subtilis DHFR, all CR atoms from the protein and all heavy atoms from the substrate and cofactor are included in the covariance correlation matrices. The covariance correlation matrices for B. subtilis DHFR were tested for convergence in the same manner as reported previously for E. coli DHFR.20 Specifically, the data set was divided in half, and a covariance matrix was calculated for each half. The covariance matrices constructed from these partial data sets were found to be in good agreement with the covariance matrix of the complete data set. The covariance correlation matrices were calculated for the reactant state (RS), transition state (TS), and product state (PS). The rank correlation analysis of the data is based on the procedure described in ref 23. This approach identifies the correlation between a thermally averaged property along the reaction coordinate and a target model. Ranks rather than numerical values are used for the calculation of these correlations. As a result, this approach is particularly robust and resistant to defects in the data. On the other hand, this approach does not provide direct information about the magnitudes of the changes in properties. In our rank correlation analysis of E. coli and B. subtilis DHFR, all thermally averaged interatomic distances among the heavy (i.e., non-hydrogen) atoms from the protein and the substrate and cofactor are correlated to a monotonically increasing linear target function from the reactant to the transition state. The resulting rank correlation matrix provides a comprehensive statistical analysis of changes in thermally averaged distances along the collective reaction coordinate. To provide information about the magnitudes of the

Coupled Motions in Dihydrofolate Reductase

J. Phys. Chem. B, Vol. 110, No. 20, 2006 10133

Figure 4. Thermally averaged structures for E. coli (red) and B. subtilis (blue) DHFR at the reactant state (RS), transition state (TS), and product state (PS).

changes, we also plot selected thermally averaged distances and angles along the collective reaction coordinate. III. Results Figure 3 depicts the free energy profiles for both the classical and the quantum mechanical treatments of the transferring hydrogen nucleus. Treatment of the transferring hydrogen nucleus quantum mechanically lowers the free energy barrier by 3.3 kcal/mol relative to the classical barrier. This reduction in the free energy barrier is ∼0.9 kcal/mol larger than the reduction observed previously in E. coli DHFR.10,11 Our previous analysis of the splittings between the lowest two hydrogen vibrational states and the energies of these states relative to the barrier along the donor-acceptor axis for the transition state configurations indicated that hydrogen tunneling along the donor-acceptor axis is significant for hydride transfer in E. coli DHFR.11 The present simulations indicate that nuclear quantum effects such as zero point energy and hydrogen tunneling are also important for hydride transfer in B. subtilis DHFR. The sequence and structural alignment scheme in PyMOL34 was used to superimpose the crystal structures of E. coli and B. subtilis DHFR. The rmsd between CR atoms was determined to be 1.53 Å. As shown in Figure 1, the tertiary structures of both enzymes are very similar. Considerable differences are observed in the turn region at residues 85-91. To a lesser extent, the βF-βG and βG-βH loops also exhibit some differences, but these loops are known to be highly mobile.3,8,19 One interesting feature that arises from comparing the structures is the almost perfect superposition of the substrate, the cofactor, and the Met20 loop, which is in the closed conformation for both crystal structures. The thermally averaged structures of E. coli and B. subtilis DHFR for the RS, TS, and PS are depicted in Figure 4. The structures were averaged using configurations selected within 2.5 kcal/mol of the RS, TS, or PS along the collective reaction coordinate. The thermally averaged structures of E. coli DHFR were obtained from previous simulations.20 The rmsd’s between CR atoms for the RS, TS, and PS of the two enzymes are 1.88, 2.03, and 1.87 Å, respectively. These rmsd values from the simulations are only slightly larger than the RMSD between the two crystal structures. These variations are caused mainly by differences in the adenosine binding loops and, to a lesser extent, differences in the βF-βG and βG-βH loops. The CR atom

rmsd between the RS and the TS for B. subtilis DHFR is 0.95 Å, while the corresponding rmsd for E. coli is 0.79 Å. On the basis of these average structures, B. subtilis DHFR does not appear to undergo any major structural changes in the backbone during hydride transfer that were not observed in E. coli DHFR. Previously, a covariance analysis was used to investigate correlated protein motions in E. coli DHFR. Radkiewicz and Brooks found that the reactant and product state ternary complexes of E. coli DHFR exhibit significant differences in the correlated and anticorrelated motions.7 Subsequently, our hybrid quantum/classical molecular dynamics simulations identified regions of the protein that undergo characteristic changes in correlated and anticorrelated motions from the reactant to the transition state.20 The motions that change from the reactant to the transition state are expected to be more directly related to catalysis. We performed the same covariance analysis for B. subtilis DHFR. Figure 5 provides the covariance correlation matrices for E. coli and B. subtilis DHFR at the RS, TS, and PS. The red regions denote highly correlated motions, while the blue regions denote highly anticorrelated motions. Several characteristics of the covariance correlation matrices are similar for E. coli and B. subtilis DHFR. The most prominent similarities are regions of correlated motion that are maintained from the RS to the TS in both enzymes. Many of these motions tend to be located near the diagonal and thus correspond to spatially close pairs of atoms. Several regions further removed from the diagonal are also highly correlated. For example, positions 130-140 (region between loops) and 105-115 (partially R-helix) are both correlated with positions 150-159 (partially β-sheet). These motions may be important for maintaining the structure. Another prominent island of correlated motions is observed between positions 1-20 (including a portion of the Met-20 loop) and 90-120 (including a portion of the βF-βG loop). All of these regions of correlated motions are qualitatively similar in E. coli and B. subtilis DHFR. Of greater interest are the similar changes that occur in both enzymes as they progress from the RS to the TS. The most striking example is the emergence of highly correlated motion between the H2F pterin ring and the NADPH nicotinamide ring at the TS. Although there are similarities between the covariance correlation matrices of E. coli and B. subtilis DHFR, there are also many important differences. The RS of E. coli DHFR exhibits many regions of both highly correlated and highly

10134 J. Phys. Chem. B, Vol. 110, No. 20, 2006

Watney and Hammes-Schiffer

Figure 5. Covariance correlation matrices for E. coli (left) and B. subtilis (right) DHFR for the (a) reactant state, (b) transition state, and (c) product state.

anticorrelated motions. In contrast, the RS of B. subtilis DHFR exhibits fewer correlated and anticorrelated motions. Furthermore, the regions of correlated motions decrease at the TS of E. coli DHFR, whereas the regions of correlated motions increase slightly at the TS of B. subtilis DHFR. Another major difference between the two enzymes is the intraligand and interligand motion. For the E. coli complex, intraligand motions of the NADPH cofactor are highly correlated at the RS and become less correlated at the TS. For the B. subtilis complex, however, the intraligand motions of NADPH are less correlated in the RS and retain a similar degree of correlation at the TS.

In addition, a small region of anticorrelated motions between the adenosine rings and the phosphate linkages of NADPH is observed in the RS of the B. subtilis complex and is absent at the TS, whereas the corresponding regions of the E. coli complex do not exhibit these types of motions. Anticorrelated motions between the substrate pterin ring and the NADPH adenosine ring diminish as the reaction proceeds from the RS to the TS in the E. coli complex, whereas the corresponding regions of the B. subtilis complex do not exhibit these types of motions. Differences are also found in the correlated and anticorrelated protein-ligand and protein-protein motions. For

Coupled Motions in Dihydrofolate Reductase

J. Phys. Chem. B, Vol. 110, No. 20, 2006 10135 TABLE 1: Thermally Averaged Donor-Acceptor Distance for the Reactant State (RS) and Transition State (TS) for E. coli and B. subtilis DHFRa property

E. coli

B. subtilis

CD-CA(RS) CD-CA(TS)

4.01 (0.38) 2.72 (0.05)

3.56 (0.20) 2.69 (0.05)

a Distances are given in angstroms, and the root-mean-square deviations are given in parentheses. The E. coli DHFR data were obtained from previous simulations.20

TABLE 2: Common Set of Residues That Play a Significant Role in the Network of Coupled Motions Leading to Configurations Conducive to Hydride Transfer for both E. coli and B. subtilis DHFRa E. coli

B. subtilis

E. coli

B. subtilis

Leu 28 Met 42 Arg 44 Ile 50 Gly 51 Arg 52 Pro 53 Ser 64 Ser 77 Glu 90

Leu 28 Met 42 Arg 44 Ile 50 Gly 51 Arg 52 Pro 53 Ser 64 Ser 78 Glu 91

Ala 19 Met 20 His 45 Trp 47 Ile 61 Leu 62 Ser 63 Gly 67 Thr 68 Val 72 Trp 74 Arg 98 Val 99 Gln 102 Gln 108 Glu 129 His 149 Ile 155

Asp 19 Leu 20 Lys 45 Phe 47 Val 61 Val 62 Tyr 63 Asp 67 Glu 69 Cys 73 Val 75 Gln 99 Leu 100 Asp 103 Asp 109 Asp 130 Tyr 150 Met 156

a The procedure used to identify this set of residues is described in the text. The residues listed in the first two columns are conserved between E. coli and B. subtilis DHFR, while the residues listed in the last two columns are not conserved.

Figure 6. Rank correlation analysis maps for (a) E. coli and (b) B. subtilis DHFR.

example, the B. subtilis complex exhibits an increase of moderate anticorrelated motions between the cofactor and positions 50-70 (including a β-sheet and a portion of the adenosine binding loop) at the TS relative to the RS, but this trend is not observed in the E. coli complex. In general, these differences in the femtosecond to picosecond fluctuations about the average structure will impact the conformational sampling of the enzymes and thus can alter the slower conformational changes required for efficient hydride transfer. We analyzed the thermally averaged conformational changes along the collective reaction coordinate with the rank correlation analysis method. This method identifies thermally averaged interatomic distances that increase or decrease virtually monotonically between the RS and the TS. A comparison of the results for E. coli and B. subtilis DHFR is provided in Figure 6. The red regions denote distances that increase monotonically, and the blue regions denote distances that decrease monotonically. Accounting for minor differences in numbering of the

residues, many of the red and blue regions are similar for the two enzymes. On the other hand, a number of features are observed in one enzyme but not in the other. These similarities and differences suggest a balance of conservation and flexibility in the thermally averaged conformational changes during hydride transfer. Although the rank correlation analysis provides useful information about the changes in thermally averaged properties along the collective reaction coordinate, it does not provide information about the magnitudes of these changes. To address this issue, we depict selected thermally averaged distances and angles along the collective reaction coordinate in Figure 7. The data for E. coli DHFR were obtained from previous simulations.10 For both E. coli and B. subtilis DHFR, the donoracceptor distance decreases as the reaction progresses from the RS to the TS and then increases as the reaction continues to the PS. As given in Table 1, the average donor-acceptor distance at the TS is 2.69 Å for B. subtilis and 2.72 Å for E. coli DHFR. In both enzymes, the donor-acceptor distance must decrease to ∼2.7 Å to enable hydride transfer. A comparison of other motions provides further insight. Figure 7 shows that the angle between the acceptor carbon and the methylene amino linkage in the substrate increases by ∼3° for E. coli DHFR and by ∼1° for B. subtilis DHFR as the reaction evolves from the RS to the TS. The smaller change in this angle for B. subtilis DHFR may be related to the differences in the Phe-31 motion. As shown in Figure 7, the distance between Cζ of Phe-31 and C11 of the substrate decreases by ∼1 Å as the reaction evolves from the RS to the TS for E. coli DHFR, but this motion is not observed for B. subtilis DHFR. Previously the motion of the Phe-31 was hypothesized to be

10136 J. Phys. Chem. B, Vol. 110, No. 20, 2006

Watney and Hammes-Schiffer

Figure 7. Selected thermally averaged geometrical properties along the collective reaction coordinate for E. coli (left) and B. subtilis (right) DHFR. From top to bottom, the properties plotted are the following: the donor-acceptor distance; the angle between the acceptor and methylene amino linkage in the substrate; the distance between Cζ of Phe-31 and C11 of the substrate; the hydrogen-bonding distance between N of Asp-122 and O of Gly-15 (red) and between O of Ile-14 and carboxamide N of NADPH (blue); the distance between Cδ of Ile-14 and the side-chain oxygen of Tyr-100. The relevant residues are labeled in Figure 8. The E. coli data were obtained from ref 10. The differences in the statistical noise are due to an improved numerical method for calculating the thermally averaged properties in B. subtilis. The B. subtilis results obtained with the numerical method used to generate the E. coli results are provided in ref 36, and the qualitative trends are the same.

coupled to the increase in the angle between the acceptor carbon and the methylene amino linkage in the substrate and the decrease in the donor-acceptor distance.10 The attenuation of this motion in B. subtilis DHFR is a possible explanation for the smaller increase in the angle and the slight decrease in the hydride transfer rate. Note that Phe-31 was found to be tightly conserved in the previous genomic analysis,10 and mutations of this residue have been found to significantly decrease the rate of hydride transfer.35 This motion was also found to be attenuated in previous simulations of the G121V mutant.15 Finally, Figure 7 illustrates that two key hydrogen bonds observed in E. coli DHFR exhibit similar trends in B. subtilis DHFR. The hydrogen bond between N of Asp-122 and O of Gly-15 weakens considerably as the reaction evolves from the reactant to the transition state, but the hydrogen bond between O of Ile-14 and the carboxamide N of NADPH remains intact. These hydrogen bonds form a link that connects the βF-βG loop to the NADPH cofactor, and the associated motions were previously hypothesized to be coupled to the decrease in the donor-acceptor distance.10 To further explore the network of coupled motions, we identified a common set of residues that play a significant role in the network of coupled motions for both E. coli and B. subtilis DHFR. A comparison of residue positions between the two

enzymes was accomplished using the positions provided by the sequence alignment. For each enzyme, we compiled a list of pairwise atomic distances for which the thermally averaged change in the distance from RS to TS was greater than 10%. Only atom pairs involving a CR, a terminal side chain atom, or an atom from the substrate or cofactor were retained. This list of atom pairs was further reduced by taking an intersection of the lists from E. coli and B. subtilis DHFR. Here we define the intersection to be atom pairs involving the same residue position and the same atom position in both enzymes. Note that atom pairs involving terminal side chain atoms from different amino acids at the same residue position were considered to be part of the intersection. In the next step of the procedure, a list of residue-residue pairs was constructed by removing all redundant atom-atom distances involving the same two residues. Finally, we compiled a list of the key residues in the network of coupled motions for both E. coli and B. subtilis DHFR by retaining only those residues that were involved in at least four different residue-residue pairs in both enzymes. The results are given in Table 2. The residues listed in Table 2 were found to move relative to other regions of the protein during the thermally averaged conformational changes conducive to hydride transfer for both E. coli and B. subtilis DHFR. We emphasize that this list of

Coupled Motions in Dihydrofolate Reductase

J. Phys. Chem. B, Vol. 110, No. 20, 2006 10137 IV. Conclusions

Figure 8. Three-dimensional structure of E. coli DHFR. As in ref 10, the residues conserved across numerous species from E. coli to human are indicated by a gradient color scheme (gray to red, where red is the most conserved). NADPH and DHF are in green and magenta, respectively. Reproduced with permission from ref 10.

Figure 9. Three-dimensional structure of E. coli DHFR, with the substrate and cofactor shown in gold. The red spheres indicate the conserved common set of residues that were found to play a significant role in the network of coupled motions leading to configurations conducive to hydride transfer for both E. coli and B. subtilis DHFR. The procedure used to identify this set of residues is described in the text.

residues represents the overlap between the two species and is not complete for either enzyme. Moreover, the residues conserved between E. coli and B. subtilis are likely to be more essential to the motions. These conserved residues are depicted as red spheres in Figure 9. Note that these residues span the region along the cofactor and the substrate. The mutation of Met-42, a conserved residue identified by this procedure, has been found experimentally to dramatically decrease the rate of hydride transfer.9

This paper presents a comparison of the coupled motions in E. coli and B. subtilis DHFR. These two enzymes were found to have 44% sequence identity, and the experimentally determined X-ray crystallographic structures and hydride transfer rates are similar. Our simulations indicate that the root-meansquare deviations among the E. coli and B. subtilis DHFR thermally averaged reactant, transition state, and product structures are ∼2 Å. In both enzymes, the donor-acceptor distance decreases to ∼2.7 Å at the transition state configurations to enable hydride transfer. These results indicate that the tertiary structures of both enzymes evolve in a similar manner during the hydride transfer reaction. Our simulations also show that the inclusion of the nuclear quantum effects for the transferring hydrogen lowers the free energy barrier by 2.4 kcal/mol in E. coli and 3.3 kcal/mol in B. subtilis DHFR. These results, along with analyses of the hydrogen vibrational states for transition state configurations, indicate that zero point energy and hydrogen tunneling effects are significant for both E. coli and B. subtilis DHFR. To elucidate the role of protein motion for hydride transfer in DHFR, we performed both a covariance analysis and a rank correlation analysis of the motions in E. coli and B. subtilis DHFR. The covariance analysis, which identifies correlations between fluctuations of pairs of atoms about an average structure, was performed for the thermally averaged reactant state, transition state, and product state. The rank correlation analysis, which identifies correlations between changes in thermally averaged interatomic distances along the reaction coordinate and a target model, was used to elucidate thermally averaged interatomic distances that increase or decrease monotonically between the reactant and the transition state. These analyses illustrate that E. coli and B. subtilis DHFR exhibit both similarities and differences in the equilibrium fluctuations and conformational changes that are correlated to hydride transfer. In addition, we identified a common set of residues that play a significant role in the network of coupled motions leading to configurations conducive to hydride transfer for both E. coli and B. subtilis DHFR. These observations suggest a balance between conservation and flexibility in the thermal motions and conformational changes during hydride transfer. Thus, the previously proposed network of coupled motions is conserved to some extent across species but is robust and versatile enough to withstand significant changes. This network of coupled motions corresponds to equilibrium, thermally averaged conformational changes along the collective reaction coordinate that facilitate hydride transfer by bringing the donor and acceptor closer together, orienting the cofactor and substrate properly, and providing an appropriate electrostatic environment. These conformational changes are caused by Brownian thermal fluctuations occurring within the confines of the protein structure. As a result, this network of coupled motions is not unique, and there are many different ways to sample configurations that are conducive to hydride transfer. Homologous protein structures, in conjunction with conformational sampling, enable enzymes with different sequences to catalyze the same hydride transfer reaction with similar catalytic efficiency. Acknowledgment. We thank Steve Benkovic for helpful discussions. We are grateful to Steve Benkovic, Alex Horswill, and Guibao Fan for providing the experimentally measured forward and backward rates of hydride transfer and to Mark Wilson and Greg Petsko for providing the coordinates from their X-ray crystallographic structure of B. subtilis DHFR. We thank

10138 J. Phys. Chem. B, Vol. 110, No. 20, 2006 Kim Wong and Pratul Agarwal for helpful discussions and providing the E. coli DHFR simulation data. We are grateful for financial support from National Institutes of Health Grant GM56207. References and Notes (1) Tousignant, A.; Pelletier, J. N. Chem. Biol. 2004, 11, 1037. (2) Benkovic, S. J.; Hammes-Schiffer, S. Science 2003, 301, 1196. (3) Falzone, C. J.; Wright, P. E.; Benkovic, S. J. Biochemistry 1994, 33, 439. (4) Huang, Z.; Wagner, C. R.; Benkovic, S. J. Biochemistry 1994, 33, 11576. (5) Wagner, C. R.; Huang, Z.; Singleton, S. F.; Benkovic, S. J. Biochemistry 1995, 34, 15671. (6) Sawaya, M. R.; Kraut, J. Biochemistry 1997, 36, 586. (7) Radkiewicz, J. L.; Brooks, C. L. J. Am. Chem. Soc. 2000, 122, 225. (8) Osborne, M. J.; Schnell, J.; Benkovic, S. J.; Dyson, H. J.; Wright, P. E. Biochemistry 2001, 40, 9846. (9) Rajagopalan, P. T. R.; Lutz, S.; Benkovic, S. J. Biochemistry 2002, 41, 12618. (10) Agarwal, P. K.; Billeter, S. R.; Rajagopalan, P. T. R.; Benkovic, S. J.; Hammes-Schiffer, S. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 2794. (11) Agarwal, P. K.; Billeter, S. R.; Hammes-Schiffer, S. J. Phys. Chem. B 2002, 106, 3283. (12) Maglia, G.; Allemann, R. K. J. Am. Chem. Soc. 2003, 125, 13372. (13) Garcia-Viloca, M.; Truhlar, D. G.; Gao, J. J. Mol. Biol. 2003, 327, 549. (14) Garcia-Viloca, M.; Truhlar, D. G.; Gao, J. Biochemistry 2003, 42, 13558. (15) Watney, J. B.; Agarwal, P. K.; Hammes-Schiffer, S. J. Am. Chem. Soc. 2003, 125, 3745. (16) Thorpe, I. F.; Brooks, C. L., III. J. Phys. Chem.B 2003, 107, 14042. (17) Rod, T. H.; Radkiewicz, J. L.; Brooks, C. L., III. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 6980. (18) Sikorski, R. S.; Wang, L.; Markham, K. A.; Rajagopalan, P. T. R.; Benkovic, S. J.; Kohen, A. J. Am. Chem. Soc. 2004, 126, 4778.

Watney and Hammes-Schiffer (19) Schnell, J. R.; Dyson, H. J.; Wright, P. E. Biochemistry 2004, 43, 374. (20) Wong, K. F.; Watney, J. B.; Hammes-Schiffer, S. J. Phys. Chem. B 2004, 108, 12231. (21) Swanwick, R. S.; Shrimpton, P. J.; Allemann, R. K. Biochemistry 2004, 43, 4119. (22) Pu, J.; Ma, S.; Gao, J.; Truhlar, D. G. J. Phys. Chem. B 2005, 109, 8551. (23) Wong, K. F.; Selzer, T.; Benkovic, S. J.; Hammes-Schiffer, S. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 6807. (24) Wilson, M. A.; Horswill, A. R.; Fan, G.; Benkovic, S. J.; Ring, D.; Petsko, G. A. Manuscript in preparation. (25) Fierke, C. A.; Johnson, K. A.; Benkovic, S. J. Biochemistry 1987, 26, 4085. (26) Billeter, S. R.; Webb, S. P.; Iordanov, T.; Agarwal, P. K.; HammesSchiffer, S. J. Chem. Phys. 2001, 114, 6925. (27) Billeter, S. R.; Webb, S. P.; Agarwal, P. K.; Iordanov, T.; HammesSchiffer, S. J. Am. Chem. Soc. 2001, 123, 11262. (28) Warshel, A. Computer Modeling of Chemical Reactions in Enzymes and Solutions; John Wiley & Sons: New York, 1991. (29) van Gunsteren, W. F.; Billeter, S. R.; Eising, A. A.; Hunenberger, P. H.; Kruger, P.; Mark, A. E.; Scott, W. R. P.; Tironi, I. G. Biomolecular simulation: The GROMOS96 manual and user guide; VdF Hochschulverlag, ETH Zurich: Zurich, 1996. (30) Scott, W. R. P.; Hunenberger, P. H.; Tironi, I. G.; Mark, A. E.; Billeter, S. R.; Fennen, J.; Torda, A. E.; Huber, T.; Kruger, P.; van Gunsteren, W. F. J. Phys. Chem. A 1999, 103, 3596. (31) Kumar, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; Rosenberg, J. M. J. Comput. Chem. 1992, 13, 1011. (32) Ichiye, T.; Karplus, M. Proteins: Struct., Funct., Genet. 1991, 11, 205. (33) Kitao, A.; Go, N. Curr. Opin. Chem. Biol. 1999, 9, 164. (34) DeLano, W. L. The Pymol Molecular Graphics System; DeLano Scientific: San Carlos, CA, 2002. (35) Chen, J. T.; Taira, K.; Tu, C. P. D.; Benkovic, S. J. Biochemistry 1987, 26, 4093. (36) Watney, J. B. The Pennsylvania State University, 2005.