New Insights into the Reaction Mechanism Catalyzed by the

Publication Date (Web): February 8, 2007. Copyright © 2007 American ... Katie L. Whalen , Kenneth B. Tussey , Steven R. Blanke , and M. Ashley Spies...
0 downloads 0 Views 563KB Size
J. Phys. Chem. B 2007, 111, 2385-2397

2385

New Insights into the Reaction Mechanism Catalyzed by the Glutamate Racemase Enzyme: pH Titration Curves and Classical Molecular Dynamics Simulations Eduard Puig,† Mireia Garcia-Viloca,*,‡ A Å ngels Gonza´ lez-Lafont,†,‡ Jose´ M. Lluch,†,‡ and Martin J. Field§ Departament de Quı´mica and Institut de Biotecnologia i de Biomedicina, UniVersitat Auto` noma de Barcelona, 08193 Bellaterra, Barcelona, Spain, and Laboratoire de Dynamique Mole´ culaire, Institut de Biologie Structurale Jean-Pierre Ebel, CEA/CNRS/UJF UMR-9075, 41 Rue Jules Horowitz, F-38027 Grenoble Cedex 1, France ReceiVed: September 27, 2006; In Final Form: December 14, 2006

The mechanism of the reactions catalyzed by the pyridoxal-phosphate-independent amino acid racemases and epimerases faces the difficult task of deprotonating a relatively low acidicity proton, the amino acid’s R-hydrogen, with a relatively poor base, a cysteine. In this work, we propose a mechanism for one of these enzymes, glutamate racemase (MurI), about which many controversies exist, and the roles that its active site residues may play. The titration curves and the pK1/2 values of all of the ionizable residues for different structures leading from reactants to products have been analyzed. From these results a concerted mechanism has been proposed in which the Cys70 residue would deprotonate the R-hydrogen of the substrate while, at the same time, being deprotonated by the Asp7 residue. To study the consistency of this mechanism classical molecular dynamics (MD) simulations have been carried out along with pK1/2 calculations on the MD-generated structures.

1. Introduction Much of the catalytic power of enzymes comes from the precise positioning of key amino acid residues that act as acids or bases in the chemical steps during which proton transfers occur. Equally important, these catalytic residues must be present in the right ionization state that the proton-transfer event requires under physiological conditions. The protonation state of acidic and basic residues of a protein along with the protonation state of any substrate will be influenced by their interactions with the environment and among themselves. Well-known examples of the perturbation of the aqueous solution pKa values of enzymatic residues by ∼5 pKa units in the active site of several enzymes involve Asp, Glu, and Lys residues but less acidic residues also have been implicated.1 As a consequence, the group of residues acting as general bases in enzymatic reactions includes not only the expected candidates, Glu, Asp, and His, but also Tyr, Cys, Lys, and even Arg, Thr, and Ser. Moreover, it has recently been demonstrated that the active site of an enzyme can be predicted by the existence of proximal residues that deviate from the aqueous solution acid/base behavior.2 Thus, understanding the catalytic function of an enzyme (especially when rate-limiting proton transfers are present) rests in part on understanding the environmental determinants of pKa values. The above considerations are particularly well illustrated by the glutamate racemase enzyme and other cofactor-independent amino acid racemases. All known amino acid racemases are found in bacteria where they convert L-amino acids into D-amino acids during cell wall biosynthesis.3 The inhibition of these * Author to whom correspondence should be addressed. E-mail: mireia@ bioinf.uab.es. † Departament de Quı´mica, Universitat Auto ` noma de Barcelona. ‡ Institut de Biotecnologia i de Biomedicina, Universitat Auto ` noma de Barcelona. § Institut de Biologie Structurale Jean-Pierre Ebel.

enzymes prevents this synthesis and leads to bacterial cell lysis and death. In Streptococcus pneumoniae,4 Escherichia coli,5 Staphylococus haemeolyticus,6 and Bacillus subtilis,7 the presence of a functional racemase enzyme has been shown to be essential for cell viability, and so these enzymes constitute a possible target for the design of new antibacterial agents.8,9 The racemase reaction requires an initial deprotonation of the substrate’s R-proton (a relatively low acidic proton), followed by a reprotonation on the opposite face of the resulting planar anionic intermediate. This is a difficult operation that nature has overcome using different chemical strategies. One of these is through the use of the pyridoxal phosphate (PLP) cofactor, as in the case of alanine racemase. However, some amino acid racemases operate in a cofactor-independent manner. In the general mechanism assumed for the reactions catalyzed by the PLP-independent amino acid racemases, two optimally located cysteines take the role of catalytic acid/base residues in the two proton transfers required to invert the amino acid’s stereochemistry.10,11 These two cysteines are able to interchange their acid/base roles depending on which enantiomer binds the active site. This is a particularly interesting case from a chemical point of view, as it implies the deprotonation of a relatively low acidic proton, the R-proton of the amino acid, by a relatively weak base such as a cysteine. Although the acid/base catalytic roles of the cysteines have been demonstrated, it is not clear how their acid/base properties are modified by their environments to make the assumed mechanism possible. Enzymes operating in this cofactor-independent fashion include glutamate racemase,10 proline racemase,12 aspartate racemase,13 and diaminopimelate epimerase.14 To our knowledge, very few crystallographic structures of the PLP-independent amino acid racemases are available. In particular, only structures of glutamate racemase (MurI) and diaminopimelate epimerase enzymes are found in the Protein

10.1021/jp066350a CCC: $37.00 © 2007 American Chemical Society Published on Web 02/08/2007

2386 J. Phys. Chem. B, Vol. 111, No. 9, 2007

Figure 1. Scheme of the active site of MurI in a given ionization state.

Data Bank (PDB). In the case of glutamate racemase, until recently, only the structure of the Aquifex pyrophillus enzyme had been resolved,15 which together with biochemical data from the Lactobacillus fermenti enzyme10 provided the first insights into the structure/function relationships of this enzyme and confirmed the general trends of the proposed mechanism. The active site of the glutamate racemase enzyme from A. pyrophillus with the D-Glu substrate is shown schematically in Figure 1. It can be seen that there are several ionizable residues in the active site, and so different active site protonation states may act during the catalytic process. The crystallographic structure did not provide a complete description of the mode of substrate binding because the enzyme was cocrystallized with an inhibitor. A recent computational study carried out by some of us,16 together with previous work of Mo¨bitz and Bruice,17 demonstrates that the substrate can adopt different binding modes that are dynamically stable. During the course of the research reported in this paper, the threedimensional structure of B. subtilis glutamate racemase (RacE) complexed with the substrate was resolved.18 Although our theoretical studies are based on the crystallographic structure of the A. pyrophillus enzyme, we will compare the results obtained with the recently resolved substrate-enzyme complex structure. In A. pyrophillus glutamate racemase, Cys70 is the catalytic base, and Cys178 is the catalytic acid for the racemization of the D-glutamate substrate. These cysteines reverse their roles when the substrate is L-glutamate. The catalytic function of these cysteines has been experimentally confirmed by several assays,19-21 but there is no definitive explanation for the ability of the enzyme to abstract the R-proton, as it remains unclear how the enzyme manages to stabilize the carbanionic intermediate. Recently, Mo¨ritz and Bruice17 have pointed out a major problem of the mechanism shown in Figure 1, which is that the cysteine thiolate, instead of abstracting the R-proton leading to the catalytic reaction, might prefer to abstract the more acidic protons of the ammonium group in the substrate main chain. However, our recent theoretical results support the proposal by Richard and co-workers22-25 that the positive charge on the ammonium group stabilizes the negative charge developed on the R-carbon of the ylide intermediate after the R-proton deprotonation.26 Information on the roles of other active site residues remains sparse. The residue Glu147 is strictly conserved among the 13

Puig et al. available amino acid sequences of MurI,21,15 but the experimental information obtained for the different species studied up to now suggests that its role in the catalytic mechanism might be different in each of them. The crystallographic structure of the A. pyrophillus enzyme15 indicates that the biologically active form of the enzyme is a dimer because a monomer’s active site contains residue Glu147′ from an adjacent monomer. (From now on the prime symbol denotes a residue belonging to this second monomer.) On the basis of a 100-fold decrease in the kcat value for the L f D direction in the E147′N mutant, Hwang et al.15 proposed that Glu147′ assists the Cys178 thiol in abstracting the R-proton of L-Glu. Their results for KM in the D f L direction also indicate that Glu147′ may participate in the binding of the D-substrate. On the contrary, the L. fermenti enzyme has been reported to be active as a monomer on the basis of size exclusion chromatography,21 and its Glu152 (Glu147′ in A. pyrophillus) residue does not seem to have a catalytic role. The recently reported crystallographic structure of B. subtilis RacE has a trimeric asymmetric unit. A dimer is identified as the biologically active species although the structure shows that the active site is not in the interface between two monomers.18 A comparison of this structure with the structure of the A. pyrophillus enzyme reveals important conformational changes that involve different positions of Glu153 (Glu147′ in A. pyrophillus) and that have the net effect of reducing the exposed solvent-accessible area of the substrate in the B. subtilis structure. In the modeled MurI (A. pyrophillus)-D-Glu structure of Hwang et al.,15 there is a hydrogen bond interaction between Glu147′ and the main chain substrate carboxylate. Our previous work, based on the A. pyrophillus crystal structure, indicates that (1) the hydrogen bond interaction between Glu147′ (modeled as protonated) and the main chain R-carboxylate group of the substrate is stable, as it is maintained throughout our molecular dynamics (MD) simulations of the Michaelis complex,16 and (2) the mechanism proposed by ref 22, in which the more abundant zwitterionic glutamate binds to the enzyme and is protonated on its R-carboxylate by Glu147′, is plausible.26 The results of our second study demonstrate that this proton transfer lowers the pKa of the R-hydrogen thereby facilitating its abstraction by the cysteine.26 Thus, from our previous results, one might speculate both catalytic and binding roles for Glu147, in agreement with the proposal of the authors of the A. pyrophillus crystallographic structure. In contrast, in the resolved RacE (B. subtilis)-D-Glu complex the side chain of this residue (Glu153) forms a long water-mediated contact with the side chain of the glutamate substrate.18 This difference led the authors to support the idea that Glu153 is more likely to have a binding role than a catalytic role in the B. subtilis variant,18 in agreement with Glavas and Tanner’s mutagenesis results for the L. fermenti enzyme.21 The structure of the A. pyrophillus enzyme (Figure 1) shows that Asp7, which belongs to the group of strictly conserved residues,21 is in the active site. Glavas and Tanner21 proposed that one of the possible roles of the Asp7 residue could be to lower Cys70’s pKa by stabilizing its thiolate form by means of hydrogen bonding (Figure 2a). However, our classical MD simulations with a neutral Asp7 residue showed that it could be important for substrate binding, by hydrogen bonding to the substrate δ-carboxylate,16 and not for its stabilization of the Cys70 thiolate. In fact, the idea that Asp10 (the corresponding residue in the L. fermenti variant) could be functioning as an acid that protonates or hydrogen-bonds to the δ-carboxylate of D-glutamate was not excluded by Glavas and Tanner. An alternative possibility, also proposed by them21 (Figure 2b), is

New Insights into the Glutamate Racemase Mechanism

J. Phys. Chem. B, Vol. 111, No. 9, 2007 2387

Figure 2. Mechanisms proposed to assist the catalytic cysteines: (a) mechanism A, stabilization of a thiolate by hydrogen bonding, or (b) mechanism B, general base assistance involving a neutral thiol.

that Asp7 could assist the deprotonation of D-glutamate by participating in general base catalysis with neutral Cys70 thiol. As mentioned above, we studied the mechanism schematized in Figure 2a in a previous work,26 but we had not explored the second mechanism (Figure 2b). Finally, the recently resolved B. subtilis RacE structure shows that the Cys74 thiol group is close to the side chain of Asp10 (3.9 Å), which led the authors to support the mechanism in which RacE utilizes a neutral cysteine thiol.18 There are also different proposals about the role that another strictly conserved residue of the active site, His180, could play in the mechanism.15,21 Glavas and Tanner’s mutagenesis studies led them to conclude that His186 (His180 in A. pyrophilus) may assist Cys184 (Cys178 in the A. pyrophilus variant). In the recently resolved structure of the B. subtilis variant18 the side chain of His187 is close to Cys185 (4.2 Å) and simultaneously forms a hydrogen bond with Glu153. This is consistent with both catalytic and binding roles of His180. In summary, all this controversy seems to indicate that, depending on the species, both Glu147 and His180 could assist Cys178 in the L f D direction, which is not surprising in view of (1) the variability of the quaternary structure between the members of this family and (2) the fact that Glu147′ belongs to a different monomer in the active site of the A. pyrophilus enzyme whereas it is part of the same monomer as the rest of the active site residues in the B. subtilis variant. In summary, several studies on MurI have shown that the protonation state of the active site residues is of crucial importance to catalysis. Although a recent computational work2 has determined the titration curves of the ionizable residues of the A. pyrophillus enzyme, that research was aimed at locating the enzyme active site by identifying clusters of residues with perturbed titration curves. This they did, thereby confirming the location of the MurI active site. However, their results are not conclusive from a mechanistic point of view, as the analysis was done for the X-ray structure of the enzyme without the presence of the substrate. The purpose of the work in this paper was to calculate and analyze the titration curves of several residues of the glutamate racemase enzyme in both the presence and the absence of the D-Glu substrate and for different points along the reaction mechanism schematized in Figure 2a. Our results serve to achieve a qualitative understanding of the catalytic role of the active site residues. In addition, classical MD calculations were carried out to confirm the feasibility of the mechanistic proposals established from the titration curves, which also involve the reaction mechanism shown in Figure 2b.

The paper is organized as follows: First, in section 2, we outline the statistical methods used to calculate and analyze the titration curves, and then we describe the specific computational details of the study in section 3. This is followed, in section 4, by a discussion of the results and, finally, a summary of the main ideas that we have obtained from this and previous studies on the glutamate racemase enzyme in section 5. 2. Theoretical Background The acid/base residues involved in an enzymatic mechanism exist in equilibrium between their protonated and unprotonated forms. The relative stabilities of the two forms can be determined as follows at a particular pH. Consider a macromolecule with N ionization sites in equilibrium with protons in the surrounding solution. Each site is capable of gaining or losing one proton. Let us assume that the reference state of the protein is one in which all sites are neutral. The ionization state of site i is specified by xi, where xi ) 1 or 0 depending on whether this site is ionized or neutral. Since a given site is classified as an acid (neutral and deprotonated forms are possible) or as a base (neutral and protonated forms are possible), the total number of possible ionization microstates of the macromolecule is 2N, and an ionization microstate, k, can be conveniently represented by a N-dimensional vector b xk ) (xk1, xk2, ..., xkN). Each ionization microstate at a given pH and temperature, T, is characterized by its free energy, Gk, and the excess of charge caused by the N ionization process, nk ) ∑i)1 γixki , where γi is +1 for bases and -1 for acids. The equilibrium ionization probability of site i averaged over all microstates is given by27 2N

∑ xki exp(-βGk)Λn

k

θi )

k)1

(1)

Z

where β ) (kBT)-1, Λ ) 10-pH, and the macrocanonical partition function Z is given by 2N

Z)

exp(-βGk)Λn ∑ k)1

k

(2)

The ionization process is represented as the gain or the loss of a proton at the titrating site, and it is assumed that the only significant contributions to the ionization free energy come from electrostatic interactions (Gk ) Gkelec). The electrostatic free energy of a given ionization microstate in the macromolecule

2388 J. Phys. Chem. B, Vol. 111, No. 9, 2007

Puig et al.

is usually calculated with respect to reference states for each titrating site. These are normally model compounds, such as single amino acid residues.28 The pKa model of a model compound in aqueous solution corresponds to the pH where θi ) 1/2. In the special case of non-interacting sites, the free energy of microstate b xk can be written as N

Gkelec ) -kBT ln 10

int γixki pK a,i ∑ i)1

(3)

where pK int a,i is known as the intrinsic pKa of site i because it only takes into account two of the three electrostatic charging terms: the Born solvation energy of the charge at site i in the macromolecule and the interaction of the charge at site i with nontitrating “background” charges (e.g., the peptide group dipoles).29,28 At this stage, the electrostatic interactions between titrating sites are not included. pK int a,i is a useful intermediate quantity as it is the pKa of the site when all other titrable sites in the macromolecule are neutralized. Equation 1 can then be rewritten as 2N

θi )

∑ k)1

N

xki exp(ln 10

int γjxkj pK a,j - nk pH ln 10) ∑ j)1

(4)

Z

where Z is 2N

Z)

∑ k)1

N

int γjxkj pK a,j - nk pH ln 10) ∑ j)1

exp(ln 10

(5)

The additive property of Gkelec in eq 3 allows one to group terms in eq 2 and rewrite the ionization probability of each site as

10γi(pK a,i -pH) int

θi )

(6)

1 + 10γi(pK a,i -pH) int

which is algebraically equivalent to the Henderson-Hasselbalch (HH) equation describing sigmoidal titration curves. However, when sites interact, the additive property of eq 3 does not hold, and θi is no longer given by the simple formula xk of eq 6. If we suppose that Gkelec of a ionization microstate b can be written in terms of intrinsic contributions and contributions arising from pairwise interactions between the titrable groups,28 then we obtain N

Gkelec ) -kBT ln 10

∑ i)1

int γixki pK a,i +

1

N

N

γiγjxki xkj Wij ∑ ∑ 2 i)1 j*i

(7)

where Wij is the interaction free energy between the ionized forms at sites i and j. The free energy in eq 7 can then be used in eqs 1 and 2 to calculate the equilibrium ionization probability of site i

θi )

1

2N



Z k)1

(

nk pH ln10 where

N

xki exp ln 10 β 2

int γjxkj pK a,j ∑ j)1 N

N

∑ ∑ j)1 l*j

)

γjγlxkj xkl Wjl

(8)

2N

Z)

∑ k)1

(

N

exp ln 10

int γjxkj pK a,j - nk pH ln 10 ∑ j)1

β

N

N

)

∑ ∑ γjγlxkj xkl Wjl 2 j)1 l*j

(9)

The summation in eqs 8 and 9 contains 2N terms and is impractical to compute when N is large. For this reason, we adopt the following common strategy. First those ionization states that make no significant contributions to the summations are selected.30 If the pH is far above (or below) that at which a particular site titrates, then that site will almost always be unprotonated (or protonated), and ionization microstates b xk having that site protonated (or unprotonated) will be negligible in the summations. At each pH then, we first precalculate the maximum and minimum possible θi values for each site. The maximum (θmax,i) and minimum (θmin,i) possible ionization probabilities of site i correspond to having the largest number of sites ionized with an opposite charge or with a charge of the same sign, respectively. If θmax,i is less than (or equal to) 0.05, then the site is taken as entirely un-ionized, and if θmin,i is greater than (or equal to) 0.95, then the site is regarded as entirely ionized.30 Once the fixed and the variable sites have been identified for a given pH, the intrinsic free energies of the variable sites are calculated taking into account the influence of the fixed charges. At each pH there can be a different set of fixed and variable sites. The ionization probability of each ionizable site i is calculated by means of the cluster method.31,32 The basic idea of this approximation consists of separating the ionizable sites into different groups (or clusters) according to the strength of their electrostatic interactions. The electrostatic free energy of the sites within a cluster is calculated with eq 7, whereas the cluster-cluster interactions are modeled by a mean field approximation, in which the ionization equilibrium of site I is influenced by site j (not belonging to the cluster I that contains site i) according to θjWij. The electrostatic free energy of cluster I (Gkelec,I) is then given by NI

Gkelec,I

)

Gkelec,I(0)

+

∑ ∑ γixki γjθjWij i)1 j∉I

(10)

where b xk represents an ionization microstate of cluster I, NI is the number of ionization sites in cluster I, and Gkelec,I(0) corresponds to the electrostatic free energy of cluster I (eq 7) at the ionization state b xk when the other clusters are not taken into account. Within the cluster approximation, the number of NC possible ionization states that have to be considered is ∑I)1 2NI, where NC is the number of clusters, instead of 2N. This reduction in the number of ionization states renders the calculation feasible for large N as long as the maximum number of sites per cluster is limited to approximately 10. The final expression for the ionization probability of each ionizable site i is

θi )

1 Z

2 NI

(

∑ xki exp -βGkelec,I(0) -

k)1

nk pH ln 10 -

β 2

NI

)

∑ ∑ γjxkj γlθlWjl (11) j)1 l∉I

New Insights into the Glutamate Racemase Mechanism where 2NI

Z)

that the distribution function corresponding to a HH site resembles a normal distribution.

(

∑ exp -βGkelec,I(0) - nk pH ln10 -

k)1

J. Phys. Chem. B, Vol. 111, No. 9, 2007 2389

β 2

NI

γjxkj γlθlWjl ∑ ∑ j)1 l∉I

)

(12)

Equations 10 and 11 are solved self-consistently for the θi values at a given pH and the titration curves of the ionizable sites generated by repeating the calculation as a function of pH. The pK1/2 value for each site may then be obtained from the corresponding titration curve (when θi ) 1/2, pK1/2 ) pH). As indicated previously, the titration curves of individual sites in biomolecules do not generally follow a HH sigmoidal shape, due to the electrostatic interactions between sites. These nonHH titration curves have been analyzed here by considering them as perturbed titration curves with respect to HH behavior. To quantify the degree of perturbation of a titration curve, a statistical analysis of the first derivative of the equilibrium ionization probability of site i versus the pH can be applied.33,2 The shape of a distribution, in our case dθi/d(pH), can be characterized by a set of n numbers called the moments of the distribution.34 These distribution functions are automatically normalized since the range of θi values necessarily implies that the area under the first derivative function is unity. Hence, we can define the nth central moment of our distribution by

µn )

(dθi/d pH)(pH - pH)n ∑ pH

(13)

where pH is the first central moment (or mean value) of the distribution defined as

pH )

pH(dθi/d pH) ∑ pH

(14)

From the general expression of eq 13 the skewness and kurtosis parameters can be defined. The skewness is an adimensional number that characterizes the degree of asymmetry of a distribution around its mean and is related to the third moment of the distribution by

skewness )

µ3 σ3

(15)

where σ is the standard deviation. For a HH site, the third moment and, hence, the skewness are equal to zero. However, it will be nonzero for a site that interacts strongly with other ionizing groups when the strength of this interaction in the range pH < pK1/2 is different from that in the range pH > pK1/2. The kurtosis is also an adimensional quantity and is related to the fourth moment of the distribution

kurtosis )

µ4 σ4

-3

(16)

It measures the relative peakedness or flatness of a distribution relative to the normal distribution, whose kurtosis equals to zero. In comparison to a normal distribution, a site with positive kurtosis values has higher population in the tails and less population in the region around the mean (where it is peaked), whereas a site with negative kurtosis has higher population around the mean (where it is flat) and less in the tails. Note

3. Computational Details 3.1. pH Titration Curves Calculations. In this paper we have determined the pH titration curves for two model systems, which were built to explore the two mechanistic proposals shown in Figure 2. The first model, which will be called model A in this paper, has been used in our two previous papers on the glutamate racemase enzyme. In brief, the aim of the first paper was to clarify the mode of binding of the substrate in the active site of MurI, starting from different models of the substrate-enzyme Michaelis complex, and running classical MD simulations.16 Simulations with one of these complexessthat built from the MurI (A. pyrophillus)-D-Glu structure which itself was modeled by Hwang et al.15 from the X-ray structure with an inhibitor (D-Gln) in the active sitesgave potentially reactive structures suitable for studying reaction mechanism A. Importantly, in this first model, Asp7 and Glu147′ are initially protonated, the substrate main chain is in its zwitterionic form, and Cys70 is deprotonated. The second paper26 used one of the reactive structures coming from our previous simulations to explore the substrate activation mechanism proposed by Richard and coworkers.22 This consists of a proton transfer from Glu147′ to the substrate carboxylate main chain to give the activated structure schematized in Figure 1. To do this, we computed minimum energy reaction paths (MEPs) for the R-proton abstraction from both the activated and the nonactivated substrates using quantum mechanical/molecular mechanical (QM/MM) potentials35-37 and the conjugate peak refinement (CPR) algorithm. In the work presented here, we employed the energetically most favorable MEP for the activated substrate and carried out pK1/2 calculations on seven structures selected from it. The computational details for the building of the model and the MEP generation are given in the previous studies16,26 and will not be repeated here. As mentioned, this first model corresponds to mechanism A shown in Figure 2a, for which we have a set of structures belonging to a MEP for the R-proton abstraction. As a complement to the work on mechanism A, we also performed pK1/2 calculations on a new model (model B) of the reactant complex corresponding to the second mechanistic proposal (mechanism B in Figure 2). The computational details of this second model are given in the next section. For both models (A and B), the protocol used to determine the pK1/2 values and to obtain the theoretical titration curves was the one outlined in the previous section and is essentially similar to that described by Gilson et al.,38 by Antosiewicz,28 and by David.32 It is assumed that the environmental effects, which determine Gkelec,I, are predominantly electrostatic in origin and are calculated by solving the linearized PoissonBoltzmann equation for the electrostatic potential at each site of a given cluster and for different protonation states. The first step of the pK1/2 calculations is to set the protein at its reference (neutral) state, in which the acids are protonated and the bases are deprotonated. This required the generation of appropriate structure files, which was done by reinitializing or deleting the Cartesian coordinates of the protein hydrogen atoms and redefining the topologies of the protein residues to their neutral state. The HBUILD facility of the CHARMM program39 was then used to add the required hydrogens. The side chain and the main chain carboxylate groups of the substrate were also modeled as neutral. The electrostatic potential calculations were carried out with the UHBD program.40

2390 J. Phys. Chem. B, Vol. 111, No. 9, 2007 The electrostatic potentials were determined at 358 K with an ionic strength of 150 mM. The solvent-accessible region was defined as the volume that could be swept out by a spherical probe of radius 1.4 Å whereas the radius accessible to ions was set to 2.0 Å. Dielectric constants of 80 and of 20 were used for the solvent and the protein, respectively. A four-step method of electrostatic focusing for increased accuracy was employed with consecutive grid sizes of 5.0, 1.0, 0.5, and 0.25 Å. All of the atomic charges and radii were taken from the UHBD program. A temperature of 358 K was also employed for the cluster algorithm, and the maximum number of sites per cluster was set to 13. The ionization probabilities for the sites were determined at pH values from -6 to +18 in increments of 0.1. The histidine residues required special treatment as they have two different protonation sites, Nδ1 and N2. The following procedure was adopted. First, the ionization probability calculations were performed with the histidine Nδ1 sites filled. If, as a result of this calculation, a histidine is protonated on its N2 atom at pH ) 8.5 then no ambiguity exists as both histidine sites are filled. To decide which sites to fill for neutral histidines, the potentials created by all of the other residues at the Nδ1 and N2 atoms are calculated. The site with the lower potential is protonated, and the ionization probability calculation is repeated again for this new set of histidine protonation states. The process is repeated until a consistent set of histidine states is obtained. 3.2. Molecular Dynamics on the Second Model. As mentioned above, a second model (model B) of the MurI-DGlu Michaelis complex was built to check the plausibility of mechanism B in Figure 2b. The model is also based on the 1B74 monomer structure of A. pyrophillus glutamate racemase.15 The dimer was created using the Build Crystallographic Symmetry Tool from the Swiss-PdbViewer. The 1B74 structure was crystallized with the inhibitor D-glutamine, but the final Cartesian coordinates released to the PDB for this ligand were modeled on the basis of the electronic density results and the proposed reaction mechanism of MurI.15 In our model B, the NH2 group of the glutamine side chain has been changed to OH to represent the substrate D-glutamate. The corresponding structure is quite similar to the one we called the D-glutamate model in a previous work.16 The differences are on the Cys70 residue, which has been modeled neutral, and the Asp7, which has been modeled ionized. Previously, we proposed that Glu147′ activates the substrate by donating a proton to its R-carboxylate.26 The results obtained for model A in a previous work support the plausability of this activation mechanism in MurI, which was suggested by Richard and co-workers for the general case of the PLP-independent amino acid racemases.22 Taking into account those previous results26 and as there cannot be a proton tranfer during a classical MD simulation, we modeled two different states to explore the two different scenariossbefore the Glu147′ residue transfers the proton (state 1) and once the proton has been transferred (state 2). Thus, in state 1, the Glu147′ residue is neutral, and the R-carboxylate of the glutamate substrate is ionized, and in state 2, the Glu147′ residue is ionized, and the R-carboxylate of the substrate is neutral. Both systems, corresponding to states 1 and 2, were protonated using the HBUILD facility in the CHARMM program.39 They were then solvated with a box of previously equilibrated water molecules with a density of 0.0334 molecules/Å3, centered at the geometric center of the protein-substrate complex. The initial dimensions of the water box were 101 × 71 × 62 Å3, which is enough to add an 8 Å shell of solvent on each side of

Puig et al. the protein. Water molecules that were within 2.2 Å of any heavy atom of the protein were deleted. All of the crystallographic waters were kept; the final number of water molecules in our model was 13 275. The systems were minimized to allow the water molecules to adjust around the protein. The steepest descent (SD) algorithm was used for 50 steps, followed by 50 more steps using the adopted basis Newton-Raphson (ABNR) algorithm. The atoms of the protein were constrained with harmonic force constants of 50.0 and 20.0 kcal/(mol Å2) during the SD and ABNR minimizations, respectively. Classical MD simulations in the isothermal-isobaric (NPT) ensemble at 1 atm were carried out with periodic boundary conditions (PBCs) by means of the CRYSTAL module of CHARMM (version 31). First, the solvated protein systems were heated gradually without constraints to 358 K in 150 ps. A cutoff of 11.5 Å was used for the nonbonded interactions along with a switch function in the region from 10.5 to 11.5 Å to smoothly reduce the interaction energy to zero. The nonbonded pair list and the image list were built on the basis of group separations and were updated every 30 and 120 ps, respectively. The relative dielectric constant was set to 1. The equations of motions were propagated using the leapfrog integrator with a time step of 1 fs and using the extended system constant pressure and temperature algorithm implemented in CHARMM. All bond lengths involving hydrogen atoms were constrained by the SHAKE algorithm.41 Once the systems were equilibrated, structures were collected every 0.5 ps for later analysis during 1 ns for both state 1 and state 2. The all-atom CHARMM22 force field42 was used to represent the protein and the substrate, and the three-point-charge TIP3P43 model was used for water. For state 2, we derived the charges of the protonated R-carboxylate of the D-Glu substrate as follows. The substrate was partitioned in a molecular mechanics moiety that included the side chain atoms Cγ, Hγ1, Hγ2, Cδ, O1, and O2 and treated with the CHARMM22 force field, and a quantum mechanical moiety that included the rest of atoms and was treated with the AM1 Hamiltonian.44 The QM/MM frontier was placed at the Cβ atom, treating it as a generalized hybrid orbital (GHO) atom.45 A series of 20 ABNR minimizations (starting from structures of the equilibrated substrateenzyme Michaelis complex) of the glutamate were carried out, each of 1000 steps, and the Mulliken charges were calculated. The charges (in a.u.) of the structure with the set of charges closest to the average Mulliken charges over all the minimized structures were taken (N, -0.04; Hτ1, 0.28; Hτ2, 0.28; Hτ3, 0.28; CR, -0.19; HR, 0.20; Cβ, -0.18; C, 0.70; Oβ, -0.47; OH, -0.36; H, 0.32). 4. Results and Discussion 4.1. pK1/2 Calculations on the Minimum Energy Path of Mechanism A. pK1/2 calculations have been performed on the seven structures obtained along the MEP calculated from two minimized structures: the reactant complex structure corresponding to the activated D-Glu substrate (with a proton on the carboxylate main chain as schematized in Figure 1) and the product complex of the racemization. The first columns of Tables 1 and 2 show the pK1/2 values along the path with and without substrate in the active site, respectively, and for selected residues of the glutamate racemase enzyme. The residues have been selected and ordered on the basis of maximum to minimum pK1/2 variation (in absolute value) along the MEP (ninth column). The last column lists the pK1/2 values for the free residues in open exposure to aqueous solution, which are the pKa model values for the residues. To start the discussion, we

New Insights into the Glutamate Racemase Mechanism

J. Phys. Chem. B, Vol. 111, No. 9, 2007 2391

TABLE 1: pK1/2 Values, pK1/2 Variation, and pKa model Values for the Substrate and Selected Residues with the Substrate (Glu) in the Active Centera D-Glu-

site

MurI

L-Glu-

2

3

4

5

6

0.0

14.0 19.5 21.5 21.6 22.8

D-Glu-

MurI ∆pK1/2 pKa model

Cys70 10.4 11.3 12.9 13.1 13.2 13.2 13.3 Cys178 12.7 12.6 12.4 12.0 11.9 11.6 10.6 δ-OOC3.3 3.3 3.6 3.9 4.1 3.9 4.1 Glu Glu147′ 5.3 5.2 5.1 4.8 4.9 5.0 5.0 Asp7 5.0 5.2 5.1 4.9 4.7 4.8 4.8 Tyr39 17.4 17.1 17.2 17.2 17.2 17.2 17.3 Tyr163 14.4 14.4 14.4 14.4 14.4 14.1 14.1 Cys139 10.2 10.2 9.9 9.9 9.9 9.9 10.0 R-OOC- -1.4 -1.5 -1.5 -1.6 -1.7 -1.6 -1.6 Glu Lys119 11.3 11.3 11.2 11.2 11.3 11.3 11.4 His180 6.0 5.9 5.9 5.9 5.9 5.9 5.8 Cys54 12.0 11.9 11.9 11.9 11.9 11.9 11.8 Asp33 2.1 2.0 1.9 1.9 1.9 1.9 1.9 Asp199 2.1 2.0 2.0 2.0 2.0 1.9 1.9 Arg36 15.8 15.8 15.7 15.8 15.8 15.7 15.6 Arg294 14.0 14.0 13.8 13.8 13.8 13.8 13.8 energy

TABLE 2: pK1/2 Values Calculated without the Substrate in the Active Center, pK1/2 Variation, and pKa model Values for Selected Residuesa

2.9 2.1 0.8

8.18 8.18 4.25

0.5 0.5 0.3 0.3 0.3 0.3

4.25 3.65 10.07 10.07 8.18 2.19

0.2 0.2 0.2 0.2 0.2 0.2 0.2

10.53 6.00 8.18 3.65 3.65 12.48 12.48

2.2

L-Glu-

site

MurI

2

3

4

5

6

MurI

Tyr39 Cys178 Tyr30 Cys70 Cys139 Lys185 Lys119 Glu99 Glu147' Glu196 Cys54 Cys430 Asp7 Tyr433 Tyr282 Tyr163 Lys98 Lys61 Lys269 His180

17.0 10.5 14.5 10.7 10.0 12.4 11.2 3.4 3.4 3.0 11.5 12.2 3.1 14.5 15.5 13.9 11.0 12.5 12.8 5.9

16.6 10.4 14.6 10.9 10.1 12.4 11.1 3.4 3.3 3.0 11.7 12.2 3.2 14.5 15.5 13.9 10.9 12.5 12.9 5.9

16.7 10.4 14.6 10.8 9.8 12.4 11.1 3.4 3.2 3.0 11.7 12.2 3.1 14.5 15.5 14.0 10.9 12.6 12.9 5.9

16.7 10.4 14.8 10.6 9.8 12.4 11.1 3.4 3.2 3.0 11.7 12.2 3.1 14.5 15.4 14.0 10.9 12.6 12.9 5.8

16.7 10.2 14.6 10.8 9.8 12.4 11.2 3.4 3.2 3.0 11.7 12.2 3.2 14.5 15.4 14.0 10.9 12.6 12.9 5.8

16.7 10.1 14.6 10.9 9.8 12.4 11.2 3.5 3.3 2.9 11.7 12.2 3.2 14.5 15.4 14.0 10.9 12.6 12.9 5.8

16.7 10.2 14.6 10.7 9.9 12.6 11.3 3.6 3.2 2.8 11.7 12.0 3.0 14.6 15.4 13.9 10.9 12.6 12.9 5.8

energy

0.0

14.0 19.5 21.5 21.6 22.8

2.2

∆pK1/2 pKa model 0.4 0.4 0.3 0.3 0.3 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.1 0.1 0.1 0.1 0.1 0.1 0.1

10.07 8.18 10.07 8.18 8.18 10.53 10.53 4.25 4.25 4.25 8.18 8.18 3.65 10.07 10.07 10.07 10.53 10.53 10.53 6.00

a

The last row indicates the relative potential energies (in kcal/mol) along the MEP of mechanism A (Figure 2). The prime symbol indicates that the residue belongs to the other monomer.

will take as a reference the results in Table 1, which correspond to the substrate-enzyme complex. Importantly, most of the active site residues presented in the Introduction as having a role in binding and/or catalysis appear in it. In general, the comparison of the calculated pK1/2 and the measured pKa model in solution indicates that the protein environment enhances the acidity differences between the residues with respect to the situation in aqueous solution. That is, the residues that are less/ more acidic in aqueous solution diminish/increase, respectively, their acidity in the protein. The exceptions are Asp7 and Glu147′, which become less acidic as will be commented upon below. As expected, the two catalytic cysteines (70 and 178) are the residues showing the largest pK1/2 variation along the reaction path. Cys70 increases its pK1/2 value as the reaction proceeds, from 10.4 to 13.3, whereas Cys178 decreases its values from 12.7 to 10.6. Thus, in reactants, Cys70 is more likely to be deprotonated, and Cys178 is more likely to be protonated as expected for the reaction going from the D- to the L-glutamate enantiomer. When these calculations are repeated for the same structures but excluding the substrate’s influence the cysteines show a quite different behavior. Table 2 shows the pK1/2 calculations for each structure without the substrate, and as can be seen in comparison with Table 1, the systematic drift has been lost, and the cysteines show similar pK1/2 values throughout the reaction path. The pK1/2 variation of each cysteine is decreased from 2.9 to 0.3 for Cys70 and from 2.1 to 0.4 for Cys178. These results show that the substrate in the active site is able to modulate the acidity of the catalytic cysteines along the reaction path depending on the acid/base role expected for each cysteine. Although the pK1/2 variation of the catalytic cysteines supports mechanism A, their absolute pK1/2 values in the D-Glu-MurI and L-Glu-MurI complexes indicate that the fraction of the thiolate species is very low, which is at odds with the proposal that the thiolate groups of Cys70 and Cys178 abstract the R-proton of the D-Glu and L-Glu substrates, respectively. After the cysteines, the δ-carboxylate of the glutamate substrate is the ionizable site showing the largest pK1/2 variation along the path in Table 1. It is important to mention that the

a

The last row indicates the relative potential energies (in kcal/mol) along the MEP of mechanism A (Figure 2). The prime symbol indicates that the residue belongs to the other monomer.

δ-carboxylate of the glutamate substrate is hydrogen-bonded throughout the MEP to the Asp7 carboxylate, which is protonated in model A (Figure 1). The interaction distance between the Hδ of Asp7 and the O2 of the δ-carboxylate group oscillates between 2.02 and 2.10 Å along the path. The values in Table 1 suggest that the proton is more likely to be bonded to the carboxylate group of Asp7. This is also one of the residues that shows a large pK1/2 variation (0.5), and its pK1/2 values are shifted by approximately 1.5 units from the pKa model of an aspartic acid residue. Table 2 shows that this shift is not observed when the substrate is absent. Thus, these results support the hypothesis that Asp7 has a role in the binding of the substrate by interacting with its δ-carboxylate group during the reaction. Glu147′ shows the same pK1/2 variation as Asp7 along the path and also has a positive shift of its pK1/2 value with respect to a glutamate in aqueous solution. In contrast, it has a negative shift when the substrate is not included in the calculations. This different result is because the δ-carboxylate of Glu147′ is strongly influenced by the substrate via a hydrogen bond to its main chain R-carboxylate. The length of this hydrogen bond is even shorter than the substrate side chain-Asp7 hydrogen bond interaction mentioned above because the distance between the hydrogen bonded to the substrate main chain carboxylate and the O2 atom of the δ-carboxylate of Glu147′ oscillates between 1.93 and 2.05 Å along the path. It has previously been reported that the strong interaction between titrable groups that are in areas of low solvent accessibility might result in pK1/2 values outside the usual pH range in protein studies (2-12).46 This could explain the negative pK1/2 values of the substrate R-carboxylate. Comparison of the pK1/2 values of the substrate R-carboxylate and the δ-carboxylate of Glu147′ shows that the thermodynamically more stable position for the hydrogen-bond proton is bound to Glu147′, which indicates that activation of substrate by proton transfer to its R-carboxylate has a free energy cost. The other residues in Table 1 show a much smaller variation of the pK1/2 values along the path. Among them, it is just worth mentioning the evolution of the pK1/2 value of His180 along

2392 J. Phys. Chem. B, Vol. 111, No. 9, 2007

Puig et al.

Figure 3. Titration curves of Cys70 (in red) and Cys178 (in green) calculated with the substrate in the active center for the structures along the racemization reaction pathway of mechanism A.

the path, which becomes more acidic in going from D-Glu to L-Glu, indicating that it might act as a proton donor to the thiolate group of Cys178 after (or while) this residue reprotonates the R-carbon (Figure 1). 4.2. Titration Curves and Statistical Parameters Corresponding to Mechanism A. The titration curves of Cys70 and Cys178 in the presence of the substrate in the active center are shown for each structure of the MEP in Figure 3. It can be seen that Cys70 becomes less acidic and Cys178 becomes more acidic as the reaction proceeds. The titration curves coincide at

approximately the third structure. Thus, during the racemization reaction the enzyme enhances the acid/base properties that each cysteine should have for the reaction to proceed from the D- to the L-glutamate enantiomer. In Table 3, the skewness and kurtosis parameters (eqs 15 and 16, respectively) are shown for the residues that have the largest pK1/2 variation along the path. The values of the skewness indicate that the titration curves of all of the residues have substantial asymmetry. It is interesting to note the systematic behavior of the cysteines’ kurtosis values and the different

New Insights into the Glutamate Racemase Mechanism

J. Phys. Chem. B, Vol. 111, No. 9, 2007 2393

TABLE 3: Skewness and Kurtosis Values for the Residues Which Show the Largest Variations of pK1/2 Values along the MEP of Mechanism A Cys70 D-Glu-MurI

2 3 4 5 6 L-Glu-MurI

Cys178

Asp7

Glu147′

skewness

kurtosis

skewness

kurtosis

skewness

kurtosis

skewness

kurtosis

0.30 0.03 -0.40 -0.53 -0.56 -0.70 -0.24

-0.84 -0.86 -0.28 0.52 0.62 0.86 1.66

0.18 0.01 0.12 0.19 0.23 0.42 0.22

2.40 1.95 1.39 0.89 0.84 1.26 -0.14

-0.75 -0.84 -0.77 -0.63 -0.59 -0.60 -0.52

-0.06 0.22 0.16 0.01 -0.02 0.10 0.01

-0.96 -0.91 -0.81 -0.72 -0.68 -0.73 -0.67

1.43 1.42 1.24 1.06 0.93 1.11 1.07

behaviors that they have depending upon whether the cysteine is Cys70 or Cys178. The kurtosis of Cys70 goes from negative (-0.84) to positive values (1.66) as the reaction proceeds. The titration curve distributions, which are just the slopes of the curves (section 2), are flat at the reactant and become peaked at the product. In general, the flatter the distribution (negative values of kurtosis), the smoother the titration curve over a wide range of pH and, as a consequence, the more readily the residue can act as a catalytic acid/base residue. In the case of Cys70 this occurs for the first three structures (Figure 3 and Table 3). Thus, at the reactant, the enzyme modifies Cys70’s titration curve in such a way that protonated and deprotonated states can coexist over a wide pH range. These results agree with the fact that Cys70 acts at an early stage of the reaction to deprotonate the R-carbon. In contrast, the kurtosis of Cys178 goes from positive values (2.40) at the reactant to a negative value at the product (-0.14). The distribution of the Cys178 titration curve becomes flat at the end of the reaction path when it reprotonates the R-carbon. The skewness and kurtosis parameters are also shown for the residues Asp7 and Glu147′ in Table 3. Both residues have large negative values of the skewness parameter along the reaction path. The skewness parameters were calculated for all of the aspartic and glutamic residues at the reactant structure of the MEP. Asp7 and Glu147′ showed the most negative skewness of all of the aspartic and glutamic residues, respectively. However, it is difficult to extract catalytic information from the kurtosis parameters obtained for Glu147′ and Asp7 as they show an erratic behavior. 4.3. Stabilization of Cys70 and an Alternative Model. The Cys70 pK1/2 values in Table 1 are shifted by up to 5 units relative to the pKa model value of a cysteine residue in aqueous solution. The cysteine acidity is thus diminished, and it becomes more capable of abstracting the R-proton. As can be seen from the titration curves (Figure 3) of the cysteines at the first structure (D-Glu-MurI), the probablity of being ionized is higher for Cys70 than that for Cys187 at the reaction’s optimum pH (8.5),47 and the situation is the opposite at the product complex (the seventh structure, L-Glu-MurI). However, the titration curves in Figure 3 indicate that for all the structures along the path there are low yields of the unprotonated thiol. To find an explanation for these results, which indicate that the glutamate racemase active site is not able to stabilize the ionized forms of Cys70 and Cys178, we have searched for residues that could hydrogen-bond to Cys70 as well as for residues that could perturb its acid/base properties by chargecharge interactions or charge-dipole interactions. At the reactant structure of the MEP there is no residue close enough to the thiolate sulfur to donate a proton in model A. The closest residue is Asp7 with a Cys70-Sγ-Asp7-Cγ distance of 4.64 Å, followed by Thr72-Oγ1 at 5.12 Å, Ser74-Oγ at 6.06 Å, and Ser8-Oγ at 6.41 Å. There are not positive charges close to the thiolate that could stabilize its negative charge. His180 is the closest at 12.32

Å, followed by Lys82 and Lys119 at 14.08 and 14.86 Å, respectively. The Cys70 sulfur is thus quite isolated in the active site of glutamate racemase, and we think that its very low acidity might come from the desolvation effect suffered by this residue in the apolar active site of glutamate racemase. The only interactions that have been found for the Cys70 residue are with the backbone NH dipoles of the adjacent residues (Asn71, Thr72, Ala73, and Ser74), but they are probably not strong enough to stabilize the sulfur thiolate. Thus, in proposed mechanism A, the stabilization of the thiolate sulfur is still unresolved. The titration curves of Cys70 and Cys178 tell us that there will not be appreciable amounts of deprotonated cysteine as there is no residue stabilizing the sulfur anion. Taking into account all of this information, we decided to explore the second mechanistic proposal of Glavas and Tanner,21 which has been mentioned in the Introduction and is schematically represented in Figure 2b. In mechanism B, Cys70 would remain mostly neutral during the reaction while acting as an assisted general base. Thus, following the suggestion of Glavas and Tanner and in view of the positive shift of the pK1/2 values reported in Table 1 for Asp7, we argue that Cys70 deprotonates the R-proton while it is in turn deprotonated by the Asp7 residue. The two deprotonations would occur in a concerted manner to minimize the amount of time that Cys70 remains deprotonated. Thus, as the first stage of our research on this new mechanistic proposal, we explore here the stability of its reactant state by means of a new model, denoted in section 3 as model B, in which Cys70 is neutral and Asp7 is ionized. 4.4. Molecular Dynamics Simulations with Neutral Cys70. Molecular dynamics simulations were carried out with a neutral Cys70 and an ionized Asp7 to confirm the feasibility of proposed mechanism B. Two different ionization states of the active site were modeled that differ in the ionization states of Glu147′ and the substrate’s R-carboxylate. In state 1, Glu147′ is protonated, and the substrate’s R-carboxylate deprotonated, whereas in state 2, Glu147′ is deprotonated, and the substrate’s carboxylate is protonated. These two different protonation states have been modeled to explore the prereactive contacts before (state 1) and after (state 2) the Glu147′ residue protonates the substrate’s R-carboxylate. As mentioned previously, this protontransfer process has been demonstrated to activate the substrate prior to the racemization reaction.22,26 4.4.1. State 1, the Zwitterionic Ligand. Figure 4 displays a trajectory snapshot from the PBC MD simulation of state 1 of model B. The active site residues and some key interactions are marked on it, the evolutions of some of which along the trajectory are shown in Figure 5. Note that in this model the glutamate substrate is in its zwitterionic form, which is thought to be the predominant state when it enters the active site because it is the predominant state in solution. During a long enough MD simulation the hydrogen bond between the substrate’s R-carboxylate and the Glu147′ side chain is not dynamically stable (pink line in Figure 5). This fact can

2394 J. Phys. Chem. B, Vol. 111, No. 9, 2007

Figure 4. Trajectory snapshot of the PBC MD simulation corresponding to state 1 of model B.

be explained by taking into account that the residue Glu147′ and the rest of the active site residues belong to different monomers within the dimer. Thus Glu147′ is not rigidly anchored to the active site and is able to move significantly inside the deep pocket formed between the two monomers. During the first 500 ps this hydrogen bond holds although some intense oscillations of the corresponding distance appear. Afterward, this hydrogen bond breaks, being substituted by interactions of the Glu147′ side chain with waters present in the pocket and a hydrogen bond between the substrate’s R-carboxylate and the NH Tyr39 backbone. We could expect that the proton transfer from Glu147′ to the substrate’s R-carboxylate occurs while the corresponding hydrogen bond exists (during the first 500 ps in our MD simulation), although our MM force field does not allow us to verify this point. Other important well-maintained interactions during the MD simulation are the ones between the strictly conserved residues Ser8 and Asp7. Two hydrogen bonds are alternatively formed between the carboxylate oxygen atoms of Asp7 and either the Ser8-HN or the Ser8-Hγ1 atoms. We propose that Ser8 is able to orient the Asp7 side chain through these interactions to a reactive position (Figure 4). The position of Asp7 is crucial because this residue in turn orients Cys70 and the glutamate substrate, as will be discussed next. The Asp7 residue interacts with the catalytic Cys70-Hγ1. The interaction is very labile, going randomnly from distances of 6.0 to 2.0 Å (red line in Figure 5). The interaction is weak as a thiol group is a poor hydrogen donor. Nevertheless, the Asp7 residue is the only residue close enough to interact with the Cys70 residue. No other interactions have been found for this cysteine during the simulation. Importantly, some of the conformations generated in the MD simulation could lead to proton transfer from the thiolate to the Asp7 residue. In this simulation, the Cys70 residue is still a bit far from the R-hydrogen of the substrate (green line in Figure 5). We will see that in state 2 (see next section), involving a deproto-

Puig et al. nated Glu147′ and a protonated substrate’s R-carboxylate, the R-hydrogen of the substrate is already close to the Cys70 residue. 4.4.2. State 2, the Carboxylate Protonated Ligand. Figure 6 displays a trajectory snapshot from the PBC MD simulation of state 2 of model B, and the evolutions of important interaction distances are shown in Figure 7. In this protonation state, the R-carboxylate of the glutamate substrate is protonated. We assume that this proton comes from the Glu147′ residue, which has been modeled deprotonated. In this model, the interaction between the Asp7-Oδ2 atom and the Cys70-Hγ1 atom is also very labile (red line in Figure 7), at least for the first 190 ps of the simulation. A comparison of the red lines of Figures 5 and 7 during the first part of the simulations indicates that the number of configurations that could lead to the proton transfer from Cys70 to Asp7 has increased in going from state 1 to state 2. It is noteworthy that in the recently resolved B. subtilis RacE structure18 (the only one with the substrate in the active site) the thiol group of Cys74 is 3.9 Å from the side chain of Asp7. Thus, our model is able to visit configurations similar to that in this X-ray structure, despite being built from an inhibitor-enzyme structure that has been proposed as a noncatalytic conformation of the enzyme.18 The results of our simulation do not support this last proposal. After the first 190 ps of the simulation of model B in state 2, the Cys70 residue changes its conformation and interacts with the backbone carbonyl oxygen of Asp7. These last conformations are not reactive, but nevertheless, we would expect the proton transfer to have already occurred by that time. The glutamate substrate is also closer to the Cys70 residue than in the state 1 simulation, and the Glu147′-OOC-D-Glu hydrogen bond no longer exists. Actually, in this model, the neutral state of the substrate main chain carboxylate enhances the formation of a strong salt bridge interaction between the positive ammonium group of the substrate and the negative Asp7, which is maintained throughout the simulation (blue line in Figure 7). This distance is much shorter than that in state 1 (blue line in Figure 5). These results suggest that this interaction orients the substrate and prepares it for proton abstraction by moving it toward the Cys70 thiol. Thus, the green line in Figure 7 shows that during the first 190 ps of the simulation there are some configurations in which the Cys70-Sγ atom is close enough to interact with the R-hydrogen of the substrate. Note that this is the same part of the trajectory during which the Cys70-Sγ atom interacts with the Asp7-Oδ2 atom. The interaction between the Cys70-Sγ and the R-hydrogen is expected to be stronger if Cys70 is deprotonated and eventually would lead to R-hydrogen abstraction. In summary, the structural results obtained here support the role of neutral Cys70 thiol as a general base assisted by Asp7. Cys70 would abstract the substrate’s R-hydrogen while at the same time having its Cys70-Hγ1 abstracted by Asp7. In addition, the already mentioned strong charge-charge interaction of the ammonium group protons with the Asp7 carboxylate prevents the Cys70-Sγ atom abstracting the more acidic protons of this positive group. Thus, the MD results presented here support the ionized state of the substrate’s amino group as a requirement for the enzymatic reaction, in contrast to the objections raised by Mo¨ritz and Bruice17 who doubted whether a D/L-Glu substrate with a protonated amino group could be inverted to give the L/D-Glu product. 4.5. pK1/2 Calculations on the New Model. The structural results presented in the last section indicate that mechanism B of Figure 2 is plausible. In this section, we present (in Table 4)

New Insights into the Glutamate Racemase Mechanism

J. Phys. Chem. B, Vol. 111, No. 9, 2007 2395

Figure 5. Evolution of selected interatomic distances from the PBC MD simulation of state 1 of model B: Asp7-Cys70 in red, Cys70-R-proton in green, NH3-D-Glu-Asp7 in blue, and Glu147′-H2-δ-OOC-D-Glu in pink. Units corresponding to the x- and y-axes are picoseconds and a˚ngstroms, respectively.

Figure 6. Trajectory snapshot of the PBC MD simulation corresponding to state 2 of model B.

the results of pK1/2 calculations carried out on one of the structures generated by the MD simulation of state 2 of model B. The most notable difference between this structure and the ones employed for generating the results in Table 1 is the disposition of the Asp7 residue. In the new structure the Asp7 residue no longer maintains a hydrogen bond with the δ-carboxylate of the substrate, but it interacts with its amino group by a salt bridge (Figure 6). This results in an even higher pK1/2 shift for this residue, up to 6.1, which is 1.1 units less acid than the MEP reactant structure (Table 1). Thus, both the structural results shown in Figure 7 and the pK1/2 value of Table 4 suggest that, in the structures generated by MD simulations of model B, the Asp7 residue could more easily act as a base to deprotonate Cys70. Taking into account our previous results obtained for the MEP structures of model A (Table 1), which suggest a binding role for Asp7, and our results for model B, which indicate the possible role of its carboxylate as an assistant base, we find it

reasonable to draw the following picture of the different roles of the Asp7 carboxylate in the racemization reaction catalyzed by MurI. First, it catalyzes proton exchange on Cys70, ensuring that the active deprotonated Cys is in rapid equilibrium with the inactive protonated form, and second, once it has been protonated, it stabilizes the carbanionic intermediate through a hydrogen bond to the side chain δ-carboxylate of the substrate. Finally, although we have not simulated the L f D racemization, it is reasonable to think that Asp7 might reprotonate Cys70 while this catalytic residue, now acting as an acid, reprotonates the carbanionic intermediate formed after the R-proton abstraction. The cysteine residues are slightly less acidic in this model and show similar pK1/2 values (13.0 for Cys70 and 13.5 for Cys178). The decreased acidity of the cysteines will make the abstraction of the R-proton of the substrate easier but make it more difficult for the Asp7 residue to deprotonate Cys70. Interestingly the titration curves of the two cysteines are now very similar (results not shown). This can be understood as the MD-generated structure is more symmetric than the reactant structure of the MEP as both cysteines are protonated. Nevertheless, we expect Cys178’s acidity to be modulated along the reaction path in a similar way as along the MEP, once (or at the same time as) Cys70 is deprotonated by the Asp7 residue. Another remarkable structural difference of this model is that the Glu147′ residue interacts with the His180 residue (Figure 6) instead of forming a hydrogen bond with the R-carboxylate of the substrate. This has important consequences for the pK1/2 values of Glu147′, His180, and the R-carboxylate of the substrate (as can be seen by comparing Tables 1 and 4) and will be commented on below. The average and standard deviation for the Glu147′-His180 interaction distance are 5.2 and 1.2 Å, respectively. The value of the standard deviation shows that the distance fluctuates significantly and for some configurations is close to the value of 2.8 Å found in the B. subtilis enzyme X-ray structure.18 Thus, this is a second indication (the first was the Cys70-Asp7 interaction of the previous section) that our model is able to visit configurations similar to the only existing X-ray structure of the enzyme-substrate complex. The Glu147′His180 interaction makes Glu147′ more acidic and His180 less acidic than in the structures of model A (Table 1). However, the behavior of the pK1/2 values of His180 along the MEP (Table 1) suggest that its high pK1/2 value in the model B Michaelis complex (Table 4) could also decrease during the D f L reaction. This would permit it to protonate Cys178 when this residue acts as an acid to reprotonate the R-carbon on the L face. As an additional consequence of the structural difference

2396 J. Phys. Chem. B, Vol. 111, No. 9, 2007

Puig et al.

Figure 7. Evolution of selected interatomic distances from the PBC MD simulation of state 2 of model B: Asp7-Cys70 in red, Cys70-R-proton in green, and NH3-D-Glu-Asp7 in blue. Units corresponding to the x- and y-axes are picoseconds and a˚ngstroms, respectively.

TABLE 4: pK1/2 Values Calculated for Mechanism B (Figure 2b) with the r-Carboxylate of the Substrate Protonated and pKa model as a Reference Cys70 Cys178 R-OOC-Glu His180 Asp7 δ-OOC-Glu Glu147'

pK1/2

pKa model

13.0 13.5 8.6 7.5 6.1 4.5 1.4

8.18 8.18 2.19 6.00 3.65 4.25 4.25

mentioned at the beginning of the previous paragraph, the pK1/2 value of the R-carboxylate of the substrate is raised from negative values in model A (Table 1) up to 8.6 in model B (Table 4). This means the the ordering of the acidities of Glu147′ and the R-carboxylate of the substrate is exactly opposite in the two models. Contrary to what has been found for model A structures, the higher pK1/2 value of the R-carboxylate of the substrate in model B indicates that the proton is more stable in this carboxylate. These results highlight how dependent ionization probabilities are on the structure and lead us to conclude that the Glu147′ and R-OOC-Glu carboxylates could interchange their acid/base roles during the racemization reaction. Thus, if the mechanism presented above for the completion of the D f L reaction occurred (via the double proton transfer of His180 to Cys178 and from Cys178 to the R-carbon), His180 would become neutral and would probably lose its interaction with Glu147′. The latter would then interact more strongly with the R-OOC-Glu carboxylate. The reformation of the hydrogen bond between R-OOC-Glu (now protonated) and Glu147′ would strongly change their acidities and eventually facilitate the proton transfer from the L-Glu substrate to Glu147′. This proposal is supported by the results shown in Table 1. 5. Conclusions In this paper we have explored the plausibility of two mechanisms proposed for the racemization reaction catalyzed by the A. pyrophillus glutamate racemase enzyme (Figure 2). We have determined the ionization probabilities (and the pK1/2 values) of the active site residues for reactant complexes of the two mechanisms, which have been modeled and equilibrated by means of MD simulations. For mechanism A, we have also carried out pK1/2 calculations along a minimum energy path determined in a previous study.26 The results obtained for mechanism A indicate that when the catalytic Cys70 and Cys178 are modeled as ionized and neutral,

respectively, in the D-Glu-MurI complex (and vice versa in the L-Glu-MurI complex) the active site is able to modulate their acidities depending on the racemization direction. However, mechanism A appears to be disfavored because of the low stability of thiolate in the active site of A. pyrophillus MurI. Mechanism B was explored by means of MD simulations of a model in which both cysteines were neutral residues and Asp7 was ionized. Both structural and pK1/2 results obtained for the equilibrated Michaelis complex of this model support the following mechanism that consists of several proton transfers: After the activation of the substrate by means of a proton transfer from Glu147′ to its R-carboxylate, the Asp7 residue deprotonates the Cys70 residue while this cysteine abstracts the R-proton of the substrate. The substrate is then reprotonated by the Cys178 residue on the opposite stereochemical face, which is in turn reprotonated by His180. Overall, the global reaction is a proton transfer from His180 to the Asp7 residue along with inversion of the substrate’s R-carbon. The pK1/2 values of these two residues are quite similar due to the high pK1/2 value of the Asp7 residue, and thus the global reaction is thermodynamically possible. However, the determination of the relative rates and the degree of synchronicity of all of these proton transfers would require a dynamical study of the whole process. Our results also provide a qualitative picture that leads us to support roles for the important residues in the active site of MurI: Asp7, Glu147′, and His180. Thus, Asp7 has three roles: (1) to assist the deprotonation of D-Glu by participating in general base catalysis with the neutral Cys70 thiol; (2) to bind and stabilize the substrate along the reaction path via a hydrogen bond to its δ-carboxylate; and (3) to reprotonate Cys70 while this residue reprotonates the R-carbon in the L f D direction. For Glu147′, we propose two roles: (1) to activate the bound substrate by donating a proton to its carboxylate main chain in a step prior to the racemization process; and (2) to stabilize the cationic form of His180. Finally, His180 is proposed to assist Cys178 in the deprotonation of L-Glu. Importantly, these conclusions are in general agreement with the proposals derived from previous experimental work on the glutamate racemase enzymes of three different species;15,18,21 in particular, they agree with the structural data arising from the only X-ray structure of the enzyme that is available with substrate in the active site,18 in spite of the fact that the theoretical model used here has not been built from this structure. As a more general conclusion, the results presented here show that the behavior of an enzyme depends critically on the ionization state of its active site residues and thus justifies the

New Insights into the Glutamate Racemase Mechanism calculation of their ionization probabilities to shed light on the enzyme catalytic mechanism. Acknowledgment. We are grateful for financial support from the Spanish Ministerio de Educacio´n y Ciencia and the Fondo Europeo de Desarrollo Regional through Project No. CTQ200507115/BQU and the Generalitat de Catalunya (Grant No. 2005SGR00400). M.G.-V. thanks the Ramon y Cajal program for financial support. E.P. also thanks all of the people in the Laboratoire de Dynamique Mole´culaire for their hospitality. Supporting Information Available: Analysis of the rootmean-square deviations of the active site residues along the MD simulations. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Schlippe, Y. V. G.; Hedstrom, L. Arch. Biochem. Biophys. 2005, 433, 266. (2) Ko, J. J.; Murga, L. F.; Andre, P.; Yang, H. Y.; Ondrechen, M. J.; Williams, R. J.; Agunwamba, A.; Budil, D. E. Proteins 2005, 59, 183. (3) Walsh, C. T. J. Biol. Chem. 1989, 264, 2393. (4) Baltz, R. H.; Hoskins, J. A.; Solenberg, P. J.; Treadway, P. J. U. S. Patent 5,981,281, 1999. (5) Doublet, P.; van Heijenoort, J.; Mengin-Lecreulx, D. Biochemistry 1994, 33, 5285. (6) Pucci, M. J.; Thanassi, J. A.; Ho, H. T.; Falk, P. J.; Dougherty, T. J. J. Bacteriol. 1995, 177, 336. (7) Kobayashi, K.; Ehrlich, S. D.; Albertini, A.; Amati, G.; Andersen, K. K.; Arnaud, M.; Asai, K.; Ashikaga, S.; Aymerich, S.; Bessieres, P.; Boland, F.; Brignell, S. C.; Bron, S.; Bunai, K.; Chapuis, J.; Christiansen, L. C.; Danchin, A.; Debarbouille, M.; Dervyn, E.; Deuerling, E.; Devine, K.; Devine, S. K.; Dreesen, O.; Errington, J.; Fillinger, S.; Foster, S. J.; Fujita, Y.; Galizzi, A.; Gardan, R.; Eschevins, C.; Fukushima, T.; Haga, K.; Harwood, C. R.; Hecker, M.; Hosoya, D.; Hullo, M. F.; Kakeshita, H.; Karamata, D.; Kasahara, Y.; Kawamura, F.; Koga, K.; Koski, P.; Kuwana, R.; Imamura, D.; Ishimaru, M.; Ishikawa, S.; Ishio, I.; Lecoq, D.; Masson, A.; Mauel, C.; Meima, R.; Mellado, R. P.; Moir, A.; Moriya, S.; Nagakawa, E.; Nanamiya, H.; Nakai, S.; Nygaard, P.; Ogura, M.; Ohanan, T.; O’Reilly, M.; O’Rourke, M.; Pragai, Z.; Pooley, H. M.; Rapoport, G.; Rawlins, J. P.; Rivas, L. A.; Rivolta, C.; Sadaie, A.; Sadaie, Y.; Sarvas, M.; Sato, T.; Saxild, H. H.; Scanlan, E.; Schumann, W.; Seegers, J. F. M. L.; Sekiguchi, J.; Sekowska, A.; Seror, S. J.; Simon, M.; Stragier, P.; Studer, R.; Takamatsu, H.; Tanaka, T.; Takeuchi, M.; Thomaides, H. B.; Vagner, V.; van Dijl, J. M.; Watabe, K.; Wipat, A.; Yamamoto, H.; Yamamoto, M.; Yamamoto, Y.; Yamane, K.; Yata, K.; Yoshida, K.; Yoshikawa, H.; Zuber, U.; Ogasawara, N. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 4678. (8) Dios, de A.; Prieto, L.; Martin, J. A.; Rubio, A.; Ezquerra, J.; Tebbe, M.; deUralde, B. L.; Martin, J.; Sanchez, A.; LeTourneau, D. L.; McGee, J. E.; Boylan, C.; Parr, T. R.; Smith, M. C. J. Med. Chem. 2002, 45, 4559. (9) Doublet, P.; van Heijenoort, J.; Mengin-Lecreulx, D. Microb. Drug Resist. 1996, 2, 43. (10) Tanner, M. E. Acc. Chem. Res. 2002, 35, 237. (11) Yoshimura, T.; Esaki, N. J. Biosci. Bioeng. 2003, 96, 103. (12) Buschiazzo, A.; Goytia, M.; Schaeffer, F.; Degrave, W.; Shepard, W.; Gregoire, C.; Chamond, N.; Cosson, A.; Berneman, A.; Coatnoan, N.; Alzari, P. M.; Minoprio, P. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 1705. (13) Abe, K.; Takahashi, S.; Muroki, Y.; Kera, Y.; Yamada, R. H. J. Biochem. 2006, 139, 235. (14) Lloyd, A. J.; Huyton, T.; Turkenburg, J.; Roper, D. I. Acta Crystallogr., Sect. D: Biol. Crystallogr. 2004, 60, 397.

J. Phys. Chem. B, Vol. 111, No. 9, 2007 2397 (15) Hwang, K. Y.; Cho, C. S.; Kim, S. S.; Sung, H. C.; Yu, Y. G.; Cho, Y. J. Nat. Struct. Biol. 1999, 6, 422. (16) Puig, E.; Garcia-Viloca, M.; Gonza´lez-Lafont, A.; Lopez, I.; Daura, X.; Lluch, J. M. J. Chem. Theory Comput. 2005, 1, 737. (17) Mobitz, H.; Bruice, T. C. Biochemistry 2004, 43, 9685. (18) Ruzheinikov, S. N.; Taal, M. A.; Sedelnikova, S. E.; Baker, P. J.; Rice, D. W. Structure 2005, 13, 1707. (19) Gallo, K. A.; Tanner, M. E.; Knowles, J. R. Biochemistry 1993, 32, 3991. (20) Glavas, S.; Tanner, M. E. Biochemistry 1999, 38, 4106. (21) Glavas, S.; Tanner, M. E. Biochemistry 2001, 40, 6199. (22) Rios, A.; Amyes, T. L.; Richard, J. P. J. Am. Chem. Soc. 2000, 122, 9373. (23) Richard, J. P.; Williams, G.; O’Donoghue, A. C.; Amyes, T. L. J. Am. Chem. Soc. 2002, 124, 2957. (24) Richard, J. P.; Amyes, T. L. Bioorg. Chem. 2004, 32, 354. (25) Richard, J. P.; Amyes, T. L. Curr. Opin. Chem. Biol. 2001, 5, 626. (26) Puig, E.; Garcia-Viloca, M.; Gonza´lez-Lafont, A.; Lluch, J. M. J. Phys. Chem. A 2006, 110, 717. (27) Onufriev, A.; Case, D. A.; Ullmann, G. M. Biochemistry 2001, 40, 3413. (28) Antosiewicz, J.; Briggs, J. M.; Elcock, A. H.; Gilson, M. K.; McCammon, J. A. J. Comput. Chem. 1996, 17, 1633. (29) Bashford, D.; Karplus, M. Biochemistry 1990, 29, 10219. (30) Bashford, D.; Karplus, M. J. Phys. Chem. 1991, 95, 9556. (31) Antosiewicz, J.; McCammon, J. A.; Gilson, M. K. J. Mol. Biol. 1994, 238, 415. (32) David, L. Ph.D. Thesis, Universite´ Joseph Fourier-Grenoble I, France, 1996. (33) Ko, J.; Murga, L. F.; Wei, Y.; Ondrechen, M. J. Bioinformatics 2005, 21, 1258. (34) Press, W. H. Numerical Recipies in Fortran77: The Art of Scientific Computing, 2nd ed.; Cambridge University Press, Cambridge, U. K., 1992, p 604. (35) Warshel, A.; Levitt, M. J. Mol. Biol. 1976, 103, 227. (36) Singh, U. C.; Kollman, P. A. J. Comput. Chem. 1986, 7, 718. (37) Field, M. J.; Bash, P. A.; Karplus, M. J. Comput. Chem. 1990, 11, 700. (38) Gilson, M. Proteins 1993, 15, 266. (39) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187. (40) Madura, J. D.; Briggs, J. M.; Wade, R. C.; Davis, M. E.; Luty, B. A.; Ilin, A.; Antosiewicz, J.; Gilson, M. K.; Bagheri, B.; Scott, L. R.; McCammon, J. A. Comput. Phys. Commun. 1995, 91, 57. (41) Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (42) MacKerell, A. D., Jr.; Bashford, D.; Bellott, M.; Jr., R. L. D.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; JosephMcCarthy, D.; Kuchnir, L.; Kuczera, K.; Lau, F. T. K.; Mattos, C.; Michnick, S.; Ngo, T.; Nguyen, D. T.; Prodhom, B.; III, W. E. R.; Roux, B.; Schlenkrich, M.; Smith, J. C.; Stote, R.; Straub, J.; Watanabe, M.; Wio´rkiewicz-Kuczera, J.; Yin, D.; Karplus, M. J. Phys. Chem. B 1998, 102, 3586. (43) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (44) Dewar, M. J. S.; Zoebisch, E. G.; Healy, E. F.; Stewart, J. J. P. J. Am. Chem. Soc. 1985, 107, 3902. (45) Gao, J.; Amara, P.; Alhambra, C.; Field, M. J. J. Phys. Chem. A 1998, 102, 4714. (46) Blachut-Okrasinska, E.; Lesyng, B.; Briggs, J. H.; McCammon, J. A.; Antosiewicz, J. M. Eur. Biophys. J. 1999, 28, 457. (47) Kim, S. S.; Choi, I. G.; Kim, S. H.; Yu, Y. G. Extremophiles 1999, 3, 175.