Computational Mechanistic Studies Addressed to the Transimination

Mar 29, 2011 - Henrique Silva Fernandes , Maria João Ramos , Nuno M. F. S. A. Cerqueira ... Carla S. Silva Teixeira , Henrique S. Fernandes , Pedro A...
0 downloads 0 Views 4MB Size
ARTICLE pubs.acs.org/JCTC

Computational Mechanistic Studies Addressed to the Transimination Reaction Present in All Pyridoxal 50 -Phosphate-Requiring Enzymes N. M. F. S. A. Cerqueira, P. A. Fernandes, and M. J. Ramos* REQUIMTE, Faculdade de Ci^encias, Universidade do Porto, Rua Campo Alegre, s/n, 4169-007 Porto, Portugal ABSTRACT: The pyridoxal-50 -phosphate-dependent enzymes (PLP enzymes) catalyze a myriad of biochemical reactions, being actively involved in the biosynthesis of amino acids and amino acid-derived metabolites as well as in the biosynthetic pathways of amino sugars and in the synthesis or catabolism of neurotransmitters. Although the scope of PLP-catalyzed reactions initially appears to be bewilderingly diverse, there is a simple unifying principle: In the resting state, the cofactor (PLP) is covalently bonded to the amino group of an active site lysine, forming an internal aldimine. Once the amino substrate interacts with the active site, a new Schiff base is generated, commonly referred to as the external aldimine. Only after this step, the mechanistic pathway for each PLPcatalyzed reaction diverges. In this paper, density functional methods have been applied to investigate this common step present in all PLP-dependent enzymes—the transimination reaction. The results indicate that the reaction involves three sequential steps: (i) formation of a tetrahedral intermediate with the active site lysine and the amino substrate bonded to the PLP cofactor; (ii) nondirect proton transfer between the amino substrate and the lysine residue; and (iii) formation of the external aldimine after the dissociation of the lysine residue. The overall reaction is exothermic (12.0 kcal/mol), and the rate-limiting step is the second one with 12.6 kcal/mol for the activation energy.

1. INTRODUCTION Pyridoxal 50 -phosphate (PLP) is a derivative of vitamin B6 and acts as a cofactor in a myriad of chemical reactions involving amino acids.1 The enzyme commission (EC) has already catalogued more than 150 distinct enzymatic activities of this type of enzyme, which includes decarboxylations, racemisations, transiminations, retroaldo cleavages, and β or γ eliminations. The study of PLP enzymes is one of the most fascinating frontiers in enzymology, not only because of their unrivaled versatility as catalysts but also because they are involved in many cellular processes.2 Their importance is further underscored by the number of receptors that have been identified as drug targets. For example, inhibitors of γ-aminobutyric acid aminotransferase (GABA ATase) are used in the treatment of epilepsy,3 serine hydroxyl methyl transferase (SHMT) has been identified as a target for cancer therapy,4 and inhibitors of ornithine decarboxylase (ODC) are employed in the treatment of African sleeping sickness.5 Functional defects in PLP enzymes have also been implicated in several pathologies, such as homocystinuria.6 Understanding the function of this important group of enzymes has thus become an important milestone in medicine and biorelated research areas to develop new molecules capable of impairing enzymatic activity and especially to design improved protein-based catalysts. A comprehensive understanding of PLP-related enzymes or even their classification in different families is not a straightforward task, because beyond the wide variety of reactions that they can catalyze, these enzymes have diverse quaternary structures. Some enzymes can be found active as monomers, others as dimers, and some of them as tetramers or even hexamers.7,8 Although the scope of the PLP-catalyzed reactions initially appears to be bewilderingly diverse, there is a simple unifying principle: In the resting state the cofactor is covalently attached to the ε-amino group of an r 2011 American Chemical Society

active site lysine, forming an internal aldimine. Once the R-amino substrate interacts with the active site, the lysine residue dissociates from PLP, and the substrate becomes covalently bonded to it, generating a new Schiff base with PLP. This intermediate is commonly called the external aldimine. Only after this step, the mechanistic pathway for each PLP-catalyzed reaction diverges as it is depicted in Scheme 1. The conversion between the internal (lysine PLP-imine) and external aldimine (substrate PLP-imine) has therefore a preponderant role in all PLP-related enzymes, and it is a crucial step in their activation. Moreover, as the inverse reaction is required for the enzymatic turnover, this reaction is critical toward the overall activity of these enzymes. Experimental studies have provided some clues regarding this step. All PLP intermediates (internal and external aldimines) have distinct absorption bands and change only by a few nanometers from system to system, which makes their identification and characterization straightforward.9,10 The postulated mechanism proposes that at physiological conditions, the imine linkage of the internal aldimine is protonated in order to form a more electrophilic and reactive iminium ion (Scheme 1). Upon amino substrate binding, it is suggested that the nitrogen N10 of the substrate attacks carbon C8 of the internal aldimine. Subsequently, the lysine residue dissociates from the PLP endorsing the formation of the external aldimine. During this process, it is proposed that the dissociation of the lysine residues and the attachment of the amino acid substrate is not accomplished in a single step. Instead, it might involve the formation of a transient but stable geminal diamine intermediate.

Received: April 26, 2010 Published: March 29, 2011 1356

dx.doi.org/10.1021/ct1002219 | J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Scheme 1. Currently Accepted Mechanism for the Transimination Reaction: A Common Step in the Catalytic Mechanism of All PLP-Dependent Enzymesa

a

P represents the phosphate group and A an R-amino substrate.

