Mutational Effects on the Catalytic Mechanism of Cellobiohydrolase I

ARTICLE pubs.acs.org/JPCB

Mutational Effects on the Catalytic Mechanism of Cellobiohydrolase I from Trichoderma reesei Shihai Yan, Tong Li, and Lishan Yao* Key Lab of Biofuels, Qingdao Institute of Bioenergy and Bioprocess Technology, Chinese Academy of Sciences, Qingdao 266101, PR China ABSTRACT: QM/MD simulations are performed to study mutational effects on the glycosylation step of the oligosaccharide hydrolysis catalyzed by Trichoderma reesei cellobiohydrolase I. The potential of mean force along the reaction pathway is determined by the umbrella sampling method. A detailed mechanism is developed to illustrate the decrease in activity of the mutants. Our calculations demonstrate that (1) the E212Q mutation increases the overall activation barrier by ∼4.0 kcal/mol, while the D214N mutation causes ∼0.4 kcal/mol increase of the barrier, and (2) there is only one transition state identified in the wild type (WT) and D214N mutant, while two transition states exist in the E212Q mutant for the glycosylation process. The results explain the experimental observation that the E212Q mutant loses most of its hydrolysis capability, while the D214N mutant only reduces it slightly compared to the WT. Further analysis suggests that the proton transfer from Glu217 to O4 and the glycosidic bond cleavage between subsites þ1 and
1 are concerted, facilitating the subsequent nucleophilic attack of Glu212 on C10 in subsite
1. Our QM/MD study illustrates the importance of the prearrangement of the active site and provides atomic details of the enzymatic catalytic mechanism.

1. INTRODUCTION Cellulose, the most abundant renewable biomass source on earth, is an insoluble homopolymer linked by unbranched β-1,4 glycosidic bonds. The biodegradation of cellulose gives birth to glucose, which has crucial functions in biological systems. Furthermore, glucose can be fermented to ethanol or other chemicals to replace nonrenewable petroleum. The growing needs for utilization and recovery of natural resources have prompted extensive scientific efforts toward a better recognition of the basic mechanisms behind cellulose degradation. Cellulases, the enzymes biodegrading the cellulose, usually consist of two domains, a large catalytic domain (CD) and a small carbohydrate-binding module (CBM), combined by a glycosylated linker peptide.1,2 Efficient enzymatic conversion of crystalline cellulose is important for an environmentally sustainable bioeconomy.3 The major bottleneck in the conversion of crystalline cellulose-based biomass to ethanol or other chemicals is the low activity of cellulases. It is now widely accepted that a better understanding of the cellulase-catalyzed hydrolysis of cellulose is highly desired and remains critical to enabling a successful bioethanol industry. The catalytic mechanisms for glycosyl hydrolases (GHs) are described in detail in several excellent reviews.2,4
12 Cellulolytic enzymes, like all GHs, hydrolyze glycosidic bonds via the mechanism of general acid/base catalysis, whereby two critical amino acid residues are needed, a proton donor and a nucleophilic assistant. Generally, the two residues playing these roles are glutamic and aspartic acid. So far, two distinct reaction r 2011 American Chemical Society

mechanisms for the GH superfamily enzymes have been proposed: inverting and retaining hydrolysis (Scheme 1). This scheme demonstrates that several covalent and noncovalent interactions are involved in the bond breaking and formation processes during the enzyme catalysis with both mechanisms. In the first mechanism (inverting hydrolysis), the RO group gets a proton from the protonated acidic residue and shifts away from the cellulose chain. Simultaneously, the negatively charged acidic residue activates water for nucleophilic attack on the carbocation center by electrostatic attraction.10,13
16 It has been suggested, based on the experimental observations and computational results, that the second mechanism is composed of two separate steps.13,14,17,18 Initially, along with proton transfer from the donor residue to the leaving RO group, a covalent glycosylenzyme intermediate is generated by the nucleophilic attack of the negatively charged residue on the glycosyl carbocation center. The original proton donor activates a water molecule by capturing one of its protons in the subsequent step. The hydrolysis product is released from the active site by nucleophilic attack of the water hydroxyl group. Accordingly, the retaining mechanism includes the separate glycosyl decomposition and hydrolysis steps. The predominant structural difference between inverting and retaining glycosidases is the distance between the proton donor and the nucleophile. The average distance between Received: January 13, 2011 Revised: March 30, 2011 Published: April 08, 2011 4982

dx.doi.org/10.1021/jp200384m | J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B

ARTICLE

Scheme 1

Figure 1. Schematic representation of T. reesei Cel7A (PDB entry 8CEL) with a cellooligomer bound. The acidic amino acid residue, Glu217, is above the cellulose chain, acting as the proton donor. Glu212, the nucleophile, is below the active site. Asp214 is behind the chain.

