A QM/MM study of the catalytic mechanism of α-1,4-glucan lyase from the red seaweed Gracilariopsis lemaneiformis

Abstract

α-1,4-Glucan lyase (GLase, EC 4.2.2.13), a unique glycoside hydrolase family member, specifically cleaves the α-1,4-glucosidic linkages in glycogen, starch and malto-oligosaccharides to produce 1,5-anhydro-D-fructose from the non-reducing end. Previous studies have proved that GLase belongs to the retaining glycoside lyase, and the catalytic reaction contains both the glycosylation and deglycosylation/elimination steps, in which a covalent glycosyl–enzyme intermediate is involved. On the basis of the newly reported crystal structure of GLase (2X2I) and the speculated mechanism, the whole catalytic cycle of GLase has been studied by using a QM/MM method. Calculation results indicate that the whole catalytic cycle contains five elementary steps. Firstly, the aspartic acid residue D665 acts as an acid to protonate the glycoside oxygen, which is concerted with the cleavage of the glycoside bond. Then, the residue D553 functions as the nucleophile to attack the anomeric carbon to form the glycosyl–enzyme intermediate. Different from the retaining α-glucosidases whose glycosylation is a typical concerted process, the glycosylation process of glycosidic lyase follows a stepwise mechanism. For the deglycosylation/elimination step, two cases with or without the maltotriose group in the active site were considered. The departure of the maltotriose can facilitate the proceeding of this process. The deprotonated aspartic acid residue D553 further acts as a catalytic base to abstract the C2-proton of the glucosyl residue. The proton abstraction in the deglycosylation/elimination step is calculated to be the rate-limiting step of the whole catalytic reaction, which corresponds to the energy barriers of 20.69 and 18.53 kcal mol⁻¹ for both of the two cases.

---

1. Introduction