The presence of a geminal diamine intermediate is assumed to be common in most PLP-dependent enzymes. However, there is some controversy about the mechanism by which the transimination reaction occurs. Snell and Jenkins were the first ones to propose the formation of a geminal diamine intermediate.11 This mechanism was supported by further experimental work.12,13 Other studies suggest the involvement of a two-fold additionelimination type of mechanism instead14,15 or even the possible involvement of the phosphate group of PLP in the proton shuttle mechanism.14,7 However, the involvement of the 50 -phosphate of PLP in the transimination process was later on discarded, because this group is buried in the protein and cannot interact directly with the region where the transimination reaction occurs. Several studies proposed that the main function of this group is to act as an anchor to hold the PLP group in the active site.16 This is rather important after the formation of the external aldimine, in which the PLP group becomes disconnected from the enzyme (it was earlier bonded through Lys69), and therefore it requires some sort of interaction that maintains it bonded to the active site and aligned in a proper orientation for effective catalysis. Yet, there are some PLP-dependent enzymes in which the 50 -phosphate group is indeed involved in the catalysis, as it is the case of the enzyme glucogen phosphorylase17 and GDP-4-keto-deoxymannose-3-dehydratase (CoID).18 However, in those cases the substrates are not amino acid-related molecules, and therefore the transimination reaction does not occur. All the other mechanisms remain as open possibilities, and several studies were conducted aiming to elucidate the most favorable pathway for the transimination reaction present in the PLP-dependent enzymes. The high rate complexity of the enzymatic transimination has made it very difficult to study the transimination reaction by experimental means. Even so several works have supported the presence of the geminal diamine intermediate in the transimination reaction.1922 However, the correct mechanism that could elucidate its formation and involvement in the transimination reaction is still poorly understood. During the past decade, several computational studies have emerged that have tried to explain the transimination reaction at the atomic level detail. The

studies performed by Mu~noz et al.23 were pioneers in this field and revealed that the conversion between the internal aldimine into the external aldmine requires the direct participation of at least one water molecule. In this proposal, it is suggested that the phenolate oxygen of the PLP should receive a proton from the amino substrate (through the water molecule) to favor the formation of the geminal diamine intermediate. Another study, performed by Zhao et al., proposes that the transimination reaction occurs through the direct proton transfer between both amino groups that are bonded to the PLP. However, prior to the formation of the geminal diamine intermediate, it is required the direct participation of the 50 -phosphate group of PLP.24 Despite the enormous interest that this reaction has received, none of these studies have explored the participation of key active site amino acids in the transimination reaction. In fact, analyzing several X-ray structures, we found that near the PLP cofactor there are two conserved residues that are capable of catalyzing this reaction, i.e., Cys360 and Tyr389 (considering the PDB code 2OO0). Both residues are pointing to the place in which the transimination reaction should occur. In addition, there is one water molecule that is pointing to the same reaction spot and can favor the transimination reaction without the direct involvement of the phenolic oxygen of the PLP, as it has been proposed by Mu~noz et al.25 In this paper, we explore all these possibilities using quantum mechanics calculations in order to enhance the knowledge about the transimination reaction and highlight the most favorable mechanism involved in this process.

2. METHODOLOGY 1. Model. The model system used in this work was based on the X-ray structure 2OO0 determined by Dufe et al.,26 which contains the human ornithine decarboxylase. This structure has a good resolution (1.9 Å) but lacks the substrate inside the active site. To acquire this information, we superimposed the X-ray structure 1F3T, which contains the ornithine descarboxylase from Trypanosoma brucei complexed with the putrescine (the ODC’s reaction product) in the binding pocket.27 Near the active site region, the two structures can be almost superimposed, 1357

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

Figure 1. Active site model taken from the PDB structure 2OO0 used. All the residues that were used in this study are shown in ball and sticks. The substrate was built based on the X-ray structure 1F3T that contains the product of the reaction. The atoms marked with F* were kept frozen during the geometry optimizations, and the symbol T* highlights the place where the truncation of the residues took place.

allowing us to model the correct position of the substrate inside the active site of the X-ray structure 2OO0. The reactant of the reaction was subsequently modeled substituting the hydrogen atom, which is attached to the carbon atom that is closer to the PLP ring, by a carboxyl group. This model was subsequently simplified, eliminating all the amino acids that do not interact directly with the PLP-imine. The final model contained the PLP, the substrate, Lys69, Cys360, and Tyr389 (Figure 1). The selected amino acids were initially truncated at the R carbon. However, in order to maintain the net of hydrogen bonds within the residues of the model, we decided to keep the main chain of each residue and protonate the carboxylate and the amino groups. The calculations that we have performed have shown that this latter approach turned out to be more satisfactory than the reverse, because it improved the robustness of the model and, in some cases, resulted in the reduction of the activation energy by 3 kcal/mol and the reaction barrier by 1 kcal/mol. In order to simplify the model, we decided to substitute the phosphate group of PLP by a methyl group. We have chosen the methyl group instead of a hydroxyl group because the first one ensures and maintains the stability of the active site, without requiring the inclusion of additional residues. Nevertheless, the differences between the barriers obtained for the first step of the studied mechanism were very similar in both cases [Ea = 5.3 and Er = 2.5 kcal/mol vs Ea = 4.5 and Er = 1.9 kcal/mol (values obtained with DFT6-31G(d))]. In spite of all simplifications the model had in total 108 atoms (Figure 1). This model was then subjected to geometry optimizations. To keep the optimized structures close to the X-ray structure, some atoms were kept frozen as depicted in Figure 1. 2. Methods. Density functional theory (DFT) calculations were performed with the Gaussian09 software package.28 All structures were fully optimized and characterized both at the B3LYP level2932 and using the new hybrid exchange correlation functional proposed by Zhao and Truhlar, M0633,34 together with the 6-31G(d) basis set. The 6-31þG(d,p) basis set was also

ARTICLE

used to optimize all geometries of the rate-limiting step and analyze the effect of polarization and diffusion functions in the hydrogen atoms. In all geometry optimizations, we first searched for the transition state starting from the reactants. This was obtained with a scan in which the reaction coordinate that we were interested in was shortened or stretched. The transition states were subsequently fully geometry optimized, starting from the structure of the higher energy point of the scans. The reactants and the products, associated with it, were determined through internal reaction coordinate (IRC) calculations. In all cases, the geometry optimizations and the stationary points were obtained with standard Gaussian convergence criteria. The transition-state structures were all verified by vibrational frequency calculations, having exactly one imaginary frequency with the correct transition vector, even using frozen atoms, which shows that the frozen atoms were almost free from steric strain. The final electronic energies were calculated using the allelectron 6-311þþG(3df,2pd) basis set, using the functionals B3LYP, M06, and M06-2X. These structures were the optimized geometries obtained with the M06/6-31G(d) level of theory. Zero-point corrections, thermal, and entropic effects (T = 310.15 K, P = 1 bar) were added to all calculated energies, with the 6-31G(d) basis set. To estimate the solvation effects of the rest of the enzyme, single point calculations on the optimized geometries were performed with IEF-PCM, as implemented in Gaussian 09,28 with the 6-311þþG(3df,2pd) basis set. This feature is of particular importance to the study of enzymatic catalysis because the use of a continuum model is normally taken as an approximation to the effect of the global enzyme environment in a reaction. A dielectric constant of ε = 4 was chosen to describe the protein environment of the active site in agreement with previous suggestions.3538 The atomic spin density distributions were calculated at the M06 level employing a Mulliken population analysis, using the 6-31G(d) basis set.

