Article pubs.acs.org/JPCB
Elucidation of the Catalytic Mechanism of a Miniature Zinc Finger Hydrolase Abir Ganguly,† Trung Quan Luong,‡ Oliver Brylski,§ Michael Dirkmann,∥ David Möller,∥ Simon Ebbinghaus,§ Frank Schulz,∥ Roland Winter,‡ Elsa Sanchez-Garcia,*,† and Walter Thiel*,† †
Max-Planck-Institut für Kohlenforschung, 45470 Mülheim an der Ruhr, Germany Fakultät für Chemie und Chemische Biologie, Physikalische Chemie I, Technische Universität Dortmund, 44227 Dortmund, Germany § Fakultät für Chemie und Biochemie, Physikalische Chemie II, Ruhr-Universität Bochum, 44780, Bochum, Germany ∥ Fakultät für Chemie und Biochemie, Organische Chemie I, Ruhr-Universität Bochum, 44780 Bochum, Germany ‡
S Supporting Information *
ABSTRACT: To improve our mechanistic understanding of zinc metalloenzymes, we report a joint computational and experimental study of a minimal carbonic anhydrase (CA) mimic, a 22-residue Zn-finger hydrolase. We combine classical molecular dynamics (MD) simulations, quantum mechanics/molecular mechanics (QM/MM) geometry optimizations, and QM/MM free energy simulations with ambient and highpressure kinetic measurements to investigate the mechanism of the hydrolysis of the substrate p-nitrophenylacetate (pNPA). The zinc center of the hydrolase prefers a pentacoordinated geometry, as found in most naturally occurring CAs and CA-like enzymes. Two possible mechanisms for the catalytic reaction are investigated. The first one is analogous to the commonly accepted mechanism for CA-like enzymes: a sequential pathway, in which a Zn2+-bound hydroxide acts as a nucleophile and the hydrolysis proceeds through a tetrahedral intermediate. The initial rate-limiting step of this reaction is the nucleophilic attack of the hydroxide on pNPA to form the tetrahedral intermediate. The computed free energy barrier of 18.5 kcal/mol is consistent with the experimental value of 20.5 kcal/mol obtained from our kinetics experiments. We also explore an alternative reverse protonation pathway for the hydrolase, in which a nearby hydroxide ion from the bulk acts as the nucleophile (instead of a zinc-bound hydroxide). According to QM/ MM MD simulations, hydrolysis occurs spontaneously along this pathway. However, this second scenario is not viable in our system, as the tertiary structure of the hydrolase lacks a suitably positioned residue that would act as a general base and generate a hydroxide ion from a nearby bulk water molecule. Hence, our combined theoretical and experimental study indicates that the investigated minimal CA mimic retains the essential mechanistic features of CA-like enzyme catalysis. The high-pressure experiments show that its catalytic efficiency can be enhanced by applying hydrostatic pressure. According to the simulations, more drastic improvements might be afforded by mutations that make the reverse protonation pathway accessible.
■
INTRODUCTION Zinc metalloenzymes are ubiquitous in biology and participate in many life processes.1 The zinc ion, which is often essential for the activity of these proteins, can play either a structural role in stabilizing a biologically relevant conformation of the protein, a catalytic role by directly participating in enzyme catalysis, or even a cocatalytic role by assisting catalysis without active involvement.2 Therefore, a fundamental understanding of catalytic zinc sites is important to elucidate the biological role of zinc. The study of catalytic zinc centers also has implications for the pharmaceutical and biotechnological industries, in the development of drug-binding enzymes, or of novel enzymes capable of carrying out biotechnologically relevant reactions.3,4 Zn-hydrolases constitute the largest family of zinc-containing enzymes. They catalyze hydrolytic reactions using an active-site water molecule or a related nucleophile. While the reaction mechanisms of Zn-hydrolases are diverse, they share common © XXXX American Chemical Society
features, extending the relevance of mechanistic studies beyond the peculiarities of a specific system. In the initial step of most Zn-hydrolase catalyzed reactions (Scheme 1), the tetrahedral catalytic zinc center formed by the Zn2+ ion, the side-chains of three amino acid residues (typically histidines), and a water molecule, becomes pentacoordinated by including the substrate. The Zn2+ ion plays a dual role of polarizing both the substrate and the catalytic water molecule. In the next step, the coordinated water molecule is activated, by the Zn2+ ion and/or a general base such as a nearby glutamate or asparagine, to generate a Zn2+-bound hydroxide ion, which acts as the active nucleophile and attacks the substrate to form a tetrahedral intermediate.5,6 The final step is a ligand exchange that yields Received: May 24, 2017 Revised: June 5, 2017
A
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B
metal center, typically giving hexacoordinated Zn2+ ions even in a protein environment.13,14 Importantly, in addition to charge transfer and polarization contributions, the coordination of Zn in biological systems may also have a sizable covalent component.6 The use of hybrid quantum mechanical/molecular mechanical (QM/MM) approaches15−19 overcomes such MMrelated deficiencies, since the Zn2+ ion and its coordination sphere are treated quantum mechanically (as accurately as needed) while the rest of the protein and the bulk solvent are described classically. Recently Srivastava and Durani20 introduced a minimal de novo zinc metalloenzyme that exhibits the catalytic activity of naturally occurring carbonic anhydrases. The 22-residue protein was designed using a synthetic αββ Zn-finger protein as a template, grafting a catalytic Zn(His)3O site at the base of the αββ fold and the substrate binding site at the mouth (Figure 1).
Scheme 1. Schematic Representation of Reaction Pathway 1a
a
In this mechanism, the Zn2+ ion simultaneously activates the carbonyl group of the substrate and the coordinated water molecule. The hydroxide formed by deprotonation of the Zn2+-bound water acts as a nucleophile.
the hydrolyzed products and regenerates the tetrahedral zinc center. In some cases, an alternate reverse protonation mechanism has been proposed,7,8 in which an activated water molecule from the bulk (instead of the Zn2+-bound hydroxide) acts as nucleophile (Scheme 2). Figure 1. De novo designed structure20 of the Zn-finger protein.
Scheme 2. Schematic Representation of Reaction Pathway 2a
The catalytic activity of the αββ enzyme was assessed with respect to the hydrolysis of the substrate p-nitrophenylacetate (pNPA), and its efficiency was found to be ∼1000-fold lower than that of naturally occurring CAs.20 While being low, this reactivity is comparable to that of the best-known artificial CA mimics.21 Naturally occurring Zn-finger proteins typically employ zinc as a structural element for stabilizing the tertiary structures.22 Owing to its simplistic design, the αββ hydrolase allows the study of a catalytic zinc center in a reduced environment. Here, we report classical and QM/MM molecular dynamics (MD) simulations, QM/MM geometry optimizations, and QM/MM free energy simulations in conjunction with ambient and high-pressure kinetic measurements, for a systematic investigation of catalysis in this de novo Zn-finger hydrolase using p-nitrophenylacetate (pNPA) as substrate. First, we determine the preferred coordination geometry of the catalytic Zn2+ ion. QM/MM MD simulations along with high-pressure kinetics measurements indicate that the catalytic Zn2+ ion prefers a pentacoordinated geometry as in most CAs and CAlike enzymes. We investigate two different pathways of pNPA hydrolysis. Our primary focus is on the generally accepted mechanism (pathway 1) for CA-like enzymes in which a Zn2+bound hydroxide/activated water acts as a nucleophile (Scheme 1). Extensive QM/MM geometry optimizations reveal a minimum energy path that is sequential, passing through a tetrahedral intermediate, again in line with the reported behavior of CAs. QM/MM free energy simulations confirm the sequential nature of this pathway and give a free energy barrier of 18.5 kcal/mol for the rate-limiting step of the reaction, which is close to the experimental value from our kinetics measurements. We also explored a less conventional, reverse protonation mechanism (pathway 2) in which a
In this “reverse protonation” mechanism, the Zn2+ ion activates only the substrate carbonyl and does not take part in the activation of the nucleophile. Upon activation by a base, a bulk water molecule acts as nucleophile. a
While zinc-assisted catalytic reactions have been studied extensively by both experimental and computational techniques,5 the active-site features that allow zinc to be catalytically active are not yet fully understood. Several attempts to design artificial functional zinc centers using approaches such as de novo design,9 introducing metal binding sites into naturally occurring proteins,10 or modifying structurally characterized zinc sites11,12 have typically resulted in metal sites with diminished affinity and little or no activity.4 On the computational side, it is challenging to properly describe the zinc coordination architecture at the force field level. Zn2+, a d10 metal ion, offers negligible ligand-field stabilization and can adopt a variety of coordination geometries depending on the environment.3 It prefers to be hexacoordinated in water, whereas in proteins tetrahedral geometries are more common. The standard nonbonded description of the Zn2+ ion in popular molecular mechanical force fields often does not provide the correct coordination geometry around the B
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B
the dRMSD calculations were Zn2+−N3 (His10), Zn2+−N3 (His15), Zn2+−N3 (His22), and Zn2+−O1 (pNPA). After the classical MD solvent equilibration, the full system was subjected to QM/MM MD simulations to generate starting structures with pentacoordinated zinc for the subsequent QM/ MM reaction path calculations (see the initial part of the Results and Discussion for a detailed description). Starting structures for QM/MM geometry optimizations were obtained from the QM/MM MD simulations after equilibration of solvent and substrate. To prepare the system for QM/MM calculations, all solvent atoms lying beyond a distance of 20 Å from the protein were deleted. Since we are interested mainly in the local changes in and around the active site, all atoms lying outside a sphere of 20 Å radius (centered at the C2 atom of the substrate) were kept fixed during the QM/MM calculations. The water molecule bound to Zn2+ was manually converted to a hydroxide ion to drive the nucleophilic attack on pNPA. The water deprotonation step was not considered explicitly in these calculations because the identity of the nucleophile that activates the Zn2+-bound water molecule is unknown. Thus, the activation of the water molecule to a hydroxide was assumed to occur prior to the nucleophilic attack of the latter on the substrate. For the quantum mechanics treatment, density functional theory (DFT) was used at the PBE0/def2-SVP level of theory. In a recent study, the PBE0 functional was shown to give the best results for a similar Zn2+ containing system.32 As in the classical MD simulations, the MM atoms were described using the CHARMM27 force field with CMAP corrections. All QM/MM calculations were performed with the ChemShell software33,34 using TURBOMOLE 6.335 for the QM calculations and DL-POLY36 for the MM calculations. The minimum free energy path (MFEP) of the catalytic reaction was calculated using an enhanced sampling technique that combines traditional umbrella sampling (US) simulations with a finite temperature string method. This method has been detailed in several recent studies;37−39 here we provide a brief summary. In this technique, the reaction path is represented as a string, essentially a curve, passing through the multidimensional space of a few chosen coordinates that are important to the reaction. The curve, initially representing a guess reaction pathway, is divided equally into a series of images, where each image corresponds to certain values of the chosen reaction coordinates. Umbrella sampling is performed on each image with restraining potentials imposed on each reaction coordinate centered on their respective current values. From the restrained simulations the average values of the reaction coordinates are calculated and the string is updated by fitting a new curve to these average values. The new string is again divided equally into a series of images to obtain a new set of centers of restraining potentials for the US simulations. This process of updating the string iteratively based on restrained dynamics is repeated until the string is converged, i.e. until the root-meansquare deviation of the updated string from the current string falls below a certain threshold. Finally, the restrained simulations from all images and iterations were put together to obtain the multidimensional free energy surface of the reaction, and the converged string on this surface corresponds to the MFEP of the reaction. The string simulations were initiated from the minimum energy path (MEP) obtained from QM/MM geometry optimizations. Representative structures along the MEP were further subjected to 0.5 ps of QM/MM MD equilibration
hydroxide or activated water molecule from the bulk acts as a nucleophile (Scheme 2). QM/MM MD simulations suggest that such reaction is extremely facile, with spontaneous hydrolysis occurring in several independent trajectories. Finally, based on our results, we suggest ways to improve the enzymatic activity of the investigated hydrolase.
■
COMPUTATIONAL METHODS As initial coordinates we used the de novo structure of a Znfinger hydrolase reported previously by Srivastava and Durani (Zn-hydrolase variant B2).20 The substrate, p-nitrophenyl acetate (pNPA), was docked into the αββ fold of the hydrolase using the program Autodock.23 The protein−substrate complex was immersed in a cubic box containing rigid TIP3P24 water molecules. The box was first neutralized by adding two Cl− ions, and then NaCl was added to the solvent to give a physiological salt concentration of 0.15 M. Following a standard equilibration procedure, in which the protein was kept fixed to its initial coordinates and only the substrate, solvent, and ions were allowed to move, production trajectories were propagated in the constant NVT ensemble, with SHAKE constraints25 imposed on the bonds involving hydrogen atoms. Periodic boundary conditions and a 1 fs MD time step were used. The smooth particle mesh Ewald algorithm26 was employed to compute the long-range electrostatic interactions and Langevin dynamics was used to maintain constant temperature during the simulations. The partial charges on the pNPA atoms were calculated using a procedure described elsewhere.27 Briefly, the electrostatic potential of the pNPA molecule was calculated at the HF/6-31G* level of theory. Thereafter, a set of partial charges, centered on the pNPA atoms, was generated to reproduce the electrostatic potential. The simulations were performed using the NAMD molecular dynamics package28 with the CHARMM27 force field29 including CMAP corrections.30 The MD trajectories were analyzed using standard tools implemented in GROMACS,31 specifically g_rms for calculating the root mean-square deviations (RMSDs) of residues along the trajectories and g_cluster for clustering the configurations in the trajectories. The trajectories were clustered using the gromos algorithm with a RMSD cutoff of 1 Å. Water distributions around the protein active site were analyzed by water isodensity plots. For this purpose, a 12 Å side cube was defined centered on the C2 (pNPA) atom. A threedimensional grid with 1 Å spacing was constructed on this cube and the occupancies of water molecules at the grid cells were calculated from the respective MD trajectories. The grid cell occupancies were normalized with respect to the total number of frames of the respective trajectories. The grid cells in the isodensity plots were colored to highlight the top 5% and 10% of water occupancies in dark and light cyan, respectively. The conformation of the hydrolase active site during the MD simulations was monitored by calculating distance RMSDs (dRMSDs) of the Zn2+ catalytic center along the MD trajectories, according to the equation: dRMSD =
∑i ∑j (dia, j − dib, j)2 N
(1)
where i and j are indices of atoms, a and b are two different configurations, and N is the number of distances considered. In this study, dRMSDs were calculated with respect to the initial de novo structure, and the active-site distances considered in C
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B during which the reaction coordinates were harmonically restrained to their initial values. During the initial iterations 26 images were considered along the string, and later the number of images was increased to 52 to ensure a better sampling along the MFEP. Restraining potentials with force constants of 100−200 kcal/mol/Å2 were used for the US simulations. The string was updated after 200 fs of MD for each image, except for the last two iterations when the string was almost converged and 1 ps of MD for each image was run to improve sampling statistics. The criterion for string convergence was chosen to be similar to that used in a previous study.38 The string was considered to be converged when the RMSD of all the reaction coordinates with respect to their respective average values from the previous 10 iterations, across all images, fell below a threshold of 0.5 Å. Further details regarding the convergence of the string are given in the Supporting Information (Figure S1). The multidimensional weighted histogram analysis method (WHAM) was used to unbias the simulations. A convergence threshold of 0.001 kcal/ mol was used for the WHAM equations.
[E]KM 1 K [E] [E] 1 1 = + = + M v0 vmax vmax [S] kcat kcat [S]
where [E] is the concentration of the enzyme, [S] is the initial substrate concentration, and νmax is the maximum achievable rate. The ratio kcat/KM is the enzymatic efficiency. The pressure effect on the rate of the reaction is obtained from the Eyring equation: ⎛ ∂ln kcat ⎞ ΔV # ⎜ ⎟ =− RT ⎝ ∂p ⎠T
EXPERIMENTAL METHODS The 22-mer peptide was obtained from Eurogentec with a purity of >90% by HPLC, confirmed by MALDI MS, and dissolved in 20 mMTris-HCl, pH 7.5 to a concentration of 2 mM. Kinetic Assays. The enzymatic activity was measured using a high-pressure stopped-flow system (HPSF-56 of Hi-Tech Scientific).40,41 The system operates in the pressure range of 0.1−200 MPa. All measurements were carried out in Tris-HCl buffer (20 mMTris, pH 7.5, with 20 μM Zn(ClO4)2) and at a peptide concentration of 20 μM (peptide: Ac-D-Pro1-Gln2Trp3-Glu4-Asp5-Ser6-Ala7-Lys8-Leu9-His10-Gly11-Gly12Gly13-Gly14-His15-Thr16-Tyr17- D -Pro18-Gly19-Arg20Asn21-His22-NH2; corresponding to variant B2 from ref 20). The pNPA concentration was varied from 0.5 to 2 mM. The kinetics of product formation (p-nitrophenol, pNP) was followed by optical absorbance at 405 nm based on experimentally determined molar extinction coefficients of pNP using UV spectroscopy (ε405 nm = 13 520 M−1 cm−1). The temperature of the samples was controlled by a thermostat and kept constant at (25 ± 0.1) °C in all experiments. The pressure was increased stepwise from 0.1 to 200 MPa. A large excess of substrate that is saturating the enzyme was used to avoid significant back reaction. Under the concentration conditions used here, the enzyme kinetics can be described by the Michaelis−Menten scheme:42,43 k1
k2
k −2
■
E + S ⎯→ES ⎯→ ⎯ E+P k−1
(4)
where p is the pressure, T is the absolute temperature, R is the ideal gas constant, and ΔV# is the activation volume of the reaction, which is the difference of the volume of the transition state and the ground state of the enzyme−substrate complex (ES), i.e., ΔV# = V# − VES. If these volumes are compressed to different extents under high pressure, a change of ΔV# with pressure will be observed. Plate Reader Experiments. Enzymatic activity was detected via absorption over time using a plate reader (Clariostar by BMG Labtech). The assay was started using the internal injection system to add substrate prior dissolved in acetonitrile and two seconds of double-orbital shaking at 200 rpm. Data points were further collected in 30 s intervals for at least 20 min to clearly distinguish between background and peptide reaction. The measurements were performed in 20 mM Tris-HCl buffer (Carl Roth) at pH 7.4, 25 °C and a final acetonitrile concentration of 2.66% (v/v). The peptide was used at a concentration of 20 μM and Zn(ClO4)2 was added in 1:1 ratio to induce Zn-finger formation. Clear 96-well polystyrene flat bottom plates (Nunc 96, Sigma-Aldrich) served as the reaction vessel for absorption measurements. Evolution of pNP during the experiment was detected measuring the absorption at 410 ± 5 nm. Determined absorption values were converted into a pNP concentration using a standard curve prepared at the given experimental conditions. The acquired data was analyzed using GraphPad Prism 6 (GraphPad Software). Linear regression of substrate concentration increase was determined over the whole time scale of the experiment. Conversion rates of each single concentration measurement including the catalytic peptide were subtracted by the corresponding background autohydrolysis of pNPA at the used concentration. Net conversion rates of each concentration measurement were then fitted using the Michaelis−Menten model implemented in the software. Standard errors were determined via Gaussian propagation of uncertainty.
■
←⎯⎯⎯
(3)
RESULTS AND DISCUSSION The active site of the Zn-finger protein was designed by Srivastava and Durani20 to mimic the active site of typical carbonic anhydrases in which the Zn2+ ion coordinates tetrahedrally to three histidine residues and the substrate. In our computational work, we first investigated the Zn2+ coordination architecture in the hydrolase, and thereafter carried out QM/MM calculations to explore the two reaction pathways considered. On the experimental side, we performed kinetic measurements to validate the findings from the computational studies. Zinc Coordination. During classical MD solvent equilibration, two water molecules from the bulk solvent entered the inner coordination sphere of Zn2+ to transform the tetrahedral
(2)
where k1, k−1, and k2 are rate constants, E is the enzyme, S is the substrate, ES is the enzyme−substrate complex, and P is the product. The initial rate of the reaction, ν0, was calculated from the slope of the linear fit of the time dependent absorbance data at 405 nm (ν0 = d[P]/dt, where [P] is the concentration of the product). As pNPA can autohydrolyze in buffer solution, we subtracted the rate of autohydrolysis from the initial reaction rates of the corresponding enzyme-catalyzed reaction to obtain the net reaction rate of the enzyme. The kinetic parameters KM = (k−1 + k2)/k1 (Michaelis constant) and kcat = k2 (rate constant of the catalysis or turnover number) were obtained by linear regression analysis of the Lineweaver−Burk plot: D
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B Zn2+ center to octahedral (Figure 2A). As discussed in the introduction, the hexacoordination of Zn2+ arises from the
Figure 3. Minimum energy paths of the reaction catalyzed by Znfinger hydrolase (pathway 1): Active-site geometries of reactant (A), intermediate (B), and product (C) obtained from QM/MM optimizations; potential energy profiles along the O(OH−)−C2(pNPA) (D) and C2(pNPA)−O2(pNPA) (E) reaction coordinates. Panels (D) and (E) correspond to the initial nucleophilic attack of the hydroxide on the pNPA carbonyl group and to the subsequent pNPA bond cleavage, respectively.
Figure 2. Coordination architecture of the Zn2+center: Water distribution in the active site obtained from QM/MM MD simulations (A); QM region considered in the QM/MM calculations (B); time evolution of the distances between Zn2+and the oxygen atoms of two initially Zn2+-bound water molecules (color coded in panel B) (C); time evolution of the distances between Zn2+ and the oxygen atoms of one water molecule and one hydroxide ion (replacing the water in blue in panel B) (D).
identical to that shown in Figure 2B, except that one of the water molecules was replaced by a hydroxide. The nucleophilic attack was driven by restraining the distance between the hydroxide oxygen atom and the C2 atom of pNPA (for atom labels, see Figure 2B). The potential energy profile for this scan (Figure 3D) indicates the formation of an intermediate, which is a true minimum on the potential energy surface. It corresponds to a structure in which the C2 atom of pNPA has a tetrahedral geometry with the C2(pNPA)−hydroxide bond fully formed and the C2(pNPA)-O2(pNPA) bond still intact (Figure 3B). To obtain the final product (Figure 3C), the pNPAC2−O2 bond of the fully optimized intermediate state was gradually increased until the bond is fully broken (Figure 3E). The geometric and energetic parameters of the optimized reactant, intermediate, and product states are provided in Table S1. The MEP implicates a sequential reaction mechanism passing through a tetrahedral intermediate in accordance with similar CA-like enzymes. Reaction Pathway 1: CA-like Mechanism, Minimum Free Energy Path (MFEP) Calculations. As discussed in Computational Methods, representative structures along the MEP were used to construct the initial pathway for the string simulations. The QM region considered in these simulations was the same as the one employed in the QM/MM MEP calculations. Seven active-site distances were chosen as reaction coordinates (see Figure 4A). As expected, only the distances that correspond to hydroxide nucleophilic attack (r1) and pNPA bond cleavage (r2) turned out to be relevant. The other five coordinates did not change significantly during the reaction (Figure S3). The 2D free energy surface of the catalytic reaction, obtained by projecting the seven-dimensional free energy surface calculated from the string simulations onto the space of reaction coordinates r1 and r2, is illustrated in Figure 4B. The solid black line on the free energy surface corresponds to the converged string or the minimum free energy path (MFEP) of the reaction. The 1D free energy profile along the
classical description of the Zn2+ ion and does not necessarily reflect its correct coordination state. To circumvent this issue we adopted the following strategy: After the classical MD solvent equilibration, the system was first subjected to 5 ps of QM/MM MD solvent equilibration. The QM region in these simulations comprised 52 atoms and consisted of the side chains of the three histidine residues (His10, His15, and His22), the side chain of a nearby serine residue (Ser6), the pNPA molecule, the Zn2+ ion, and the two water molecules coordinating to the Zn2+ ion (Figure 2B). During the QM/MM MD simulations we noted that, of the two water molecules coordinated to zinc at the MM level, the less tightly bound one drifts away from the inner coordination sphere of the metal ion (Figure 2C). The remaining Zn2+-bound water molecule was thereafter replaced by a hydroxide ion and the system was subjected to an additional 20 ps QM/MM MD simulation, during which the hydroxide ion remained strongly bound to Zn2+ and the other water molecule settled in the outer coordination shell of Zn2+(Figure 2D). We repeated this protocol for three different initial configurations and in each case observed a similar behavior (Figure S2). This finding strongly suggests that, in the αββ hydrolase, the Zn2+ ion prefers a pentacoordinated geometry, with a single water molecule or a hydroxide ion in its first coordination sphere. Reaction Pathway 1: CA-like Mechanism, Minimum Energy Path (MEP) Calculations. The reactant state was obtained from a QM/MM geometry optimization starting from a QM/MM MD solvent-equilibrated configuration in which Zn2+ is pentacoordinated. The first coordination sphere of Zn2+ thus consists of the three histidines, the substrate, and the hydroxide ion (Figure 3A). Starting from the optimized reactant state, the MEP of the reaction was computed by performing extensive QM/MM optimizations along suitable reaction coordinates. The QM region in these calculations was E
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B
4D: at the first stage of the reaction, the r1 coordinate decreases sharply and the r2 coordinate remains practically unchanged. In the intermediate region, with r1 ≈ 1.4 Å, the hydroxide is fully bound to C2(pNPA) and the latter assumes a tetrahedral geometry. Thereafter, C2−O2 bond cleavage occurs to form the final product: in the final stage of the MFEP, the r2 coordinate increases and the r1 coordinate remains unchanged. The free energy barriers for the first (OH− nucleophilic attack) and second (pNPA cleavage) steps of the reaction are 18.5 and 4.2 kcal/mol, respectively. As reported below, the experimentally determined rate constant of the Zn-finger hydrolase catalytic reaction is (8.5 ± 2.1) × 10−4 s−1. Using transition state theory, the reaction rate constant can be expressed in terms of the free energy barrier according to the equation k = A e−ΔG/RT, where A is the Arrhenius prefactor, R is the universal gas constant, and T is the absolute temperature. Assuming a prefactor of 1 ps−1, the experimentally measured rate constant corresponds to a free energy barrier of 20.5 ± 0.3 kcal/mol at 25 °C, in close agreement with the value obtained from the string simulations. Reaction Pathway 2: Reverse Protonation Mechanism. As discussed in the Introduction, for some zinccontaining metalloenzymes such as thermolysin and carboxypeptidase A, an alternative mechanism has been proposed to explain their experimentally observed kinetic behavior.7,8 In this “reverse protonation mechanism”, the Zn2+-bound hydroxide is coordinated too strongly to the metal ion to act as a nucleophile. Instead, an activated water molecule from the bulk attacks the substrate to affect the hydrolysis, and the only role of Zn2+ is to activate the carbonyl group of the substrate. We decided to assess the viability of this pathway in the minimal Zn-finger hydrolase. For this purpose, we first scrutinized the classical MD simulations to check if bulk water molecules move into the vicinity of the pNPA carbonyl so that they can potentially act as
Figure 4. QM/MM finite temperature string simulations indicate a sequential mechanism for Zn-finger catalyzed hydrolysis (reaction pathway 1). Active-site distances considered as reaction coordinates (A). 2D free energy surface of the catalytic reaction obtained by projecting the seven-dimensional free energy calculated from the simulations onto the space of reaction coordinates r1 and r2 (B). The inset zooms into the region of the tetrahedral intermediate (B). The dotted and solid black lines correspond to the initial and final (converged) strings, respectively (B). The converged string represents the minimum free energy path (MFEP) of the reaction (B). 1D free energy profile along the MFEP (C). Evolution of coordinates r1 and r2 along the MFEP (D).
MFEP is shown in Figure 4C, and the evolution of the coordinates r1 and r2 is shown in Figure 4D. The free energy profiles are consistent with the MEP calculations. They also indicate that the rate-limiting step of pathway 1 is the nucleophilic attack of the hydroxide to pNPA and that the pathway is sequential, passing through a shallow tetrahedral intermediate (Figure 4B−D). Along the MFEP, the hydroxide initially attacks pNPA. This can be seen in Figure
Figure 5. Water isodensity plots obtained from five independent MD trajectories of 25 ns each (A−E, dark cyan−top 5%, light cyan−top 10%); spontaneous pNPA cleavage observed in QM/MM MD simulations (F). These latter simulations were started from configurations with a water molecule close to C2(pPNPA) that was replaced by hydroxide (see text). The solid and dashed lines refer to the O(OH−)−C2(pNPA) and C2(pNPA)−O2(pNPA) distances, respectively, while the red, blue, and green colors correspond to independent simulations with different initial configurations. F
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B
evaluation of the data using the Michaelis−Menten equation compared to the double reciprocal plotting in the Lineweaver− Burke analysis.44 Catalytic Activity under High Hydrostatic Pressure. The pressure effect on the enzymatic activity was investigated in the range of 0.1−200 MPa. Figure 7 shows the pressure
nucleophiles upon activation. The distributions of water around the pNPA carbonyl in the various MD simulations are illustrated in Figure 5A−E; evidently, water molecules seldom come close to C2(pNPA). Nonetheless, to study the most favorable mechanistic scenario, we conducted the following computational experiment. We chose four configurations in which a water molecule approached the C2(pNPA) atom to within 3.5 Å, replaced this water molecule with a hydroxide ion, and initiated independent QM/MM MD simulations. The O(OH−)−C2(pNPA) and C2(pNPA)−O2(pPNPA) distances along the various QM/MM MD trajectories are depicted in Figure 5F. Obviously, we observe a spontaneous hydrolytic reaction in each trajectory, suggesting that this mode of nucleophilic attack is extremely facile, mainly because the manually inserted hydroxide anion is a powerful nucleophile, and also because the Zn2+ ion polarizes the pNPA carbonyl group and makes it more susceptible toward nucleophilic attack. Hence, the alternative reverse protonation pathway would become viable if a bulk water molecule could come close to the substrate and be activated enough to behave like a hydroxide anion. The currently studied αββ hydrolase does not possess a suitably positioned residue that would act as a general base and assist in such activation of a water molecule: the only possible candidate, a serine residue near the substrate, does not act in this manner and remains a spectator in all QM/MM calculations. Experimental Studies of the Enzyme under Ambient Conditions. The catalytic activity of the enzyme was monitored through the release of pNP. At 2.66% (v/v) acetonitrile, the turnover number kcat was determined as (8.5 ± 2.1) × 10−4 s−1, the Michaelis constant KM as (5.8 ± 2.1) mM and the catalytic efficiency kcat/KM as (0.148 ± 0.064)/M·s (see Figure 6 and Table S2). As discussed above, the kcat value is in
Figure 7. Pressure dependence of the kinetic constants for the hydrolysis of pNPA: turnover number, kcat (blue); Michaelis constant, KM (red); and catalytic efficiency, kcat/KM (black).
Figure 6. Michaelis−Menten kinetics of the Zn-finger hydrolase (20 μM) including 20 μM Zn(ClO4)2 and 2.66% MeCN in 20 mMTrisHCl buffer. Kinetics were derived from four independent measurements at each concentration. The data represents the average of each background corrected concentration. Error bars correspond to standard deviations.
dependence of kcat, KM, and kcat/KM, for the hydrolysis of pNPA catalyzed by the enzyme. We find that increasing the pressure does not affect KM, while it increases kcat as well as kcat/KM. The deduced value for the activation volume ΔV# is negative, indicating that the turnover number of the hydrolysis reaction is enhanced upon compression. A negative ΔV# value implies a smaller compression of the enzyme−substrate complex (ES) relative to the transition state (ES#). We find ΔV# to be −18 ± 2 cm3 mol−1, which is of the order of the volume of a single water molecule (18 cm3 mol−1). This value is similar to ΔV# of amyloid fibril esterase enzymes with the same substrate42 and of α-chymotrypsin in the hydrolysis of N-succinyl-L-phenylalanine-p-nitroanilide.45 Our results show that pressure effectively enhances the enzymatic activity of the de novo designed enzyme. These experimental findings are in line with our computational results, which indicate that the transition state for reaction pathway 1 is more compact than the reactant state. We performed a crude estimate of the change of volume of the catalytic center during the reaction (see the Supporting Information), assuming that the catalytic center consists of the substrate, the Zn2+ ion, and either a Zn2+ bound water (in the reactant state) or hydroxide ion (in the transition state). On the basis of the QM/MM optimized geometries of the reactant and transition, we find a negative activation volume ΔV# of −39 Å3 (corresponding to −23 cm3 mol−1), in qualitative agreement with experiment, which can be taken as further support for the proposed carbonic anhydrase-like mechanism.
excellent agreement with the computed barrier, which lends support to the proposed mechanism. Compared to previous literature data (variant B2),20 our measured kcat and KM values are larger by about one and 2 orders of magnitude, while the measured catalytic efficiency is smaller by about 1 order of magnitude. These differences could result from different substrate concentration ranges used to determine the kinetic constants. While the previous experiments covered a narrow concentration range (10−100 μM) at undefined acetonitrile portion, our current experiments yield reliable constants over a broad range of concentration (0.2−3.2 mM). This allows direct
CONCLUSIONS We conducted computational and experimental studies to investigate catalysis in an in silico designed Zn-finger protein, a minimal system capable of carrying out CA-like hydrolysis. Using classical and QM/MM MD simulations, QM/MM geometry optimizations, and advanced QM/MM free energy simulations, we established that the catalytic Zn2+ center in the Zn-finger hydrolase prefers a pentacoordinated geometry, like the Zn2+ centers found in naturally occurring CAs and CA-like enzymes. We investigated two pathways of substrate hydrolysis. In analogy to the commonly accepted mechanism for CA-like
■
G
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
The Journal of Physical Chemistry B
■
enzymes, we found a sequential pathway, in which a Zn2+bound hydroxide acts as a nucleophile and the hydrolysis proceeds through a tetrahedral intermediate. According to the QM/MM free energy simulations, the rate-limiting step of this reaction is the initial nucleophilic attack of the hydroxide on the substrate, which has a free energy barrier of 18.5 kcal/mol. These computational results are fully consistent with kinetic measurements under ambient and high-pressure conditions. We also explored an alternative reverse protonation mechanism in which a hydroxide (activated water molecule) from the bulk acts as the nucleophile (instead of the Zn2+-bound hydroxide). When such a hydroxide from the bulk is positioned favorably near the substrate carbonyl that is polarized by the Zn2+ ion, hydrolysis occurs spontaneously in our QM/MM MD simulations. However, this pathway does not appear to be viable in the current enzyme architecture, as the tertiary structure of the hydrolase is rather flexible and lacks a suitably positioned residue that could act as a general base and activate a bulk water molecule. Hence, it is extremely unlikely to find a hydroxide from the bulk near the substrate carbonyl in our system. It is conceivable, however, that the enzymatic activity could be increased by targeted mutagenesis of residues in the vicinity of the substrate carbonyl to favor the activation of a noncoordinated water molecule and thereby promote the reverse protonation pathway. Furthermore, the experimental results show that the catalytic efficiency can also be enhanced by increasing the hydrostatic pressure.
■
REFERENCES
(1) Lipscomb, W. N.; Sträter, N. Recent Advances in Zinc Enzymology. Chem. Rev. 1996, 96, 2375−2434. (2) McCall, K. A.; Huang, C.-c.; Fierke, C. A. Function and Mechanism of Zinc Metalloenzymes. J. Nutr. 2000, 130, 1437S− 1446S. (3) Vallee, B. L.; Auld, D. S. Zinc Coordination, Function, and Structure of Zinc Enzymes and other Proteins. Biochemistry 1990, 29, 5647−5659. (4) Zastrow, M. L.; Pecoraro, V. L. Designing Hydrolytic Zinc Metalloenzymes. Biochemistry 2014, 53, 957−978. (5) Bertini, I.; Luchinat, C. The Reaction Pathways of Zinc Enzymes and Related Biological Catalysts. In Bioinorganic Chemistry; Bertini, I., Gray, H. B., Lippard, S. J., Valentine, J. S., Eds.; University Science Books: Mill Valley, CA, 1994. (6) Xu, D.; Cui, Q.; Guo, H. Quantum Mechanical/Molecular Mechanical Studies of Zinc Hydrolases. Int. Rev. Phys. Chem. 2014, 33, 1−41. (7) Mock, W. L.; Stanford, D. J. Arazoformyl Dipeptide Substrates for Thermolysin. Confirmation of a Reverse Protonation Catalytic Mechanism. Biochemistry 1996, 35, 7369−7377. (8) Parkin, G. Synthetic Analogues Relevant to the Structure and Function of Zinc Enzymes. Chem. Rev. 2004, 104, 699−768. (9) DeGrado, W. F.; Summa, C. M.; Pavone, V.; Nastri, F.; Lombardi, A. De Novo Design and Structural Characterization of Proteins and Metalloproteins. Annu. Rev. Biochem. 1999, 68, 779−819. (10) Lu, Y.; Berry, S. M.; Pfister, T. D. Engineering Novel Metalloproteins: Design of Metal-Binding Sites into Native Protein Scaffolds. Chem. Rev. 2001, 101, 3047−3080. (11) Dhanasekaran, M.; Negi, S.; Sugiura, Y. Designer Zinc Finger Proteins: Tools for Creating Artificial DNA-Binding Functional Proteins. Acc. Chem. Res. 2006, 39, 45−52. (12) Nomura, A.; Sugiura, Y. Hydrolytic Reaction by Zinc Finger Mutant Peptides: Successful Redesign of Structural Zinc Sites into Catalytic Zinc Sites. Inorg. Chem. 2004, 43, 1708−1713. (13) Obst, S.; Bradaczek, H. Molecular Dynamics Simulations of Zinc Ions in Water Using CHARMM. J. Mol. Model. 1997, 3, 224−232. (14) Mujika, J. I.; Mulholland, A. J.; Harvey, J. N. Computational Methods: Modeling of Reactivity in Zn-Containing Enzymes; John Wiley & Sons, Ltd, 2006. (15) Estiu, G.; Suárez, D.; Merz, K. M. Quantum Mechanical and Molecular Dynamics Simulations of Ureases and Zn β-Lactamases. J. Comput. Chem. 2006, 27, 1240−1262. (16) Ackerman, S. H.; Gatti, D. L. Biapenem Inactivation by B2Metallo β-Lactamases: Energy Landscape of the Hydrolysis Reaction. PLoS One 2013, 8, e55136. (17) Dal Peraro, M.; Llarrull, L. I.; Rothlisberger, U.; Vila, A. J.; Carloni, P. Water-Assisted Reaction Mechanism of Monozinc βLactamases. J. Am. Chem. Soc. 2004, 126, 12661−12668. (18) Dal Peraro, M.; Vila, A. J.; Carloni, P.; Klein, M. L. Role of Zinc Content on the Catalytic Efficiency of B1Metallo β-Lactamases. J. Am. Chem. Soc. 2007, 129, 2808−2816. (19) Nechay, M. R.; Valdez, C. E.; Alexandrova, A. N. Computational Treatment of Metalloproteins. J. Phys. Chem. B 2015, 119, 5945−5956. (20) Srivastava, K. R.; Durani, S. Design of a Zinc-Finger Hydrolase with a Synthetic αββ Protein. PLoS One 2014, 9, e96234. (21) Zastrow, M. L.; Peacock, A. F. A.; Stuckey, J. A.; Pecoraro, V. L. Hydrolytic Catalysis and Structural Stabilization in a Designed Metalloprotein. Nat. Chem. 2012, 4, 118−123. (22) Laity, J. H.; Lee, B. M.; Wright, P. E. Zinc Finger Proteins: New Insights into Structural and Functional Diversity. Curr. Opin. Struct. Biol. 2001, 11, 39−46. (23) Morris, G. M.; Huey, R.; Lindstrom, W.; Sanner, M. F.; Belew, R. K.; Goodsell, D. S.; Olson, A. J. AutoDock4 and AutoDockTools4: Automated Docking with Selective Receptor Flexibility. J. Comput. Chem. 2009, 30, 2785−2791. (24) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926−935.
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.7b05027. Supplementary figures and tables with further computational results, Cartesian coordinates of representative QM regions, and supplementary experimental data (PDF)
■
Article
AUTHOR INFORMATION
Corresponding Authors
*E-mail:
[email protected]. *E-mail:
[email protected]. ORCID
Roland Winter: 0000-0002-3512-6928 Elsa Sanchez-Garcia: 0000-0002-9211-5803 Walter Thiel: 0000-0001-6780-0350 Author Contributions
The manuscript was written through contributions of all authors. All authors have given approval to the final version of the manuscript. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS We thank Dr. Kinshuk Raj Srivastava for providing the initial coordinates of the Zn-finger hydrolase. The authors acknowledge the German Research Foundation for support within the Cluster of Excellence RESOLV EXC1069 (S.E., F.S., R.W., E.S.G., W.T.) and the Collaborative Research Center CRC1093 (E.S.-G.). E.S.-G. acknowledges the support from the Boehringer Ingelheim Foundation (Plus-3 Program). H
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B
(45) Mozhaev, V. V.; Lange, R.; Kudryashova, E. V.; Balny, C. Application of High Hydrostatic Pressure for Increasing Activity and Stability of Enzymes. Biotechnol. Bioeng. 1996, 52, 320−331.
(25) Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H. J. C. Numerical Integration of the Cartesian Equations of Motion of a System with Constraints: Molecular Dynamics of N-Alkans. J. Comput. Phys. 1977, 23, 327−341. (26) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N· log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (27) Luo, G.; Andricioaei, I.; Xie, X. S.; Karplus, M. Dynamic Distance Disorder in Proteins Is Caused by Trapping. J. Phys. Chem. B 2006, 110, 9363−9367. (28) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kalé, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781− 1802. (29) MacKerell, A. D.; Bashford, D.; Bellott, M.; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; et al. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586−3616. (30) Mackerell, A. D.; Feig, M.; Brooks, C. L. Extending the Treatment of Backbone Energetics in Protein Force Fields: Limitations of Gas-Phase Quantum Mechanics in Reproducing Protein Conformational Distributions in Molecular Dynamics Simulations. J. Comput. Chem. 2004, 25, 1400−1415. (31) Van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: Fast, Flexible, and Free. J. Comput. Chem. 2005, 26, 1701−1718. (32) Vasilevskaya, T.; Khrenova, M. G.; Nemukhin, A. V.; Thiel, W. Mechanism of Proteolysis in Matrix Metalloproteinase-2 Revealed by QM/MM Modeling. J. Comput. Chem. 2015, 36, 1621−1630. (33) ChemShell, a Computational Chemistry Shell; http://www. chemshell.org/ (June 3, 2017). (34) Metz, S.; Kästner, J.; Sokol, A. A.; Keal, T. W.; Sherwood, P. ChemShellA Modular Software Package for QM/MM Simulations. WIREs Comput. Mol. Sci. 2014, 4, 101−110. (35) Ahlrichs, R.; Bär, M.; Häser, M.; Horn, H.; Kölmel, C. Electronic Structure Calculations on Workstation Computers: The Program System Turbomole. Chem. Phys. Lett. 1989, 162, 165−169. (36) Forrester, T. R.; Smith, W. DL-POLY Program; Daresbury Laboratory: Daresbury, Warrington, England, 1996. (37) Rosta, E.; Nowotny, M.; Yang, W.; Hummer, G. Catalytic Mechanism of RNA Backbone Cleavage by Ribonuclease H from QM/MM Simulations. J. Am. Chem. Soc. 2011, 133, 8934−8941. (38) Ganguly, A.; Thaplyal, P.; Rosta, E.; Bevilacqua, P. C.; HammesSchiffer, S. Quantum Mechanical/Molecular Mechanical Free Energy Simulations of the Self-Cleavage Reaction in the Hepatitis Delta Virus Ribozyme. J. Am. Chem. Soc. 2014, 136, 1483−1496. (39) Zhang, S. X.; Ganguly, A.; Goyal, P.; Bingaman, J. L.; Bevilacqua, P. C.; Hammes-Schiffer, S. Role of the Active Site Guanine in the glmS Ribozyme Self-Cleavage Mechanism: Quantum Mechanical/Molecular Mechanical Free Energy Simulations. J. Am. Chem. Soc. 2015, 137, 784−798. (40) Rosin, C.; Estel, K.; Hälker, J.; Winter, R. Combined Effects of Temperature, Pressure, and Co-Solvents on the Polymerization Kinetics of Actin. ChemPhysChem 2015, 16, 1379−1385. (41) Bugnon, P.; Laurenczy, G.; Ducommun, Y.; Sauvageat, P.-Y.; Merbach, A. E.; Ith, R.; Tschanz, R.; Doludda, M.; Bergbauer, R.; Grell, E. High-Pressure Stopped-Flow Spectrometer for Kinetic Studies of Fast Reactions by Absorbance and Fluorescence Detection. Anal. Chem. 1996, 68, 3045−3049. (42) Luong, T. Q.; Kapoor, S.; Winter, R. PressureA Gateway to Fundamental Insights into Protein Solvation, Dynamics, and Function. ChemPhysChem 2015, 16, 3555−3571. (43) Eisenmenger, M. J.; Reyes-De-Corcuera, J. I. High Pressure Enhancement of Enzymes: A Review. Enzyme Microb. Technol. 2009, 45, 331−347. (44) Dowd, J. E.; Riggs, D. S. A Comparison of Estimates of Michaelis-Menten Kinetic Constants from Various Linear Transformations. J. Biol. Chem. 1965, 240, 863−869. I
DOI: 10.1021/acs.jpcb.7b05027 J. Phys. Chem. B XXXX, XXX, XXX−XXX