Why Do Divalent Metal Ions Either Promote or Inhibit Enzymatic Reactions?

Divalent metal ions are essential to many enzymatic reactions involving nucleic acids, but their critical and specific role still needs to be uncovered. Restriction endonucleases are a prominent group of such metal-requiring enzymes. Large scale accurate simulations of Mg- and Ca-BamHI elucidate the mechanism of the catalytic reaction leading to DNA cleavage and show that it involves the concerted action of two metal ions and water molecules. It is also established that what is decisive for the dramatically different behavior of magnesium (a cocatalyst) and calcium (an inhibitor) are kinetic factors and not the properties of the prereactive states of the enzymes. A new perspective is opened for the understanding of the functional role of metal ions in biological processes.

Many enzymatic reactions involving nucleic acids are known to critically rely on the cooperative action of divalent metal cofactors, such as magnesium, that are incorporated into the protein as positively charged ions (1). Although quite well established, their functional role has not yet received a satisfactory explanation and is still a matter of debate and speculation. Among these reactions, the ones more often investigated experimentally pertain to the family of enzymes called type II restriction endonucleases. These enzymes are highly site-specific; their biological function is that of cleaving DNA at both strands (2), and they all require metal cofactors. This is a process of primary biochemical importance. For example, this process provides protection to bacterial cells against bacteriophage infection, occurs naturally during programmed cell death (apoptosis), and has been exploited in genetic engineering.
The experimental information accumulated over the years (see e.g. Refs. 1-10) mostly concerns the structural properties of the protein-DNA complexes and the effect of mutations on their activity. More is needed to unravel the precise mechanism driving the cleavage of the phosphodiester bond at the DNA strands (11). Although there is little doubt that the reaction takes place via a concerted nucleophilic substitution, it is still unclear whether the attacking nucleophile is a water molecule or a hydroxide ion. More importantly the activity of these enzymes crucially depends on the presence and identity of the metal ions, 1 but simple features, such as ionic size, electronegativity, or average coordination number, do not simply correlate with it. Other issues concern the extent to which the characteristics of the initial state of the reaction or the presence of water molecules can affect its development. Here we make definitive progress toward answering all these questions with the aid of accurate computer simulations that provide the so far missing insight into the mechanism of the catalytic reaction and a reliable evaluation of its energetics.
The validity of any computational approach to enzymatic processes critically relies on a quantum mechanical description of the active site, but the modeling often suffers from the size limitation typical of ab initio calculations. To overcome this restriction, combined quantum-classical approaches have been proposed for a long time (12)(13)(14)) that partition the system into two subsystems: the close environment of the reaction site, treated quantum mechanically (QM), 2 and the outer region described with a molecular mechanics model (MM). However, only recently has this method started to provide the required level of accuracy (15)(16)(17). Our QM/MM approach combines two well established methodologies, namely density functional theory (DFT) (18) and the GROMOS force field (19), and includes a new careful treatment of the electrostatic interactions (20).
We apply this method to the study of the restriction endonuclease BamHI (1, 4 -6) and focus on the comparison of the activities of magnesium and calcium, which are experimentally known to act as cocatalyst and inhibitor, respectively, of the DNA cleavage reaction. This drastically different behavior of magnesium and calcium is quite common, and the reasons for it are largely unknown. The advantage that this specific system provides for a thorough investigation consists in the availability of x-ray structures that could be used as reference for our simulations: these structures refer to the complexes obtained by cocrystallizing the protein-DNA complex with calcium and with manganese 3 and can be considered as representative of the pre-and postreactive states, respectively. In particular, by comparing the two states, the structure of the protein appeared to be largely insensitive to the reaction as evidenced by the overall root mean square displacement of only 0.33 Å between the C␣ positions of the pre-and postreactive states. Two metal atoms (A and B) were observed at the active site, corroborating earlier proposals for a two-metal mechanism (4), the nature of which will emerge clearly from our results.