3. RESULTS AND DISCUSSION This section will be divided in two main parts. In the first part, we will discuss the most favorable pathway for the transimination process at the M06-2X/6-311þþG(3df,2pd)//M06/6-31G(d) level. In the second part, we will discuss a small benchmarking study that was performed to compare the results that were obtained with the very popular B3LYP and two hybrid meta exchange correlation functionals of Zhao and Truhlar, M06 and M06-2X. 3.1. Transimination Mechanism. Our first task in this study was to review all the available experimental data concerning this subject. In this analysis, we found three interesting X-ray structures of the enzyme ornithine decarboxylase (enzyme that catalyzes the decarboxylation of ornithine to putrescine). Each PDB reveals the PLP in three different states that correspond to snapshots of the active site during the initial steps of the catalytic process (Figure 2). The first PDB structure (PDB code: 7ODC)31 clearly shows a lysine residue bonded to the PLP cofactor. This should correspond to the initial state of the enzyme (internal aldimine) and is a common characteristic of all PLP-dependent enzymes. The PDB structure 1F3T20 shows the lysine residue dissociated from the PLP cofactor. This structure resembles what should be found in the external aldimine but with the putrescine in place of the ornithine. The active site region of the mutant PDB structure 1SZR13 shows a structure where the PLP cofactor is bound simultaneously to the lysine residue and to the amino substrate. This 1358

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Figure 2. Active site topology of three PDB files of ornithine decarboxylase available in the protein databank (* mutant PDB structure of ODC).

Scheme 2. First Step of the Transimination Reaction of PLP-Dependent Enzymesa

a

P stands for a phosphate group.

structure suggests that the transimination reaction might not occur in a single step but rather in three subsequent steps, i.e.,: (i) formation of the tetrahedral intermediate with the lysine and the ornithine bonded to the PLP cofactor; (ii) proton transfer between both amino substrates; and (iii) subsequent dissociation of the lysine residue with the concomitant formation of the external aldimine. The enzyme ornithine descarboxylase (ODC) was used in the subsequent sections as a model of all PLP-dependent enzymes. Since the transimination reaction is a common feature in all PLPdependent enzymes, the mechanism described in the following sections can be transferred unequivocally to any PLP-dependent enzyme. Step 1: Formation of the Tetrahedral GeminalDiamine Intermediate. Taking into account the mutant X-ray structure 1SZR we have tested in this step whether the formation of the geminaldiamine intermediate occurs, and if such an intermediate is a stable compound (Scheme 2). The computational results have shown that once the ornithine enters the active site (the amino substrate), it interacts with the PLP nearby carbon C8 and oxygen O1 (2.60 Å). The optimized structure of the reactant indicates that in the initial state, Lys69 remains tightly bonded to carbon C8 of the PLP through the NH group (1.33 Å). The bond length between atoms C8 and N9 is characteristic of secondary amines (1.31 Å), which means that the lysine residue is attached to PLP by a single bond

(bond lengthCdN ∼ 1.28 Å). This is also emphasized by the covalent nature of the double bond between atoms carbon C2 (0.40 au) and oxygen O1 (0.75 au), characteristic of the carbonyl group (1.25 vs 1.23 Å). All these results allow us to conclude that resonance structure b from scheme 2 is the one that describes better what must be found in the reactants of this reaction. The positive charge of the system is mainly distributed along the extended π systems of the PLP ring (0.86 au). It is interesting to note that the charge of carbon C8 is slightly more positive (0.36 au) than the other atoms of the PLP ring (on average ∼ 0.28 au). This effect turns carbon C8 more prone to accept the nucleophilic attack of the amino substrate, therefore, favoring the reaction. This result is in agreement with previous suggestions11 that point to the fact that atom C8 of PLP becomes bonded to the R-amino group of the substrate. The transition-state structure (Figure 3) of this reaction is characterized by an imaginary frequency of 165i cm1. The optimized structure indicates that the ε-amino group of Lys69 remains tightly bonded to carbon C8 of the PLP (1.35 Å). The R-amino group of ornithine comes closer to the same center (1.92 Å), and carbon C8 changes its hybridization from sp2 to sp3, remaining slightly positively charged (0.39 au). The formation of such a tetrahedral intermediate is favored by the pyridine ring that acts as an electron sink, thus allowing to concentrate the excess of negative charge around oxygen O1 (0.64 au) and to maintain the electropositive nature of carbon C8 (0.40 au). 1359

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation In the optimized geometry of the products, both amino groups become covalently bound to carbon C8 of the PLP, generating a stable geminaldiamine intermediate, as it was observed in the mutant protein 1SZR14. The bond length between the nitrogen N9 of Lys69 and PLP slightly elongates to 1.41 Å (1.33 Å before), while the distance between the nitrogen N10 from the substrate and PLP gets shortened to 1.56 Å. Due to the formation of the tetrahedral intermediate, the bond length between carbons C7 and C8 becomes slightly elongated to 1.53 Å (1.44 Å before), and the same is also true for the bond length between carbon C2 and oxygen O1 (1.25 Å vs 1.23 Å in the reactants). Oxygen atom O1 starts to interact very closely with the proton of the R-amino group of ornithine through a hydrogen bond (1.59 Å), and this creates a pseudoring as depicted in Scheme 2. This rearrangement stabilizes this region, allowing the positive charge of the Ramino group of the substrate and the negative charge of oxygen O1 to spread along the pseudoring (0.14 au). This reaction is characterized by an activation energy of 10.8 kcal/mol and is exothermic by 5.7 kcal/mol. An overall evaluation of the reaction allows us to conclude that the PLP group has a central role in this reaction by delocalizing the negative charge through its π system. The charge becomes mainly concentrated at the carbonyl group that simultaneously allows guiding and aligning the amino substrate with carbon C8.

