Modeling and Simulation of Batch Sugarcane Alcoholic Fermentation Using the Metabolic Model

Abstract

The present work sought to implement a model different from the more traditional ones for the fermentation process of ethanol production by the action of the fungus Saccharomyces cerevisiae, using a relevant metabolic network based on the glycolytic Embden–Meyerhof–Parnas route, also called “EMP”. We developed two models to represent this phenomenon. In the first model, we used the simple and unbranched EMP route, with a constant concentration of microorganisms throughout the process and glucose as the whole substrate. We called this first model “SR”, regarding the Portuguese name “sem ramificações”, which means “no branches”. We developed the second model by simply adding some branches to the SR model. We called this model “CR”, regarding the Portuguese name “com ramificações”, which means “with branches”. Both models were implemented in MATLAB^TM software considering a constant temperature equal to 32 °C, similar to that practiced in sugar and ethanol plants, and a wide range of substrate concentrations, ranging from 30 to 100 g/L, and all the enzymes necessary for fermentation were already expressed in the cells so all the enzymes showed a constant concentration throughout the fermentation. The addition of common branches to the EMP route resulted in a considerable improvement in the results, especially predicting ethanol production closer to what we saw experimentally. Therefore, the results obtained are promising, making adjustments consistent with experimental data, meaning that all the models proposed here are a good basis for the development of future metabolic models of discontinuous fermentative processes.

1. Introduction

A wide range of fermented products are still produced by batch or sequential batch processes using easily accessible raw materials. In the case of ethanol production, biofuel is of interest to compete with fossil fuels, since its use brings a reduction in the emission of harmful gases, such as SO_x and NO_x. This is similar to the use of one of the most widely used processes in the production of this alcohol, the Melle-Boinot process, a form of fed-batch fermentation developed by Fermin Boinot, in the town of Melle, located in the Nouvelle-Aquitaine region of France, in the first half of the 20th century. Furthermore, in the production of ethanol, the yeast Saccharomyces cerevisiae is the most used microorganism, even if other microorganisms have a better potential in laboratory situations, such as the bacterium Zymomonas mobilis [1,2,3,4].

Due to the importance of this type of process, from the second half of the 20th century, a series of works seeking to understand these interactions began to be developed, especially after Gaden [5], who developed an empirical analysis (black box) of batch fermentation processes, verifying whether it was possible to associate the production of a metabolite with cell growth or not. This relatively simple approach had good results, despite ignoring a large part of the interactions that occur inside the cells. Because of this, numerous works were developed based on the work of Gaden and the kinetic model of Monod, making additions that sought to insert any metabolic knowledge.

Among the models that began seeking to insert these metabolic details, we can mention the model developed by Sonnleitner and Käppeli [6], developed based on a quasi-ideal Monod kinetic model. In such a model, Sonnleitner and Käppeli simplified the relevant metabolism into three stoichiometric relationships: two containing purely oxidative relationships, and one reductive. Furthermore, this model already considered the possibility of alcoholic fermentation not being fueled by a single substrate: it considered the possibility of multiple substrates, such as the possibility of ethanol being used by Saccharomyces cerevisiae as a carbon source. Despite the relevant results achieved in such a model, it still considered the great metabolic complexity in an extremely simplified way.

In addition to the model by Sonnleitner and Käppeli [6], several other relevant models for alcoholic fermentation were developed in a black box or with the addition of some metabolic details. Among such models, we can cite the work of Birol et al. [7], Kostov et al. [8], Sainz et al. [9], and Freitas et al. [10]. The works by Birol et al. [7] and Kostov et al. [8] have the merit of analyzing a large number of black box models for alcoholic fermentation, highlighting the Monod, Andrews, Naock, and Hinshelwood models. Fermentation processes conducted in batch and fed-batch form generate a considerably complex dynamic interaction between microorganisms and the environment around them, as studied by Freitas et al. [10] using genetic algorithms, differential evolution, and real-time dynamic optimization. These interactions have been generally studied and well analyzed in the form of black box analyses for some decades, as highlighted by Sainz et al. [9], but they fail to consider that each of the cells of a microorganism is capable of carrying out hundreds of reactions simultaneously. It was fair to the model proposed by Sainz et al. [9] that such interactions began to be considered. Each of these reactions is extremely well controlled by the action of substances that generate the expression of or inhibit very specific biocatalysts—enzymes [11,12].