the catalytic residues is 5 ( 0.5 Å for the retaining and 10 ( 0.5 Å for the inverting enzymes.10,11,19 The large separation is necessary in inverting enzymes to allow the substrate and water to bind simultaneously between the carboxyl groups. Cel7A from Trichoderma reesei (T. reesei), also known as cellobiohydrolase I (CBHI), hydrolyzes the β-1,4 linkages of a cellulose chain from its reducing end via a retaining mechanism liberating the product, which consists of roughly 63% β-cellobiose and 37% R-cellobiose due to the anomeric equilibrium.1,20
26 It has been confirmed that Cel7A provides most hydrolytic power among the T. reesei enzymes contributing to cellulose degradation.27 The structural determination of the CD has paved the way for detailed investigations of the catalytic mechanism.17,28,29 Taking the position of the cleaved glycosidic bond as the reference point, the glucosylbinding subsites are numbered according to the generally accepted subsite-naming convention.30 The active site of CBHI contains three carboxylate residues (Glu212, Asp214, and Glu217) that appear strategically positioned for participation in catalysis (Figure 1). The proximity of Glu217 to the O4 atom in subsite þ1 suggests that this residue may act as a general acid by donating a proton, while Glu212 plays the role of a nucleophile during a double displacement at the anomeric carbon atom C10 . The third residue, Asp214, between the two glutamic acid residues, may be responsible for the correct positioning and protonation of Glu212. To assess the importance of these three carboxylate residues, several studies have been carried out by means of site-directed mutagenesis, where each carboxylate residue was replaced by its isosteric amide.28,29 The most dramatic effects were observed for the E217Q and E212Q mutations, while the smallest effect was shown for the D214N mutant, which still displays some activity. Although the structures have been determined in atomic detail, the intrinsic mechanism for the distinct activity reduction of the mutants remains unknown. To understand the mutational effect on enzyme catalysis, in principle one has to

study the whole process, including substrate binding, the catalytic reaction (glycosylation and deglycosylation), and product release. However, if the step controlling the reaction rate is known (assuming it is the same for the wild type and mutants), one can greatly reduce the amount of work by studying only this step. It has been proposed that glycosylation is rate-limiting for the substrate 2-chloro-4-nitrophenyl-β-lactoside.28 Therefore, illustrating the mechanism of the glycosylation step is critical to understand the mutational effects. Here, we apply a quantum mechanics molecular dynamics (QM/MD) simulation method to study this process for the wild type and the E212Q and D214N mutants. E217Q was not studied here since it has a similar effect as the E212Q mutant. To better mimic the hydrolysis of cellulose, we use a cellulose nanomer as the substrate. Computational modeling approaches can provide a powerful tool for understanding enzyme catalysis.17,31 QM/MD methods that combine the accuracy of quantum mechanics (QM) and the efficiency of molecular mechanics (MM) are an effective tool to understand characteristics of macromolecular systems.32 It allows for chemical bond breaking/formation and dynamics of the active center while including the effects of the fluctuating protein environment. In this work, QM/MD simulations for the above-described enzyme
substrate complexes were performed. The umbrella sampling method was used to move along a reaction coordinate appropriate to simulate the glycosylation step, and the weighted histogram analysis method (WHAM)33
35 was utilized to calculate the potential of mean force profile along this reaction coordinate. The simulations were carried out on a modeled structure, 8CEL, which is based on the protein X-ray crystal structures of the E217Q mutant complexed with cellohexaose and the E212Q mutant complexed with cellotetraose and cellopentaose, respectively. Two mutants, E212Q and D214N, were generated as described in Section 2. The computational details, including the preparation of the WT and mutant QM/MD structural models and the description of the hybrid QM/MD methodology, are presented in Section 2. Section 3 collects the results and discusses the mechanistic details. Our conclusions are summarized in Section 4.

2. COMPUTATIONAL DETAILS Preparation of the Simulation System. The initial coordinates of the WT were prepared according to a theoretical model, PDB entry 8CEL (Figure 1).29 The E212Q and D214N mutants were generated computationally, Glu212 Oε2 f Gln212 Nε2 for E212Q and Asp214 Oδ2 f Asn214 Nδ2 for D214N, respectively. 4983

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B E217 and D214 were assumed to be neutral based on the X-ray structure.28,29 All prepared systems were neutralized by adding Naþ ions with the Amber tool and were solvated in a rectangular box with a 12 Å buffer distance between the solvent box wall and the nearest solute atoms. For water molecules, the TIP3P model is utilized. Each of these three systems consists of ∼65 000 atoms. Throughout this work, all of the simulations were performed with the AMBER11 molecular dynamics package.36 The system was first minimized and equilibrated and then simulated for 1 ns with periodic boundary conditions, constant temperature at 300 K, and pressure of 1 atm. The FF99SB Amber37 and GLYCAM_06 force fields38 are employed to model the enzyme and the substrate, respectively. Hybrid QM/MD Approach. The active site of CBHI contains three carboxylate residues, Glu212, Asp214, and Glu217. The cellulose hydrolysis at the active site is initialized by the proton transfer from the donor residue Glu217 carboxyl group to O4, which holds the glucosyl residue to the reducing side of the glycosidic bond. The system was partitioned into QM and MM regions. The QM region, including three carboxylate residues as well as the glycosyl rings in subsites
1 and þ1, was treated by the PM3 semiempirical method.39,40 The MM region contains all other protein residues, solvent water molecules, as well as the glycosyl rings except the subsites
1 and þ1. The last snapshot from the classical MD simulation was employed for the subsequent QM/MD simulations (with the step size of 0.001 ps). The umbrella sampling and weighted histogram analysis method (WHAM)33
35 were utilized to determine the free-energy profile for the proton transfer and the nucleophilic attack. The Glu217 Hε2 proton transfer to the O4 atom was divided into ∼30 windows (with the Hε2
O4 distance as the reaction coordinate), each of which was simulated for 50 ps. Totally, ∼1.5 ns QM/MD simulations were carried out to characterize the proton transfer process in the WT and mutants, E212Q and D214N. The nucleophilic attacking process, corresponding to the bond formation between C10 in subsite
1 and Oε1 of Glu212 (or Gln212), was simulated with about 35 windows (with the C10
Oε1 distance as the reaction coordinate) and 50 ps for each window. Harmonic restraints with force constants of 400 kcal/mol Å2 were employed for proton transfer reaction and the subsequent nucleophilic attacking process by shortening the Hε2
O4 and C10
Oε distances, respectively. The configurations after 10 ps in each window were collected for data analysis. The probability distribution along the reaction coordinate was determined for each window and pieced together with the WHAM approach,15,16,41
44 to calculate the potential of mean force (PMF) profile. The error analysis was performed by dividing the gathered data into two equal blocks.

3. RESULTS AND DISCUSSION QM/MD simulations were carried out after classical molecular dynamics equilibration, which was performed for the fully solvated WT and the mutants at the NPT condition (300 K and 1 atm). The primary geometric parameters of the enzyme substrate complex (ES), the transition state of proton transfer process (TSP), the intermediate after proton transfer (IM), the transition state of nucleophilic attacking reaction (TSN), and the glycosylenzyme intermediate (GI) are all collected in Tables 1
3, as extracted from the QM/MD simulations for the WT and mutants. The PMFs for the WT and the mutants determined by the QM/MD simulations and umbrella sampling method are

ARTICLE

Table 1. Primary Geometric Parameters (Distances in Å) of the Enzyme Substrate Complex (ES), the Transition State of the Proton Transfer Process (TSP), the Intermediate after Proton Transfer (IM), the Transition State of the Glycosylation Reaction (TSN), and the Enzyme Glycosyl Intermediate (GI) for the WTa η1

107

N (Arg

a

ES

TSP

IM

TSN

GI 2.90

)
O3(
2)

2.90

2.89

2.89

2.90

Nη2(Arg107)
O6(
3)

3.05

3.02

3.02

4.72

4.72

Oη(Tyr145)
O2(
2)

2.89

2.85

2.85

2.79

2.79

Oδ2(Asp179)
O6(
3) Oη(Tyr247)
O6(
2)

2.55 3.87

2.56 4.18

2.56 4.18

2.54 3.91

2.54 3.91

O(Asp368)
O2(
3)

3.19

3.23

3.23

3.10

3.10

Oε2
Hε2(Glu217)

0.96

1.04

2.08

2.14

1.98

Hε2(Glu217)-O4

2.28

1.24

0.97

0.97

0.97

Oε2(Glu217)
O60

3.97

3.62

2.98

2.83

2.94

O4
C10

1.42

1.47

3.41

3.26

3.89

C10
Oε1(Glu212)

5.02

4.96

3.51

3.68

2.82

C10
Oε2(Glu212) Oε2(Glu212)
O20

3.60 2.74

3.48 2.74

2.15 3.06

2.34 3.02

1.44 3.22

Oε1(Glu212)
O20

3.16

3.18

2.70

2.70

2.64

Oδ2(Asp214)
Oε2(Glu212)

2.68

2.68

2.69

2.69

3.43

Nε2(His228)
O3

2.86

3.27

3.61

3.39

3.38

Nε2(His228)
Oδ1(Asp214)

3.09

3.17

3.03

3.06

3.08

Oδ1(Asp214)
O3

3.24

2.66

3.16

3.12

3.18

N(Ser174)
Oε1(Glu212)

2.90

2.89

3.07

3.08

3.30

Oγ(Ser174)-Nδ1(His228) Oγ(Ser174)
Nε2(Gln175)

3.50 3.06

3.57 2.99

2.96 3.03

3.24 3.11

2.88 3.21

Nε2(Gln175)
O20

3.23

3.50

3.22

3.21

3.14

Oε1(Gln175)
Nη1(Arg251)

3.11

3.28

3.07

3.10

3.10

Oε1(Gln175)
Nη2(Arg251)

2.88

2.89

2.88

2.87

2.85

Oδ1(Asp259)
Nε(Arg251)

2.84

2.85

2.84

2.86

2.84

Nη1(Arg394)
O1(þ2)

3.37

3.40

3.73

3.89

3.67

Nη1(Arg394)
O5(þ2)

2.95

2.90

2.93

3.24

2.95

Nη1(Arg394)
O6(þ2) Nη2(Arg394)
O1(þ2)

4.21 3.10

4.05 2.95

4.05 3.06

4.23 3.12

3.94 3.01

Nη2(Arg394)
O5(þ2)

3.96

3.84

3.65

3.52

3.64

Nη2(Arg394)
O6(þ2)

5.41

5.34

5.19

4.89

5.00

The data are averaged over the snapshots.

shown in Figure 2. The variation of the distance of O4 and C10 , O4
C10 , in the WT and mutants during the glycosylation process is presented in Figure 3. Figure 4 illustrates the enzyme
substrate interactions in the tunnel of CBHI. 3.1. Potential of Mean Force. To reveal the mechanism behind the mutants’ activity loss, the PMFs of the glycosylation process, which is composed of the proton transfer and the nucleophilic attack, for the WT and the mutants are determined employing the QM/MD simulations and umbrella sampling. The obtained PMFs are shown in Figure 2, which was drawn by setting the free energies of substrate complexes at zero. It should be emphasized that this origin is used only to compare the PMF barriers more easily. The binding free energies of the three substrate complexes are not determined here and most likely are different. As a result, the direct comparison of the absolute PMF of the three complexes is not meaningful. For the WT, a high free energy activation barrier of 32.6 kcal/mol needs to be overcome during the proton transfer process. This barrier is higher than 4984

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B

ARTICLE

Table 2. Primary Geometry Parameters (Distances in Å) of the Enzyme Substrate Complex (ES), the Transition State of the Proton Transfer Process (TSP), the Intermediate after Proton Transfer (IM), the Transition State of the Glycosylation Reaction (TsN), and the Enzyme Glycosyl Intermediate (GI) for the E212Q Mutanta ES η1

107

N (Arg

IM

TSN

GI η1

107

reactant

TSP

IM

TSN

GI 2.90

)
O3(
2)

2.93

2.90

2.98

2.89

2.87

N (Arg

)
O3(
2)

2.88

2.95

2.92

2.90

Nη2(Arg107)
O6(
3)

3.23

3.04

3.17

3.14

3.09

Nη2(Arg107)
O6(
3)

2.95

3.28

3.03

3.14

3.14

Oη(Tyr145)
O2(
2)

2.79

2.82

2.81

2.84

2.84

Oη(Tyr145)
O2(
2)

2.85

2.86

2.85

2.98

2.93

2.64 390

2.68 4.12

2.62 4.52

2.62 4.38

2.61 4.58

Oδ2(Asp179)
O6(
3) Oη(Tyr247)
O6(
2)

3.51 3.28

2.59 2.77

2.66 2.75

2.54 2.97

2.55 2.75

2.92

2.75

2.82

2.71

2.69

O(Asp368)
O2(
3)

3.85

3.54

3.59

3.56

3.66

0.98

1.78

1.85

3.38

3.10

Oε2
Hε2(Glu217)

0.96

1.66

1.79

1.78

1.81

Hε2(Glu217)
O4

1.82

1.02

0.96

0.95

0.95

Hε2(Glu217)
O4

1.97

1.04

0.98

0.97

0.96

Oε2(Glu217)
O60

4.63

3.11

3.06

4.56

4.58

Oε2(Glu217)
O60

3.83

3.15

3.27

3.46

3.05

O4
C10

1.44

4.05

3.08

5.13

5.11

O4
C10

1.43

2.88

3.47

4.92

4.91

C10
Oε1(Gln212)

3.46

3.33

2.80

2.09

1.47

C10
Oε1(Glu212)

4.34

3.51

3.39

2.29

1.39

C10
Nε2(Gln212) Nε2(Gln212)
O20

4.80 3.86

3.67 3.70

4.22 4.04

3.34 3.04

2.99 2.81

C10
Oε2(Glu212) Oε2(Glu212)
O20

3.42 3.19

2.92 3.47

2.76 3.43

4.26 4.50

3.36 4.17

3.52

4.33

3.55

3.29

3.08

Oε1(Glu212)
O20

2.84

2.76

2.71

3.34

3.07

2.69

2.69

2.75

2.66

2.78

Nδ2(Asn214)
Oε2(Glu212)

4.82

4.92

4.80

4.75

4.70

Nε2(His228)
O3

3.17

3.29

3.39

3.32

3.32

Nε2(His228)
O3

3.16

3.23

3.21

3.18

3.08

Nε2(His228)
Oδ1(Asp214)

3.04

3.14

2.99

2.97

2.95

Nε2(His228)
Oδ1(Asn214)

5.59

5.45

5.21

4.73

5.15

3.08

3.18

3.33

3.34

3.30

Oδ1(Asn214)
O3

3.23

3.23

2.93

2.82

3.24

5.88

5.68

5.75

5.93

5.99

N(Ser174)
Oε1(Glu212)

3.12

3.08

3.14

6.17

6.07

O (Ser )
N (His ) Oγ(Ser174)
Nε2(Gln175)

2.89 3.45

3.05 3.72

2.99 3.26

3.01 3.58

3.15 5.20

N(Ser174)
Nδ1(His228) N(Ser174)
Nε2(Gln175)

3.61 2.95

3.57 2.99

3.58 2.98

2.91 3.24

3.00 3.42

Nε2(Gln175)
O20

3.16

3.40

3.33

3.54

3.66

Nε2(Gln175)
O20

3.29

3.47

3.66

3.62

3.60

)

4.85

4.57

4.73

4.79

4.89

Oε1(Gln175)
Nη1(Arg251)

3.24

3.41

3.59

3.62

3.49

Oε1(Gln175)
Nη2(Arg251)

4.70

4.60

4.60

4.66

5.14

Oε1(Gln175)
Nη2(Arg251)

2.89

2.92

2.91

3.07

2.94

Oδ1(Asp259)
Nε(Arg251)

3.21

3.19

3.10

4.44

4.51

Oδ1(Asp259)
Nε(Arg251)

3.42

3.64

3.48

3.60

3.59

Nη1(Arg394)
O1(þ2)

3.88

3.99

3.91

3.85

3.95

Nη1(Arg394)
O1(þ2)

4.01

4.33

4.37

4.13

4.22

)
O5(þ2)