Figure 3. Optimized structure of the transition state of step 1.

ARTICLE

These two effects have an important role in the reaction: (i) The first favors the ammoniumcation nature of the carbon located at position 8 (Scheme 2); and (ii) the second enhances the formation of a strong hydrogen bond between oxygen O1 and the amino group of ornithine that helps to align the substrate inside the active site. In the product of this reaction, this strong hydrogen bond is preserved, and it favors the creation of a pseudoring between PLP and the amino substrate that stabilizes the formation of the geminaldiamine intermediate. It may be interesting to note that the same type of interaction between the lysine residue and PLP was never observed in our calculations. This occurs because at the beginning of the reaction, the lysine is located in a perpendicular plane to that of the PLP and in the opposite direction of the carbonyl group. This behavior is in agreement with the available X-ray structures as depicted in Figure 4. Comparing the optimized structure of the tetrahedral geminal diamine obtained by computational means and the one obtained experimentally (X-ray structure 1SZR), we can find many similarities (Figure 4). Both amino groups of Lys69 and ornithine become covalently attached to carbon C8 of the PLP at similar distances (1.41 and 1.56 Å vs 1.31 and 1.50 Å). The same is also true across all the bond lengths from the PLP ring, which emphasizes the equivalence between the theoretical model and the X-ray structure. It must be said that there is a very good match between both structures presented in Figure 4. The only difference lies on the conformations adopted by Lys69 and ornithine for this study. The model contains all the residues that are directly involved in the reaction (in total they account to 108 atoms) but lacks all the others (much less important) that sometimes help to keep the conformation and orientation observed in the active site of the PDB structure 1SZR. However, this is not a crucial aspect for this reaction, because no active site residue is involved in it. Therefore taking into account the similarity observed between all the bond lengths between carbon C8 of PLP and Lys69 and ornithine, we can conclude that the structures presented in Figure 4 are very similar. Step 2: Proton Transfer between the Amino Groups of the Tetrahedral GeminalDiamine Intermediate. To complete the transimination process, Lys69 must dissociate from the PLP cofactor. Because of the formation of the tetrahedral intermediate, this reaction is not straightforward, and a proton must be first transferred from the R-amino group of the substrate to the ε-amino group of Lys69 (Scheme 3).

Scheme 3. Possible Proton-Transfer Pathways (a and b) between the Amino Acid Substrate and the Lysine Residuea

a

Pathway a consists of one step. Pathway b consists of two steps (dashed line): b1 and b2 with X = Tyr, Cys, or H2O. 1360

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Figure 4. Left: X-ray structure of the PDB entry 1SZR. Right: Optimized structure of the products of the first step (the geminaldiamine intermediate).

Figure 5. Active site topology of the X-ray structure 1SZR (1SZRa:PLP from chain D; 1SZRb:PLP from chain C) and 1F3T (PLP from chain A).

Multiple pathways can be drawn for this step, as is depicted in Scheme 3. The easiest one involves the direct proton transfer from the R-amino group of ornithine (N10) to the ε-amino group of Lys69 (N9). This reaction corresponds to pathway a of Scheme 3. Pathway b of the same scheme, involves the presence of a neighbor and proton donor/acceptor active site residue (residue X in Scheme 3) capable of catalyzing the proton transfer. Two residues can play such role in the mechanism, considering the PDB structures 1SZR and 1F3T of ornithine descarboxylase, and they are Cys360 and Tyr389 (Figure 5). The first residue is at 3.01 Å from carbon C8 of PLP. Tyr389 is not so close to this center (3.58 Å), but the flexibility of the side chain still allows it to interact with both amino groups (figure 5). Another hypothesis is the involvement of one water molecule in the reaction (Figure 5), similar to what was proposed by Mu~noz et al.23 A closer inspection of the crystal structure 1SZRb reveals the presence of the water molecule 1477 (1SZR numbering) that is very close to the site where the reaction happens (but different to what was modeled by Mu~ noz et al.). Moreover, this molecule is stabilized by several residues, such as Asp361, Asp332, and Tyr 389, suggesting that it might be conserved in the active site. In order to model these reactions and understand which is, from a kinetic and thermodynamic point of view, the most favorable pathway, we have created four scenarios of the active

site. Each scenario contains the residues/water molecule capable of catalyzing the proton transfer. The conformation of each residue was carefully chosen beforehand using a set of rotamer libraries and considering the protein surroundings. Scenario 1 corresponds to pathway a of Scheme 3 and does not contain any residues between the two amino groups (direct pathway). Scenarios 24 contain a tyrosine (X = OH), a cysteine (X = SH), and a water molecule (X = H2O), respectively. These three scenarios will be used to study pathway b of Scheme 3, and the acquired information will allow us to understand which residue is more likely to catalyze the proton transfer. The optimized structures of the reactants are very similar between each scenario. Both amino groups remain covalently bonded to carbon C8 of the PLP, but while the ε-amino group of Lys69 (0.25 au) is at 1.40 Å from carbon C8, the R-amino group of ornithine (0.15 au) is at 1.61 Å from the same atom. The discrepancies that are observed in the bond lengths are a consequence of the protonation state of each amino group. The R-amino group from ornithine has two protons, and when it is bonded to carbon C8, it becomes more electropositive (0.12 au), a situation that results in a weaker chemical bond. The ε-amino group of Lys69 has only one proton, and therefore the nitrogen is more electronegative (0.54 au), a situation that results in a stronger chemical bond between the nitrogen atom and carbon C8 stronger. 1361

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Figure 6. Optimized transition-state structures of step 2 of the transimination reaction of scenario 1 (direct pathway), 2 (catalyzed by Tyr389), 3 (catalyzed by Cys360), and 4 (catalyzed by water1477).