Among the most important metabolic routes for numerous fermentative processes, we can mention the Embden–Meyerhof–Parnas route, or the EMP route. Completely described in the 1940s, the EMP pathway is present in a large number of organisms and describes a series of spontaneous reactions that convert glucose to pyruvate. The pyruvate generated in this sequence can be converted into a wide range of other substances, including ethanol, which can be produced under anaerobic conditions by reducing pyruvate [11,13,14,15].

Considering the importance of ethanol in the current economic and environmental scenario, as well as the development of models that are more faithful to reality, the present work sought to model a batch alcoholic fermentation process. We used kinetic equations derived from the EMP route of Saccharomyces cerevisiae, considering substrate concentrations lower than those used industrially to avoid the appearance of significant inhibitory effects. The production of ethanol still needs further clarification at the metabolic level. Thus, we propose two distinct models: the SR model, built only with the reactions present in the EMP route followed by alcoholic fermentation, and the CR model, which consists of using the SR model as a basis, adding the ramifications of trehalose, glycerol, succinate, and acetate. Furthermore, the development of models based on as much metabolic information as possible can benefit countless other industrial sectors besides sugar and alcohol: an understanding of the behavior of possible ramifications is something that would greatly benefit the food and pharmaceutical industries, for example. In the case of the food industry, this knowledge would make better control of the processes and the selection and development of new strains of microorganisms possible, making it possible to obtain larger amounts of specific substances, especially those involved in aromas and flavors.

2. Method and Experimental Procedures

Figure 1 below shows the sequence of reactions selected for the composition of the model. In short, there are the transport of the substrate from the extracellular environment to the intracellular environment, the EMP route, and alcoholic fermentation. Additionally, highlighted are the enzymes responsible for each of the reactions, in blue; products or resulting secondary reagents, in green; and the intracellular energy molecules ADP, ATP, NAD⁺, and NADH, in red. A list of symbols and meanings is available in the end of this paper.

In Figure 1, there is the existence of two substrates, glucose (GLC_ex) and fructose (FRU_ex), where the “ex” subscript indicates that these substances are found in an extracellular environment. Even though the yeast metabolizes these two sugars and both are originally present in the raw material, used at the ratio of 52.5% glucose and 47.5% fructose, in this work, we assumed that the only substrate available is glucose for reasons of simplification. Furthermore, there is only one metabolic product of interest, ethanol, EtOH.

In addition to these substances, there is a series of intracellular metabolites: ${\mathrm{GLC}}_{\mathrm{in}}$, $\mathrm{ADP}$, $\mathrm{ATP}$, $\mathrm{TRE}$, ${\mathrm{FRU}}_{\mathrm{in}}$, $\mathrm{G}6\mathrm{P}$, $\mathrm{F}6\mathrm{P}$, $\mathrm{F}16\mathrm{P}$, $\mathrm{DHAP}$,$\mathrm{GAP}$, $\mathrm{GLY}$, ${\mathrm{P}}_{\mathrm{i}}$, $\mathrm{BPG}$, $\mathrm{NAD}$, $\mathrm{NADH}$, $\mathrm{P}3\mathrm{G}$, $\mathrm{P}2\mathrm{G}$, $\mathrm{PEP}$, $\mathrm{PYR}$, $\mathrm{ACA}$, $\mathrm{SUC}$, $\mathrm{ACE}$, $\mathrm{T}6\mathrm{P}$, and $\mathrm{F}26\mathrm{P}$. Thus, starting from the sequence of metabolic reactions shown in Figure 1, we propose the following set of differential equations for such a model:

$$\frac{dGL{C}_{ex}}{dt}=-v_{HXT,G}$$

(1)

$$\frac{dEtOH}{dt}=v_{ADH}$$

(2)

$$\frac{dGL{C}_{in}}{dt}=v_{HXT,G}-v_{HK}-2v_{TRE1}$$

(3)

$$\frac{dATP}{dt}=v_{PGK}+v_{PYK}-v_{HK}-v_{PFK}-2v_{TRE2}$$

(4)

$$\frac{dTRE}{dt}=v_{TRE2}-v_{TRE1}$$

(5)