3.05

3.17

3.06

3.05

3.10

Nη1(Arg394)
O5(þ2)

2.96

3.13

3.15

3.00

3.06

Nη1(Arg394)
O6(þ2) Nη2(Arg394)
O1(þ2)

3.10 2.96

2.97 2.96

3.01 2.94

2.94 2.91

2.93 2.98

Nη1(Arg394)
O6(þ2) Nη2(Arg394)
O1(þ2)

3.18 3.12

2.95 3.19

2.93 3.11

2.97 3.10

2.92 3.11

Nη2(Arg394)
O5(þ2)

3.36

3.30

3.28

3.26

3.19

Nη2(Arg394)
O5(þ2)

3.49

3.37

3.27

3.41

3.37

4.79

4.57

4.63

4.56

4.47

Nη2(Arg394)
O6(þ2)

4.51

4.24

4.24

4.42

4.31

Oδ2(Asp179)
O6(
3) Oη(Tyr247)
O6(
2) O(Asp368)
O2(
3) ε2

ε2

O
H (Glu

217

)

Oε1(Gln212)
O20 δ2

214

O (Asp

δ1

ε1

212

)
O (Gln

)

214

O (Asp

)
O3

N(Ser174)
Oε1(Gln212) γ

ε1

η1

228

η1

175