The distance between carbon C7 and carbon C8 of PLP remains stabilized around 1.52 Å, and the positive charge delocalized around the PLP ring (0.85 au), as before. The interaction between oxygen O1 and the R-amino group of ornithine is maintained similarly to what is observed in the products of the first step of the mechanism. In all scenarios, hydrogen HN (Scheme 3) of the R-amino group of ornithine is on average 2.25 Å away from the nitrogen atom of the ε-amino group of Lys69. In scenarios 24, this hydrogen atom makes an additional hydrogen bond with the sulfur atom of Cys360 (2.69 Å), with the oxygen atom of Tyr389 (2.15 Å), and with the oxygen atom of water molecule 1477 (1.88 Å), respectively. Hydrogen HX of the cysteine and tyrosine residues do not establish a hydrogen bond with the nitrogen atom from Lys69, as it would be expected. Such types of interaction only occur when the water molecule is placed between the two amino groups (3.3 Å). The transition state of scenario 1 (Figure 6) is characterized by an imaginary frequency of 1609i cm1. In this scenario, proton HN is found halfway between both amino groups (∼1.35 Å), and the

charge distribution kept the same trend that was observed beforehand. In scenarios 2 and 3 (Figure 6), the transition states are characterized by imaginary frequencies of 1350 and 1114 cm1, respectively. Proton HX gets closer to the nitrogen atom of Lys69 (1.28 Å at scenario 2 and 1.08 Å at scenario 3), whereas hydrogen HN is halfway between the nitrogen atom of ornithine and the oxygen atom of tyrosine (1.28 and 1.22 Å, respectively) in scenario 2 and the sulfur atom of cysteine (1.29 and 1.70 Å, respectively) in scenario 3. Accordingly, the sulfur atom of Cys360 (scenario 3) becomes more electronegative (0.27 vs 0.04 au in the reactants), while the oxygen atom of Tyr389 (scenario 2) becomes more electropositive (0.16 vs 0.27 au in the reactants). The charge distribution around the PLP cofactor and the pseudoring remains unchanged in both cases. The transition state of scenario 4 (Figure 6) is characterized by an imaginary frequency of 1476i cm1. Both protons, HN (0.64 au) and HX (0.60 au) are equally shared between both amino groups and the oxygen atom of the water molecule (1.25 Å on 1362

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Figure 7. Energies involved in all of the four scenarios used to study the second step of the transimination reaction. Scenario 1: direct proton shuttle between Lys69 and the amino substrate. Scenarios 24: model reactions in which Tyr360, Cys360, or Water1477 are involved in the proton transfer between the amino substrate and Lys69, respectively (activation energies are colored in black and the reaction energies in light gray).

average), forming a tightly bound complex. This means that in the transition state of scenario 4, we have a H3Oþ ion (0.55 au) within both amino groups (0.70 au). Similarly to what happens in scenarios 2 and 3, the charge distribution around the PLP group and the pseudoring remains unchanged (þ0.86 au), In the product of each reaction, the proton that was previously bound to the R-amino group of ornithine (N10) transfers to the ε-amino group of Lys69 (N9) (HX in scenarios 24 and HN in scenario 1) (on average 1.03 Å bond length). Both amino groups remain bound to carbon C8 of PLP, but the bond length between the R-amino group of ornithine and carbon C8 of PLP is now shorter (1.42 vs 1.56 Å), whereas the bond length between the ε-amino group of Lys69 and carbon C8 of PLP elongated to 1.6 Å (1.41 Å beforehand). In scenarios 24, HN becomes covalently bonded to the oxygen atom of Tyr389 (0.98 Å), to the sulfur atom of Cys360 (1.35 Å), and to the hydroxyl group of the water molecule (1.00 Å). Similar to what happens in the reactants and in the transition state, the charge distribution around the PLP ring remains practically unchanged in all scenarios (0.85 au). The charge of the pseudoring in the products of the reaction becomes more electronegative (0.22 au) than was observed in the reactant and the transition state of this reaction (0.19 au). The graphic shown in Figure 7 resumes the energetic profile obtained from all the studied reactions for step 2 of the transimination reaction. The results show that under physiological conditions, all of the studied pathways can occur and efficiently catalyze the proton transfer between both amino groups of the tetrahedral geminaldiamine intermediate. All the studied pathways have one characteristic in common: The reactions are almost thermoneutral, which indicates the feasibility of the reaction in both directions, as it is predicted experimentally.12 However, some reactions are more favorable than others. Accordingly, the involvement of Tyr389 (scenario 2)

or of the water molecule (scenario 4) in the proton transfer tends to decrease considerably the activation energy. The direct proton transfer of pathway a (scenario 1) or the involvement of Cys360 (scenario 3) in the reaction is less favorable, accounting to activations energies higher than 20 kcal/mol. The energies involved in scenarios 2 and 4 are comparable (Eascenario2 = 14.9 vs Eascenario4 = 12.6 kcal/mol and Erscenario2 = 4.4 vs Erscenario4 = 2.0 kcal/mol), but it is evident that the reaction is more favored when the water molecule is directly involved in the reaction. This can be explained considering the volume and the flexibility of the water molecule that, when compared to Tyr398, promote a better and closer interaction of both amino groups (∼1.2 Å). As a consequence, the difference between the energies of the reactants and the transition state tends to decrease, favoring the reaction. However, the differences observed in these energies are not sufficient to disclose which one should be preferred. Therefore, both of them should be equally capable of catalyzing the proton transfer. This means that if water is available in the active site, then it can catalyze the reaction. Otherwise, Tyr389 catalyzes the reaction without requiring a significant energetic cost. Looking at the pKa values of the residues involved in the proton transfer, we can understand why Tyr389 and the water molecule are better proton shuttles than Cys360. Under physiological conditions, the pKa values of a cysteine, a tyrosine, and the water molecule are respectively 8.09, 10.07, and 15.74. This means that Cys360 has a greater tendency to be found in the dissociated form, rather than Tyr389 or the water molecule. This behavior makes Cys360 a good proton donor but not a good proton acceptor. Accordingly, Cys360 is a worst proton shuttle than Tyrs389 and the water molecule. Looking at the charge distribution and the bond lengths in the optimized models of the transition-state structures, we can conclude basically the same. 1363

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Scheme 4. Charge Distribution and Bond Lengths in the Transition State of Scenarios 4 and 3

Scheme 5. Third Step of the Transimination Reactiona

a

P stands for phosphate group.