COMPUTATIONAL METHOD
The QM subsystem was treated in the framework of DFT (18) using Becke's approximation (21) for the exchange energy functional and the * The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked "advertisement" in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.
Lee-Yang-Parr approximation (22) for the correlation energy functional. The QM calculations were performed with the Car-Parrinello molecular dynamics (CPMD) code 4 using a plane-wave basis set for the description of the valence electron wave functions (energy cut-off 70 Rydberg) and norm-conserving angular momentum-dependent pseudopotentials for the description of the core-valence interaction. The MM subsystem was treated within the GROMOS force field (19), which is well established (see e.g. Ref. 23) for the determination of both the structure and the dynamics of proteins in aqueous solutions (24). The parametrization used corresponds to the 43A1 version of the force field. The interface between the two subsystems was modeled using the link-atom method (14). The QM subsystem was enclosed in a box (of edges 43.36, 48.55, and 39.21 atomic units) and treated as an "embedded cluster" namely using the Hockney method (25) and accounting for the influence of the MM subsystem via explicit calculations of van der Waals and electrostatic interactions. The electrostatic potential of the MM atoms was calculated in the QM region exactly within a cut-off radius of 15 Å and including a Poisson-Boltzmann reaction field correction (26) outside. The influence of the QM atoms on the MM atoms was accounted for via the electrostatic interaction. The QM region was represented by a distribution of effective charges, which we obtained from fitting to the ab initio electrostatic potential using a novel and suitable procedure that avoids common errors due to inconsistency between quantum mechanically derived parameters and the rest of the force field (20).
The initial structure is taken from the x-ray structure of the calcium complex (6) (Protein Data Bank entry 2BAM), and hydrogen atoms are added as required by the force field (19). The system was solvated in a rectangular box of 14,819 simple point charge water molecules (24). Sodium ions were added in proximity of the DNA chain to provide charge neutralization.
The optimization procedure (of the electronic and geometric variables) consisted of three consecutive steps. First the structure of the QM system was optimized keeping the positions of the MM atoms fixed, consistent with the experimental observation of a minor difference between pre-and postreactive states. Next the MM geometry was relaxed with molecular dynamics runs and optimized. The interaction with the QM system (fixed) was modeled as discussed above. Finally the structure of the QM subsystem was optimized again. Only minor changes were observed in these last two steps.

RESULTS AND DISCUSSION
In our simulations, the system is immersed in water to provide a realistic model (see Fig. 1) for the study of the enzy-matic reaction. As the starting configuration, we used the crystallographic structure of the BamHI-DNA complex with calcium and optimized its structure (Fig. 2b) in water; we then repeated this procedure after replacing calcium with magnesium (Fig. 2a). As attacking agent, we considered a hydroxide generated either by a water molecule (4, 6) deprotonated in situ (intrinsic mechanism) or by an external source (extrinsic mechanism) (27). Fig. 2 refers only to the modeling of the latter reaction. In this case, only minor geometrical differences emerged in the initial states of the two enzymes apart from the expected decrease of the metal-oxygen distances (on average from 2.6 to 2.2 Å) on passing from calcium to magnesium. Water molecules are part of the metal coordination shell and are crucial to guarantee the preferred coordination number for both the A and B ions. If the hydroxide is replaced by a water molecule as attacking species, more significant differences between the A and B ion environments in the initial state as well as between those of calcium and magnesium are observed. For example, the magnesium B ion loses the phosphate O-2 oxygen as binding partner in a 5-fold coordination. However, the structural discrepancies we have detected in the prereactive states of the two metal-requiring proteins are minor and do not reflect a significantly different propensity of calcium and magnesium to donate electrons or polarize the surrounding bonding pattern. This is best illustrated by the spatial distribution of the charges obtained from the topological analysis of the electron density along the lines of Bader's theory (28) (Fig. 3). The atomic charges associated with the A and B ions are the same, and only a small difference is found between calcium and magnesium. Both A and B ions help stabilize a cluster distribution of negative charges localized at the active site. This is in agreement with the experimental observation of an enhancement of the stability of the protein-DNA complex in the presence of divalent ions (5). In the two-metal atom mechanism suggested in Ref. 6, the role of the A ion is that of stabilizing hydroxyl ions generated by the Glu-113-driven deprotonation of water, and the failure of Ca 2ϩ to act as cocatalyst was interpreted as failure to provide the concentration of metalhydroxyl ions necessary for the nucleophilic attack because of an acidity lower than that of either Mg 2ϩ or Mn 2ϩ . The stabilizing role of the A ion is confirmed by our calculations for magnesium, but the same is true for calcium. From our anal- ysis, it has become clear that the active site does not contain any relevant feature in the initial state of the binding that could provide a simple explanation for the dramatically different behavior of magnesium and calcium.
We proceeded to investigate the reaction mechanisms and evaluate the energy profiles. For the phosphodiester cleavage by the hydroxide, be it extrinsic or intrinsic, we chose its approach distance to the phosphate as reaction coordinate. For the deprotonation of the water molecule, we selected its distance from the protonating amino acid (Glu-113 in Fig. 2). In the case of magnesium, the mere bond breaking induced by the hydroxide costs 14 kcal/mol in the extrinsic and 29 kcal/mol in the intrinsic process. Substantially higher energy barriers are obtained when magnesium is replaced by calcium, namely 31 and 52 kcal/mol, respectively. The presence of calcium instead of magnesium also increases, although to a lesser extent, the energy barrier relative to deprotonation of water in situ, namely from 20 to 26 kcal/mol. On the other hand, the cost to transfer a hydroxide from the bulk solution to the metal coordination shell, which is not explicitly calculated here, is expected to be only slightly different in the two cases. An estimate of the free energy change involved in this process can be obtained (29) from the difference between the pK a of the metal ion-bound water and the experimental external pH (30). This  (a and c) and calcium (b and d) complexes. All distances are in Å. 3. Prereactive states of both magnesium (a) and calcium (b) complexes. The atomic charges obtained using Bader's approach (28) are shown in the range Ϫ2.0 (deep blue) to ϩ2.0 (deep red). Values outside this range were condensed in the extremes for sake of clarity, and aliphatic hydrogen atoms were omitted. Also displayed is the charge distribution projected onto a plane containing the scissile bond. Note the strong similarity all over the region of the active site. Bader's atomic charges may differ only slightly from the nominal ionic charges, e.g. they are ϩ2.0 for magnesium and ϩ1.9 for calcium. amounts to about 5 kcal/mol for magnesium and 7 kcal/mol for calcium. In conclusion, our results unambiguously show that the kinetics of the reaction is decisive for the dramatically different behavior of the two divalent ions. Our calculations also show that the extrinsic mechanism is the most plausible one, which is consistent with recent model calculations of the free energy gain in the presence of the enzyme for the formation of the hydroxide in the two cases (27). An activation free energy of ϳ18 kcal/mol was derived from the measured value of the rate constant of the phosphodiester bond cleavage reaction (30) catalyzed by the magnesium enzyme. If we subtract the amount due to the inclusion of the hydroxide in the system from it, we obtain a value of ϳ13 kcal/mol. Comparison with our calculated value of ϳ14 kcal/mol for the energy barrier indicates that the entropy contribution is minor.