Carbohydrates have numerous significant functions in organisms. Besides serving as structural components and providing the vital energy necessary for living bodies, carbohydrates also play vital roles in immune responses via trimmed peptidoglycans.^1–3 In general, the utilization of carbohydrates in organisms depends largely on the breakdown of the stable glycosidic bond. Glycoside hydrolases (GHs) are extremely common enzymes, which can specifically catalyze the hydrolysis of glycosidic bonds such as O-, N- and S-linked glycosides with high efficiency, and many GHs can enhance the catalytic rate on the order of 10¹⁷-fold than the nonenzymatic (natural) reaction.⁴ GHs have been detected in nearly all domains of life, from bacteria and fungi to plants and mammals. These GHs so far discovered have been classified into more than 130 glycoside hydrolase families (http://www.cazy.org) based on their conserved amino acid sequence motifs.⁵ According to the stereochemical difference at the anomeric center of the catalytic reaction, these enzymes can be classified into two types: inverting and retaining.⁶ Both mechanisms require the participation of oxocarbenium ion-like transition states and a pair of enzymatic residues, typically carboxylate residues at the active sites.⁷ In inverting enzymes, they employ a single-displacement mechanism, and two residues, with a distance of about 10.5 Å to allow the substrate and a water nucleophile gripped between them, act as acid and base respectively.^7–9 Whereas in retaining enzymes, they utilize a double-displacement mechanism via a covalent glycosyl–enzyme intermediate, in which one carboxylate residue functions as the catalytic acid/base and another one in a distance of about 5.5 Å acts as the nucleophile.⁷ The formation of the covalent glycosyl–enzyme intermediate has been verified by X-ray crystal structures and theoretical calculations.^10–14 Generally most of the enzymes require the participation of an extra water molecule, but some non-hydrolytic enzymes also exist, such as the α-1,4-glucan lyases, which cleave the glycosidic bond with no hydration to form anhydrofructoses in their enol forms.

The α-1,4-glucan lyases (EC 4.2.2.13, GLases) was early purified and characterized from red seaweeds, and now the activity of lyases has been widely founded in eukaryotes.^15,16 Owing to the similarities in cleaving C1–O glycosidic bonds and sequential structure, GLases have been classified into the family 31 of glucan hydrolases. Since no water molecule takes part in the reaction and the product 1,5-anhydrofructose is clearly different from those of other glycosidases, GLases form the subgroup 2 of GH31.¹⁶

GLases are retaining glycoside hydrolases,¹⁷ which catalyze the degradation of starch, maltose and glycogen via a nonhydrolytic pathway to produce 1,5-anhydro-D-fructose (1,5AnFru) from the non-reducing end, and release a free glucose, as shown in Scheme 1. Previous experimental studies by kinetic isotope effects have suggested that the GLases-catalyzed reaction follows a two-step mechanism which involves both the glycosylation and deglycosylation steps.^17–19 In the glycosylation step, a glycosyl–enzyme intermediate (the glycosylated GLase) is formed via a transition state with significant oxocarbenium ion characteristic, which is identical to the glycosylation step of the retaining glycoside hydrolases of this family (Scheme 2A).^17–20 The formation of the covalent glycosyl–enzyme intermediate has been observed by trapping with 5-fluoro-β-L-idopyranosyl fluoride, and subsequent detection and sequencing of the labeled peptide by mass spectrometry.²¹ But for the deglycosylation step, it was proposed that the enzyme–substrate intermediate undergoes a syn-elimination reaction of an adjacent proton and the nucleophile to product the un-saturated product 1,5AnFru,¹⁹ as shown in Scheme 2B. Since no water molecule participates in the deglycosylation step, the mechanism of GLases is obviously different from that of α-glycosidases (Scheme 2B′), which release a glucose as the final product.

Scheme 1 Catalytic reaction of glucan lyase.

Scheme 2 Comparison of the catalytic mechanisms of alpha-1,4-glucan lyase (A and B) and alpha-glycosidase (A and B′).

In 2013, H. J. Rozeboom et al. resolved the crystal structure of α-1,4-glucan lyase from the Red Seaweed Gracilariopsis lemaneiformis, in its native form and in complex with four inhibitors acarbose (Acr), 1-dexoynojirimycin, 5FβIdoF and 5FαGlcF.²⁰ Based on the analyses of active site structures, they further proposed that the catalytic nucleophile (D553) which forms the covalent enzyme–glycosyl intermediate in the glycosylation step also acts as the catalytic base in the elimination step. Their results resolved the long-time ambiguity on the nature of the proton abstracting residue.

Although useful information about the mechanism of α-1,4-glucan lyases has been obtained, open questions still remain unresolved. For example, by using substrates p-nitrophenyl α-D-glucopyranoside (PNPαGlc) and α-D-glucopyranosyl fluoride (αGlcF), it was concluded that the deglycosylation step can not be rate-limiting.^17–19 But if the substrate 5-fluoro-α-D-glucosyl fluoride was used, the deglycosylation step was suggested to be rate-limiting.^17,18 As for the deglycosylation/elimination step, in 2003 Lee et al. have suggested an E1-like E2 concerted mechanism.¹⁸ But a base responsible for the H-2 abstraction is required, which was not clear at that time. The most plausible candidate is the catalytic nucleophile itself (D553) in the glycosylation step. Until 2013 it was again suggested by the crystal structure that the catalytic nucleophile (D553) is in a correct position to act also as a base to abstracts the proton of the C2 carbon atom.²⁰ Since D553 plays a dual role in the catalytic cycle and no other residue is available to abstract the C-2 proton, and the C-2 proton should be acidified by the largely broken of C–O bond allowing its facile removal, the E1-like E2 concerted mechanism should be reconsidered. Besides, how the catalytic reaction proceeds are still unknown. To answer these questions, in this work, the combined quantum mechanics and molecular mechanics (QM/MM) method has been employed to explore the catalytic mechanism of α-1,4-glucan lyases. The QM/MM method has been widely applied in elucidating the catalytic mechanism of extended systems, including enzymes.^22–34 The overall view of the QM/MM approach is that the whole system was divided into two related regions, the QM region which contains the reactive part and the MM region which contains the environment part. The QM region is described by quantum mechanical while the MM region is treated by a classic force field.^35–40 In this methodology, the bond breaking and forming, and the effect of surrounding enzymatic residues could be reasonably described. Based on our calculations, the structures of intermediates and transition states, the energetic details as well as the roles of key residues were determined.

2. Computational details

2.1 Computational model

The initial atomic coordinates were derived from a recently issued X-ray crystal structure of α-1,4-glucan lyase of GH31 family from red seaweed Gracilariopsis lemaneiformis (EC 4.2.2.13, PDB code 2X2I), which was a complex of α-1,4-glucan lyase with a potent inhibitor acarbose (Acr).²⁰ The crystal structure and the active pocket are shown in Fig. 1a and b. Due to the high similarity between the inhibitor and the natural substrate glycogen, the residues in the active pocket are supposed to situate in their right positions to catalyze the reaction. Specifically, the valienamine residue of the inhibitor is bound to the subsite −1 of the enzyme. The nitrogen atom of valienamine moiety, which mimics the glucosidic oxygen atom at the non-reducing end of maltotetraose, makes hydrogen bond to residue D665 with a distance of 2.88 Å. The carboxyl group of D553, which is 3.46 Å away from the C1 atom at the non-reducing end, is correctly positioned for attacking C1. For subsite +1, residue R649 forms one hydrogen bond with the OH-2 group with a distance of 2.88 Å. D239 makes hydrogen bonds to the hydroxyl groups of the 6-deoxy-glucose moiety of acarbose. Besides, the side chain of residue D412 forms two hydrogen bonds with the OH-4 and OH-6 of valienamine residue. Overall, these residues are believed to be crucial for the biological catalysis of GLase, which will be described by QM method.

Fig. 1 (a) The crystal structure of the GLase (PDB code: 2X2I); (b) the structure of the proposed active pocket of the GLase crystal. Possible hydrogen bonds were shown in black dash lines.

To get the enzyme–substrate Michaelis complex, the acarbose in the active site was firstly deleted manually. Then, the maltotetraose, a kind of maltooligosaccharides, was chosen as the substrate and was placed into the active site by the help of the Autodock program.⁴¹ Before docking, the maltotetraose molecule was optimized with Gaussian 09 package at the level of B3LYP/6-31G(d,p).^42–46 When docking, the framework of the protein was kept rigid while all the flexible bonds of the maltotetraose were kept free. The grid box scale was set to 50 Å in all three dimensions by the use of Grid module, and the grid spacing was 0.375 Å. 100 independent docking runs were simulated and the representative conformation which performed the same conformation as the inhibitor was chosen as the start point of the following calculations. The protonation states of all titratable residues were determined based on the experimental condition, and the pK_a values were predicted by the PROPKA 3.1 program.⁴⁷ The final structure was checked carefully by the VMD program.⁴⁸ The missing hydrogen atoms were generated by using the HBUILD module in the CHARMM package.⁴⁹ Subsequently, the obtained enzyme–substrate model was surrounded by a water sphere of 45 Å radius which contains 9192 water molecules, while the crystallographic water molecules were kept in their original positions. Then the system was neutralized using 27 sodium cations at random positions. Finally, a neutral system was formed, which contains 42 150 atoms. The resulting system was equilibrated with a succession of minimizations, and a 10 ns MD simulations were performed with the CHARMM22/CMAP all-atom force field⁵⁰ for the protein. The water molecules were described by TIP3P model program.⁵¹ The root-mean-squared deviation (RMSD) of the protein during the MD simulation was derived and shown in Fig. S1.† The dynamics trajectory was basically stable after 6 ns with the RMSD value of 1.3 Å, which reveals that the backbone of the protein only changes slightly.

2.2 QM/MM calculations

During the QM/MM calculations, the QM region includes two sugar rings, the side chains of Asp239, Asp412, Asp553 and protonated Asp665, Asn459 and Arg649, as shown in Fig. 2. The total charge of the QM region is −2, and the remaining atoms were defined as MM subsystem.

Fig. 2 The selected QM region in our QM/MM calculations for Model 1.

ChemShell package⁵² integrating TURBOMOLE⁵³ and DL-POLY programs⁵⁴ were employed throughout the QM/MM calculations. The atoms in QM region and the MM atoms within a distance of 16 Å from the glucosidic oxygen atom of maltotetraose were kept free during the geometrical optimization, while the remaining MM atoms were kept frozen. The boundaries between QM region and MM region were handled by the link atom approach with a charge shift scheme.^55,56 In order to avoid the hyperpolarization of the QM wave function, the electronic embedding scheme⁵⁷ was used to incorporate the MM atomic point charges into the one-electron Hamiltonian of the QM calculation.

Geometry optimizations were implemented by the hybrid delocalized internal coordinates (HDLC) optimiser⁴⁹ in Chemshell. To optimize the geometry, the QM part was described by the B3LYP functional and 6-31G(d,p) basis set with Turbomole module,⁵³ and the MM region was described by the CHARMM22 force field⁵⁰ with DL-POLY program.⁵⁴ During the QM/MM calculations, the limited memory Broyden–Fletcher–Goldfarb–Shanno (L-BFGS) algorithm^58,59 was used for minima search and the partitioned rational function optimization (P-RFO) algorithm⁶⁰ was used for transition state search. The L-BFGS algorithm is one of Quasi-Newton methods and approximates to the inverse Hessian matrix using a limited memory variation of the BFGS update. It is an excellent method for large scale optimization and well suited for optimization problems with many variables. After the geometry optimization, the single-point calculations were performed at B3LYP/6-311++G(2d,2p) level to get more accurate energies. To correct the B3LYP energies for dispersion, the DFT-D3 program was used to calculate the empirical dispersion correction (B3LYP-D).⁶¹ Thus, all the reported energies in this work are single-point energies at B3LYP/6-311++G(2d,2p) including empirical dispersion correction. Based on our calculations, all the energies of dispersion correction are very small, as shown in Table S3,† and therefore the contribution of the dispersion correction to the calculated energy profile is neglectable.

3. Results and discussion

3.1 Structure of the GLase in complex with maltotetraose

It is well known that the initial structure may greatly influence the calculation results.^62–64 To guarantee our selection of the initial computational model to be representative, a series of snapshots were taken from the MD simulation at an interval of 200 ps from 8 ns to 10 ns. These eleven snapshots were firstly optimized by using QM/MM method, as shown in the superposition of the active sites of optimized geometries (Fig. S2†). To find a representative structure, an average structure was firstly derived from the 11 optimized geometries, and then the RMSDs of the 11 structures relative to the average one was calculated, which are shown in Fig. S2.† Since the RMSD at 9.6 ns corresponds to the smallest value, therefore, the structure at 9.6 ns was thought to be representative and was used as the computational model for the following studies. The structure of enzyme–maltotetraose complex (denoted as R) used for the following calculations is shown in Fig. 3.

Fig. 3 The optimized active site structure of GLase in complex with the substrate, which is denoted as the reactant (R) in our calculations. The atoms in QM region are displayed in ball and stick models. Some key distances are shown in angstrom. Hydrogen bonds were shown in black dash lines.

Compared with the crystal structure of enzyme–inhibitor complex, there are some slight deviations on the positions of the residues, which may be mainly caused by the structural difference between the substrate maltotetraose and inhibitor acarbose. After optimization, the proposed catalytic nucleophile D553 is also located at the same position as that of the crystal structure. In the active center, the non-reducing end (glucose unit) is gripped by the two catalytic residues D553 and D665. The distance between D553 and D665 is 6.4 Å (not shown in Fig. 3), which is slightly larger than that (5.8 Å) of the crystal structure of α-1,4-glucan lyase.²⁰ It should be noted that the distance between these two conserved residues in the typical retaining GHs enzyme is approximately 5.5 Å.⁷ The protonated residue D665 makes a strong hydrogen bond to the glucosidic oxygen atom with a distance of 1.90 Å, which is supposed to facilitate the protonation of glucosidic oxygen atom. One amino group of residue R649 forms two hydrogen bonds with the C2 and C3′ hydroxyl groups. Besides, the C3′ hydroxyl is within a hydrogen bonding distance from D239. D412 forms two hydrogen bonds with Glc1. These hydrogen bond interactions are in accordance with those of the crystal structure of enzyme–inhibitor complex, which make the substrate tightly bounded to the enzyme with specific orientation, and are vital to enzymatic reaction.²⁰

3.2 Reaction paths

Based on the optimized structure of the maltotetraose complex (R), the catalytic mechanism of GLase was systematically studied by using QM/MM method. Since the whole catalytic cycle contains two successive steps involving glycosylation and deglycosylation/elimination, two parts of the mechanism will be respectively discussed in the following sections.

3.2.1 Glycosylation process. According to our calculations and the proposed mechanism,²⁰ the glycosylation contains two elementary steps (Scheme 3), including the protonation of glycosidic oxygen and the formation of GLase–substrate intermediate.

Scheme 3 The catalytic reaction of GLase based on our calculations.

The optimized structures of the transition states (TS1, TS2) and intermediates (IM1, IM2) are displayed in Fig. 4. As shown in Fig. 4, the glycosidic oxygen is firstly protonated by the carboxyl proton of D665. In TS1, the distances of H_A–O_A and H_A–O₁ (Fig. 4) change to 1.60 Å and 1.00 Å, respectively, and the length of α-1,4-glycosidic bond (C1–O₁) increases from 1.43 Å to 1.76 Å, which means in TS1 the proton (H_A) of D665 is basically shifted to the glycosidic oxygen, and the glycosidic bond has been greatly loosened. Thus, we can conclude that the protonation of glycosidic oxygen greatly facilitates the broken of glycosidic bond. Besides, the distance between the C1 of glycose and O_B of carboxyl of D553 decreases from 3.63 Å to 3.31 Å.

Fig. 4 Optimized structures of transition states and intermediates for the glycosylation step. Distances are given in angstrom. Hydrogen bonds were shown in black dash lines.

Clear changes were also found in the sugar ring. For comparison, some key structural parameters and charges are listed in Table S1 of ESI.† In TS1, the bond length of C1–O5 decreases from 1.38 Å to 1.30 Å, implying that the C1–O5 bond is changing from single to double bond. In addition, the dihedral angle defined by H1–C1–O5–C2 (θ₁) changes from −117.4° to −136.9°, indicating a tendency of in-plane configuration of H1, C1, O5 and C2, which is in consistent with the hybridization of anomeric carbon changing from sp³ to sp². Owing to protonation of glycosidic oxygen, the net charge of C1 changes from 0.54e to 0.67e, and that of O5 changes from −0.67e to −0.63e. All these changes indicate the formation of an oxocarbenium-ion-like transition state. In IM1, the glycosidic bond has been completely broken with the distance of C1–O1 increased to 2.70 Å. In this concerted process, the proton transfer occurs prior to the cleavage of the glycosidic bond. Besides, the bond length of C1–O5 decreases to 1.25 Å, the dihedral angle of θ₁ increases to 174.1°, and the charges of the C1 and O5 change to 0.77 and −0.57e, respectively, which indicate the formation of the oxocarbenium ion transient intermediate. Interestingly, the largest positive charge on the anomeric carbon does not present on the TS1 but in IM1 in this mechanism, which was also discovered by Biarn X. and his coworkers in the study of 1,3-1,4-β-glucanase.⁶⁵

The following step is the formation of glycosyl–enzyme intermediate (IM2), in which the C1 atom is attacked by the carboxyl O_B atom of D553. In TS2, the distance between the O_B of D553 and the C1 of sugar ring decreases to 2.13 Å. In IM2, this distance further decreases to 1.49 Å to form a new glycosidic bond, and the anomeric carbon inverts its configuration.

The energy profile of the glycosylation step is shown in Fig. 5. It is shown that the concerted protonation of glycosidic oxygen and cleavage of glucosidic bond is calculated to be endothermic by 7.96 kcal mol⁻¹ with an energy barrier of 15.78 kcal mol⁻¹. The following step, the formation of glycosyl–enzyme intermediate, is very easy to occur. It is not surprising for a ligation of a carboxyl and an oxocarbenium ion intermediate.

Fig. 5 Energy profiles of glycosylation and deglycosylation/elimination steps. The energies are calculated at the level of 6-311++G(2d,2p) and corrected by DFT-D3.

For this process, the possibility that the nucleophilic attack of D553 on the anomeric carbon (C1) occurs firstly to trigger the reaction was also explored. Starting from the initial optimized structure R, the reaction coordinate (RC), defined as d = r₁(O1⋯C1) − r₂(C1⋯O_B), was scanned in a step size of 0.06 Å. It was found that, with the reaction coordinate gradually changing from −2.17 Å to −1.04 Å (the peak of the curve in Fig. S3†), r₁ increases from 1.44 Å to 1.91 Å, which means the glycoslic bond is nearly broken. But after this point, the system suddenly reaches to a state, in which the O_B of carboxyl of D553 departs from the C1 of the Glc1, and the proton of D665 completely attaches to the glucosidic O1. In other words, accompanying the gradually cleavage of glycosidic bond, the protonation of the glycosidic oxygen occurs automatically. Besides, the highest energy in the scanning curve is as high as 27.44 kcal mol⁻¹ at the level of B3LYP/6-31G(d,p). All these results suggest the pathway that the nucleophilic attack occurs firstly to be unlikely.

3.2.2 Deglycosylation/elimination process. After the glycosylation process, the α-1,4-glycosidic bond of maltotetraose is completely broken and the non-reducing end forms a new covalent bond with D553. But the remaining moiety (maltotriose) still forms hydrogen bonds with D239 and D665. Thus, for the deglycosylation process, two cases should be considered. In the first case, the remaining maltotriose group still presents in the active site, and in another case the maltotriose group has diffused into the solution. To fully elucidate this process, two models (Model 1 and Model 2) of the glycosyl–enzyme complex have been constructed.

Model 1 was constructed from the previously formed glycosyl–enzyme intermediate (IM2), in which the maltotriose group still exists in the reactive center. Based on our calculations, three elementary steps was determined in the deglycosylation process. The optimized structures of transition states and intermediates are shown in Fig. 6, and the energy profile is displayed in Fig. 5.

Fig. 6 The optimized structures of transition states, intermediates and product for deglycosylation step. Distances are given in angstrom.

In the first step, the scissile glycosidic bond of the glycosyl–enzyme intermediate is broken, generating the oxocarbenium ion intermediate IM3 via transition state TS3. This step is considered to be indispensable for the following elimination reaction. It is the oxocarbenium ion character of this species that leads to substantial acidification of the C-2 proton allowing its relatively facile removal. In TS3, the glycosidic bond loosens gradually with the bond length changing from 1.49 Å (IM2) to 1.98 Å. Besides, the C1–O5 bond length decreases from 1.35 Å to 1.27 Å (Table S1†), and the dihedral angle defined by H1–C1–O5–C2 (θ₁) increases to 156.5°. All these changes support the formation of the oxocarbenium ion character transition state. In IM3, the length of glycosidic bond further increases to 2.27 Å, implying the glycosidic bond has already broken. The bond length of C1–O5 is 1.35 Å and the dihedral angle of θ₁ is 172.8°, implies that the C1 atom has been changed to sp² hybridized. The calculated energy barrier of this process is 8.79 kcal mol⁻¹.

The following step is the rotation of the carboxyl group of residue D553. To abstract the C2-proton of the glycosyl group, the carboxyl group of residue D553 undergoes a rotation of 73.1° around the C–C bond to adjust the orientation of O_C. In IM3, the distance between O_C and C2-proton is 2.47 Å. After rotation, this distance decreases to 2.06 Å (IM4). More importantly, in IM4 the carboxyl group of D553 is almost co-planar with the C1, C2 and C2-proton of the glycosyl group, which is favorable for the proton abstraction. The calculated energy barrier of this rotation is only 3.49 kcal mol⁻¹, which is very easy to occur.

The last step is the abstraction of the C2-proton of glycosyl group by residue D553 to generate the final product 1,5-anhydrofructose (1,5AnFru) via TS5. In TS5, the distance between the C2-proton and the O_C atom of D553 decreases from 2.06 Å to 1.27 Å. Accompanying this abstraction, the bond length of C1–C2 of glycosyl group changes from 1.50 Å in IM4 to 1.39 Å in TS5, and finally changes to 1.32 Å in the final product. In P, the dihedral angles θ₁ and θ₂ change to −176.9° and 171.6°, respectively, further suggesting the formation of a double bond between C1 and C2. Besides, the hydrogen bond between D665 and the glycosyl group in +1 site is disappeared owing to conformational change of D665. This elementary step is calculated to be rate-determining with an energy barrier of 20.69 kcal mol⁻¹. And the whole deglycosylation/elimination process for this model corresponds to an overall energy barrier of 28.92 kcal mol⁻¹. The relative energy of the product (P) is 1.49 kcal mol⁻¹, indicating that the whole reaction cycle is slightly endothermic.

Besides, after the glycosylation step, the remaining maltotriose may diffuse to the solution, and only the covalent glucosyl residue remains in the active site. Thus, on the basis of intermediate IM2, a new computational model (Model 2) was constructed. Firstly, the remaining maltotriose group was removed from the active pocket. Then, an 8 ns MD simulation was performed to equilibrate the whole system. The time dependence of RMSD is shown in Fig. S4.† The final representative snapshot was employed for the following QM/MM calculations.

The obtained computational model (denoted as IM2-1) from QM/MM calculations is shown in Fig. 7. One can see that the distance between O_C atom of D553 and the C2-proton of glucoside group is 2.61 Å, which is slightly shorter than that of the intermediate structure (IM2). The glucosidic bond of the intermediate is almost unchanged (1.49 Å vs. 1.50 Å). The selected QM region for this model is shown in Fig. 8. Since there is no direct interaction between D239 and the glucoside group, D239 was not included in the QM region.

Fig. 7 The optimized active site structure of GLase–substrate intermediate (IM2-1), in which the maltotriose has been removed.

Fig. 8 The selected QM region in our QM/MM calculations for Model 2.

Starting from IM2-1 (Model 2), the structures of transition states (TS3-1, TS4-1 and TS5-1), intermediates (IM3-1, IM4-1) and product (P-1) for deglycosylation/elimination step were further optimized, which are shown in Fig. 9. The energy profile is shown in Fig. 10.

Fig. 9 The optimized structures of transition states, intermediates and product for deglycosylation step by using Model 2. Distances are given in angstrom.

Fig. 10 The energy profile of the deglycosylation step for Model 2.

One can see that the deglycosylation/elimination step also contains three elementary steps. Firstly, the glucosidic bond is cleaved to generate an oxocarbenium-ion-like intermediate (IM3-1) via TS3-1. As shown in Fig. 9, in TS3-1 the length of O_B–C1 increases from 1.50 Å to 1.77 Å, and the distance between the carboxyl O_C of residue D553 and the C2-proton of glycosyl decreases from 2.61 Å to 1.95 Å. Other bond lengths also display minor changes. For example, the length of C1–O5 decreases form 1.34 Å to 1.29 Å (see Table S2†), and length of C1–C2 decreases from 1.54 Å to 1.52 Å. Moreover, the dihedral angle of θ₁ increases from 127.3° to 141.1°, implying the nearly in-plane configuration of H1, C1, O5 and C2 atoms and sp² hybridization character of C1. The charge of C1 increasing from 0.55e to 0.64e further reveals the formation of an oxocarbenium-ion-like transition state. In IM3-1, the glucosidic bond is further loosened to 2.32 Å, the bond length of C1–O5 decreases to 1.24 Å and θ₁ increases to 170.7°. After the cleavage of glucosidic bond, the carboxyl group of D553 undergoes a self-rotation around C–C bond, making its O_C more close to the C2-proton of glycosyl group, which corresponds to an energy barrier of 12.06 kcal mol⁻¹. Then, the C2-proton of glycosyl group was abstracted by the carboxyl O_C of D553, generating the final product (P-1) via a six-member ring transition state (TS5-1). In TS5-1, the bond length of C1–C2 decreases from 1.50 Å in IM4-1 to 1.39 Å (Table S2†), and dihedral angles θ₁ and θ₂ change to −173.1° and −141.7° respectively, indicating the formation of a partial double bond between these two atoms. In P-1, the C2-proton of glycosyl group has been completely abstracted by the carboxyl O_C of D553. These four atoms (C1, C2, O5, and H1) are coplanar with θ₁ and θ₂ being 179.9° and 178.9°, respectively. Based on the energy profile, the abstraction of C2-proton by D553 is calculated to be the rate-limiting step with an energy barrier of 18.53 kcal mol⁻¹. But the whole deglycosylation/elimination step corresponds to an overall energy barrier of 23.93 kcal mol⁻¹. By comparing the calculations of Model 1 and Model 2, we can conclude that the deglycosylation/elimination step is easier to occur if the cleaved maltotriose group departs from the binding pocket.

As is well known, the partition of the QM region and MM region is significant to QM/MM calculation. To verify our selection of QM region is rational, a larger QM region model was also constructed, in which residue W551 and the amino of the peptide bond of M554 were also included, as shown in Fig. S5.† The QM/MM optimized structures are shown in Fig. S6–S8,† and the calculated energy profiles are shown in Fig. S9.† One can see that there are no obvious changes between the two QM region models, which indicate that our selection for QM region is reasonable.

4. Conclusions

In this work, the catalytic mechanism of α-1,4-glucan lyase from the red seaweed Gracilariopsis lemaneiformis has been investigated by using QM/MM method. On the basis of our calculations, the most significant conclusions can be summarized as follows:

(1) The glycosylation process starts from the protonation of glycosidic oxygen by the nearby residue D665, which is followed by the formation of a GLase–substrate intermediate. The protonation of glycosidic oxygen and cleavage of glycosidic bond follows a concerted but asynchronous mechanism, corresponding to an energy barrier of 15.78 kcal mol⁻¹. The formation of glycosidic bond between the substrate and D553 is quite easy to occur. Different from the retaining α-glucosidases whose glycosylation are typical concerted processes, the glycosylation process of glycosidic lyase follows a stepwise mechanism. The different distances (6.4 Å vs. 5.5 Å) between the two catalytic residues D553 and D665 in glycosidic lyase and retaining GHs enzymes may be the main reason to result in the different mechanism.

(2) The deglycosylation/elimination process is calculated to be the rate-limiting step of the whole catalytic reaction, in which the abstraction of C2-proton of the sugar ring by D553 corresponds to the highest energy barrier. D553 plays a dual role, which acts as nucleophile in the glycosylation step and base in the elimination step, respectively.

Acknowledgements

This work was supported by the Natural Science Foundation of China (21173129, 21373125) and Project of Shandong Province Higher Educational Science and Technology Program (no. J12LD10).

References

1. A. Kobata, Glycobiology: An expanding research area in carbohydrate chemistry, Acc. Chem. Res., 1993, 26, 319–324.
2. C.-H. Wong, R. L. Halcomb, Y. Ichikawa and T. Kajimoto, Enzymes in Organic Synthesis: Application to the Problems of Carbohydrate Recognition (Part 1), Angew. Chem., Int. Ed. Engl., 1995, 34, 412–432.
3. R. A. Dwek, Glycobiology: Toward Understanding the Function of Sugars, Chem. Rev., 1996, 96, 683–720.
4. R. Wolfenden, Benchmark reaction rates, the stability of biological molecules in water, and the evolution of catalytic power in enzymes, Annu. Rev. Biochem., 2011, 80, 645–667.
5. B. L. Cantarel, P. M. Coutinho, C. Rancurel, T. Bernard, V. Lombard and B. Henrissat, The carbohydrate active enzymes database (CAZy): an expert resource for glycogenomics, Nucleic Acids Res., 2009, 37, D233–D238.
6. D. E. Koshland, Stereochemistry and the mechanism of enzymatic reactions, Biol. Rev. Cambridge Philos. Soc., 1953, 28, 416–436.
7. D. L. Zechel and S. G. Withers, Glycosidase mechanisms: anatomy of a finely tuned catalyst, Acc. Chem. Res., 2000, 33, 11–18.
8. J. McCarter and S. G. Withers, Mechanisms of enzymatic glycoside hydrolysis, Curr. Opin. Struct. Biol., 1994, 4, 885–892.
9. Q. Wang, R. W. Graham, D. Trimbur, R. A. J. Warren and S. G. Withers, Changing enzymatic reaction mechanisms by mutagenesis: conversion of a retaining glucosidase to an inverting enzyme, J. Am. Chem. Soc., 1994, 116, 11594–11595.
10. G. J. Davies, V. M. Ducros, A. Varrot and D. L. Zechel, Mapping the conformational itinerary of beta-glycosidases by X-ray crystallography, Biochem. Soc. Trans., 2003, 31, 523–527.
11. J. Li, L. Du and L. Wang, Glycosidic-Bond Hydrolysis Mechanism Catalyzed by Cellulase Cel7A fromTrichodermareesei: A Comprehensive Theoretical Study by Performing MD, QM, and QM/MM Calculations, J. Phys. Chem. B, 2010, 114, 15261–15268.
12. B. C. Knott, M. Haddad Momeni, M. F. Crowley, L. F. Mackenzie, A. W. Götz, M. Sandgren, S. G. Withers, J. Ståhlberg and G. T. Beckham, The mechanism of cellulose hydrolysis by a two-step, retaining cellobiohydrolase elucidated by structural and transition path sampling studies, J. Am. Chem. Soc., 2014, 136, 321–329.
13. J. Wang, Q. Hou, L. Dong, Y. Liu and C. Liu, QM/MM studies on the glycosylation mechanism of rice BGlu1 β-glucosidase, J. Mol. Graphics Modell., 2011, 30, 148–152.
14. S. Badieyan, D. R. Bevan and C. Zhang, Probing the Active Site Chemistry of β-Glucosidases along the Hydrolysis Reaction Pathway, Biochemistry, 2012, 51, 8907–8918.
15. S. Yu, L. Kenne and M. Pedersen, α-1,4-Glucan lyase, a new class of starch/glycogen degrading enzyme. I. Efficient purification and characterization from red seaweeds, Biochim. Biophys. Acta, Gen. Subj., 1993, 1156, 313–320.
16. S. Yu, K. Bojsen, B. Svensson and J. Marcussen, α-1,4-Glucan lyases producing 1,5-anhydro-D-fructose from starch and glycogen have sequence similarity to α-glucosidases, Biochim. Biophys. Acta, Protein Struct. Mol. Enzymol., 1999, 1433, 1–15.
17. S. S. Lee, S. Yu and S. G. Withers, Mechanism of action of exo-acting α-1,4-glucan lyase: a glycoside hydrolase family 31 enzyme, BiologiaBratislava., 2005, 16, 137–148.
18. S. S. Lee, S. Yu and S. G. Withers, Detailed dissection of a new mechanism for glycoside cleavage: α-1,4-glucan lyase, Biochemistry, 2003, 42, 13081–13090.
19. S. S. Lee, S. Yu and S. G. Withers, α-1,4-Glucan lyase performs a trans-elimination via a nucleophilic displacement followed by a syn-elimination, J. Am. Chem. Soc., 2002, 124, 4948–4949.
20. H. J. Rozeboom, S. Yu, S. Madrid, K. H. Kalk, R. Zhang and B. W. Dijkstra, Crystal structure of α-1,4-glucan lyase, a unique glycoside hydrolase family member with a novel catalytic mechanism, J. Biol. Chem., 2013, 288, 26764–26774.
21. Y.-M. Kim, W. Saburi, S. Yu, H. Nakai, J. Maneesan, M.-S. Kang, S. Chiba, D. Kim, M. Okuyama, H. Mori and A. Kimura, A novel metabolic pathway for glucose production mediated by R-glucosidase-catalyzed conversion of 1,5-anhydrofructose, J. Biol. Chem., 2012, 287, 22441–22444.
22. M. Zheng and D. Xu, Catalytic mechanism of hyaluronate lyase from spectrococcus pneumonia: Quantum Mechanical/Molecular Mechanical and Density Functional Theory studies, J. Phys. Chem. B, 2013, 117, 10161–10172.
23. B. C. Knott, M. H. Momeni, M. F. Crowley, L. F. Mackenzie, A. W. Götz, M. Sandgren, S. G. Withers, J. Ståhlberg and G. T. Beckham, The mechanism of cellulose hydrolysis by a two-step, retaining cellobiohydrolase elucidated by structural and transition path sampling studies, J. Am. Chem. Soc., 2014, 136, 321–329.
24. X. Sheng and Y. Liu, A QM/MM study of the catalytic mechanism of nicotinamidase, Org. Biomol. Chem., 2014, 12, 1265–1277.
25. C. Li and S. Shaik, How do perfluorinated alkanoic acids elicit cytochrome P450 to catalyze methane hydroxylation? An MD andQM/MM study, RSC Adv., 2013, 3, 2995–3005.
26. B. Kang, H. Shi, S. Yan and J. Y. Lee, Solvent effect on electron and proton transfer in the excited state of a hydrogen bonded phenol–imidazole complex, RSC Adv., 2014, 4, 38551–38557.
27. A. Bobovská, I. Tvaroška and J. Kóňa, A theoretical study on the catalytic mechanism of the retaining α-1,2-mannosyltransferase Kre2p/Mnt1p: the impact of different metal ions on catalysis, Org. Biomol. Chem., 2014, 12, 4201–4210.
28. G. Ma, L. Dong and Y. Liu, Insights into the catalytic mechanism of dTDP-glucose 4,6-dehydratase from quantum mechanics/molecular mechanics simulations, RSC Adv., 2014, 4, 35449–35458.
29. K. Gkionis, S. T. Mutter and J. A. Platts, QM/MM description of platinum–DNA interactions: comparison of binding and DNA distortion of five drugs, RSC Adv., 2013, 3, 4066–4073.
30. X. Sheng, Y. Liu and R. Zhang, A theoretical study of the catalytic mechanism of oxalyl-CoA decarboxylase, an enzyme for treating urolithiasis, RSC Adv., 2014, 4, 35777–35788.
31. X. Sheng and Y. Liu, QM/MM Study of the Reaction Mechanism of the Carboxyl Transferase Domain of Pyruvate Carboxylase from Staphylococcus aureus, Biochemistry, 2014, 53, 4455–4466.
32. X. Chen, G. Ma, W. Sun, H. Dai, D. Xiao, Y. Zhang, X. Qin, Y. Liu and Y. Bu, Water Promoting Electron Hole Transport between Tyrosine and Cysteine in Proteins via a Special Mechanism: Double Proton Coupled Electron Transfer, J. Am. Chem. Soc., 2014, 136, 4515–4524.
33. Q. Hou, X. Sheng and Y. Liu, QM/MM studies of the mechanism of unusual bifunctional fructose-1,6-bisphosphate aldolase/phosphatase, Phys. Chem. Chem. Phys., 2014, 16, 11366–11373.
34. Q. Hou, L. Du, J. Gao, Y. Liu and C. Liu, QM/MM study on the reaction mechanism of O6-alkylguanine-DNA alkyltransferase, J. Phys. Chem. B, 2010, 114, 15296–15300.
35. H. M. Senn and W. Thiel, QM/MM methods for biomolecular systems, Angew. Chem., Int. Ed., 2009, 48, 1198–1229.
36. J. L. Gao, S. H. Ma, D. T. Major, K. Nam, J. Z. Pu and D. G. Truhlar, Mechanisms and free energies of enzymatic reactions, Chem. Rev., 2006, 106, 3188–3209.
37. Q. Cui and M. Karplus, Catalysis and specificity in enzymes: a study of triosephosphat isomerase (TIM) and comparison with methylglyoxal synthase, Adv. Protein Chem., 2003, 66, 15–372.
38. J. L. Gao and D. G. Truhlar, Quantum mechanical methods for enzyme kinetics, Annu. Rev. Phys. Chem., 2002, 53, 467–505.
39. M. J. Field, P. A. Bash and M. Karplus, A combined quantum-mechanical and molecular mechanical potential for molecular dynamics simulations, J. Comput. Chem., 1990, 11, 700–733.
40. A. Warshel and M. Levitt, Theoretical studies of enzymic reactions: dielectric, electrostatic and steric stabilization of the carbonium ion in the reaction of lysozyme, J. Mol. Biol., 1976, 103, 227–249.
41. G. M. Morris, D. S. Goodsell, R. S. Halliday, R. Huey, W. E. Hart, R. K. Belew and A. J. Olson, Automated docking using a Lamarckian genetic algorithm and an empirical binding free energy function, J. Comput. Chem., 1998, 19, 1639–1662.
42. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. Montgomery, Jr, J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, O. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski and D. J. Fox, Gaussian 09, Revision A.1, Gaussian, Inc., Wallingford, CT, 2009.
43. S. H. Vosko, L. Wilk and M. Nusair, Accurate spin-dependent electron liquid correlation energies for local spin density calculations: a critical analysis, Can. J. Phys., 1980, 58, 1200–1211.
44. C. Lee, W. Yang and R. G. Parr, Development of the ColleSalvetti correlation-energy formula into a functional of the electron density, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785–789.
45. A. D. Becke, Density-functional thermochemistry. III. The role of exact exchange, J. Chem. Phys., 1993, 98, 5648–5652.
46. P. J. Stephens, F. J. Devlin, C. F. Chabalowski and M. J. Frisch, Ab initio calculation of vibrational absorption and circular dichroism spectra using density functional force fields, J. Phys. Chem., 1994, 98, 11623–11627.
47. M. H. M. Olsson, C. R. Søndergard, M. Rostkowski and J. H. Jensen, PROPKA3: Consistent Treatment of Internal and Surface Residues in Empirical pK_a Predictions, J. Chem. Theory Comput., 2011, 7, 525–537.
48. W. Humphrey, A. Dalke and K. Schulten, VMD: visual molecular dynamics, J. Mol. Graphics Modell., 1996, 14, 33–38.
49. B. R. Brooks, R. E. Bruccoleri, B. D. Olafson, D. J. States, S. Swaminathan and M. Karplus, CHARMM: A program for macromolecular energy, minimization, and dynamics calculations, J. Comput. Chem., 1983, 4, 187–217.
50. A. D. MacKerell Jr, D. Bashford, M. Bellott, R. L. Dunbrack Jr, J. D. Evanseck, M. J. Field, S. Fischer, J. Gao, H. Guo, S. Ha, D. Joseph-McCarthy, L. Kuchnir, K. Kuczera, F. T. K. Lau, C. Mattos, S. Michnick, T. Ngo, D. T. Nguyen, B. Prodhom, W. E. Reiher III, B. Roux, M. Schlenkrich, J. C. Smith, R. Stote, J. Straub, M. Watanabe, J. Wiórkiewicz-Kuczera, D. Yin and M. Karplus, All-atom empirical potential for molecular modeling and dynamics studies of proteins, J. Phys. Chem. B, 1998, 102, 3586–3616.
51. W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, Comparison of simple potential functions for simulating liquid water, J. Chem. Phys., 1983, 79, 926–935.
52. P. Sherwood, A. H. de Vries, M. F. Guest, G. Schreckenbach, C. R. A. Catlow, S. A. French, A. A. Sokol, S. T. Bromley, W. Thiel, A. J. Turner, S. Billeter, F. Terstegen, S. Thiel, J. Kendrick, S. C. Rogers, J. Casci, M. Watson, F. King, E. Karlsen, M. Sjøvoll, A. Fahmi, A. Schäfer and C. Lennartz, QUASI: a general purpose implementation of the QM/MM approach and its application to problems in catalysis, J. Mol. Struct.: THEOCHEM, 2003, 632, 1–28.
53. R. Ahlrichs, M. Bär, M. Häser, H. Horn and C. Kälmel, Electronic structure calculations on workstation computers: The program system turbomole, Chem. Phys. Lett., 1989, 162, 165–169.
54. W. Smith and T. R. Forester, DL_POLY_2.0: a general-purpose parallel molecular dynamics simulation package, J. Mol. Graphics Modell., 1996, 14, 136–141.
55. M. J. Field, P. A. Bash and M. Karplus, A combined quantum mechanical and molecular mechanical potential for molecular dynamics simulations, J. Comput. Chem., 1990, 11, 700–733.
56. A. H. de Vries, P. Sherwood, S. J. Collins, A. M. Rigby, M. Rigutto and G. J. Kramer, Zeolite structure and reactivity by combined quantum-chemical-classical calculations, J. Phys. Chem. B, 1999, 103, 6133–6141.
57. S. R. Billeter, A. J. Turner and W. Thiel, Linear scaling geometry optimisation and transition state search in hybrid delocalized internal coordinates, Phys. Chem. Chem. Phys., 2000, 2, 2177–2186.
58. J. Nocedal, Updating quasi-Newton matrices with limited storage, Math. Comput., 1980, 35, 773–782.
59. D. C. Liu and J. Nocedal, On the limited memory BFGS method for large scale optimization, Math. Program., 1989, 45, 503–528.
60. A. Banerjee, N. Adams, J. Simons and R. Shepard, Search for stationary points on surfaces, J. Phys. Chem., 1985, 89, 52–57.
61. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu, J. Chem. Phys., 2010, 132, 154104.
62. I. F. Thorpe and C. L. Brooks III, Conformational substates modulate hydride transfer in dihydrofolate reductase, J. Am. Chem. Soc., 2005, 127, 12997–13006.
63. M. Klähn, S. Braun-Sand, E. Rosta and A. Warshel, On possible pitfalls in ab initio quantum mechanics/molecular mechanics minimization approaches for studies of enzymatic reactions, J. Phys. Chem. B, 2005, 109, 15645–15650.
64. P. Hu and Y. K. Zhang, Catalytic mechanism and product specificity of the histone lysine methyltransferase SET7/9: An ab initio QM/MM-FE study with multiple initial structures, J. Am. Chem. Soc., 2006, 128, 1272–1278.
65. X. Biarnés, A. Ardèvol, J. Iglesias-Fernández, A. Planas and C. Rovira, Catalytic itinerary in 1,3-1,4-β-glucanase unraveled by QM/MM metadynamics. Charge is not yet fully developed at the oxocarbenium ion-like transition state, J. Am. Chem. Soc., 2011, 133, 20301–20309.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c4ra09758k

---