251

)
N (Arg

394

N (Arg

η2

δ1

174

O (Gln

394

N (Arg a

TSP

Table 3. Primary Geometry Parameters (Distances in Å) of the Enzyme Substrate Complex (ES), the Transition State of the Proton Transfer Process (TSP), the Intermediate after Proton Transfer (IM), the Transition State of the Glycosylation Reaction (TSN), and the Enzyme Glycosyl Intermediate (GI) for the D214N Mutanta

)
O6(þ2)

The data are averaged over the snapshots.

that obtained by Li et al.18 for the WT CBHI-catalyzed glycosylation of cellobiose. In Li et al.’s work, the energy barrier of 14.1 kcal/mol is obtained with a simplified model system including three residues E212, D214, and E217 and cellobiose employing the density functional theory (DFT) method (B3LYP/6-31G (d,p)45). Several factors can contribute to this discrepancy. The three residues in the DFT calculations were allowed to move freely in the search for the transition state, which is artificial, since the position of these residues and cellobiose is restrained by the enzyme. In fact, the barrier increases to 23.9 kcal/mol, when the whole protein is included, closer to that obtained in this work. Furthermore, the barrier from Li et al.’s work is in energy while in our work it is free energy, and different QM methods used may also contribute to the barrier difference. Compared to the barrier of proton transfer, the barrier of the nucleophilic attacking step is small, only ∼0.4 kcal/mol. Thus, the whole process from the ES to GI can be taken as a single transition state (TSP) reaction, and the proton transfer process is

a

The data are averaged over the snapshots.

rate-limiting for the WT. The error analysis indicates that the free energy profile converges well, suggesting the sampling is sufficient. Similar to the WT, the glycosylation of D214N essentially has one transition state (TSP) as seen in Figure 2. The free energy activation barriers of the proton transfer process and the nucleophilic attack reaction are higher for the D214N mutant by ∼0.4 kcal/mol than those of WT. The slightly higher energy barrier indicates a somewhat lower catalytic activity. In addition, the GI of D214N is lower in energy by about 4.0 kcal/mol than its ES, as is favorable for the glycosylation step. The small 0.4 kcal/mol barrier increase corresponds to the kcat loss of 50% based on the Arrhenius equation, close to the experimental data where after 100 h incubation of bacterial microcrystalline cellulose (BMCC) with the D214N mutant the amount of reducing sugar is ∼35% less than that from incubation with the WT. It is worth noting that the comparison is very qualitative since the product quantity difference in experiment also depends on the 4985

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B

ARTICLE

Figure 2. Complete free energy profile for the proton transfer process (a) and nucleophilic attacking reaction (b) determined by the QM/MD simulations and umbrella sampling. The statistical error is estimated by averaging the free energy difference between 10 and 30 ps and 30
50 ps. The differences of the GI PMF values in the three enzyme complexes reflect different relative product binding affinity from that of the reactant.

binding kinetics and thermodynamics of the substrate/product to the WT and D214N enzyme, which are unknown. For the mutant E212Q, the free energy activation barrier of the proton transfer process is 31.2 kcal/mol, slightly lower than that of the WT. However, there is a higher barrier (8.4 kcal/mol) in the nucleophilic attacking step. The energy of the nucleophilic attacking transition state (TSN) is higher by more than 4.0 kcal/mol than the proton transfer transition state, TSP. Therefore, two transition states (TSP and TSN) exist during the process explored here, and the TSN determines the reaction rate. Another point should be mentioned that the GI of the E212Q mutant is thermodynamically unfavorable because its energy is 13 kcal/mol higher than its ES. Therefore, the hydrolysis rate should be slower than that catalyzed by the WT. The experimental study showed that after 100 h incubation of BMCC with the E212Q mutant the amount of reducing sugar is only ∼6% of that yielded from the incubation with the wild type enzyme,28 which corresponds to a barrier increase of 1.7 kcal/mol if the mutation only impacts kcat. This value is comparable to the barrier increase obtained in our work. The above analyses illustrate qualitatively the loss of the enzyme activity in the mutants, which confirms the previous experimental results.28,29 The small variations in PMF during the proton transfer process illustrate the slight influence brought by mutations, while distinct effects can be found for the nucleophilic attacking reaction and the GI free energy relative to the ES. From the structural point of view, it is reasonable since residues 212 and 214 are close to the nucleophilic attack reaction center and are somewhat removed from the proton transfer site. The second high free energy barrier of the E212Q mutant is conceivably caused by the net charge change of this residue from
1 to 0. Since the target of the nucleophilic attack C10 is positively charged, the weakening of the electrostatics between Oε (Q212) and C10 can no longer overcome the repulsion between the two atoms and corresponding fragments as they approach each other. Therefore, a new energy barrier is observed for E212Q but not for the WT and D214N. A more detailed geometric analysis is provided below to illustrate the mutational effects on the enzyme activity. 3.2. Active Site Geometry. Various proteins carry out different biological functions owing to their special geometric arrangements. The biological function of a protein is an inherent property of the structure, which is the basis of and tightly related to the function of the molecule. For the hydrolysis of cellulose by Cel7A, the reaction mechanism is determined by the coupling

mode and the conformation of the enzyme
substrate complex with the active site located between subsites
7 and þ2. The two critical residues, a proton donor and a nucleophilic assistant, lie between subsites
1 and þ1. Therefore, the interactions approximate to these two subsites are vital and play key roles during the hydrolysis process. Here, the couplings in subsites
3,
2, þ2, as well as
1 and þ1 are considered in detail. 3.2.1. Binding Subsites
3 and
2. The hydrogen bond interactions are summarized in Table 1 for the WT. The Oδ2(Asp179)
HO6(
3) hydrogen bond, between Oδ of the Asp179 side chain carboxyl group and the O6H hydroxyl of glucosyl in subsite
3, which was observed experimentally by Tones et al.,29 is a short, strong hydrogen bond (SSHB) in the ES state. The SSHB characteristic of this interaction remains during the whole proton transfer and the nucleophilic attack processes, as reflected by the short distance in Table 1. O6 in subsite
3 participates in another hydrogen bond with Nη2 of the Arg107 residue, Nη2H(Arg107)
O6(
3), which weakens in the GI state. The third hydrogen bond, found in subsite
3 between the O2 hydroxyl group and backbone carbonyl O of residue Asp368, O(Asp368)
HO2(
3), changes slightly during the reaction. In subsite
2, O3 of glucosyl forms a hydrogen bond with the Nη1H of Arg107. This hydrogen bond can be observed from the ES to the GI state. During the proton transfer and nucleophilic attack processes, two tyrosine residues, Tyr145 and Tyr247, located on the opposite sites of subsite
2, take part in the hydrogen bonds with O2 and O6, respectively, with the former stronger than the latter. Besides these two hydrogen bonds, glycosyl in subsites
2 also aligns approximately parallel to the tryptophan ring of Trp367. Such an arrangement minimizes the interaction energy and should be responsible for the twist of the glycosyl in this subsite to be almost perpendicular to that in subsite
4. All the hydrogen bonds and van der Waals interactions discussed here can be found in the E212Q and D214N mutants (Tables 2 and 3). It can be inferred that the influence of mutations to subsites
3 and
2 as well as their ambient interactions is small. 3.2.2. Binding Subsites
1 and þ1. One of the most prominent phenomena that occurs between subsites
1 and þ1 is the glycosidic bond O4
C10 cleavage as a result of the protonation of O4 by Hε2(Glu217). The variations of the glycosidic bond, O4
C10 , in the WT and mutants along the proton transfer coordinate and the nucleophilic attack reaction are presented in Figure 3. Plots (a)
(c) illustrate the variations 4986

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B

ARTICLE

Figure 3. Variations of the O4
C10 distance with the simulation time. (a)
(c) represent the distance changes during the proton transfer process for the wild-type enzyme
substrate complex and the mutants, respectively. (d)
(f) denotes the alterations during the nucleophilic attack process, respectively.

during the proton transfer process, and the changes in the nucleophilic attack step are collected in (d)
(f). When the length of Hε2(Glu217)
O4 approaches 1.0 Å, the O4
C10 bond lengthens significantly. In other words, the covalent O4
C10 bond cleavage is concerted with the Hε2(Glu217)
O4 bond formation, after which the cellobiose dissociates from the cellulose. Accompanying the O4
C10 bond cleavage, the Cδ(Glu212) 3 3 3 C10 distance decreases distinctly, favoring the subsequent nucleophilic attack. This phenomenon is observed in the WT and mutants. Therefore, the mutations maintain the property that the proton transfer and the glycosidic bond O4
C10 cleavage occur synchronously, followed by the nucleophilic attack. This is consistent with the earlier QM/MM results.17,18,32,46 In WT, the O4
C10 distance lengthens along with the decrease of the Cδ(Glu212)
C10 distance during the

nucleophilic attack of Glu212, which indicates that the cellobiose moves away from the active site. A similar increase of the distance is observed in E212Q, while in the D214N mutant it fluctuates around 4.8 Å during the nucleophilic attack process. At the end of the glycosylation reaction, a novel covalent bond is generated between C10 and Oε2 of the Glu212 carboxyl group (Tables 1, 2, and 3). The tables also indicate that this carboxyl group interacts with the O2H group in subsite
1 during the nucleophilic attack. In the WT, when the distance between Hε2 (Glu217) and O4 decreases from 2.28 to 0.97 Å, the Oε2
Hε2 distance of residue Glu217 increases continuously. This proton covalently links to O4 in the IM state (Table 1). The variation of the newly formed Hε2(Glu217)
O4 bond is very small in the subsequent nucleophilic attack process. The same phenomenon can also be observed in the E212Q and D214N mutants. The electrostatic 4987

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B

Figure 4. Snapshots of the enzyme
substrate intermediate. The configurations of residues 212, 214, and 217 as well as the ring conformations of subunit
1 are shown.

interaction between Oε2(Glu217) and O60 H in the ES state changes into the hydrogen bond in the IM state. This hydrogen bond still can be found after glycosylation, as is supported by the experimental X-ray structures.29 The same phenomenon occurs in D214N, while the hydrogen bond generated in IM disappears in the GI state of the E212Q mutant (Tables 1, 2, and 3). These variations are coincidental with the changes of the PMFs, but whether this interaction is related to the free energy barrier remains unclear. Another point should be mentioned that, in agreement with the previous report,28 a strong hydrogen bond exists between Asp214 (protonated) and Glu212 carboxyl groups in the ES of WT. This hydrogen bond is maintained until the formation of the covalent C10
Oε2(Glu212) bond in the GI state where it is weakened significantly (Table 1). This observation suggests that Oδ2 of Asp214 positions Oε2 of Glu212 to facilitate the nucleophilic attack process. The Oδ2H(Asp214)
Oε1(Gln212) hydrogen bond is maintained throughout the glycosylation step in the E212Q mutant. However, its counterpart in the D214N mutant, the Oδ2H(Asn214)
Oε1(Glu212) hydrogen bond, cannot be observed. Besides the hydrogen bonds mentioned above, a complex hydrogen bond network is also observed in the subsites
1, þ1, and their surroundings in the WT (Table 1). The O3H is hydrogen bonded to Oδ1(Asp214) and remains throughout the whole glycosylation process. Nε2 of His228, located around subsite þ1, acts as a hydrogen bond donor to the Asp214 carboxyl group, which forms another hydrogen bond with the O3H group in subsite þ1. These two hydrogen bonds, unchanged during the proton transfer and nucleophilic attack processes, effectively restrain the position of His228 with the assist of a third strengthening hydrogen bond, OγHγ

ARTICLE

(Ser174) 3 3 3 Nδ1(His228). Oγ(Ser174) is a hydrogen bond acceptor to the Nε2Hε22 group of residue Gln175 which is also hydrogen bonded to O20 as a donor, while NH(Ser174) forms a hydrogen bond to Oε1(Glu212) in the ES of WT (2.90 Å). The NH(Ser174) 3 3 3 Oε1(Glu212) hydrogen bond weakens during the reaction, especially along the covalent C10
Oε2(Glu212) bond formation in the nucleophilic attack process. All these H-bonds can be observed in mutants except that the Nε2H(His228)
Oδ1(Asp214) hydrogen bond disappears in the D214N mutant, and NH(Ser174) 3 3 3 Oε1(Glu212) disappears during the nucleophilic attack reaction in D214N and does not exist at all in E212Q (Tables 1, 2, and 3). In addition, Oε1(Gln175) couples bidentately with the Nη1
Hη12 and Nη2
Hη22 groups of Arg251 through a hydrogen bond which only exists in the WT and D214N mutant. In the WT, the distance between two carboxylate groups of the donor residue (Glu217) and the nucleophile (Glu212) is about 6.0 ( 0.5 Å, in good agreement with the commonly observed results in previous experimental and theoretical investigations.4,10,11,18,19 This distance, in the E212Q mutant, lengthens from 7.0 to 9.0 Å following the proton transfer and then is reduced back to ∼7.0 Å during the nucleophilic attack process. The distance change of these two carboxylate groups may be related to the high free energy barrier during the nucleophilic attack process. In D214N, the distance between two carboxylate groups of Glu212 and Glu217 is around 7.5 ( 0.5 Å, longer than that in WT, but shorter than that in E212Q on average. At the beginning of the proton transfer, the ring in subsite
1 (the ring labeled in Scheme 1) is in chair conformation and tilts by 60° with the plane defined by the other three parallel neighboring rings (subsites
2, þ1, and þ2). As the reaction proceeds, the tilt angle decreases rapidly. The glycosyl ring in subsite
1 varies from chair (4C1) to skew-boat (1S3) conformation, through a half-chair (4H3) transition state during the proton transfer reaction. These phenomena are shared by the WT and the mutants. The corresponding IM configurations are shown in Figure 4, and a common characteristic is observed: the C1, C2, C5, and O5 atoms form a plane. This is in good agreement with the previous reports.7,17 This ring conformation changes back to 4 C1 after the nucleophilic attack for WT.18 However, a 2,5B conformation comes into being in the GI state for the E212Q and D214N mutants. 3.2.3. Binding Subsite þ2. The Nε(Arg251) is hydrogen bonded with the carboxylate group of the Asp259 residue in three ES forms. During the whole reaction, the variation of this hydrogen bond is small in the WT and D214N mutant. However, the change is substantial during the nucleophilic attack step in E212Q. Arg394, located at the cellobiose exit, interacts with O1, O5, and O6 of subsite þ2. In the ES, O1 locates between two NH2 groups of Arg394, while the O5 and O6 are positioned inside the channel. During the proton transfer and the nucleophilic attack processes, cellobiose approaches toward Arg394 and moves out of the active site. Further study will be needed to prove whether the interaction between the two facilitates the departure of cellobiose.

4. CONCLUSIONS On the basis of the QM/MD simulation results, the configurations and the coupling details around the enzyme active site are analyzed. Along with the transfer of a proton from the donor residue, Glu217, to the glycosidic O4, the covalent CO bond, O4
C10 , is broken, and the cellobiose departs from the active site. At the same time, the distance between C10 and the 4988

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989

The Journal of Physical Chemistry B nucleophilic group Glu212 (or Gln212) decreases, which facilitates the subsequent nucleophilic attack reaction. A detailed mechanism is presented to illustrate the decrease of activity of the mutants during the hydrolysis of cellulose. The dramatic loss of enzyme catalytic capability of the E212Q mutant is due to the overall activation barrier increase of ∼4.0 kcal/mol. In contrast, the D214N mutation only leads to an increase of the activation energy by ∼0.4 kcal/mol, which accounts for the slight decrease of the catalytic efficiency of this mutant. The detailed geometric analysis suggests that a hydrogen bond network contributes to the active site interaction and enzymatic catalytic process. The distance of the two carboxylate groups of Glu212 and Glu217 is lengthened for the mutations.

’ AUTHOR INFORMATION Corresponding Author

*Tel.: þ86 532 80662792. Fax: þ86 532 80662778. E-mail: [email protected].

’ ACKNOWLEDGMENT We are thankful to Supercomputing Center of Chinese Academy of Sciences (CAS) for providing the computer resources and time. This work was supported by 100 Talent Project, the Knowledge Innovation Program of the CAS (Grant No. KSCX2-EW-J-10), the Director Innovation Foundation of Qingdao Institute of Biomass Energy and Bioprocess Technology of CAS, and the Foundation for Outstanding Young Scientist in Shandong Province (No. BS2010NJ020). ’ REFERENCES (1) Divne, C.; Stahlberg, J.; Reinikainen, T.; Ruohonen, L.; Pettersson, G.; Knowles, J. K. C.; Teeri, T. T.; Jones, T. A. Science 1994, 265, 524. (2) Mosier, N. S.; Hall, P.; Ladisch, C. M.; Ladisch, M. R. Adv. Biochem. Eng. Biotechnol. 1999, 65, 23. (3) Vaaje-Kolstad, G.; Westereng, B.; Horn, S. J.; Liu, Z.; Zhai, H.; Sorlie, M.; Eijsink, V. G. H. Science 2010, 330, 219. (4) Zechel, D. L.; Withers, S. G. Acc. Chem. Res. 2000, 33, 11. (5) Zechel, D. L.; Withers, S. G. Curr. Opin. Chem. Biol. 2001, 5, 643. (6) Vasella, A.; Davies, G. J.; Bohm, M. Curr. Opin. Chem. Biol. 2002, 6, 619. (7) Vocadlo, D. J.; Davies, G. J. Curr. Opin. Chem. Biol. 2008, 12, 539. (8) Rye, C. S.; Withers, S. G. Curr. Opin. Chem. Biol. 2000, 4, 573. (9) Crich, D. Acc. Chem. Res. 2010, 43, 1144. (10) McCarter, J. D.; Withers, S. G. Curr. Opin. Struct. Biol. 1994, 4, 885. (11) Davies, G.; Henrissat, B. Structure 1995, 3, 853. (12) Sinnott, M. L. Chem. Rev. 1990, 90, 1171. (13) Knowles, J. K. C.; Lentovaara, P.; Murray, M.; Sinnott, M. L. J. Chem. Soc., Chem. Commun. 1988, 1401. (14) Claeyssens, M.; Tomme, P.; Brewer, C. F.; Hehre, E. J. Febs Lett. 1990, 263, 89. (15) Petersen, L.; Ardevol, A.; Rovira, C.; Reilly, P. J. J. Phys. Chem. B 2009, 113, 7331. (16) Koivula, A.; Ruohonen, L.; Wohlfahrt, G.; Reinikainen, T.; Teeri, T. T.; Piens, K.; Claeyssens, M.; Weber, M.; Vasella, A.; Becker, D.; Sinnott, M. L.; Zou, J. Y.; Kleywegt, G. J.; Szardenings, M.; Stahlberg, J.; Jones, T. A. J. Am. Chem. Soc. 2002, 124, 10015. (17) Liu, J.; Wang, X.; Xu, D. J. Phys. Chem. B 2010, 114, 1462. (18) Li, J.; Du, L.; Wang, L. J. Phys. Chem. B 2010, 114, 15261. (19) Wang, Q. P.; Graham, R. W.; Trimbur, D.; Warren, R. A. J.; Withers, S. G. J. Am. Chem. Soc. 1994, 116, 11594. (20) Vrsanska, M.; Biely, P. Carbohydr. Res. 1992, 227, 19.

ARTICLE

(21) Barr, B. K.; Hsieh, Y. L.; Ganem, B.; Wilson, D. B. Biochemistry 1996, 35, 586. (22) Bu, L. T.; Beckham, G. T.; Crowley, M. F.; Chang, C. H.; Matthews, J. F.; Bomble, Y. J.; Adney, W. S.; Himmel, M. E.; Nimlos, M. R. J. Phys. Chem. B 2009, 113, 10994. (23) von Ossowski, I.; Stahlberg, J.; Koivula, A.; Piens, K.; Becker, D.; Boer, H.; Harle, R.; Harris, M.; Divne, C.; Mahdi, S.; Zhao, Y. X.; Driguez, H.; Claeyssens, M.; Sinnott, M. L.; Teeri, T. T. J. Mol. Biol. 2003, 333, 817. (24) Boer, H.; Teeri, T. T.; Koivula, A. Biotechnol. Bioeng. 2000, 69, 486. (25) Kipper, K.; Valjamae, P.; Johansson, G. Biochem. J. 2005, 385, 527. (26) Igarashi, K.; Koivula, A.; Wada, M.; Kimura, S.; Penttila, M.; Samejima, M. J. Biol. Chem. 2009, 284, 36186. (27) Baker, J. O.; Ehrman, C. I.; Adney, W. S.; Thomas, S. R.; Himmel, M. E. Appl. Biochem. Biotechnol. 1998, 70
2, 395. (28) Stahlberg, J.; Divne, C.; Koivula, A.; Piens, K.; Claeyssens, M.; Teeri, T. T.; Jones, T. A. J. Mol. Biol. 1996, 264, 337. (29) Divne, C.; Stahlberg, J.; Teeri, T. T.; Jones, T. A. J. Mol. Biol. 1998, 275, 309. (30) Biely, P.; Kratky, Z.; Vrsanska, M. Eur. J. Biochem. 1981, 119, 559. (31) Warshel, A. Annu. Rev. Biophys. Biomol. 2003, 32, 425. (32) Saharay, M.; Guo, H.; Smith, J. C. PLoS One 2010, 5, e12947. (33) Kumar, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; Rosenberg, J. M. J. Comput. Chem. 1992, 13, 1011. (34) Souaille, M.; Roux, B. Comput. Phys. Commun. 2001, 135, 40. (35) Ferrenberg, A. M.; Swendsen, R. H. Phys. Rev. Lett. 1988, 61, 2635. (36) Case, D.A.; D., T. A.; Cheatham, T. E., III; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R. C. W.; Zhang, W.; Merz, K. M.; Roberts, B.; Wang, B.; Hayik, S.; Roitberg, A.; Seabra, G.; K., I.; Wong, K. F.; Paesani, F.; Vanicek, J.; Wu, X.; Brozell, S. R.; Steinbrecher, T.; G., H.; Cai, Q.; Ye, X.; Wang, J.; Hsieh, M.-J.; Cui, G.; Roe, D. R.; Mathews, D. H.; S., M. G.; Sagui, C.; Babin, V.; Luchko, T.; Gusarov, S.; Kovalenko, A.; Kollman, P. A. AMBER 11; University of California: San Francisco, 2010. (37) Hornak, V.; Simmerling, C. Proteins 2003, 51, 577. (38) Kirschner, K. N.; Yongye, A. B.; Tschampel, S. M.; GonzalezOuteirino, J.; Daniels, C. R.; Foley, B. L.; Woods, R. J. J. Comput. Chem. 2008, 29, 622. (39) Aliev, A. E.; Courtier-Murias, D. J. Phys. Chem. B 2007, 111, 14034. (40) Fong, P.; McNamara, J. P.; Hillier, I. H.; Bryce, R. A. J. Chem. Inf. Model. 2009, 49, 913. (41) Wang, S. L.; Hu, P.; Zhang, Y. K. J. Phys. Chem. B 2007, 111, 3758. (42) Lu, Z. Y.; Zhang, Y. K. J. Chem. Theory Comput. 2008, 4, 1237. (43) Wu, R. B.; Wang, S. L.; Zhou, N. J.; Cao, Z. X.; Zhang, Y. K. J. Am. Chem. Soc. 2010, 132, 9471. (44) Zheng, H.; Wang, S. L.; Zhang, Y. K. J. Comput. Chem. 2009, 30, 2706. (45) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. B: Condens. Matter 1988, 37, 785. (46) Liu, D. J.; Nimlos, M. R.; Johnson, D. K.; Himmel, M. E.; Qian, X. H. J. Phys. Chem. A 2010, 114, 12936.

4989

dx.doi.org/10.1021/jp200384m |J. Phys. Chem. B 2011, 115, 4982–4989