$$\frac{dG6P}{dt}=v_{HK}-v_{PHI}-2v_{TRE2}$$

(6)

$$\frac{dF6P}{dt}=v_{PHI}-v_{PFK}$$

(7)

$$\frac{dF16P}{dt}=v_{PFK}-v_{ALD}$$

(8)

$$\frac{dDHAP}{dt}=v_{ALD}-v_{GLY}-v_{TPI}$$

(9)

$$\frac{dGAP}{dt}=v_{ALD}-v_{GAPDH}+v_{TPI}$$

(10)

$$\frac{dGLY}{dt}=v_{GLY}$$

(11)

$$\frac{d{P}_i}{dt}=-v_{GAPDH}$$

(12)

$$\frac{dBPG}{dt}=v_{GAPDH}-v_{PGK}$$

(13)

$$\frac{dNADH}{dt}=v_{ACE}+v_{GAPDH}+3v_{SUC}-v_{GLY}-v_{ADH}$$

(14)

$$\frac{dP3G}{dt}=v_{PKG}-v_{PGM}$$

(15)

$$\frac{dP2G}{dt}=v_{PGM}-v_{ENO}$$

(16)

$$\frac{dPEP}{dt}=v_{ENO}-v_{PYK}$$

(17)

$$\frac{dPYR}{dt}=v_{PYK}-v_{PDC}-2v_{SUC}$$

(18)

$$\frac{dACA}{dt}=v_{PDC}-v_{ADH}-v_{ACE}$$

(19)

$$\frac{dSUC}{dt}=v_{SUC}$$

(20)

$$\frac{dACE}{dt}=v_{ACE}$$

(21)

The model generated from Equations (1)–(21) constitutes what will be called the CR model, that is, the most complete model that considers the existence of ramifications in the EMP route followed by alcoholic fermentation. On the other hand, the simplest model, called the SR model, is the model for which the velocities of Equations (5), (11), (20), and (21) will be equal to zero, that is, the SR model is the model where the EMP route and alcoholic fermentation occur directly and without deviations. Furthermore, as the focus of the present work is the development of a model for alcoholic fermentation, and cell reproduction in alcoholic fermentation in an anaerobic medium is low, we opted for the simplified use of kinetic relationships for cell growth, as will be discussed later.

A problem that may occur when using all these reactions in a metabolic model relates to the concentrations of energetic molecules, such as AMP, ADP, ATP, T6P, and NAD. These substances are used in numerous other reactions that occur in cell metabolism, being generated and consumed, thus maintaining such substances at approximately constant levels in healthy cells. Therefore, as chosen in the model proposed by van Eunen et al. [16], we decided to maintain the levels of these substances as constant in the present work. In addition to these substances, F26P is also another substance of importance in the metabolic network in question, but it is difficult to measure and ended up having its concentration considered as constant as well. In

Table 1. Substances with constant concentrations in the model application. Source: van Eunen et al. [16].

Substance	Concentration (mM)	Substance	Concentration (mM)
$ATP$	3	$T6P$	0.2
$ADP$	1	$F26P$	0.014
$AMP$	0.3	$NAD$	1.59

, we list the fixed concentrations for such substances.

In addition to constant concentrations, we need to employ a series of initial concentrations for the relevant metabolites in order to solve the model. The use of good initial values is essential for a good model response; however, finding practical values for the initial concentration of metabolites is complex. Teusink [17] and van Eunen [16] used a similar set of initial concentrations of substances, but some concentrations used by them were considerably high compared to the data presented in other works, such as Sato et al. [18], Ruoff et al. [19], Casei et al. [20], and Peeters et al. [21]. Thus, based on the data presented in such studies, we list the initial concentrations of metabolites in

Table 2. Initial concentrations of the substrate, ${\left[\mathrm{i}\right]}_0$, product, and internal metabolites from alcoholic fermentation. Source: van Eunen et al. [16]; Sato et al. [18]; Teusink et al. [17]; Ruoff et al. [19]; Casei et al. [20]; Peeters et al. [21].

Substance (mM)	[i]0 (mM)	Substance (mM)	[i]0 (mM)
GLCex	*	P3G	1.09
GLCin	0.1213	P23G	0.15
G6P	0.50	PEP	0.11
F6P	0.20	PYR	0.25
F16P	1.80	ACA	0.04
TRIO	0.50	EtOH	0.00
BPG	0.05	NADH	0.29