From Scheme 4, we can see that when the water molecule is used for the proton transfer, the bonds that are created/cleaved are placed within a pseudoring in which the charge is well delocalized (a similar behavior is observed with Typr389). The same does not occur with Cys360. The ring is largely deformed, and the charge is mainly located on the sulfur atom (0.27 au), while the ε-amino of Lys69 remains slightly positively charged. These results show that the proton shuttle seems to be only favored when the atom that is involved in the proton transfer can fit between both amino groups and that it can behave as a donor and acceptor of electrons. Only when these conditions are ensured, it is observed the formation of a pseudoring that enhances the delocalization of electrons in order to enhance the proton transfer. Step 3: Formation of the External Aldimine. The third step of the transimination reaction involves the formation of the external aldimine. This reaction is the reverse of the first step and involves the dissociation of Lys69 from PLP (Scheme 5). In the reactants, both amino groups remain covalently bonded to carbon C8 of PLP. However, due to the differences in the protonation state of both nitrogen atoms, the bond length between carbon C8 of PLP and the ε-amino group of Lys69 is more stretched (1.60 Å) than the bond length between carbon C8 and the amino group of ornithine (1.42 Å). The full complex

retains the tetrahedral geometry around carbon C8 of PLP, and oxygen O1 continues to interact very closely with the proton of the R-amino group through a hydrogen bond (1.92 Å). This rearrangement continues to stabilize this region, allowing the positive charge of the amino group and the negative charge of oxygen O1 to spread along with the atoms of the ring, which in total accounts for ∼0.03 au The transition state of this reaction is characterized by an imaginary frequency of 1618i cm-1 (Figure 8). In this structure the bond length between the R-amino group of ornithine (N10) and carbon C8 of PLP decreases to 1.35 Å, whereas the bond length between the ε-amino group of Lys69 (N9) and carbon C8 elongates to 1.98 Å. The tetrahedral geometry around carbon C8 of PLP breaks up, and ornithine adopts a conformation parallel to the PLP ring and interacts with oxygen atom O1 (1.84 vs 1.91 Å in the reactants). Consequently, the nitrogen atom N9 becomes more electronegative (0.01 au) than in the reactants (0.19 au). The charge around the pyridine ring and the pseudoring remains unchanged. In the products of this step, the external aldimine is obtained. The tridimensional structure resembles what is observed in the PDB structure 1F3T that contains a similar intermediate but with putrescine instead of ornithine (Figure 5). Lys69 is now disconnected 1364

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

Figure 8. Optimized transition state of the third step of the transimination reaction.

from the PLP cofactor (4.40 Å), while ornithine becomes tightly bound to carbon C8 of PLP (1.31 Å). The hydrogen bond between oxygen atom O1 and the R-amino group is still present (1.74 Å). This rearrangement continues to stabilize this region, allowing the positive charge of the amino group and the negative charge of oxygen O1 to spread along the atoms that compose the pseudoring (∼0.03 au). This means that the resonance structure b of Scheme 5 is the one that describes better the product of this reaction. The formation of the external diamine requires a very small activation energy (2.2 kcal/mol), and the reaction is exothermic in 4.3 kcal/mol. Comparing this step with the first step of the mechanism, this reaction is very similar but has a lower activation barrier of about 8.6 kcal/mol. Such behavior can be explained taking into account the type of bond that is cleaved/formed in each reaction. While in the first step, the bond that is formed/cleaved involves the ornithine residue (the substrate), in the last step the chemical bond that is cleaved/formed involves a lysine residue (which is part of the enzyme). Although these two residues are very similar (ornithine lacks only one CH2 group in the side chain, when compared to lysine), the way they bind to the PLP group is quite different. Lysine binds to PLP through the NH2 group located in the side chain (ε-amino group), whereas ornithine binds to PLP with the amino group located in the main chain (R-amino group). This means that the binding of ornithine to PLP is from a steric point of view less favorable, as the neighboring carboxylic group hinders the approach of the amino group to PLP. This explains why the activation energy of the first step is higher than that of the last step. Looking at the reaction energies, we also see that both reactions are exothermic. In the first step, this means that in spite of the steric effect that results from the binding of ornithine to PLP, the formation of the geminal diamine intermediate is very favorable. In the last step, we have the opposite situation, i.e., the formation of the external aldimine is more stable than the geminal diamine intermediate. This happens because once ornithine binds to PLP, its amino group makes a strong hydrogen bond with the carbonyl group of PLP, and this interaction overcomes the steric penalty arising from the approach of ornithine. In the products, the same type of interaction exists, which stabilizes the formation of the external aldimine and favors the dissociation of the lysine residue. The stabilization of the

ARTICLE

ornithinePLP complex arises from the formation of a pseudoring that seems to favor the delocalization of the charge around it. It must be noted that the same type of interaction between the lysine residue and PLP was never observed. This occurs because at the beginning of the reaction, the lysine is located in a perpendicular plane to that of the PLP and in the opposite direction of the carbonyl group. This is in agreement with the available X-ray structures and the model used in this study kept the same orientation, which underscores its robustness. 3.2. Functional Benchmarking: B3LYP vs M06 Family. In order to understand if the energetic profile that was obtained in this study could be influenced by the functional that was used, we performed a small benchmarking exercise comprising B3LYP and the two hybrid meta exchangecorrelation functionals M06 and M06-2X. The latter have shown to be very accurate for thermodynamics and kinetics. In addition, we have also compared the influence of the functional and the basis set in the geometries. To evaluate the influence of the functional in the optimized geometries in this process, all the transition states were recalculated, and the products and reactants of each reaction were obtained through IRC calculations. The obtained results have shown that the differences between the optimized structures obtained with M06/6-31G(d) and with B3LYP/6-31G(d) amount to less than 0.7 Å. In addition, a visual inspection of each minimum revealed that the atoms that are involved in the formation or cleavage of chemical bonds have a very small root-mean-square deviation (rmsd, below 0.12 Å). These results show, that in this type of system and when studying this type of chemistry, the optimized geometries that are obtained with B3LYP are similar to those that are obtained with the M06 functionals. We have also tested the effect of the basis set in the geometry optimizations. For this purpose, we have reoptimized all the species of the second step of the mechanism (rate limiting step) with the functional M06 and the 6-31þG(d,p) basis set (a doubleζ basis set augmented with polarized functions for heavy atoms and hydrogens). The results revealed that there is not a significant difference in the optimized geometries. When compared with those obtained with the 6-31G(d) basis set, the rmsd of the reactants (rmsd = 0.15 Å), TS (rmsd = 0.08 Å) and products (rmsd = 0.07 Å) are indeed very small and on average below 0.1 Å. In order to evaluate if the energetic profile of the transimination reaction did differ from the functional that was used, the very popular B3LYP functional was used instead of the hybrid meta exchangecorrelation M06 functional to recalculate the final energies of the optimized models. The M06-2X functional was also used to check for the effect of doubling the HF exchange, which is known to affect the barriers. For this purpose we used the M06 geometries and recalculated the energy using B3LYP and M06-2X, with the 6-311þþG(3df,2pd) basis set. The results are presented in Figure 9. From all the employed functionals, the barriers obtained with M06-2X functional were always within the values obtained with B3LYP and M06. In general there are no significant differences between them. Exceptions are limited to the activation energy of the second step, scenario 2, and the reaction energy of the third step. Comparing the values obtained with M06-2X and M06, we observe that M06 tends to result in higher barriers. Comparing the values obtained with the M06-2X and the B3LYP functionals, we observe that B3LYP underestimates the barriers. As the experimental values for the barriers are unknown, we cannot pinpoint exactly which functional is giving us the most exact result. 1365

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Figure 9. Energetic profile for the transimination reaction calculate with the B3LYP, M06, and M06-2X functionals with the 6-311þþG(3df,2pd) basis set.

