Dynamic Connection between Enzymatic Catalysis and Collective Protein Motions.

Enzymes employ a wide range of protein motions to achieve efficient catalysis of chemical reactions. While the role of collective protein motions in substrate binding, product release, and regulation of enzymatic activity is generally understood, their roles in catalytic steps per se remain uncertain. Here, molecular dynamics simulations, enzyme kinetics, X-ray crystallography, and nuclear magnetic resonance spectroscopy are combined to elucidate the catalytic mechanism of adenylate kinase and to delineate the roles of catalytic residues in catalysis and the conformational change in the enzyme. This study reveals that the motions in the active site, which occur on a time scale of picoseconds to nanoseconds, link the catalytic reaction to the slow conformational dynamics of the enzyme by modulating the free energy landscapes of subdomain motions. In particular, substantial conformational rearrangement occurs in the active site following the catalytic reaction. This rearrangement not only affects the reaction barrier but also promotes a more open conformation of the enzyme after the reaction, which then results in an accelerated opening of the enzyme compared to that of the reactant state. The results illustrate a linkage between enzymatic catalysis and collective protein motions, whereby the disparate time scales between the two processes are bridged by a cascade of intermediate-scale motion of catalytic residues modulating the free energy landscapes of the catalytic and conformational change processes.

Protein motions involved in enzymatic catalysis span a wide range of spatial and temporal scales, from thermal fluctuations that equilibrate rapidly along the catalytic reaction coordinate on one hand to the slow collective motions that accompany enzymatic catalysis, ligand binding, and allosteric regulation[4−7] on the other hand. Thermal fluctuations typically take place locally at the active site of an enzyme in response to changes in reacting species and have been studied using various kinetic, spectroscopic, and simulation approaches.[9,10] While the role of collective protein motions in substrate binding, product release, and regulation of enzymatic activity is generally understood, their role in catalysis has been challenging to understand due to the difficulty in attributing a direct functional role to such motions or quantifying their effects.[11−14] For this reason, many studies have focused on identifying a phenomenological link between slow protein motions and enzymatic catalysis.[15−18] Furthermore, elucidating the detailed mechanism and sequence of events among slow conformational change, ligand binding and/or release, and catalysis remains a theoretical and experimental challenge. Overcoming this challenge with a combined computational and experimental approach, here, we examine conformational plasticity and dynamics in the active site of Escherichia coli adenylate kinase (AdK) and explore the mechanisms that link enzymatic catalysis to large-scale collective motion of the enzyme, achieved through the modulation of their respective free energy landscape.

In the cell, AdK regulates cellular adenosine triphosphophate (ATP)/adenosine diphosphate (ADP) homeostasis by catalyzing the reversible conversion of ATP and adenosine monophosphate (AMP) to two molecules of ADP.[19] X-ray crystallography and nuclear magnetic resonance (NMR) spectroscopy have established that during the catalytic cycle, the enzyme undergoes an open-to-closed conformational change of two substrate binding subdomains [ATPlid and AMPlid (Figure )]. A third subdomain, called CORE, is static but crucial for the thermal stability of the enzyme.[20] Subdomain opening occurs on the same time scale as the catalytic turnover[3] and is thus rate-limiting to the overall catalytic turnover. On this basis, extensive studies, both experimental[8,21−24] and theoretical,[25−27] were conducted to elucidate the mechanisms of the open-to-closed conformational change of the enzyme. However, despite various site-directed mutagenesis studies,[8,28−31] the catalytic mechanism of the enzyme remains not fully understood at the atomic level, including the roles of active site residues[32] and the metal ion in catalysis,[33] and the nature of the reaction transition state.[34] For example, Kerns et al. recently studied the role of the catalytic Mg2+ ion as an electrostatic pivot in phosphoryl transfer.[8] Nevertheless, the precise role of Mg2+ in the stabilization of the phosphoryl transfer transition state remains unknown, and its study requires a quantum mechanical description of the catalytic reaction to fully characterize the nature of the transition state (TS). In addition, the coupling between enzymatic catalysis and intrinsic protein motions, which occur at different time scales,[4,35,36] remains a fundamental and open question in enzymology.

Figure 1

Open-to-closed conformational change of E. coli adenylate kinase (AdK). (a) Ribbon diagram of E. coli AdK in the open (left, Protein Data Bank entry 4AKE(1)) and closed conformations (right, Protein Data Bank entry 7APU). (b) Architecture of the closed conformation active site with bound ATP and AMP. The five active site arginine residues, one lysine residue, and the bound ATP and AMP molecules are shown in a ball and stick model. In the figure, the two ADPs in the 7APU structure are modified to ATP and AMP.

In this study, we present a mechanistic connection between the large-scale opening motion of the enzyme adenylate kinase (i.e., opening of the ATP and AMPlids) and its catalytic phosphoryl transfer reaction, discovered by a combined application of NMR, X-ray crystallography, site-directed mutagenesis, and free energy simulations. This connection is achieved through a change in their respective free energy landscapes occurring in response to the change in the other. For example, the free energy barrier of the opening of the enzyme differed between the reactant and product states of the catalytic phosphoryl transfer reaction. Likewise, the free energy of the catalytic reaction was affected by the relative opening of the enzyme. In addition, arginine residues in the active site of the enzyme changed their orientation and interactions along the catalytic reaction. These alterations led to the change in the orientations of the two lids relative to the CORE subdomain, orienting them for faster opening in the product state than in the reactant state. In this way, the catalytic reaction in the active site of the enzyme affects the slow closed-to-open protein motion, without the events on two different time scales being directly coupled. Therefore, the coupling mechanism presented in this work is different from the direct coupling mechanism in which certain rapid protein motions, such as the rate-promoting vibrations proposed by Schwartz and co-workers,[9,37] occur concurrently with the crossing of the chemical barrier. In addition, this type of coupling mechanism is not exclusive to this enzyme but has been proposed in motor proteins, most notably in F1-ATPase, where the catalytic reaction modulates the free energy landscape of the rotation of its central rotor subunit, without the events taking place concurrently.[38] In the following, we present the Materials and Methods, Results, and Discussion.