* The initial value to be chosen.

below.

The importance of choosing the initial values used for solving the systems of differential equations found in the present work should be highlighted. The system of differential equations proposed to solve the problem treated in the present work consists of a set of stiff differential equations, which generates a highly nonlinear system of differential equations. That said, the initial conditions used in solving the problem are extremely important, since a set of initial values can easily lead to non-convergence of values or even to mathematically correct results that are unrealistic in practice. Thus, a search was carried out in other works in the literature for initial concentration values for the various substances found in the system proposed here, looking for those that produced more stable and biologically viable results. Therefore, a series of initial values for the substances was tested in the model proposed here, noticing three types of behaviors, used as groups. The first group of values included values of the initial concentration of certain metabolites, such as fructose-6-phosphate and pyruvate, which presented a very abrupt drop in concentrations in the initial periods of fermentation, generating instability in the resolution, including negative concentrations. The second group of values included initial values that presented the opposite behavior, with an abrupt growth, generating concentrations biologically impossible to see inside a cell—in the case of glucose-6-phosphate, fructose-1,6-bisphosphate, and trioses. The third group of values generated a more stable response consistent with data found in the literature.

The initial substrate concentration in the present work varied between 30 g/L, 75 g/L, and 100 g/L. Another point we considered in this type of modeling is related to the cell concentration. Under anaerobic conditions with low substrate concentrations, Saccharomyces cerevisiae does not reproduce as much as under aerobic conditions. This consideration makes it common for the cell concentration to be set as constant in such fermentation assays; however, as the cell concentration relates directly to the enzyme concentration, the consideration of a variable amount of cells certainly enriches the model. Therefore, in the present model, we considered the cell concentration to follow the Andrews and Noack model, a simple model derived from the Monod model that showed a good predictive capacity according to the work of Kostov et al. [8].

$$\mu =\frac{1}{X}\frac{dX}{dt}={\mu }_0\frac{\left[S\right]}{1+\frac{K_{SX}}{S}+\frac{S}{K_{iS}}}$$

(22)

This Andrews and Noack model is a black box model, that is, it does not consider peculiarities regarding the fungus metabolism. The growth metabolism is extremely complex and not addressed in the present work. The purpose of inserting this growth model is just to add dynamic behavior to the cell concentration. Furthermore, the Andrews and Noack model takes into account a possible inhibitory effect on cell growth due to excess substrate in the culture medium, something important to take into account, since inhibitory effects on growth tend to appear even at lower levels of the substrate concentration.

With this, it is possible to start solving the system of equations. Both the kinetic equations and the parameters used are in the Appendix A and Appendix B. The values of the kinetic parameters for the action of the enzymes hexokinase (HK), phosphohexoisomerase (PHI), phosphofructokinase (PFK), fructose-bisphosphate aldolase (ALD), triose-phosphate isomerase (TPI), glyceraldehyde-3-phosphate dehydrogenase (GAPDH), bisphosphoglycerate mutase (PGK), phosphoglycerate mutase (PGM), enolase (ENO), pyruvatokinase (PYK), pyruvate decarboxylase (PDC), and alcohol dehydrogenase (ADH) were taken from Teusink et al. [22], Berthels et al. [23], van Eunen et al. [16], and Smallbone et al. [24].

For HXT, there is the limiting step of the fermentation process: that is, it is expected that this step presents a greater variety of possible values than the other steps. Thus, we carried out a manual search within a range of possible values for the kinetic constants, coming from the same works mentioned above. These constants are called apparent, as they vary with the initial conditions of the fermentation process, and the model proposed here does not aim to incorporate repression or expression terms into the behavior of proteins.

The models were validated using data from Acorsi et al. [25]. In the batch fermentation carried out in such works, the culture media used were prepared based on diluted final sugarcane honey and inverted at 50 °C with the addition of 1 mL of invertase solution for each 100 mL of the desired medium, where the enzyme solution contained 4 mg/mL of invertase. The raw material was provided by “Usina Santa Terezinha Iguatemi Unit” and collected directly from the output of the continuous honey centrifuge B, being the same raw material that supplies the plant’s fermentation vats.