FIG.
A complete understanding of the difference in the calculated energy barriers can be obtained by investigating the active site in the region of the "transition state" (Fig. 2, c and d). These structures show the presence of a trigonal bipyramid for the reacting moiety with both the attacking and leaving groups in axial positions. This is consistent with the geometry of the transition state involved in the hydrolysis of phosphate esters (31). Here the structural differences between the calcium and magnesium environments are much more pronounced than in the initial state, and their role in increasing the barrier for the transition to happen can now be detected. From our simulations of the active magnesium enzyme, it becomes clear that the role of the B ion is crucial in stabilizing the pentacoordinated transition state of the phosphate. In doing this, the B ion reduces its coordination with oxygen from six to four. This is what allows a water molecule to exit the coordination sphere of the metal and position itself to stabilize, via hydrogen bonding, another water molecule, which drives the final step of the reaction, i.e. the protonation of the O-3Ј leaving group of the phosphate. An oxygen coordination shell of four is not so uncommon for magnesium. This is not the case for calcium, which tends to overcoordinate (with respect to a 6-fold coordination) rather than undercoordinate. In fact, when the system is forced to accommodate the pentacoordinated phosphate, resulting in a loosening of the bonding between the phosphate and the metal ion, the calcium B ion does not readjust its coordination but responds to keep all possible oxygen partners. Calcium thus fails to stabilize a situation that would be favorable for the transition to happen. As Fig. 2d shows, the A ion also tends to increase its coordination, resulting in a strained rearrangement of the atoms around it in contrast to the case of the magnesium enzyme in which the A ion keeps a quasisymmetric octahedral structure during the entire reaction (see Fig. 2, a and c). All this explains the much higher energy required by the calcium complex with respect to the magnesium complex in allowing the final reaction to take place.
These simulations have elucidated the mechanism by which magnesium ions coadjuvate the BamHI enzyme in cleaving the phosphodiester bond in DNA strands and finally unraveled why calcium ions inhibit the catalytic reaction. These findings have confirmed earlier suggestions (6) that a crucial role in the reaction is played by water molecules as well as by a second metal ion, beyond the one binding the attacking agent, and elucidated their concerted mode of action. This leads to a transition state with a configuration that perturbs the coordination shell of the metal ions in a way that is compatible with the binding mode of magnesium but not with that of calcium. On this basis, barriers lower than those for the calcium complex can reasonably be expected for complexes incorporating manganese, cobalt, zinc, and cadmium, which are more similar to magnesium and also act as cocatalysts (6). Calculations from first principles were necessary to account for the nontrivial difference of the electronic structure of the two enzymes that drastically affects the kinetics of the reaction rather than the static properties of the prereactive states. The success of this investigation is largely due to our computational method, which relies on a careful description of all relevant interactions and on a sufficiently large size of the QM subsystem. In fact, as we have verified, reducing the size of the QM subsystem (50 -100 atoms) drastically changes the nature of the frontier orbitals and thus the chemical behavior of the system, and ignoring dominant interactions with the MM region makes it unstable. The large discrepancy found in the values of the energy barriers of the enzyme with the two different metals (ϳ20 kcal/mol) is such that future enhancements of the accuracy of these calculations (e.g. an improvement of the DFT functional or a dynamic search for the transition state) will not alter the explanation we have provided here for their different chemical behavior. These results also reveal a new perspective on the large variety of enzymatic reactions relying on metal ion cofactors (1) and offer an efficient and accurate method that can now be confidently extended to verify several mechanistic hypotheses and to calculate specific physical quantities such as reaction activities.