A coarse-grained model for synergistic action of multiple enzymes on cellulose
© Asztalos et al.; licensee BioMed Central Ltd. 2012
Received: 29 March 2012
Accepted: 21 June 2012
Published: 1 August 2012
Skip to main content
© Asztalos et al.; licensee BioMed Central Ltd. 2012
Received: 29 March 2012
Accepted: 21 June 2012
Published: 1 August 2012
Degradation of cellulose to glucose requires the cooperative action of three classes of enzymes, collectively known as cellulases. Endoglucanases randomly bind to cellulose surfaces and generate new chain ends by hydrolyzing β-1,4-D-glycosidic bonds. Exoglucanases bind to free chain ends and hydrolyze glycosidic bonds in a processive manner releasing cellobiose units. Then, β-glucosidases hydrolyze soluble cellobiose to glucose. Optimal synergistic action of these enzymes is essential for efficient digestion of cellulose. Experiments show that as hydrolysis proceeds and the cellulose substrate becomes more heterogeneous, the overall degradation slows down. As catalysis occurs on the surface of crystalline cellulose, several factors affect the overall hydrolysis. Therefore, spatial models of cellulose degradation must capture effects such as enzyme crowding and surface heterogeneity, which have been shown to lead to a reduction in hydrolysis rates.
We present a coarse-grained stochastic model for capturing the key events associated with the enzymatic degradation of cellulose at the mesoscopic level. This functional model accounts for the mobility and action of a single cellulase enzyme as well as the synergy of multiple endo- and exo-cellulases on a cellulose surface. The quantitative description of cellulose degradation is calculated on a spatial model by including free and bound states of both endo- and exo-cellulases with explicit reactive surface terms (e.g., hydrogen bond breaking, covalent bond cleavages) and corresponding reaction rates. The dynamical evolution of the system is simulated by including physical interactions between cellulases and cellulose.
Our coarse-grained model reproduces the qualitative behavior of endoglucanases and exoglucanases by accounting for the spatial heterogeneity of the cellulose surface as well as other spatial factors such as enzyme crowding. Importantly, it captures the endo-exo synergism of cellulase enzyme cocktails. This model constitutes a critical step towards testing hypotheses and understanding approaches for maximizing synergy and substrate properties with a goal of cost effective enzymatic hydrolysis.
Biofuel production from lignocellulosic materials is considered to be a promising option to substantially reduce the dependence on petroleum [1–3]. The conversion of lignocellulosic biomass (agronomic residues, paper wastes, energy crops) into ethanol consists of the extraction and pretreatment of cellulose from the biomass, hydrolysis (the enzymatic breakdown of crystalline cellulose fibers into monomer glucose) and finally the fermentation of glucose to ethanol. Current approaches mainly differ from one another in the method of pretreatment. Cost-competitive production of ethanol is currently prevented by the low efficiency of converting cellulose into glucose . Greater efficiency may be achievable through improvements in hydrolysis.
Enzymatic hydrolysis of cellulose is a complex reaction. In the classical model, the heterogeneous catalytic cleavage of the glycosidic bond takes place on a crystalline cellulose surface and requires the cooperative action of three classes of aqueous enzymes, collectively known as cellulases. These are (i) endoglucanases, (ii) exoglucanases or cellobiohydrolases and (iii) β-glucosidases. Recently, it has been proposed that oxidative enzymes (monoxygenases) may also play a role in cleaving glycosidic bonds, although this new mechanism may be restricted to certain types of microbes . It is widely accepted that endoglucanases cleave β-1,4-D-glycosidic bonds at random sites within both amorphous and crystalline polysaccharide chains, creating new chain ends on the cellulose surface [6–11]. Exoglucanases prefer to hydrolyze crystalline cellulose chains by acting on the free chain ends and releasing cellobiose units in a processive manner. Soluble cellobiose units are then converted into glucose by β-glucosidases. Consequently, these enzymes display strong synergy [12–18].
The classical chemical kinetics assumption of uniformly mixed systems does not hold in the case of enzymatic hydrolysis of cellulose fiber, as it is heterogeneous in nature. Such reactions are rather characterized by time-dependent rate constants and non-uniform concentration variation of reacting species. Although kinetic models [11, 19–22] have been used to explain various features of the enzymatic hydrolysis of cellulose, they fail to account for spatial details of the cellulose substrate as well as the specificity of binding sites. Recently, Zhou and colleagues [23–25] proposed a “morphology-plus-kinetics rate equation approach” that explicitly captures the hydrolytic evolution of cellulose substrate. In addition, a kinetic model was developed based on population-balance equations in which a distribution of chains with different chain-lengths was explored [26, 27]. Yet, these models give little insight regarding the action of cellulases at the molecular level.
It is imperative to develop spatial models of cellulose degradation because spatial effects such as enzyme crowding on the cellulose surface have been shown to lead to a reduction in hydrolysis rates. In order to account for the spatial heterogeneity of the system during cellulose hydrolysis, a cellular automata model  was developed to study the effect of different parameters such as enzyme binding and hydrolysis on the overall kinetics of cellulose by the cellulases. Alternatively, all-atom molecular dynamics (MD) simulations can provide details of molecular level events at high precision. Recent MD simulation studies [29–32] have proven effective for understanding enzyme-substrate binding, processivity and activity. However, because of length and time scale limitations, it is not currently possible to simulate the entire crystalline cellulose degradation process using all-atom MD simulations.
We have developed a coarse-grained stochastic model that captures the interaction of endo- and exo-cellulases with crystalline cellulose at a mesoscopic level. This model was specifically designed to improve our understanding of the molecular-level details of the enzymatic hydrolysis of crystalline cellulose. This paper introduces the basic framework and demonstrates how this model can be an effective and easily modifiable testing platform for new hypotheses based on experimental data on various cellulase components and substrate characteristics. By capturing the reactive nature of the cellulose substrate and the activities of non-complexed cellulases at the molecular level, this method forms a bridge between all-atom MD studies and deterministic reaction-rate approaches. To the best of our knowledge, it is the first model that is able to relate the synergetic action of multiple enzymes to molecular level details such as the hydrogen bond network of a cellulose substrate.
The overall efficiency of the heterogeneous catalysis that occurs in the enzymatic hydrolysis of crystalline cellulose depends on factors such as adsorption, desorption, diffusion rates on the insoluble cellulose substrate, and processivity. In our model, catalysis is broken down into distinct parts related to different kinetic events (chemical reactions) performed by individual particles (enzymes). Specifically, we include the following reactive events: adsorption of cellulases on the solid cellulose substrate, inter-chain hydrogen bond breaking, hydrolysis of glycosidic bonds, and desorption of cellulases from cellulose. These reactions constitute the main elements of this model, and their realization is achieved by following and updating the state (based on certain predefined rules discussed below) of each individual particle in the system as it evolves in time. The actions of cellulases are modeled based on the most abundant endoglucanase (EG I) and the two cellobiohydrolases (CBH I and CBH II) secreted by the filamentous fungus Trichoderma reesei, as these three enzymes have been widely studied [10, 11, 34–36] and have been the target for improvement/design for efficient biodegradation [37, 38].
State variables of a glucose unit a
covered by endo
covered by exo-R
covered by exo-N
locked by exo-R
locked by exo-N
The first parameter, P 1 informs whether the glucose unit belongs to the cellulose substrate (‘nonsoluble’) or is in the aqueous phase (‘soluble’). Initially, all glucose molecules are part of the cellulose surface (P 1 = 0), and as the simulation progresses, they become soluble (P 1 = 1) in the form of simple sugars: glucose, cellobiose or cellotriose. The length of the soluble oligomers can be easily modified.
The second parameter, P 2 specifies whether the glucose unit constitutes the nonreducing end (NE) or the reducing end (RE) of a glucan chain. Also, when a glycosidic bond in the middle of a chain is hydrolyzed, two new ends are created, one new NE and one new RE.
Parameter P 3 informs whether an endo-cellulase covers the glucose unit. Similarly, the values of parameters P 4 and P 5 indicate whether an exo-cellulase hydrolyzing from the reducing end (exo-R) or an exo-cellulase hydrolyzing from the nonreducing end (exo-N) covers the glucose unit. At any time, only one cellulase is allowed to cover a specific glucose, during which glycosidic and hydrogen bonds may be cleaved or broken. A glucose unit, which is not covered by any cellulases, may become locked by a processive cellulase, exo-R or exo-N. This is specified by parameters P 6 and P 7 , respectively. A locked glucose unit only facilitates the binding of the locking exo-R (exo-N); it does not constitute an available binding site for any other cellulases, nor can its glycosidic or hydrogen bonds be cleaved or broken, until the locking cellulase binds directly to it.
Rate constants characterizing endo-cellulases
adsorption rate constant
desorption rate constant
a higher desorption rate constant than k off
rate constant for breaking a single hydrogen bond
rate constant for hydrolyzing a glycosidic bond
Endo-cellulases, once adsorbed to the cellulose surface, can break inter-chain hydrogen bonds, hydrolyze glycosidic bonds and desorb from the substrate into solution. Each of these chemical reactions is essentially a Poisson process that takes place at a specific, constant rate defined by the propensity function of that reaction .
Rate constants and relevant parameters for characterizing and exo-cellulase
adsorption rate constant
desorption rate constant
a higher desorption rate constant than k off
probability for the exo-cellulase to desorb from cellulose
time for an exo-cellulase to break a single inter-chain hydrogen bond
time for an exo-cellulase to hydrolyze a glycosidic bond and slide one cellobiose unit along the glucan chain
time for an exo-cellulase to remain at a certain location
Exo-cellulases, once adsorbed to the cellulose surface, can break hydrogen bonds, hydrolyze glycosidic bonds, slide along a chain, and desorb from the cellulose into solution. For simplicity, in the following description we restrict ourselves to actions carried out by an exo-R cellulase. The ones carried out by an exo-N cellulase are essentially the same.
During adsorption, the inter-chain hydrogen bonds of the glucose units beneath the cellulase are instantly broken (Figure 5a). Similar to an endo-cellulase, desorption of an exo-cellulase from the cellulose sheet can take place at any time. If any of the glycosidic bonds between the glucose units covered by the cellulase in the middle glucan chain is already hydrolyzed, the cellulase desorbs from the surface within an exponentially distributed time interval with rate parameter k offfast . This assumption was built into the model in order to account for the continuity of the chain entering the tunnel of CBH I. The frequency of cellulase dissociations from the substrate is set by a probability α. Results of high speed AFM measurements  indicate that once the cellulase is adsorbed to a free end, it continues to process the chain until it reaches the end of the chain. In this light, α is usually set to a relatively low value. The ‘Spatial Conditions’ branching point in Figure 4 checks whether the cellulase can further hydrolyze the glucan chain. When another cellulase obstructs its diffusion, it may stay in place for a time of length t stay . It may also desorb within an exponentially distributed time with parameter k off , if glucose units are missing from the middle chain or if it has reached the nonreducing end of the chain (and similarly for the exo-N).
The processive movement of exo-cellulases is modeled as a two-step process: (i) the cellulase hydrolyzes the glycosidic bond positioned at the active site, followed by (ii) its sliding along the processed chain by one cellobiose unit, instantly breaking inter-chain hydrogen bonds  (Figure 5b). The active site of an exo-R is considered to be the glycosidic bond between the second and third glucose unit from the reducing end of the middle, processed glucan chain. If the two glucose units in front of the cellulase belong to the substrate, and none of the six glucose units in front of the cellulase (two consecutive ones in three chains) are occupied by any other cellulases, these six monomers become locked, implying that they could be covered by the exo-cellulase. Locked glucose units are not considered available binding sites for endo– or exo–cellulases. The processing time of exo–cellulases is calculated as t process = t move + t hbbreak *nrHb, where t move specifies the time during which the cellulase hydrolyzes the glycosidic bond at the active site and slides along the middle glucan chain by one cellobiose unit. Here nrHb is the number of present hydrogen bonds between the six locked monomers. The same strategy is employed for exo-N except that the processivity is along the opposite direction towards the left-hand side. Furthermore, an exo-cellulase is never allowed to productively bind to chains in the middle of the cellulose surface. The exo-R cellulases are only allowed to productively bind to a free reducing end, while the exo-N cellulases are only allowed to productively bind to a free non-reducing end.
Chemical reactions that involve endo-cellulases (adsorption, desorption) or are catalyzed by endo-cellulases (inter-chain hydrogen bond breaking, hydrolysis of glycosidic bond) are modeled as Poisson processes. This choice has been motivated by the knowledge of the specific reaction rate constants.
A different modeling approach is taken in the case of exo-cellulases. While adsorption of exo-cellulases to the cellulose substrate is modeled as a Poisson process, each adsorbing event is combined with an instantaneous structural transformation of the surface as hydrogen bonds below the cellulase are assumed to break at the very moment of adsorption. The processivity of exo-cellulases is defined through specific rules: when certain spatial conditions are met, the cellulase slides along the glucan chain by one cellobiose unit leaving behind a soluble cellobiose and instantly changing the cellulose surface (breaking intact hydrogen bonds) beneath itself. We assume that each exo-cellulase has the same, constant processing velocity so the only difference in processing times of one exo-cellulase from another originates in the number of intact hydrogen bonds covered by respective cellulases. Finally, desorption of exo-cellulases is also modeled as a Poisson process with a constant rate parameter.
Each chemical reaction is modeled as a discrete event occurring instantaneously while the state of the system remains unchanged between two consecutive events. The events associated with each of the adsorbed cellulases are stored in a priority queue, here referred to as the master queue of events (Figure 6). They are sorted by the simulated time at which they should occur. The simulation runs until the surface degrades to a specific degradation threshold, which is usually set to 100%, or the point at which the substrate is not able to adsorb more enzyme particles.
Input parameters used in simulations unless otherwise stated
Adsorption rate constant
Desorption rate constant
Higher desorption rate constant
Glycosidic bond hydrolysisc
Hydrogen bond breakingd
Table 4 also includes the molecular weights of EG I and CBH I used in the simulation. The molecular weights of glucose (μ = 180.15588 g/mol) and anhydroglucose (μ = 162.1406 g/mol) molecules are essential in the calculation of the soluble sugar concentration. The cellulose substrate is composed of five glucan chains, each of them having 4000–5000 glucose monomers. The area of a cellobiose unit  was set to A G2 = 5.512 × 10 -19 m 2 . The initial (molar) enzyme concentration was calculated from the number of cellulases present in the system and the volume of the system, V = Sd, where S is the cellulose surface area in m2 and d is of the order of μm. In most of the runs the initial enzyme concentration was set to be [E] 0 = 2 μM. In the following section, unless stated otherwise, we used the parameter values listed in Table 4 and the above enumerated initial conditions.
The output consists of the time evolution of the concentrations of glucose, cellobiose and cellotriose present in the aqueous phase, the adsorption density, and the number of available binding sites per gram of cellulose. Results were averaged over 10 replicas of the system. The conversion is quantified by counting the number of all the soluble glucose molecules including those in cellobiose and cellotriose and dividing that by the initial number of glucose molecules.
Even though we obtain reasonable qualitative agreement, it cannot be simply improved by just using a few experimentally determined kinetic parameters. Since our model is a detailed one, there are several other parameters that need to be optimized as well to get quantitative correspondence. For example, when we use the experimentally determined rates for adsorption (kon = 4.2 * 104 s-1- based on association constant Ka = 1.4*106 M-1 s-1 and koff) and desorption (koff = 0.03 M-1 s-1) for EG-I from cellulose crystals, we observe that the enzyme hydrolyzes the cellulose crystals completely in under 15 hours. Other physical reasons may also contribute. One of the reasons for the fast processing rate observed in the model is because we only process two-dimensional crystals of cellulose, while cellulose crystals are three-dimensional in nature. In three-dimensional crystals, there are multiple layers of cellulose chains in the crystal, and not all of the glycosidic bonds are available as substrate for the enzyme to process from the beginning. Rather the inner layers of the crystal are available as substrate only after the layers above them are partially processed. Another reason could be that the rate constant for decrystallization (or chain separation) of cellulose crystals by the enzymes is based on the theoretical timescales for breaking of hydrogen bonds in an aqueous environment, which is quite rapid.
The relatively constant processing time of the exo-cellulases is the result of using a simple coarse-grained description of processing time for exo-cellaloses since no rates have been measured for specific events associated with processivity. Values in these calculations are set such that kon, koff >> 1/tmove for exo-cellulases. Thus the binding of the enzyme to the substrate is at equilibrium. These parameters can be easily adjusted to match up with any forthcoming experimental observations. In addition, the total concentration of the substrate (reducing ends for exo-R or non-reducing ends for exo-N) does not change with time until the whole chain is processed. This ensures that the concentration of enzyme-bound substrate remains nearly constant with time for the exo-cellulases resulting in a nearly constant rate of processing until the end. We expect these effects to be reduced in three-dimensional crystals of cellulose in which multiple layers of cellulose chains have to be processed by the exo–cellulases as the chains get hydrolysed in a staggered fashion in these crystals. It has also been reported  that cellulose-binding modules bind to insoluble non-crystalline cellulose with a 10-20-fold greater affinity than to cello-oligosaccharides and/or soluble polysaccharides. Future expansion of this model will incorporate a non-constant adsorption rate of enzymes that would depend on the length of the cello-oligosaccharides; this will bring further complexity to the model. In addition, incorporation of stochasticity in the processing of the cellulose chain by the exo-cellulases and better estimates for rates of decrystallization of the cellulose crystal could lead to better agreement with experimental hydrolysis rates.
The increase of free chain ends produced by endo-cellulases is the primary source of the endo-exo synergy observed in the model. Figure 12b illustrates how this effect contributes to a large increase in the percentage of adsorbed exo-cellulases. This percentage is constant when exo-cellulases degrade the substrate alone, while in the presence of endo-cellulases it grows to higher values, contributing to a fast and efficient degradation of the substrate.
In an effort to complement both all-atom molecular dynamics and coarse-grained simulation tools, we have developed a an agent-based the dynamical, functional model capturing the surface chemical reaction of cellulose hydrolysis by enzymes at the molecular level. This model accounts for heterogeneous enzymatic hydrolysis reactions occurring on the substrate surface (a reaction taking place in dimensions less than three), and incorporates key factors controlling it that are different from those in an aqueous environment. The catalysis process is broken down into distinct parts related to different kinetic events carried out by individual particles. These events are essentially chemical reactions taking place on the surface of cellulose (adsorption, breaking inter-chain hydrogen bonds, cleaving glycosidic bonds, desorption) and constitute the main elements of this model. Reactions are monitored by following and updating the state (based on a set of predefined rules) of each individual particle in the system. Simulation results showed good qualitative agreement with experimental data.
The agreement with experimental data can be improved by obtaining better experimental estimates of the parameters in Table 4 and by extending the current model to three dimensions. Initial experiments that can greatly benefit the model are those that probe the kinetics of different steps for the individual domains of the cellulases (Carbohydrate Binding Domain and Catalytic Domain) separately. These experiments need to quantify the binding affinities, k on and k off . Then similar measurements on the entire cellulases for binding of the same substrate will help to verify the role of individual domains and provide a measure of productive and non-productive binding. These measurements need to be carried out for pure cellulose substrates of different shapes morphology, degree of polymerization (DP) and partially digested states.
Finally, our model only simulates the degradation of a single cellulose crystal layer, a feature that should be extended to capture the degradation of a whole cellulose crystal. The major effect from a three-dimensional model is expected to be that the substrate would not be completely accessible at the same time for the cellulases to digest. Also, such a three-dimensional model can capture the possibility that floating sheets of detached substrate may slow cellulases from reaching a larger surface where more efficient digestion is possible. The current model does not account for surface diffusion, which is likely to be important based on results reported by Jervis and colleagues , who showed that diffusion does not limit enzyme activity. Fortunately, these deficiencies are not of a fundamental nature because our model is easily extendable and can incorporate them as well as additional properties of various cellulase systems on different types of cellulose surfaces. Importantly, this approach could be broadened to other classes of cellulases and even to cellulosomes as additional experimental data becomes available. For this reason we believe that this model constitutes a significant contribution to the ability to simulate the complicated reactions involved in cellulose degradation.
aEach row corresponds to one of the seven parameters characterizing one monomer while the columns represent numerical values the parameters can take. Each entry of the table denotes a distinct condition of a glucose unit. For a detailed explanation, please see text. bSee reference . cSee reference . dSee reference .
Exocellulase that processes from nonreducing end
Exocellulase that processes from reducing end
The authors gratefully acknowledge the financial support from LANL LDRD program and CNLS. Also we appreciate computational support from LANL Institutional Computing. We thank Jennifer Macke for critical reading of the manuscript.
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.