Figure 10. Most favorable pathway for the transimination reaction (P stands for phosphate group). The results were obtained withy the M06-2X/6311þþ(3df,p2d) level of theory.

However, these results show that B3LYP tends to give lower values for the chemical barrier, when compared to the values obtained with the M06 family. This is in agreement with other previous studies, and therefore we may conclude that the barriers obtained with M06 and M06-2X should be more exact than those that are obtained with B3LYP. As the M06-2X functional includes twice the amount of exchange compared with M06 functional (which usually is favorable to describe activation energies) and resulted in energy values that were between the B3LYP and M06 extremes, we choose to use these values in the description of the reaction steps on the previous sections. The most important aspect of this study is that the chosen pathway (kinetically most favorable) is always the same, independently of the functional. The functionals affect the accuracy of the activation and the reaction energies but

do not affect the discrimination between different hypotheses for the catalytic mechanism.

4. CONCLUSIONS The external aldimine is the common central intermediate for all enzymatic and nonenzymatic reactions that are catalyzed by the PLP cofactor. Divergence in reaction specificity occurs from this point, which means that the formation of external aldimines from the internal PLP-aldimines represent the first level of catalysis in all PLP-dependent reactions. The results obtained in this work show that the transimination reaction is very favorable, but it is not accomplished by a single step as it is generally accepted involving instead three subsequent 1366

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation

ARTICLE

Figure 11. Conformational rearrangements adopted by the amino substrate (ornithine) and Lys69 around the PLP ring. Left: Conformational rearrangement observed in the external and internal aldimines. Right: Conformational rearrangement adopted in all the transition-state structures studied in this paper.

steps, as depicted in Figure 10. This is in agreement with previous NMR studies performed by Chan-Huot.22 Multiple pathways are possible for the conversion of the geminaldiamine intermediate into the external aldimine. In this article, we explored all reasonable pathways including the direct migration of the proton between both amino groups, the involvement of several active site residues, or even the participation of a water molecule. The results have shown that the rate of the reaction is lower if a water molecule or Tyr389 are involved in this process. The direct proton transfer or the involvement of the catalytic active site residue Cys360 were shown to be less favorable. The most favorable pathway occurs when the water molecule is directly involved in the reaction. The overall reaction accounts for 10.8 kcal/mol and is exothermic by 12.0 kcal/mol. The second step of the reaction is the rate-limiting one amounting to 12.6 kcal/mol. In Figure 10, it is displayed the most favorable pathway obtained for the transimination reaction (only scenario 4 of step 2 is displayed). These results also point out the importance of the PLP in the course of the reaction. It allows the interchange of carbon C8 between sp2 and sp3 hybridizations, without requiring a significant energetic cost. This key feature favors the formation of the tetrahedral intermediate that is the key driving force behind the conversion between the internal and external aldimines. During this process, the excess of charge in the system becomes lodged at oxygen O1 of the PLP ring. This effect is very important during the transimination reaction because it not only allows the attraction of the substrates to PLP but it also serves as a guide during the binding/dissociation of the amino substrate/lys69 residue to carbon C8 of the PLP cofactor. Moreover, this atom makes a strong hydrogen bond with the amino group of the amino substrates, favoring the alignment of both parts of the molecule on the same plane. This rearrangement, improves the stereoelectronic effect of the system, ensuring the maximum overlap of the extended π system and the stabilization of the full system. This is in agreement with Dunathan’s hypothesis39 postulated almost 50

years ago in which he predicted that the most active form of the external aldimine had a cisoide conformation underneath the same plane in order to favor the occurrence of the subsequent reactions (Figure 11, Structure a). Furthermore, Dunathan’s proposed that the bonds that are formed/broken in the PLP system should adopt a perpendicular plane to that of the pyridoxal imine system. All the transition-state structures presented in this article show exactly the same sort of conclusion, which confirms once again the accuracy of the early proposals made by Dunathan (Figure 11, Structure b). In the last 10 years, two different proposals for the transimination reaction have been suggested. One of those proposals was suggested by Salva et al.,23 in which he suggested that the transimination process occurs in seven steps, requiring the participation of one/two water molecules. That mechanism differs substantially from the one presented in this work since it required the direct participation of oxygen O1 as the key intermediate that shuttles the proton transfer between the Lys69 and the substrate. The mechanism presented here is accomplished in fewer steps (only three instead of seven) and is from an energetic point of view more favorable (Ea ∼ 10 vs Ea ∼ 20 kcal/ mol). In addition, it shows that in the absence of water molecules, nearby the active site, Tyr389 can catalyze this reaction. Another study was performed by Zhao et al.,24 in which he proposed that the transimination reaction requires the direct participation of the phosphate group of PLP, adopting a similar role that is played by oxygen O1 in the mechanism proposed by Salva. The energies involved in that process are comparable to those presented in this paper, but they do not involve the formation of the geminaldiamine intermediate that is observed experimentally. From all the analyzed data, we can conclude that the mechanism presented here is, both from thermodynamic and kinetic points of view, more favorable than the previous suggestions and includes the formation of all intermediates that are observed experimentally. We believe therefore that this mechanism should be general for all PLP-dependent enzymes, corresponding to the 1367

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368