We carried out fermentation tests from 30 to 100 g/L in a BIOSTAT^® B reactor, a bioreactor with a 5 L-capacity vessel with built-in mechanical agitation and temperature control, with the temperature controlled by the water distribution. In addition, this bioreactor features an automatic sampler, which makes sample collection possible, considerably reducing the risk of contamination during sampling. The other fermentations used smaller Kitasato-type reactors. We carried out the fermentation tests at a temperature of 32 °C [25].

For the quantification of sugars present in the fermentation medium at a given moment, as the wort used presented its sugars in the form of glucose and fruit, a DNS methodology modified for a wavelength of 600 nm was used [26,27,28]. For the quantification of ethanol, a VARIAN 330 with a Porapak Q column was used. The inlet and detector were kept at 120 °C, and the column was kept at 100 °C. Helium was used as a carrier gas, at a flow rate of 18.75 mL/min. Finally, 1 µL samples were injected to obtain the amount of ethanol present.

3. Results and Discussion

In Figure 2, Figure 3 and Figure 4, we present the results obtained by the SR model with the conditions described above.

In the three cases analyzed, we obtained a good fit for both the cell concentration and substrate consumption. However, for product formation, the results found were very poor: in none of the cases was product formation minimally close to the experimental values.

When we added the ramifications to the model, that is, when we run the model with branches, we obtained the results shown in Figure 5, Figure 6 and Figure 7.

Note that the prediction of product formation improves considerably with the addition of branches, even using extremely simple expressions such as those assumed in the present work. Thus, there is an indication that the addition of well-described ramifications to the metabolic model can further improve the results.

Considering the difference obtained between the models, we can compare the kinetic constants that most interfere in the response of the models, starting with the kinetic constants of the limiting step of alcoholic fermentation, i.e., the transport of glucose from the extracellular medium to the intracellular medium. Starting with the saturation constant K_GLC, we noticed that the increase in the concentration of the substrate in the medium decreased the microorganism’s affinity for glucose: in the analyzed fermentations, medium- to low-affinity transporters were found, with this constant varying from 6.95 to 54 mM. In both models, we found the highest affinity in the situation with the lowest substrate concentration, matching the need for yeast cells to facilitate substrate entry under nutrient-scarce conditions. Furthermore, we noted that the CR model needed slightly higher saturation constants, except in the fermentation with a higher initial substrate concentration. Furthermore, in the more concentrated fermentation, the saturation constant for the two cases with a greater availability of microorganisms in the CR model was repeated, something not seen in the SR model. This point may also indicate the need to investigate whether yeast can express different transport proteins at different stages of fermentation, according to the availability of the substrate in the medium: in the present work, the concentrations of relevant enzymes and proteins were constants, since including kinetic expressions for their expression in different phases of the fermentation process is a future step for the development of metabolic models.

As for the maximum substrate conversion speed, there was an increase in this parameter with increasing concentration, and in the case of the CR model, the maximum glucose transport speed was reached at a lower concentration of substrate available in the medium. Better results for substrate consumption can eventually be obtained by working without considering a very simple substrate composed only of glucose.

4. Conclusions

The present work aimed to propose a model for a batch alcoholic fermentation process using a series of kinetic equations derived from the metabolism of Saccharomyces cerevisiae. Considering the ease and efficiency of the fungus Saccharomyces cerevisiae in carrying out alcoholic fermentation, in addition to the fact that this fungus is a quasi-model organism, we successfully assembled a model capable of predicting sugar consumption and the formation of ethanol using the relevant metabolic network. This setup was certainly simple because the fungus considered in the present study is a quasi-model microorganism, and the central metabolic route, the EMP route, is well studied. The biggest difficulties arose with the ramifications of this route, which have not been sufficiently studied and lack information. However, even the simplified use of these ramifications showed an improvement in the results of the proposed models, especially for the formation of ethanol, which was closer to reality with the CR model. This work can still be improved by giving more detail to the kinetics of the relevant branches. Andrews and Noack’s model sufficiently represents cell growth, but it is a black box model; if we add the kinetic relationships from metabolism, its contribution would be improved. Thus, we have in our hands a basic metabolic model that shows promise for the representation and study of batch alcoholic fermentation, which can be improved and optimized in the future.