Materials and Methods

Protein Expression and Purification

E. coli adenylate kinase (AdK) was overexpressed using a previously published protocol[39] based on a self-inducing plasmid (pEAK91). Plasmids carrying the R2A, R36A, R36K, R88A, R88K, R123A, R123K, R131A, R156K, R167A, and R167K mutations were acquired from Genescript, in which each mutation, except R156K, was designed to examine the impact of side chain removal and NMR assignment. The R156K mutant was designed on the basis of ref (8), and in this study, it was prepared to provide a consistent set of kinetic parameters. The expressed proteins were purified by affinity chromatography on a Blue sepharose column, followed by gel filtration. To obtain a uniformly 15N-labeled enzyme for the 15N-edited NMR experiments, expression cells were grown on minimal medium (M9) with 15NH4Cl as the sole nitrogen source.

X-ray Crystallography

AdK at a concentration of 18.3 mg/mL (∼776 μM) was mixed with 5 mM AMP and 5 mM guanosine triphosphophate (GTP) in a 30 mM MOPS buffer at pH 7 containing 50 mM NaCl. Diffraction quality crystals were grown by the hanging drop method in 24-well Linbro plates. Mixed on a coverslip were 2 μL of AdK, preincubated with AMP and GTP, which was a part of a screen of different combinations and concentrations of AMP, ADP, and ATP with GMP, GDP, and GTP, and 2 μL of the buffer containing 0.2 M NH4CH3CO2, buffered with 100 mM CH3CO2Na adjusted to pH 4.6 and 30% (v/v) PEG 4000 as the precipitant. The coverslip with the crystallization drop was inverted and sealed over a well containing 0.5 mL of the precipitant buffer. No Mg2+ ions were added at any stage, to slow catalytic reactions in the crystallization setup. Crystals grew at 291.2 K within 5 days to a size of 0.5 mm × 0.1 mm × 0.1 mm. Crystals were mounted free floating suspended in a nylon loop and vitrified by chock cooling to 100 K in a nitrogen gas stream (Oxford CryoSystems Ltd.). X-ray diffraction data were collected at a wavelength of 0.9750 Å onto a Pilatus 6M detector (Dectris) at beamline ID23-1 of the European Synchrotron Facility (ESRF) in Grenoble, France. The diffraction intensities were indexed, integrated with XDS,[40] and prepared for scaling with pointless.[41] Aimless[42] was used to scale and merge the intensities that were converted to structure factors with cTruncate.[43] Approximately 5% of the data was used for Rfree calculation.

Phase determination and structure solution were carried out by molecular replacement (MR) using E. coli AdK [Protein Data Bank (PDB) entry 1AKE(44)] as the search model in Phaser.[45] AdK structure refinement against data extending to 1.36 Å was carried out with phenix.refine[46,47] (Phenix program package version 1.17.1-3660). Manual model building was carried out with Coot (version 0.8.9.2 EL).[48] Refinement and model building were iterated until Rfree and Rfactor converged. The data processing and refinement statistics are listed in Table S1.

1Hε–15Nε Correlation and 15N Relaxation NMR Spectroscopy