Journal of Chemical Theory and Computation first chemical transformation that is catalyzed by all PLP requiring enzymes that have amino acids as substrates.

’ AUTHOR INFORMATION Corresponding Author

*E-mail: [email protected].

’ REFERENCES (1) Toney, M. D. Arch. Biochem. Biophys. 2005, 433, 279–287. (2) Eliot, A. C.; Kirsch, J. F. Annu. Rev. Biochem. 2004, 73, 383–415. (3) Kleppner, S. R.; Tobin, A. J. Expert Opin. Ther. Targets 2001, 5, 219–239. (4) Thorndike, J.; Pelliniemi, T.; Beck, W. Cancer Res. 1979, 39, 3435–3440. (5) Wang, C. C. Annu. Rev. Pharmacol. Toxicol. 1995, 35, 93–127. (6) Mudd, S. H.; Finkelstein, J. D.; Irreverre, F.; Laster, L. Science 1964, 143, 1443–1445. (7) Jansonius, J. N. Curr. Opin. Struc. Biol. 1998, 8, 759–769. (8) Grishin, N.; Phillips, M.; Goldsmith, E. Protein Sci. 1995, 4, 1291–1304. (9) Korpela, T.; M€akel€a, M.; L€onnberg, H. Arch. Biochem. Biophys. 1981, 212, 581–8. (10) Robitaille, P.; Scott, R.; Wang, J.; Metzler, D. J. Am. Chem. Soc. 1989, 111, 3034–3040. (11) Lehtokari, M.; Puisto, J.; Raunio, R.; Korpela, T. Arch. Biochem. Biophys. 1980, 202, 533–539. (12) Snell, E. E.; Jenkins, W. T. J. Cell Comp. Physiol. 1959, 54, 161–77. (13) Tobias, P.; Kallen, R. J. Am. Chem. Soc. 1975, 97, 6530–6539. (14) Vazquez, M.; Munoz, F.; Donoso, J.; Blanco, F. J. Phys. Org. Chem. 1992, 5, 142–154. (15) Hershey, S.; Leussing, D. J. Am. Chem. Soc. 1977, 99, 1992–1993. (16) Chang, Y. C.; McCalmont, T.; Graves, D. J. Biochemistry 1983, 22, 4987–4993. (17) Fletterick, R. J.; Sprang, S. R. Acc. Chem. Res. 1982, 15, 361–369. (18) Cook, P. D.; Holden, H. M. J. Biol. Chem. 2008, 283, 4295–4303. (19) Counts, K. G.; Wong, I.; Oliveira, M. A. Biochemistry 2007, 46, 379–386. (20) Phillips, R. S.; Miles, E. W.; Cohen, L. A. Biochemistry 1984, 23, 6228–6234. (21) Schirch, L. J. Biol. Chem. 1975, 250, 1939–1945. (22) Chan-Huot, M.; Sharif, S.; Tolstoy, P. M.; Toney, M. D.; Limbach, H. Biochemistry 2010, 49 (51), 1081810830. (23) Salva, A.; Donoso, J.; Frau, J.; Mu~noz, F. Int. J. Quantum Chem. 2002, 89, 48–56. (24) Zhao, Z.; Liu, H. J. Phys. Chem. B 2008, 112, 13091–13100. (25) Ortega-Castro, J.; Adrover, M.; Frau, J.; Salva, A.; Donoso, J.; Mu~noz, F. J. Phys. Chem. A 2010, 114, 4634–4640. (26) Dufe, V. T.; Ingner, D.; Heby, O.; Khomutov, A. R.; Persson, L.; Al-Karadagi, S. Biochem. J. 2007, 405, 261–268. (27) Jackson, L. K.; Brooks, H. B.; Osterman, A. L.; GoldSmith, E.; Philips, M. Biochemistry 2000, 39, 11247–11257. (28) Frisch, M. J.; Trucks, G. W.; Cheeseman, J. R.; Scalmani, G.; Caricato, M.; Hratchian, H. P.; Li, X.; Barone, V.; Bloino, J.; Zheng, G.; Vreven, T.; Montgomery, J. A.; Petersson, G. A.; Scuseria, G. E.; Schlegel, H. B.; Nakatsuji, H.; Izmaylov, A. F.; Martin, R. L.; Sonnenberg, J. L.; Peralta, J. E.; Heyd, J. J.; Brothers, E.; Ogliaro, F.; Bearpark, M.; Robb, M. A.; Mennucci, B.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Rendell, A.; Gomperts, R.; Zakrzewski, V. G.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H. Gaussian 09; Gaussian, Inc.: Wallingford, CT, 2009. (29) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. J. Phy. Chem. 1994, 98, 11623–11627.

ARTICLE

(30) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200–1211. (31) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (32) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. B 1988, 37, 785–789. (33) Zhao, Y.; Truhlar, D. G. Theor. Chem. Acc. 2008, 119, 525–525. (34) Zhao, Y.; Truhlar, D. G. Acc. Chem. Res. 2008, 41, 157–167. (35) Chen, S.; Fang, W.; Himo, F. Theor. Chem. Acc. 2008, 120, 515–522. (36) Cerqueira, N. M. F. S. A.; Fernandes, P. A.; Eriksson, L. A.; Ramos, M. J. Biophys. J. 2006, 90, 2109–2119. (37) Amata, O.; Marino, T.; Russo, N.; Toscano, M. J. Am. Chem. Soc. 2009, 131, 14804–14811. (38) Ramos, M. J.; Fernandes, P. A. Computational enzymatic catalysis. Acc. Chem. Res. 2008, 41 (6), 689–698. (39) Dunathan, H. Proc. Natl. Acad. Sci. U.S.A. 1966, 55, 712.

1368

dx.doi.org/10.1021/ct1002219 |J. Chem. Theory Comput. 2011, 7, 1356–1368