Stochastic molecular model of enzymatic hydrolysis of cellulose for ethanol production
© Kumar and Murthy; licensee BioMed Central Ltd. 2013
Received: 17 October 2012
Accepted: 29 April 2013
Published: 2 May 2013
Skip to main content
© Kumar and Murthy; licensee BioMed Central Ltd. 2013
Received: 17 October 2012
Accepted: 29 April 2013
Published: 2 May 2013
During cellulosic ethanol production, cellulose hydrolysis is achieved by synergistic action of cellulase enzyme complex consisting of multiple enzymes with different mode of actions. Enzymatic hydrolysis of cellulose is one of the bottlenecks in the commercialization of the process due to low hydrolysis rates and high cost of enzymes. A robust hydrolysis model that can predict hydrolysis profile under various scenarios can act as an important forecasting tool to improve the hydrolysis process. However, multiple factors affecting hydrolysis: cellulose structure and complex enzyme-substrate interactions during hydrolysis make it diffucult to develop mathematical kinetic models that can simulate hydrolysis in presence of multiple enzymes with high fidelity. In this study, a comprehensive hydrolysis model based on stochastic molecular modeling approch in which each hydrolysis event is translated into a discrete event is presented. The model captures the structural features of cellulose, enzyme properties (mode of actions, synergism, inhibition), and most importantly dynamic morphological changes in the substrate that directly affect the enzyme-substrate interactions during hydrolysis.
Cellulose was modeled as a group of microfibrils consisting of elementary fibrils bundles, where each elementary fibril was represented as a three dimensional matrix of glucose molecules. Hydrolysis of cellulose was simulated based on Monte Carlo simulation technique. Cellulose hydrolysis results predicted by model simulations agree well with the experimental data from literature. Coefficients of determination for model predictions and experimental values were in the range of 0.75 to 0.96 for Avicel hydrolysis by CBH I action. Model was able to simulate the synergistic action of multiple enzymes during hydrolysis. The model simulations captured the important experimental observations: effect of structural properties, enzyme inhibition and enzyme loadings on the hydrolysis and degree of synergism among enzymes.
The model was effective in capturing the dynamic behavior of cellulose hydrolysis during action of individual as well as multiple cellulases. Simulations were in qualitative and quantitative agreement with experimental data. Several experimentally observed phenomena were simulated without the need for any additional assumptions or parameter changes and confirmed the validity of using the stochastic molecular modeling approach to quantitatively and qualitatively describe the cellulose hydrolysis.
Lignocellulosic biomass is a complex matrix of three biopolymers: cellulose (20-50%), hemicellulose (15-35%) and lignin (5-30%) [1–3]. Conversion of renewable lignocellulosic biomass into biofuels is at the heart of advanced biofuels production. Two important biochemical approaches to accomplish this involve enzymes and/or cellulolytic microorganisms. First approach involves pretreatment of lignocellulosic feedstocks [4–6], followed by hydrolysis and subsequent/simultaneous fermentation by yeasts [7, 8]. Second approach involves us of cellulolytic microbes with consolidated bioprocessing (CBP) capabilities that hydrolyze cellulose without external enzyme addition . The first approach represents state of the art in cellulosic ethanol and is the system of choice for near term (5–10 years) commercialization that is being tested at pilot and industrial scales.
Hydrolysis of cellulose and hemicellulose into sugar monomers is a critical step in the biochemical conversion of cellulose into ethanol. Cellulose is the most abundant biopolymer on earth with about 100 billion tons produced by terrestrial plants every year . Cellulose is made of glucose units linked together by β-1,4 glycosidic bonds. Due to extensive hydrogen bonding, cellulose chains form a recalcitrant crystalline structure that has a half-life of several million years at neutral pH and ambient temperatures in the absence of enzymes .
However, abundance of cellulose in diverse natural substrates has led to evolution of diverse cellulolytic organisms capable of degrading cellulose [11–14]. There are three main mechanisms of cellulose degradation found among most cellulolytic organisms . First mechanism which is used by many aerobic organisms is the extracellular secretion of free cellulases. Trichoderma reesei (T. reesei) is a well-studied fungus that employs this mechanism. Second mechanism employed by most anaerobic bacteria such as Clostridium thermocellum is the use of large multienzyme complexes called cellulosomes. Third strategy, comparable in effectiveness to above strategies, employed by organisms such as F. succinogens (anaerobic rumen bacteria) and Cytophaga hutchinsonii (aerobic soil bacterium) is presently not well understood. These organisms do not produce cellulosomes or secrete extracellular enzymes and remain tightly bound to the cellulose substrate .
Enzymatic hydrolysis process can be divided into two stages: In the first stage, long chains are hydrolyzed to form soluble oligomers, and soluble oligomers are in turn hydrolyzed to sugar monomers during the second stage of hydrolysis. First stage of hydrolysis is considered as a rate limiting step in the hydrolysis process [1, 20]. Several experimental studies have concluded that physical properties of cellulose such as crystallinity, degree of polymerization (DP) and accessible surface area, are some of the major factors responsible for controlling hydrolysis rate due to the effect on enzyme binding and substrate accessibility to the cellulase enzymes [23–25]. Crystallinity is a key factor affecting the hydrolysis of cellulose as the glycosidic bonds in crystalline regions are difficult to hydrolyze compared to those in the amorphous regions [23, 25].
Enzymatic cellulose hydrolysis process is one of the bottlenecks and key cost center (enzymes costs up to $1/gal ethanol) in the commercialization of process of cellulosic ethanol production [22, 26–30]. However, there is potential for cost reduction by improving the understanding of the process, testing enzymes and various substrates under different conditions to determine optimum hydrolysis conditions. Since conducting hydrolysis experiments is time consuming and labor intensive, a comprehensive hydrolysis model that can predict accurate hydrolysis profile of cellulose under various scenarios (substrate, enzyme, hydrolysis conditions) can be used as an important forecasting tool to understand and improve the hydrolysis process [22, 26].
There are three principal approaches to model cellulose hydrolysis by non-complexed systems [20, 27]. First approach is to fit experimental data to linear/nonlinear regression models which can be simple to construct but require large sets of experimental data.
Second approach involves formulation of mechanistic models which attempt to model some of the underlying phenomenon with simplifying assumptions [20, 22, 27, 31–37]. Rate expressions are generally described using Michaelis-Menten type enzyme kinetics with/without incorporating the effects of enzyme adsorption, temperature, pH, substrate and product inhibition. The model structure results in a set of ordinary differential equations (ODE) and model parameters are often determined by fitting model predictions to the experimental data. Substrate and product inhibition, lignin inhibition effect, cellulose accessibility and reactivity effects have been considered in some of these models. These models are reasonably accurate in predicting the experimental trends of sugar production and are currently most widely used in literature to predict cellulose hydrolysis [20, 27]. While a more general set of conditions can be simulated by incorporating additional model terms this often results in the loss of physical significance of the model terms and leads to overparameterization issues. Due to these limitations, a new parameter set must be identified whenever the substrate, pretreatment conditions, enzymes and/or process conditions are changed. Another limitation of this modeling approach is that synergistic interaction of various enzymes is difficult to determine and model accurately. Some of the models in literature do account for endo/exo and exo/exo synergistic interactions . However due to the reliance on lumped terms and inadequate detail in describing cellulase-cellulose interactions, their predictive abilities are limited .
Most of these kinetic models can not consider dynamic change in properties of substrate (accessibility of enzymes) that directly affect the enzyme-substrate interactions and rate of hydrolysis in detail. At most some of the models have considered changes in accessibility based on ratio of active/less active cellulose, amorphous/crystalline cellulose and surface area. However it is very difficult to consider accessibility changes for individual enzymes depending upon their mode of action.
Changes in substrate property, enzyme characteristics can be incorporated without the need for additional experimentation.
Concentrations of many oligomers can be tracked without increasing the complexity of model.
Chain distribution (number of chains with different chain lengths) during hydrolysis can be easily determined using SMM approach, which can provide better understanding of action mechanism and behavior of enzymes during hydrolysis.
Changes in structural properties such as number of chain ends, average degree of polymerization and crystallinity over time can be followed during hydrolysis.
The basic requirement of this approach is detailed and accurate description of substrate properties and enzymes action in the model. Although SMM approach is computationally intensive, it is more realistic and can easily incorporate the changes in conditions (substrate, enzyme or reaction conditions). However, potential of this technique has not been explored in detail for cellulose hydrolysis.
The first reported model for hydrolysis of insoluble polysaccharides using stochastic modeling approach was developed by Fenske et al. . This was a limited SMM model as it did not capture the actual structural properties of cellulose (e.g. crystallinity, DP, fibril structure) and multi enzyme dynamics. Hydrolysis was performed on a two dimensional matrix representing a single surface of cellulose with short chain length of 20. It was a theoretical study and results were not validated with experimental data. More recently, similar approach was used by Asztalos et al.  to model cellulose hydrolysis which had reasonable accuracy in predicting the hydrolysis trends for endoglucanse and CBH enzymes. Dynamic enzyme-substrate interactions were captured to some extent in the model. However the model did not include some important structural features of cellulose and had a limited usability. It was a two dimensional model in which all glucose chains were accessible to all enzymes, which in not the case in actual process. Degree of crystallinity was not considered in the cellulose structure, and consequently, activity difference of enzymes in amorphous and crystalline regions found to be significant by many researchers was not incorporated into the model [1, 20, 27]. Inhibition by cellobiose or glucose was not accounted into the model, which is very important parameter during cellulose hydrolysis. Therefore, despite these relatively recent advances, there is no SMM model to date that considers complex enzyme-substrate interactions. The objective of this study was to develop a detailed SMM model that can predict hydrolysis profiles of cellulose with high accuracy by capturing the complexities of cellulose structure and hydrolysis mechanism. The aim was to develop a general hydrolysis model that can be used for various conditions (different substrates, enzymes, hydrolysis conditions) considering structural properties of feedstock (crystallinity, degree of polymerization, accessibility), enzyme properties (mode of actions, synergism and inhibition) and most importantly dynamic changes in these properties during hydrolysis.
Correlation between experimental and model simulation data for Avicel hydrolysis by CBH I action at different conditions
25 g/L Avicel
50 g/L Avicel
Enzyme loading (mg/g glucans)
Enzyme loading (mg/g glucans)
There are several possible reasons for this divergence in the later stages of hydrolysis. Commercial Avicel substrates obtained from different sources vary in their physical properties, whereas we used average values of DP and crystallinity index for structure simulations. Another important contributor to the difference could be differences in enzyme activities. A great variability has been found among of activities of pure enzymes, depending upon degree of purification . Activity used in the model simulations were used from other literature studies, which might be different from the activity of CBH I purified in these experiments.
For all three enzymes the hydrolysis rates decrease rapidly after few initial hours of hydrolysis and then become nearly constant (Figure 6). The decrease in hydrolysis rate with cellulose conversion is a widely observed phenomenon that is expected due to morphological changes in the structure of cellulose such as decrease in glucose chains on microfibril surface. These morphological changes limit the accessibility of cellulase enzymes to glucose chains and results in rapid decline in hydrolysis rate [20, 36].
During hydrolysis of cellulose by CBH I and CBH II, cellobiose is the only major product formed, which is a strong inhibitor of these enzymes [1, 20, 27] and limits the hydrolysis rate. Concentration of cellotriose, cellobiose and glucose were higher than that of cellotetrose, cellopentose and cellohexose after 48 h hydrolysis of Avicel by EGI (Figure 7). Concentrations of cellotetrose, cellopentose and cellohexose increase during initial few hours (3–4 h) of hydrolysis and then decrease in the latter stages of hydrolysis. This trend is expected since during initial phase of hydrolysis surface glucose chains are easily accessible and endoglucanases acts randomly on them to producing short chains. However, with the progress of hydrolysis accessibility of enzymes to the glucose chains decreases. The enzymes start acting on soluble sugars and reduce their overall concentrations after an initial increase. Concentration of cellobiose and cellotriose do not decrease as EG I was assumed to act only on oligomers with DP>4. In case of hydrolysis by CBH I and CBH II, all soluble sugars except cellobiose were produced in negligible amounts (less than 0.01 g/L, data not reported).
Key Properties of simulated model cellulosic substrates
Degree of polymerization (DP)
Crystallinity index (CrI)
Number of microfibrils simulated (Total glucose molecules simulated)
Bacterial Cellulose (BC)
It can be observed from Figure 8 that model simulations agreed with the expected results and captured the inverse relationship between substrate DP and hydrolysis of cellulose. Cellulose conversion after 48 h hydrolysis of high DP substrates: bacterial cellulose (DP 1800–2000) and cotton (DP 2500–3000), which have low availability of chain ends for CBH attack, was only 10.8 and 4.4% relative to that from low DP substrate, Avicel (DP 250–300) hydrolysis. Similar results have been reported in literature from both experimental as well as modeling studies [33, 49]. Although it appears that cellulose conversion was almost similar for Avicel and filter paper, rate of hydrolysis was lower in case of filter paper. After few hours of hydrolysis, cellobiose inhibition becomes dominant and limits the hydrolysis.
Endoglucanases activity is affected by degree of crystallinity. Regions that are highly crystalline are less susceptible to hydrolysis compared to amorphous regions because of low accessibility of enzymes in these regions . This behavior was also captured by model simulations as cellulose conversion by action of endoglucanases on cotton (highly crystalline cellulose, CrI 0.85-0.95) was found 42.8% lower than that of Avicel (semi-crystalline cellulose, CrI 0.5-0.6) (Figure A6, section A6 in the Additional file 1).
Cooperative action of different enzymes, known as synergism, is one of the most important phenomenon observed in cellulose degradation [1, 18, 20, 50]. Synergistic action among different enzymes (such as endo-exo synergism, exo-exo synergism, exo-β-glucosidase synergism) have been studied; of these end-exo synergism is a highly effective synergism that has been reported in many studies [18, 20, 48, 51].
Where, ∆Cmixed is cellulose conversion obtained from mixture of ‘n’ enzymes; ∆Ci is cellulose conversion obtained from individual action of ‘ith’ enzyme.
For additional confirmation of synergistic action, two additional substrates filter paper and bacterial cellulose (substrate properties in Table 2) were simulated under similar conditions (100 g/L cellulose, 10 mg/g individual enzyme). Degrees of synergism from combined action of EGI and CBH I were calculated as 1.8, 1.9, 4.14 and 4.99 for Avicel, filter paper, bacterial cellulose and cotton respectively. The values of DS obtained from model simulations are consistent with the reported values in literature [20, 33, 37, 48, 51]. The inverse relationship between DS and substrate DP was expected and has been reported in literature [18, 20, 52]. Zhang et al.  compiled DS values from various studies (Table V of ) and reported DS values of 1.3 to 2.2 for Avicel and 4.1 to 10 for cotton and bacterial cellulose from synergism of T.reesei enzymes. When CBH I acts alone, its accessibility to chain ends is very limited and cellulose conversion is very less. During simultaneous action of EG and CBH I and/or CBH II, EG action results in creation of additional chain ends for action of CBH I and/or CBH II, which results in more effective hydrolysis (high synergism).
Glucose was the main product of hydrolysis in presence of excess β-glucosidase which converts cellobiose to glucose (Figure 11). Although concentrations of all oligomers remained low throughout the hydrolysis, cellobiose concentration was highest among all oligomers. This trend was expected since cellobiose is the main product produced by the action of CBH I and CBH II while other oligomers are mainly produced by endoglucanases only. Similar results were obtained in terms of relative amount of oligomers from Avicel hydrolysis by same enzyme mixture in the absence of β-glucosidase. Cellobiose concentration was highest followed by cellotriose, which is also assumed to be not hydrolysable by either EG, CBH I or II. However cellulose conversion after 48 h of hydrolysis in absence of β-glucosidase was 33.4%, which was 57.5% less than that in presence of β-glucosidase (78.6%). This behavior was expected because as described earlier cellobiose produced by the action of CBH I and CBH II inhibits the activity of these enzymes and reduces overall cellulose conversion. Therefore for effective cellulose hydrolysis β-glucosidase enzyme is added to the cellulase preparation T. reesei. as the enzyme mixture produced by it has low levels of β-glucosidase . Hydrolysis simulations were also performed under same hydrolysis conditions (T. reesei. enzyme mixture in excess of β-glucosidase) for cotton cellulose and the results of oligomer formation was found similar as in the case of Avicel hydrolysis. However cellulose conversion at end of 48 h hydrolysis was lower by 30.8% compared to Avicel. This observation can be attributed to the high crystallinity and higher degree of polymerization of the cotton cellulose compared to Avicel.
As expected, hydrolysis rates were found to be increasing with increase in enzyme loadings (Figure 12). For enzyme loading at 10 mg/g glucans hydrolysis rate remained constant during initial period followed by a substantial decrease and nearly constant rates thereafter. This behavior can be attributed to rapid decrease in the number of available surface glucose molecules with time. Number of accessible sites for endoglucanases decreases during hydrolysis. Initial constant rate phase was not observed at high enzyme loading (30 mg/g glucans) as all the accessible surface bonds were hydrolyzed at a very fast rate.
Cellulose was modeled as group of microfibrils made of bundles of elementary fibrils. Model incorporated all important factors of cellulose hydrolysis such as structural features, morphological changes in cellulose during hydrolysis and effect of morphological changes on the hydrolysis. Model was accurate in predicting the cellulose hydrolysis profiles from experimental studies as well as followed other trends reported in literature studies. The major advantage of this model is that any change in substrate property, enzyme characteristics can be incorporated into model after performing independent characterization experiments. Other than predicting hydrolysis profiles, the model can be used to gain insight in the dynamic formation and breakdown of cello-oligomers during hydrolysis.
A SMM model that considers structural features of cellulose and complex enzyme-substrate interactions was developed. The SMM of cellulose hydrolysis process was systematically organized into three steps, consisting 1) In-silico construction of a representative cellulose polymer model, 2) characterization of the celluolytic enzymes in the mixture and 3) modeling the enzymatic hydrolysis.
Parameters of glucose molecules in an elementary fibril (corresponds to Figure 14 )
These parameters were used to determine accessibility of enzymes depending upon their action pattern and directly affect the hydrolysis process. For example it is important to consider bond location (amorphous or crystalline region) since activity of endoglucanases enzymes is relatively low in crystalline region [18, 20]. Similarly not all endoglucanases produce glucose and the parameter “Distance_NR” helps modeling this case. EF surfaces are less accessible to enzymes compared to MF surfaces. The properties “EF Surface” (Table 3) and “MF Surface” help in identifying the bond location and simulating enzymatic action appropriately. While assigning these parameters, it was ensured that cellulose model captures all major structural properties of actual cellulose. For instance, the crystalline and amorphous regions were considered to occur as bands in elementary fibril. While the length of the crystalline region was set to 200 glucose molecule, the location was randomly chosen along the length of the elementary fibril such that the amorphous regions do not occur at the ends of the elementary fibril (Figure 1).
Mode of action for enzymatic hydrolysis of cellulose by enzymes from various classes
Mode of action
Endo-cellulase (EG, Non-processive with CBM)
• Carbohydrate binding module (CBM) binds to surface chains randomly and breaks bonds at catalytic domain (CD)
Endo-cellulase (EG, Non-processive without CBM)
• Breaks bonds on surface chains in a random pattern
Endo-cellulase (EG, Processive with CBM)
• CBM binds to surface chains randomly and break bond at CD
• Enzyme moves along the chain (towards non-reducing end or reducing end) and cuts every alternate bond releasing cellobiose until a minimum chain length is achieved
Exo-cellulase (CBH I, Processive)
• CBM attaches from reducing end only on surface chains and pulls the chain towards CD
• Chain passes through CD (tunnel like shape) and every alternate bond is broken to produce cellobiose
• Enzyme moves along the chain (towards non-reducing end) and cut every alternate bond until a minimum chain length is reached
Exo-cellulase (CBH I, Non-Processive)
• Attack from reducing end on surface chains and cuts every alternate bond to produce cellobiose
Exo-cellulase (CBH II, Processive)
• CBM attaches from non-reducing end on surface chains and pulls the chain towards CD
• Chain passes through the tunnel shaped CD and every alternate bond is broken to produce cellobiose
• Enzyme moves along the chain (towards reducing end) and cuts every alternate bond until a minimum chain length is reached
Exo-cellulase (CBH II, Non-Processive)
• Attack from non-reducing end on surface chains and breaks alternate bonds to produce cellobiose
• Acts on cellobiose and soluble oligomers (DP ≤ 6) and produce glucose by breaking bond
There can be multiple enzymes in each class identified from different organisms; however they will differ only in their hydrolytic efficiency (activity). Therefore modeling the action of these eight classes of enzymes will sufficiently capture the dynamics individual cellulases. The important aspect of this classification is that any new enzyme can be used in model after characterizing its action pattern and placing into one/multiple groups with its specific activity /activities, without making any change to the current model.
Factors affecting enzyme action can be classified as intrinsic and extrinsic factors. Intrinsic factors do not depend on substrate characteristics while extrinsic factors depend solely on substrate characteristics. Extrinsic factors (crystallinity, accessibility and average DP) are accounted in the cellulose model described in Step 1. Intrinsic characteristics such as activity and stability of the enzyme are dependent on the reactor pH and temperature. The enzyme loading (amount of enzyme/g substrate) and activity information was transformed into number of bonds hydrolyzed per unit time for each enzyme.
Where, 'Ei' is amount of 'ith' enzyme used per unit of cellulose during hydrolysis, mg/g cellulose. 'Ui' is activity of 'ith' enzyme, IU/mg enzyme. One international unit (IU) of enzyme can liberate one micromole (6.023 * 1017 molecules) of product (i.e. one micromole of bonds are broken) per minute under standard conditions . 'Gsim' is the number of glucose molecules simulated in the model. “162” is the average molecular weight of anhydrous glucose molecule in a cellulose molecule (joined by β-1-4 linkages in long chains with DP from 100 to 20000). 'Si' is stability of 'ith' enzyme under current hydrolysis conditions (temperature and pH) that will affect the activity. Value of “Si” can be calculated for any enzyme using empirical equations developed based on experiments to incorporate the effect of temperature and pH on enzyme activity. For example, Arrhenius rate relationships are commonly used to consider the effect of temperature . Value of 'Si' is a real number between 0 and 1. In the current model, value of “Si” was used as 1, indicating that the enzyme was considered to be operating at its optimum pH and temperature.
Where, fj is fraction of hydrolysable bonds for 'jth' microfibril and is dependent upon mode of action of enzymes.
Where, Nsj is number of glucose molecules on the surface of elementary fibrils in jth microfibril and 'n' is number of microfibrils simulated.
Where, NRj is number of reducing ends in jth microfibril.
Where, NNRj is number of non-reducing end in jth microfibril.
Where, NSLj is number of soluble oligomers produced from jth microfibril.
Important conditions simulated in the model for enzymatic hydrolysis of cellulose by enzymes from endo, exo and β-glucosidase classes
Endo-cellulase (EG, Non-processive action pattern by an enzyme that contains a CBM)
• Lower probability of breaking bond in crystalline region compared to those in amorphous regions
• Probability of binding to the glucose chain and hydrolysis of bonds located on EF surface but not on MF surface is 25% lower compared to the glucose chains on the both EF and MF surfaces
• A minimum number of glucose molecules, depending upon size of enzyme, are required to be unblocked (not blocked by other enzymes) and on the EF surface to enable binding (please see section A3 of Additional file 1)
• Enzymes do not act on oligomers of chain length less than five. The enzymes can hydrolyze bonds to produce glucose with a defined 90% probability. Inhibition occurs due to cellobiose and glucose [1, 54]
Exo-cellulase (CBH I, Processive action pattern by an enzyme that may/may not have a CBM.)
• Probability of binding to the glucose chain and hydrolysis of bonds located on EF surface but not on MF surface is 25% lower compared to the glucose chains on the both EF and MF surfaces
• These enzymes have higher probability of hydrolysis with shorter glucose chains compared to longer chains . A preferred chain length of 300 glucose molecules was assumed
• Enzyme may get desorbed from the glucose chains at any time. Probability of desorption in amorphous region (5.0%) was set higher than that in crystalline region (2.5%)
• Enzymes have lower probability of bond hydrolysis in the case of soluble oligomers and do not act on oligomers with DP<5
• Acts only on soluble chains (DP<6)
• Probability of action on cellobiose is higher than other soluble oligosaccharides
• Inhibited by glucose
In addition to the mode of specific action pattern of enzymes (Endo/Exo-enzymes; processive or non-processive) probability of hydrolysis was set to be dependent on enzyme characteristics such as the bond location (crystalline or amorphous region) and soluble oligomer chain length. For example, during the action by EG enzymes, the probability of β-1,4 bond hydrolysis located in crystalline and amorphous regions are different. Similarly, the probability of bond hydrolysis in soluble chains with lower DP was set lower compared to the probability of hydrolysis for chains with higher DP. These probabilities were compared with a number randomly selected from a uniform distribution. The bond at a particular location was hydrolyzed only if the random number (generated from a uniform pseudo-random number distribution 0.0-1.0) was greater than the probability of hydrolysis set using the constraints described above. For example a bond located in the crystalline region could be broken (hydrolysis event occurred) only if the probability of the hydrolysis (as indicated by the pseudo-random number from the uniform distribution 0.0-1.0) was higher than the set value for probability of crystalline bond hydrolysis by EG enzymes. Similarly EG enzyme will hydrolyze a bond next to the reducing end or the bond prior to the non-reducing end (hydrolysis of this bond will produce glucose) only if its associated probability is higher than pseudo-random number.
Inhibition effect of sugars produced was captured in the model using the same approach. Activities of both CBH and endoglucanases are inhibited by cellobiose and glucose. Inhibition properties for cellobiose (Ninhib _G2) and glucose (Ninhib _G) were assigned to each class of enzyme. During hydrolysis simulations, whenever bond location determined using random number generated was corresponding to glucose or cellobiose, counter was incremented by corresponding “Ninhib” for that enzyme. As the hydrolysis progresses, probability of encountering glucose or cellobiose increases, which results in less number of bonds hydrolysis (as counter is incremented by “Ninhib”) and inhibition effect is captured. Extent of inhibition can be controlled by adjusting these properties. For example, cellobiose is strong inhibitor compared to glucose for CBH [1, 19, 27, 50, 54]; therefore value of “Ninhib _G2” was set higher compared to “Ninhib _G” for cellobiose in the model.
Properties of all glucose molecules (mentioned in Table 3 and section A1 of Additional file 1) in the chain, where bond is broken are reassigned after productive attack (hydrolysis event occurs). Initially, only the chains on the elementary fibril surface are accessible to enzyme attack. During hydrolysis, as the bonds are broken, glucose chains having DP less than seven become soluble (i.e. soluble oligomers are removed from the surface of the cellulose) and part of chain just beneath the soluble chain is exposed and becomes accessible to enzymes.
Where, 'Wsample' is weight of total solution during hydrolysis in grams. 'S' is fraction of solids in the solution (biomass or solid loading), dimensionless. 'C' is cellulose fraction of the solid, dimensionless. Explanation about these calculations and of other oligomers is provided in detail in the supplementary material (Section A4 of Additional file 1).
The hydrolysis rates were estimated by calculating the changes in cellulose concentration (DP 6+) over a specific time interval and was expressed as g cellulose/h.
Chain distribution provides the number of chains with different chain lengths in total cellulose, which can help in understanding of hydrolysis process. Therefore chain distribution profile was generated at specific time intervals. Additionally, crystallinity index and percentage of soluble molecules were also calculated at specific time intervals to monitor the changes in these properties.
The computer algorithms for implementation of above described hydrolysis model were written in C++ language. Random number generators were used in simulation of cellulose structure and hydrolysis process .
Cellulose structure was simulated using user defined range of structural parameters: DP, CrI, number of EFs in one microfibril. This approach allows simulation of any cellulose substrate such as Avicel, filter paper (with low DP and crystallinity), bacterial cellulose, cotton fibers (long DP and highly crystalline). Hydrolysis simulation was initialized based on the weight of solution (scale of hydrolysis), solid loading, cellulose content, total enzyme loading (mg protein/g cellulose), ratio of enzymes present (out of eight classes), temperature, pH and simulation time (hydrolysis duration). Enzyme activities can be obtained from manufacturer’s data or measured on different substrates (cellulose model substrates) using standard laboratory methods . Maximum value of enzyme activity is used in the model. Data provided by user on enzyme ratio was input to model through data files. The output from model included glucose concentrations, oligosaccharide concentrations (DP 2–6), chain distribution profile, crystallinity index profile, solubility profile and data sheets for each microfibril (illustrating major properties associated with glucose molecules) at various times during hydrolysis.
The approach of using various probabilities (such as endo enzyme acting on crystalline region, cutting the bond just before non-reducing end to produce glucose) provides flexibility to simulate various enzymes and to incorporate new enzymes in the model. For example, in general EG acts preferentially in amorphous region, however in future an enzyme may be discovered/developed with high activity in crystalline region. That enzyme can be easily incorporated into current model by changing the probability for bond cleavage in the crystalline regions. Similarly other probabilities can be altered depending upon enzyme characteristics.
Cellulose structure was simulated for various model cellulose substrates such as Avicel, bacterial cellulose (BC), cotton. These substrates differ in cellulose structure in terms of DP and degree of crystallinity. For instance, Avicel is low DP cellulose (about 300), whereas cotton has very high DP value (about 2000) . Values of properties used to simulate these model cellulosic substrates are listed in Table 2[20, 33]. Multiple microfibrils were simulated for each model substrate using these properties and assuming 4–6 EF in each row and column of microfibril. A minimum of 1,000,000 total glucose molecules were simulated to represent any cellulose substrate during hydrolysis simulations. These model substrates can vary somewhat in terms of DP, crystallinity index depending upon method of preparation (e.g. Avicel will not always have DP of exact 300), therefore using a range of properties and simulating multiple microfibrils captures this variability.
All model simulations were performed for four main classes of cellulase enzymes: Non-processive endoglucanases with CBM, processive CBH I with CBM and processive CBH II with CBM and β-glucosidase. These enzymes are similar to enzyme system of Trichoderma reesei (earlier known as T. viride), comprising endoglucanases EGI, two exoglucanases CBHI and CBHII and β-glucosidase. Cellulase system of T. reesei, a soft-rot fungi, has been extensively studied [1, 3, 18, 20]. Specific activities of enzymes from T. reesei were assumed as 0.4, 0.8, and 1.6 IU/mg of EG I, CBH I and CBH II respectively . Zhang and Lynd  calculated these specific activities values from various literature studies [66–68].
Data from model simulations were compared with two sets of experimental results from literature to validate the model.
Experimental details for the experimental data are detailed in Bezerra et al. [45–47]. Briefly, Avicel was hydrolyzed using CBHI (Cel7A) purified from Celluclast 1,5 L (a commercial cellulase preparation from T. reesei, provided by Novo Nordisk A/S, Denmark), at various enzyme: substrate loadings at 40 °C in 50 mM citrate buffer (pH 4.8) for 47 h.
For model validation, cellulose conversion results were extracted from Figure 3 of Medve et al. . Briefly, Avicel (M2331from Merck, Darmstadt, Germany) was hydrolyzed using CBH I and EG II (8 mg EGII/g Avicel and 8 mg CBHI/g Avicel), purified from a commercial cellulase, Celluclast from T. reesei (Novo Nordisk, Bagsvaerd, Denmark). For model simulations, activity of EG II enzymes used in experiments were assumed to be same as EG I enzyme (0.4 IU/mg protein).
1Graduate Student, Biological and Ecological Engineering, 116 Gilmore Hall, Oregon State University, Corvallis, OR-97331. 2Assistant Professor, Biological and Ecological Engineering, 116 Gilmore Hall, Oregon State University, Corvallis, OR-97331.
Cellulose conversion by action of mixture of ‘n’ enzymes
Cellulose conversion by individual action of ‘ith’ enzyme
Cellulose fraction of the solid
Carbohydrate binding module
Concentration of glucose (g/L)
Concentration of cellobiose (g/L)
Concentration of high DP molecules (g/L)
Degree of polymerization
Degree of synergism
Amount of 'ith' enzyme used per unit of cellulose during hydrolysis (mg/g cellulose)
Fraction of hydrolysable bonds for 'jth' microfibril
Number of glucose molecules in actual sample
Number of glucose molecules simulated in the model
International unit of enzyme activity
Number of microfibrils simulated
Number of maximum bonds broken by ith enzyme on jth microfibril
Maximum number bonds broken by an enzyme per minute
Increment in counter due to inhibition
Increment in counter due to cellobiose encounter/inhibition
Increment in counter due to glucose encounter/inhibition
Increment in counter when no binding occurs
Increment in counter in case of non-productive binding
Number of non-reducing end in jth microfibril
Number of reducing ends in jth microfibril
Number of glucose molecules on the surface of elementary fibrils in jth microfibril
Number of soluble oligomers produced from jth microfibril
Ordinary differential equations
Square of the Pearson product–moment correlation coefficient
Fraction of solids in the solution
Stability of 'ith' enzyme under current hydrolysis conditions
Stochastic molecular modeling
Activity of 'ith' enzyme (IU/mg enzyme)
Volume of solution (mL)
Weight of total solution during hydrolysis (g)
Authors gratefully acknowledge the support by National Science Foundation through NSF grant No. 1236349 from Energy for Sustainability program, CBET Division.
The authors are grateful to Dr. Rui M.F. Bezerra, Associate Professor, Department of Biology and Environment, University of Trás-os-Montes and Alto Douro, Apartado 1013, 5001–801 Vila Real, Portugal for providing the experimental data of Avicel hydrolysis by CBH I.
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.