Two-dimensional (2D) 1Hε–15Nε correlation spectra for side chain arginines and backbone relaxation spectra were recorded on a Bruker Avance III HD spectrometer at 850 MHz using a triple-resonance (TXI 5 mm) cryoprobe equipped with pulsed field gradients along the x, y, and z axes. Side chain arginine spectra were measured in a 30 mM MES, 50 mM NaCl buffer at pH 5.5. Backbone relaxation experiments were performed in a 30 mM MOPS, 50 mM NaCl buffer at pH 7.0. D2O [10% (v/v)] was added to all NMR samples for the field-frequency lock. NMR spectra were processed using NMRPipe and NMRDraw[49] and analyzed using Topspin 3.6 (Bruker) and SPARKY (version 3.113; https://www.cgl.ucsf.edu/home/sparky; University of California, San Francisco, CA). The side chain correlation spectra were recorded by placing the 15N carrier frequency at 90 ppm, and the sweep width was set to 86 ppm for HSQC measurements.

The 15N relaxation experiments, T1, T2, and 15N–{1H} heteronuclear nuclear Overhauser effects (NOEs), were performed in an interleaved manner. Interleaved 2D NMR experiments based on 1H–15N correlation spectra were performed. Delays of 50, 100 (duplicate), 200, 400, 500, 600, 800, 1000 (duplicate), 1200, and 1500 ms were used for T1 measurements, and delays of 16.96, 33.92, 50.88, 84.80, 101.76, 118.72, 135.68, 152.64, 169.60 (duplicate), 186.56, and 203.52 ms were used for T2 measurements.[50] Relaxation delays of 3 and 2 s were used for the T1 and T2 measurements, respectively. Steady-state 15N–{1H} heteronuclear NOE was measured with either 5 s delays between each free induction decay or 2 s delays, followed by a 3 s series of 120° nonselective 1H pulses.[51]T1, T2, and 15N–{1H} heteronuclear NOE experiments were performed with time-domain sizes of 256 × 2048 complex points and sweep widths of 10204.08 and 2584.62 Hz along the 1H and 15N dimensions, respectively, with eight scans for T1 and T2 and 32 scans for the 15N–{1H} heteronuclear NOE experiment. Arginine side chain 15N–{1H} heteronuclear NOE experiments were performed as described with the 15N carrier frequency and sweep width set at 86 and 60 ppm, respectively. The details of the analysis of 15N relxation data are provided in the Supporting Information. Assignments of the 1Hε–15Nε correlations were accomplished by comparing NMR spectra of apo and P1,P5-di(adenosine-5′) pentaphosphate (Ap5A)-bound states with those of several AdK variants (R2A, R36A, R36K, R88A, R88K, R119A, R123A, R123K, R131A, R156K, and R167K).

System Preparation and Molecular Dynamics (MD) Simulations

Three systems were prepared on the basis of the 1AKE wild-type (WT) E. coli AdK structure and the 2CDNMycobacterium tuberculosis AdK structure[33] for ADPs and Mg2+: (1) the reactant state (RS) with ATP and AMP, (2) the transition-state mimic (TSM), and (3) the product state (PS) with two ADPs. The ATP in RS and one of the two ADPs in PS were protonated. The TSM state was prepared to mimic the geometries of the transition state (TS) of the catalytic reaction (obtained from the quantum mechanical and molecular mechanical simulations), and one of the transferred phosphoryl oxygens was protonated. In all systems, we included the Mg2+ ion, its coordinating waters, and all crystal waters of 1AKE. CHARMM22[52] and the CMAP[53] potentials were used to represent the protein and ions. For ATP, ADP, and AMP, the CHARMM27[54,55] force field was used. After solvation with a 72 Å cubic box of TIP3P[56] waters and 150 mM NaCl, each system was minimized and equilibrated for 200 ps at 298 K. Production MD simulations were performed for 300 ns at 298 K and 1 bar with the leapfrog Verlet for 2 fs time integration and SHAKE[57] applied to all bonds involving hydrogen atoms. All MD simulations were performed using the CHARMM[58]/OpenMM[59] program. The particle mesh Ewald (PME) summation method[60] was used, and van der Waals interactions were evaluated up to 12 Å. All MD simulations were performed in triplicate to achieve convergence of simulation results. See the Supporting Information for details of system preparation, force field parameters for TSM, and MD simulations.

Hybrid Quantum Mechanical and Molecular Mechanical (QM/MM) Simulation

The system for the QM/MM simulation was prepared following the procedure described above for the product system (i.e., with two ADPs), except that it was solvated with a 69 Å rhombic dodecahedron (RHDO) box of TIP3P waters. After minimization and equilibration, the QM/MM MD simulation was performed at 298 K and 1 bar, in which the two ADP molecules, one Mg2+ ion, and its coordinated waters were described by the AM1/d-PhoT QM model,[61] which was developed for the catalytic phosphoryl transfer reactions, and the rest by the CHARMM force field. Long-range electrostatic interactions were evaluated with the QM/MM-PME algorithm[62,63] and van der Waals interactions with the 10.0 Å cutoff. All QM/MM MD simulations were performed using the CHARMM program (version c43a1),[58,64] with 0.5 fs time integration and SHAKE[57] applied to non-QM atoms. The system temperature and pressure were controlled using the Nose-Hoover thermostat[65] and the extended system pressure algorithm,[66] respectively.

Starting from the path generated by the umbrella sampling (US)[67] simulations, the string method in collective variables (SMCV) simulation[68,69] was performed using the distances of the forming and cleaving P–O bonds (i.e., P–ONuc and P–OLG) as CVs. The entire path was described using 48 SM replicas (i.e., MD images). For each string replica, 100 ps QM/MM MD simulation was performed to optimize the path, followed by 400 ps QM/MM MD to calculate the free energy along the path. Thus, a total of 24 ns QM/MM MD simulations were performed for the reaction.

In the simulation, we protonated one of the γ-phosphoryl oxygens of ATP. To determine the site of protonation, the transfer of a proton from ATP to AMP was simulated using the US technique. The simulation yielded 2.4 kcal/mol as the free energy of reaction (Figure S1), suggesting ATP as the preferred site for protonation for the wild-type enzyme. (This is different in the R167A mutant, in which the protonation of AMP is favored by 3.3 kcal/mol.) In addition, the phosphoryl transfer reaction with fully deprotonated ATP and AMP was simulated (using SMCV) to test the impact of ATP protonation on the catalytic barrier and reaction free energy. The barrier and the reaction free energy determined were 10.8 ± 0.2 and 0.4 ± 0.3 kcal/mol, respectively.

SMCV Free Energy (FE) Simulation of AdK Opening

Starting from the closed-to-open conformational change path obtained from the MD simulation in PS, SMCV path optimization was performed in an iterative manner. In each iteration, a 4 ps MD simulation was performed, during which each SM image was harmonically restrained with a force constant of 100.0 kcal mol–1 Å–2 to each CV and the average free energy gradients were evaluated for the evolution of the path. In the study presented here, the CORE–ATPlid, CORE–AMPlid, and ATPlid–AMPlid distances were used as CVs to describe the path (Figure S2) and the entire path was described with 72 SM images. After the evolution of the path, the CV positions were reparametrized to have an equal distance between neighboring images, completing one path optimization cycle. Path optimization was performed for a total of 2 ns for each image (i.e., 500 iterations). Then, additional 20 ns MD simulations were performed to determine the potential of mean force along the path. Thus, a total of 1.58 μs MD simulations were performed for the conformational change FE profile.

For the opening in the reactant state, two ADPs were modified to ATP and AMP for each SM image obtained from the product-state SMCV simulations and briefly energy minimized. The rest of the path optimization and the potential mean force simulation followed the same procedure that was used for the product-state SMCV simulation. To test the convergence of the results, additional SMCV FE simulations were performed from different initial paths with essentially identical results.

Results

Phosphoryl Transfer Mechanism and Impacts of Catalytic Residue Variants

The catalytic mechanism of AdK was studied by applying the hybrid QM/MM method,[64,70] using the AM1/d-PhoT QM method[61] for the catalytically reactive region and the CHARMM force field[52] for the rest of the enzyme and water solvent. The barrier and free energy of the phosphoryl transfer reaction between ATP + AMP and two ADPs were calculated to be 10.6 ± 0.2 and −1.6 ± 0.5 kcal/mol, respectively (Figure a). These values are comparable to the <12.4 kcal/mol catalytic barrier and the −1.4 kcal/mol reaction free energy estimated from kinetic measurements,[8] suggesting the reliability of the employed QM/MM method for the study of the enzyme’s catalytic mechanism.

Figure 2

Catalytic mechanism of AdK. (a) Free energy profile of the catalytic phosphoryl transfer between ATP + AMP and two ADPs, where α denotes the normalized reaction coordinate between 0 and 1. (b) Representative structures of the reactant state (RS), transition state (TS), and product state (PS), obtained from the QM/MM simulations. (c) Kinetic parameters (kcat and KM for ATP) determined for wild-type AdK and five active site arginine variants (Table S2). Error bars indicate the calculated errors of the kcat and KM values.

The reaction begins with the nucleophilic attack of an α-phosphoryl oxygen of AMP at the γ-phosphorus of ATP. This is followed by the formation of a transition state and then the reaction product (Figure b). Details of the catalytic mechanism and the character of the transition sate are provided in the Supporting Information and Figure S3. Importantly, the active site of E. coli AdK harbors five conserved arginine residues, i.e., R36, R88, R123, R156, and R167 (Figure b), which are well positioned for efficient phosphoryl transfer and, with a bound Mg2+ ion, balance the charges of ATP and AMP. In accordance with this, their mutations to alanine, and in one case to lysine,[28] substantially reduced the catalytic activity of the enzyme (Figure c and Table S2): R88A, R123A, R156K, and R167A exhibited <1% of the wild-type activity, and R36A 17% of the wild-type activity. On the contrary, the KM values were less affected by the mutations and the values determined for ATP ranged from 52 to 120 μM (Table S2). Previously, Wolf-Watz and co-workers established that the KM value for the wild-type AdK (71 μM for ATP) can serve as a valid proxy for binding affinity, based on its similarity to the 53 μM Kd value for ATP.[71] Assuming the same for the arginine mutants, the small changes in KM indicate that the thermodynamic contribution to binding exerted by the positive charge of the arginines is small. This was unexpected, because if the effects of the arginine mutation on kcat were primarily electrostatic in nature, equally large impacts on substrate binding would be expected as seen in the mutation of another catalytic residue, lysine 13, to glutamine (K13Q), which affected both kcat and KM.[28] Therefore, there are contributions of the enzyme, other than electrostatic, to the decrease in catalytic activity in the mutants.

Likewise, significant impacts on catalytic activity have been reported from the mutation of the active site arginines in other adenylate kinases.[8,29−31] For example, in human adenylate kinase 1 (AK1), the R44A and R97A mutations (equivalent to R36A and R88A of E. coli AdK, respectively) retained 21% and 19% of the wild-type activity, respectively, and the R132A, R138A, and R149A mutations (i.e., R123A, R156A, and R167A of E. coli AdK, respectively) <0.1% of the wild-type activity.[29] Similar results were also reported for chicken AK1.[30,31] However, due to the difference in the ATPlid length, their precise location and interactions are different between AK1 and E. coli AdK, and their mutations therefore have different impacts on the catalytic activity. To this end, Table S2 provides a consistent set of kinetic data on mutations for a single and widely investigated enzyme system, i.e., E. coli AdK.

Changes in the Conformation of Active Site Residues during the Catalytic Reaction

In Figure b, we show representative structures of the reactant state (RS), transition state (TS), and product state (PS) obtained from the QM/MM simulation. These structures illustrate the reorganization of the active site residues along the catalytic reaction trajectory (Movie S1 and shown schematically in Figure S4). On the ATP side, R123, which interacts with ATP in RS, follows the transfer of the γ-phosphoryl group to AMP. It therefore maintains its interaction with both ADP molecules in PS. A conserved lysine (i.e., K13) also interacts with ATP in RS via its β- and γ-phosphoryl groups. In PS, its interaction with the γ-phosphoryl group is lost, while the interaction with the β-phosphoryl group is retained. On the AMP side, R156 interacts with ATP and AMP via their γ- and α-phosphoryl groups, respectively, in RS. With the transfer of the γ-phosphoryl group, the α-phosphoryl group of AMP turns away from R156, but the residue retains the interaction with the transferred γ-phosphoryl group. R167 undergoes a change similar to that of R156. Namely, R167 interacts with both ATP and AMP in RS, and with the transfer of the γ-phosphoryl group, its interaction with AMP weakens. In PS, the residue disengages from the α-phosphoryl group of ADP and interacts with only the transferred phosphoryl group in the AMP binding site. In Table S3, we present the hydrogen bond distances formed between the different arginine residues and ATP and AMP in the RS, TS, and PS determined from the QM/MM simulations.

In addition, helix α6 changes its orientation relative to helix α7 during the catalytic reaction, in which the two helices anchor ATPlid to the reaction substrates (Figure b and Movie S1). This change in orientation can be attributed to the change in R123, R156, and R167 interactions during the reaction; note that R123 and R156 are located in the immediate vicinity of the two helices and R167 is located in the middle of helix α7. For example, helix α6 tilts its orientation toward helix α7 following the rearrangement of R123 with the transfer of the γ-phosphoryl group. Similarly, the orientation of helix α7 can be affected by R156 and R167. As a result, the interaction between the two helices becomes stronger, which then releases the loop connecting helices α3 and α4 of AMPlid from helix α7 (Movie S1). Finally and as discussed in Dual Role of the Active Site Arginines, these changes in orientations and interactions ultimately lead to the opening of ATP and AMPlids after the reaction.

Effects of the Conformational Plasticity of R167 on Catalysis

To understand the 137-fold reduction in kcat observed for the R167A variant (Table S2), the QM/MM and (100 ns) MD simulations were carried out for the R167A mutant. Interestingly, the mutation not only increased the reaction barrier by 4.0 kcal/mol (Figure S5a,b) but also enhanced the conformational motion of the two substrate binding subdomains compared to that of the wild-type enzyme (compare Figure S5c to Figure a). In particular, while the wild-type enzyme remained closed during the entire (300 ns) MD simulations, the R167A mutant exhibited an opening, if not entirely, even during the 100 ns MD simulations. For example, in RS, the opening of AMPlid was accompanied by the opening of ATPlid. In PS, AMPlid fluctuated between the closed and (partially) open states, while ATPlid opened slowly. In both states, helix α7, harboring R167, was displaced slightly from the loop between helices α3 and α4 of AMPlid. Consistent with the wild-type QM/MM result, the displacement of helix α7 in the R167A mutant can trigger the opening of AMPlid by releasing it from helix α7. This result suggests that in the wild type, R167 contributes to the deceleration of the opening of the enzyme through the stabilization of the helix α7 orientation, while it accelerates the catalytic reaction.

Figure 3

ATPlid dynamics along the catalytic reaction trajectory. (a) Distribution of CORE–ATPlid and CORE–AMPlid distances obtained from MD simulations between the reactant state (RS), transition-state mimic (TSM), and product state (PS). The color varies from blue to red with the increasing frequency of the observed distances. (b) Comparison of ATPlid and AMPlid orientations for X-ray structures in different catalytic states. The ATPlid is colored white for the structure with a bound ATP analogue and AMP (PDB entry 1ANK(2)), light blue for the transition-state mimic structure (i.e., ADP–AlF4––AMP, PDB entry 3SR0(8)), and blue for the two-ADP structure (PDB entry 4JKY(8)). In all structures, the CORE subdomain is colored white and 2ADP and Mg2+ are shown as only bound ligands for the sake of clarity. The thickness of the tube represents the temperature factor (B-factor) of the protein backbone Cα atoms. (c) Lowest-frequency protein motion (PC1) identified from principal component analysis (PCA), whose tube thickness and color represent the amplitude of the PC1 motion. (d) Normalized projection of MD simulation trajectories onto the PC1 and PC2 vectors.

Next, to verify the change in the R167 orientation in PS, we analyzed the X-ray crystal structure of E. coli AdK with two bound ADPs in the absence of Mg2+ determined at a resolution of 1.36 Å (Figure a, right, and Table S1). The structure is essentially identical to other closed AdK structures, for example, with a 0.2 Å root-mean-square deviation (RMSD) from the 1AKE structure.[44] In addition, despite the absence of the bound Mg2+ ion, the orientation of the active site ADPs and surrounding residues was very close to the structure reported with Mg2+ by Kerns et al.,[8] thus providing meaningful details of ADP interaction at both binding sites. In our structure, the AMP binding site is partially occupied by ADP (∼90%) and AMP (∼10%) on the basis of the occupancy of the ADP β-phosphoryl group, while the ATP binding site is fully occupied by ADP. Previously, ATP and AMP were shown to bind independently to AdK, with Kd values of 53 μM for ATP and 210 μM for AMP.[71] Here, binding of ADP to AdK is quantified via isothermal titration calorimetry (Figure S6 and Table S4). The data can be fitted well with a single Kd value of 9 μM, suggesting that the two ADP molecules bind with the same binding affinity for the two binding sites.[72] The 6- and 23-fold changes in the binding affinity for ADP compared to those of ATP and AMP, respectively, suggest some level of cooperativity between the two ADP binding events.

In our ADP-bound co-crystal structure, the electron density covering R167 suggests two conformations, designated here α and β (Figure ). This heterogeneity is also observed in other X-ray structures of AdK[2,8,44,73] (Figure S7). To test the effects of the different orientations of R167 on the catalytic reaction, we performed the QM/MM simulations while restraining the Cγ–Cδ–Nε–Cζ dihedral angle (i.e., δ angle) of R167 at −66° for the β conformation. The unrestrained QM/MM simulation was close to the α conformation with an average δ angle of −171° at RS. Restraining the R167 side chain orientation increased the catalytic barrier by 4.6 kcal/mol and the reaction free energy by 1.6 kcal/mol. This suggests the importance of R167 flexibility for the low catalytic reaction barrier, possibly by stabilizing the transition state more effectively than the reactant state. Moreover, the MD simulations, which were performed in the transition-state mimic (TSM) and product states (PS) to access the behavior of the active site residues over a longer period of time, showed an enhanced fluctuation of the δ angle toward the β conformation (Figure S8). Likewise, R156 exhibited conformational heterogeneity based on the superposition of multiple X-ray structures.[8] This finding, with the results presented here, suggests the importance of conformational heterogeneity of the two arginine residues (i.e., R156 and R167) in the catalytic turnover of the enzyme.

Figure 4

Different orientations of the R167 side chain in the X-ray structure of AdK in complex with two ADPs (PDB entry 7APU): (a) the α conformation and (b) the β conformation compatible with bound ADP in the AMP binding pocket. The color scheme of the protein cartoon representation is the same as in Figure .

Role of the Mg2+ Ion in Catalysis

The Mg2+ ion bound in the active site plays an important role in achieving high catalytic efficiency, as demonstrated by the 4.4 × 104-fold reduction in kcat in the absence of Mg2+ (from 330 ± 11 to 7.5 × 10–3 ± 7.5 × 10–4 s–1) as determined from a real-time 31P NMR assay[74,75] (Figure S9). Kerns et al. also reported a similar level of reduction in kcat with a removal of Mg2+.[8] To better understand the reduction in kcat, we carried out the QM/MM and MD simulations in the absence of Mg2+. The QM/MM simulation yielded a 17.4 ± 1.3 kcal/mol reaction barrier, which is slightly lower than the 18.74 kcal/mol barrier estimated using the forward rate constant of Kerns et al.[8] and 3.1 kcal/mol lower than that of 20.5 kcal/mol estimated on the basis of the kcat value determined in this work. Compared to that of the reaction in the presence of Mg2+, the barrier is higher by 6.8 kcal/mol. This is understandable because both β- and γ-phosphoryl oxygens coordinate with the Mg2+ ion in the reactant sate (Figure c), and therefore, the ion can stabilize the negative charge developed in the leaving group oxygen at the TS.

In addition, the MD simulations performed in the absence of Mg2+ showed more suppressed ATP and AMPlid motions and thus a more closed conformation (Figure S10c). Consistent with this, NMR relaxation dispersion and pre-steady-state kinetic experiments have shown much slower lid opening in the absence of Mg2+ (by at least 3800-fold).[8] Together, this study and the study by Kerns et al. suggest that Mg2+ accelerates both the lid opening and the catalytic reaction. This can be compared to the impact of R167, which slows the opening of the enzyme, while it accelerates the catalytic reaction.

NMR Characterization of the Structural Heterogeneity of the AdK Active Site

To probe the structural heterogeneity of the active site arginine residues experimentally, we turned to the 1H–15N correlation NMR experiment. We have previously shown that backbone 1H–15N heteronuclear single-quantum coherence (HSQC) experiments are very similar between Ap5A- and ADP-saturated states, demonstrating that the overall closed conformation is the main conformation in solution for both ligands.[71] In this study, the experiment was used for the Hε proton at the Nε atom of the arginine side chain (Figure a). In AdK, five of the 13 arginines are in the active site, and assignment of the catalytic and a few additional arginines was carried out by analyzing the NMR spectra of AdK variants in which the arginines were replaced with alanine and/or lysine (Figures S11 and S12). With this approach, we could assign all catalytic arginines with the exception of R156 in the apo state and R167 in the apo and Ap5A-bound states.

Figure 5

NMR characterization of the active site plasticity of E. coli AdK. (a) Arginine side chain and notation of its heavy atoms. (b) Overlay of arginine side chain 1Hε–15Nε HSQC NMR spectra of apo (black) and Ap5A-bound (red) states of AdK showing the assignment of the catalytic arginine residues. (c) Overlay of arginine side chain 1Hε–15Nε HSQC NMR spectra of apo (black), Ap5A-bound (red), and ADP-bound (green) states of AdK. For panel c, each of the spectra is shown in Figure S13a–c. (d and e) Plots of the 15N–{1H}-heteronuclear NOEs of catalytic arginine backbone amide (maroon) and arginine side chain 15Nε–1Hε (blue) groups of (d) apo and (e) Ap5A-bound AdK, respectively. The heteronuclear NOE (Isat/Iunsat ratio of peak intensities) of each catalytic residue is plotted on the y-axis using the same values (−1.0 to 1.0) to compare the apo and Ap5A-bound AdK states. Errors associated with each NOE value from triplicate measurements are shown atop the bar graphs.

In the apo state, the catalytic arginines exhibit a relatively narrow chemical shift dispersion (black in Figure b), indicating that they are predominantly solvent exposed. Then, in the presence of the inhibitor Ap5A, which is known to stabilize the closed conformation of the enzyme,[44] the chemical shift dispersion of the arginine side chains has increased considerably (Figure b and Figure S13a). Analysis of the heteronuclear steady-state NOEs for the side chains shows that molecular motion on the picosecond to nanosecond time scale[50] is significantly restricted in the Ap5A-bound state relative to the apo state (blue bars in Figure d,e). In the apo state, the NOE values range from −0.63 for R36 to 0.31 for R123, while the assigned side chains have an average NOE of 0.70 in the Ap5A-bound state. This restriction of picosecond to nanosecond side chain motions is in line with the chemical shift perturbation (CSP) data and together suggests that the arginine side chains have well-defined conformations in the Ap5A complex. This is also consistent with a related analysis of a uridine monophosphate (UMP)/cytidine monophosphate (CMP) kinase.[32]

To benchmark the changes in the NOE values, we quantified the corresponding values of the arginine backbone 1H–15N correlations (maroon bars in panels d and e of Figure ). The average values for the apo and Ap5A-bound states were found to be 0.77 and 0.85, respectively. The most dramatic difference is seen in the apo state, where the NOE values indicate that the arginine side chains are dynamically disordered, while their backbones are well ordered. This observation is consistent with the idea that in the apo state, the side chains of all catalytic residues are projected into solvent and that the molecular motion on the picosecond to nanosecond time scale is, therefore, almost unrestricted. The situation is different in the Ap5A-bound state; i.e., the average NOE value for the side chains is well above 0.5, indicative of their order. However, it is still below that of the backbone atoms in the apo state, suggesting that although the side chains are relatively well ordered in the active site in complex with Ap5A, some additional motions on the picosecond to nanosecond time scale are present relative to the backbone. It is possible that these “additional” motions are an intrinsic property of the enzyme and contribute to its catalytic function as probed by QM/MM simulations in this work.

In contrast to the well-defined spectrum of AdK in complex with Ap5A, the chemical shift dispersion in the ADP-bound state has collapsed relative to the Ap5A state (green in Figure c and Figure S13d). However, we were unable to assign the side chain resonances in complex with ADP, and the same difficulty was previously reported for the UMP/CMP kinase.[32] Nevertheless, the narrow chemical shift dispersion is consistent with a significant dynamic disorder of the arginine side chains in the ADP complex. Although these NMR data in the absence of quantitative analysis cannot provide a molecular-level model for this dynamic disorder, it is still worth correlating this observation with the relatively insignificant changes in the KM values for the arginine variants in comparison with the wild type (Figure c and Table S2). Taken together, it is possible to envision a scenario in which binding of substrates to AdK is multivalent in nature, allowing the active site arginines to interact with the phosphates (of ATP and AMP vs two ADPs) with different binding poses. More specifically, because the catalytic reaction proceeds naturally with ADP during NMR experiments, the observed dynamic disorder may be due to the change in the arginine orientations following the catalytic reaction; that is, they relax to a new orientation on the picosecond to nanosecond time scale after the phosphoryl transfer reaction. However, more work is required to solidify such a model.

Collective Protein Motions in Catalysis

MD simulations performed at different reaction states (i.e., RS, TSM, and PS) showed that while the enzyme stayed in the closed state in all simulations, the product-state simulation exhibited a pronounced ATPlid fluctuation toward a relatively more open conformation than the reactant state (Figure a). Inspection of X-ray structures of AdK in complex with different ligands revealed a similar variation in ATPlid conformation and dynamics. For example, the structure of AdK with an ATP analogue and AMP (PDB entry 1ANK(2)) is the most closed (white in Figure b), followed by the transition-state analogue (TSA) structure with ADP, AlF4–, and AMP (PDB entry 3SR0,[8] light blue), and the structure with two ADPs is the most open (PDB entry 4JKY,[8] blue). In addition, B-factors are largest for amino acids around the ATPlid in the co-crystal structure with two ADPs, as visualized by the thickness of the tube in Figure b, which indicates a large fluctuation of ATPlid residues.

To characterize the dynamics of the ATPlid along the catalytic reaction coordinate, principal component analysis (PCA)[76] was performed using the coordinates saved in the QM/MM simulations. In agreement with Figure b, the two lowest-frequency motions identified from the analysis, PC1 and PC2, occurred along the opening direction of the enzyme (Figure c and Figure S14). In addition, the projection of the MD simulation trajectories onto the two PCA modes showed a wider conformational space exploration in TSM and PS than in RS (Figure d). Interestingly, the projection of PS (and TSM) shows a bimodal distribution, in which the larger peak on the left corresponds to the more open ATPlid and the smaller one to the right to the more closed ATPlid conformation (Figure a,d). The QM/MM simulations performed with a restraint on the distance between the ATPlid and CORE subdomains further showed a pseudolinear dependence of the reaction free energy on the CORE–ATPlid distance, while the barrier was relatively insensitive (Figure S15). These results suggest that ATPlid dynamics differs between the reactant and product states. This difference in dynamics results in the more open conformation in PS (Figure a), as a consequence of the formation of different interactions of the arginine residues in the active site, and at the same time influences the free energy landscape of the catalytic reaction.

The results presented above in turn suggest the possibility that the rearrangements of the active site residues and the resultant changes in their interactions trigger the facilitated opening of the enzyme in PS. To test this possibility, the MD simulations performed in RS and PS were extended to 1 μs. In the two MD simulations in PS, a spontaneous opening occurred within 1 μs (Movie S2), while no opening was observed in the RS MD simulations. Figure a shows the distributions of the CORE–ATPlid and CORE–AMPlid distances in different colors between the different time segments. In addition, we performed the string method free energy (FE) simulations[69] based on the conformational change path obtained in PS and determined the free energy change along the path. In agreement with the MD simulations, the FE profiles show that the opening of the enzyme has a lower FE barrier in the product state than that in the reactant state [7.4 ± 1.2 kcal/mol in PS vs 13.3 ± 1.9 kcal/mol in RS (Figure b)]. This confirms that the catalytic reaction facilitates the opening of the enzyme in the product state. This is achieved by the rearrangements of the active site residues and the resultant changes in their interactions with the reaction substrates, affecting the free energy landscape of the closed-to-open conformational transition. However, we must be careful when comparing the free energy of opening shown in Figure b to that determined by single-molecule studies.[22,24] Typically, in single-molecule studies, the free energy difference between two conformational states is estimated on the basis of their relative populations, while the FE profile shown in Figure b has a contribution of the ligand–enzyme interactions. This difference can be better understood when the free energy value based on the ligand dissociation constant (Kd) is compared with the value based on relative populations. In particular, due to the difference in the definition of the two equilibrium constants, the resulting free energy values differ in the contribution of ligand concentration. For the same reason, the population-based free energy value depends on the ligand concentration and the ligand binding affinity for the enzyme.

Figure 6

Closed-to-open conformational change. (a) Changes in CORE–ATPlid and CORE–AMPlid distances in the two 1 μs MD simulations performed in the product state. The two distances are shown in different colors with an interval of 200 ns. (b) Free energy (FE) profile of the closed-to-open conformational change described by the normalized reaction coordinate α between 0 (closed state) and 1 (open state). The FE profile shown in black is for opening in the reactant state (ATP and AMP) and in red for that in the product state (two ADPs).

Discussion

Dual Role of the Active Site Arginines

This study reveals how conformational heterogeneity around the five arginine residues in the active site of the enzyme adenylate kinase (AdK) facilitates the catalytic reaction. The arginine side chains reorganize rapidly along the catalytic reaction coordinate and, with their ability to form multiple hydrogen bonds, effectively stabilize the transition state of the reaction. For this reason, these residues are highly conserved among NMP kinases, including adenylate kinases,[77] and any mutation thereof causes a substantial reduction in catalytic activity (Table S2). Importantly, AdK’s R123 occupies a position very similar to that of the so-called “arginine finger”, R789GAP, of the Ras–RasGAP complex[78] (Figure S16). On the basis of this finding and also given that AdK inserts its catalytic arginine residues (in particular, R123, R156, and R167) into the active site through conformational change, these active site arginine residues are reminiscent of the arginine finger of RasGAP[78,79] and other ATP/GTPases.[80−83]

This study also suggests that the active site arginines control the closed-to-open conformational dynamics of AdK. As shown in Figure a, AdK explores a wider conformational space in the product state (PS) along the opening direction of the enzyme, which leads to a more favorable opening of the enzyme compared to the reactant state (RS). Our analysis of the QM/MM simulations further suggests that the enhanced motion of ATPlid in PS can be caused by the change in hydrogen bonding interactions of the arginine residues in the active site of the enzyme along the catalytic reaction trajectory (Figure b and Figure S4). For example, notable in Movie S1 and Figure b, which show the rearrangement of the R123, R156, and R167 residues during the catalytic reaction, is the change in the relative orientation of helices α6 and α7. Due to the proximity of the three arginine residues to the two α-helices, in the case of R167 in the middle of helix α7, their rearrangements during the reaction lead to the change in the orientation of the two helices. In PS, in particular, these changes trigger the release of the loop connecting helices α3 and α4 of AMPlid relative to helix α7. This release then decouples AMPlid from ATPlid and finally allows the opening of the enzyme. This sequence of events is consistent with that observed from the MD simulations in PS; i.e., a slight opening of AMPlid occurs first, followed by the opening of ATPlid (Figure a and Movie S2). If this is the case, the removal of the side chain of these residues, such as, by mutation, would result in an increased level of motion of ATP and AMPlids and eventual opening of the enzyme. This was indeed observed from the MD simulations of the R167A mutant in both RS and PS (Figure S5c). Combined with the impact on the phosphoryl transfer barrier (Figure S5a,b), this explains the >99% reduction in catalytic activity. We expect similar results with the R123A mutant, due to the orientation and interactions formed by the R123 side chain in both RS and PS and the impact of their change on the orientation of helix α6 (Figure b).

Finally, we note that the change in the time scale of conformational motion in response to the catalytic reaction has been reported in other enzymes, such as dihydrofolate reductase (DHFR)[84] and triosephosphate isomerase (TIM).[85,86] In DHFR, in addition, a variation in the conformation of active site residues was observed across the catalytic cycle of the enzyme,[87] suggesting the possible generality of the mechanism presented in this study. In these enzymes, however, the region of conformation change is rather localized around the active site of the enzyme compared to that of adenylate kinase.

Implications for the Rate-Limiting Step

The difference in the opening barrier between RS and PS has an important implication for the mechanism of the enzyme. That is, the slowest step of the reaction can be different between the forward and backward reactions (Figure ). In the forward reaction (i.e., ATP + AMP → 2ADP), the chemical step takes 9.47 μs based on the transition-state theory (TST) with a ΔF⧧ of 10.6 kcal/mol (and <0.2 ms based on the experimental ΔF⧧ of <12.4 kcal/mol[8]), while the opening step takes 3.0 μs based on the mean first passage time (MFPT) calculation.[88] Therefore, the slowest step in the forward direction is the chemical step. In contrast, it is the opening of the enzyme in the reverse reaction (450 ms MFPT). In this comparison, however, the following is worth noting. First, we used TST to estimate the time of the catalytic reaction (i.e., the rate of the catalytic reaction). In TST, because the reaction is assumed to take place only when the reactants have reached the so-called high-energy transition state, the estimated reaction time includes the time that the enzyme spends waiting for the transition-state conformation, which is the majority, while the actual reaction occurs on the time scale of the scissile/forming bond vibration at TS. Second, the MFPT rate in RS is 2 orders of magnitude larger than the opening rate measured by NMR (i.e., 3.5 ms).[3] This severe overestimation may be due to errors in the free energy profiles (Figure b) and/or limitations of the MFPT method and approximation.[88] For this reason, the precise determination of the rate of slow conformational change has remained a theoretical challenge. Nevertheless, the MFPT-based opening time in PS is comparable to the opening time scale (0.5–1.0 μs) directly observed from the MD simulations (Figure a).

Figure 7

Free energy landscape of AdK along the entire catalytic cycle and schematic representation of corresponding states. Indicated in the figure is the time estimated to cross each barrier: the mean first passage time in the conformational change and the time based on the transition-state theory in the chemical reaction step. For the opening in the reactant sate (i.e., EcRS → EoRS), the time determined on the basis of the NMR experiment[3] is also provided.

At first glance, the difference in the slowest steps seems incompatible with the rate-limiting conformational change mechanism proposed on the basis of the NMR experiments.[3] The NMR technique applied, however, probed only millisecond protein motions and thus provided no information about protein motions on a faster time scale. Furthermore, recent single-molecule experiments[24] and simulations[89,90] have shown that the opening can occur in the range of 1.8–77 μs. In contrast, in other single-molecule study,[22] 6.3 ms opening was observed in the presence of a substrate analogue (i.e., AMP-PNP and AMP). While this result is consistent with the NMR result, it can be limited by the time resolution of the single-molecule experiment to detect a faster time-scale event. Taken together, the two vastly different opening time scales we discovered in this study are consistent with the previously published experimental results and allow reconciliation between the different sets of experiments.

Conclusion

In this study, we have applied a combination of multiscale quantum mechanical and classical mechanical simulations, enzyme kinetics, NMR, and X-ray crystallography to investigate how the enzyme adenylate kinase achieves a mechanistic coupling between enzymatic catalysis and microsecond to millisecond time-scale collective protein motions. To achieve the coupling, the enzyme employs picosecond to nanosecond side chain motions, which occur locally at the active site, to influence both the catalytic reaction and the slow conformational change of the enzyme. In this way, the two events are connected, thus overcoming the disparity in their time scales. In motor proteins, such as F1-ATPase, a similar mechanism is adopted to regulate the timing of the catalytic reaction and large-scale conformational changes.