The sensitivity and significance analysis of parameters in the model of pH regulation on lactic acid production by Lactobacillus bulgaricus.

BACKGROUND: The excessive production of lactic acid by L. bulgaricus during yogurt storage is a phenomenon we are always tried to prevent. The methods used in industry either control the post-acidification inefficiently or kill the probiotics in yogurt. Genetic methods of changing the activity of one enzyme related to lactic acid metabolism make the bacteria short of energy to growth, although they are efficient ways in controlling lactic acid production.
RESULTS: A model of pH-induced promoter regulation on the production of lactic acid by L. bulgaricus was built. The modelled lactic acid metabolism without pH-induced promoter regulation fitted well with wild type L. bulgaricus (R2LAC = 0.943, R2LA = 0.942). Both the local sensitivity analysis and Sobol sensitivity analysis indicated parameters Tmax, GR, KLR, S, V0, V1 and dLR were sensitive. In order to guide the future biology experiments, three adjustable parameters, KLR, V0 and V1, were chosen for further simulations. V0 had little effect on lactic acid production if the pH-induced promoter could be well induced when pH decreased to its threshold. KLR and V1 both exhibited great influence on the producing of lactic acid.
CONCLUSIONS: The proposed method of introducing a pH-induced promoter to regulate a repressor gene could restrain the synthesis of lactic acid if an appropriate strength of promoter and/or an appropriate strength of ribosome binding sequence (RBS) in lacR gene has been designed.

Background

Lactobacillus delbrueckii subsp. bulgaricus has been widely applied in diary industry, especially as a starter for yogurt production. This microorganism is facultative anaerobic and can produce lactic acid. During the later period of fermentation and storage of yogurt, the production of lactic acid is dominated by L. bulgaricus [1]. The excessive production of lactic acid during storage makes the yogurt taste too sour and this phenomenon is called post-acidification.

At industrial level, several attempts have been made to slow down post-acidification, such as using devoid of L. bulgaricus, selecting the weak post-acidification L. bulgaricus, pasteurization after fermentation and so on [2-4]. However, these methods either control the post-acidification inefficiently or kill the probiotics in yogurt. Using genetic technology to modify the bacteria has become another option to decline post-acidification. Due to lacking molecular tools, mainly caused by the absence of a reliable transformation procedure, our understanding of the physiology and genetics of L. bulgaricus is still limited [5].

In L. bulgaricus, the absorption of lactose is processed by lactose/galactose antiport transport system [6]. Firstly, lactose is transported into cells by lactose permease encoded by lacS gene. Then the lactose is decomposed into glucose and galactose by β-galactosidase encoded by lacZ gene. The free galactose is pumped out of cells or stored in cells in the form of macromolecule i.e. carbohydrate gum. Glucose turns into pyruvate through glycolysis, and the glycolysis is transformed into lactic acid by the catalysis of lactate dehydrogenase [7-9]. The lacS gene and the lacZ gene constitute a lacSZ operon. In the downstream of lacZ gene, there is a lacR gene encoding a repressor which makes the lacSZ operon induced by lactose. However, in L. bulgaricus, the lacR gene has lost regulatory function due to the insertion of some gene fragments, resulting in constitutive expression of lacSZ operon [10].

Due to our limited knowledge, most attempts at the genetic level focused on changing the activity of one significant enzyme at a time in the metabolism of lactic acid. Druesne altered the Histidine codon to Alanine codon at the 552th locus of lactose permease, resulting in the reduction of the enzyme's activity [11]. Early in 1990, Mollet obtained the mutant L. bulgaricus without β-galactosidase activity by spontaneous deletion [12]. The new aforementioned organisms did produce less lactic acid but they had trouble in growing in milk independently due to the lack of energy. Adams selected two cold-sensitive mutants of the β-galactosidase from L. bulgaricus by using the expression system of E. coli [13]. However, their later research turned back to directly mutate and screen L. bulgaricus since it was difficult to transform with L. bulgaricus.

So far, Lactobacillus delbrueckii subsp. bulgaricus ATCC 11842 has been reported successfully transformed by electroporation, although with very low reproducibility and efficiency [5,14]. Several pH-induced promoters from Lactococcus lactis have been demonstrated, such as rcfB [15], P1 and P170 [16]. Combining the two points mentioned above, we came up with a new method to build a pH-induced promoter with a repressor gene controlling the production of lactic acid. Thus, we could turn on the switch at appropriate time. Here we analyse the parameters in the kinetic model to investigate the regulation effect of pH-induced promoter on lactic acid.

Methods

The reduced metabolic network

Here we present a model of the production of lactic acid, which includes the important enzymes, lactose permease, β-galactosidase, lactate dehydrogenase and related regulation genes lacS, lacZ, lacR, as shown in Figure 1. The process of glycolysis is reflected in one reaction, adopting the assumption that there is a constant flux from glucose to pyruvate. Pyruvate is partly catalysed to lactic acid, and the other two products, acetyl-CoA and butanedione, are not taken into account.

Figure 1

A schematic of metabolic and gene regulation network model of the production of lactic acid. The LacS, LacZ and LDH represent lactose permease, β-galactosidase, and lactate dehydrogenase respectively. When the extracellular pH reduces to 5.5, the pH-induced promoter (we chose rcfB promoter as an example) will trigger the express of LacR. The accumulation of LacR will combine with the operator sequence upstream of lacS gene, inhibiting the expression of LacS and LacZ [17].

The mathematical model

The mathematical model describes how the extracellular pH value influences the production of lactic acid. The ordinary differential equations in the model are based on Michaelis-Menten equations [18] and fundamental kinetic principles. The variables are described in Table 1. The model assumptions are as follows.

Table 1

Variables used in the model

Variable	Description	Units
LACout	Concentration of lactose out of the cells	mol/L
LACin	Concentration of lactose in the cells	mol/L
GLC	Concentration of glucose in the cells	mol/L
PYR	Concentration of pyruvate in the cells	mol/L
LA	Concentration of lactic acid out of the cells	mol/L
LR, LS, LZ	Concentration of LacR protein, lactose permease, β-galactosidase in the cells	mol/L
M	Concentration of the mRNA of lacR gene	mol/L
N	Concentration of the mRNA of lacS and lacZ gene	mol/L
X	Concentration of cell mass	g/L

(1) Each entity, described in Table 1 is given as total amount in the population.

(2) The model does not take energy metabolism into consideration.

(3) All reactions are modelled by mass action principles, except for enzyme reactions which obey Michaelis-Menten equation, and transcription which obeys saturation kinetics.

(4) The level of translation product is assumed to be proportional to mRNA levels according to Dockery and Keener [19].

(5) There is no delay in synthesis of either component or delay because of protein transportations.

(6) All substances are subject to degradation.

The metabolism of lactose

As mentioned before, the metabolism of lactose is described by Michaelis-Menten equation. Since lactose transport is reversible, a term was included to count for lactose efflux dependent on the internal lactose concentration [20]. The Vfor lactose permease and β-galactosidase is associated with the concentration of enzymes, so we use Kand the concentration of enzymes to estimate the V. Since the lac operon is usually induced by lactose, we assume that the repressor LacR protein could react with lactose. Thus we need to take the reduction into consideration when evaluating the concentration of intracellular lactose and LacR protein. As for the utilization of sugar, L. bulgaricus belongs to the type of homofermentation. Theoretically, this microorganism could totally convert glucose to lactic acid, but the conversion of glucose is about 80%-90% in practice [21]. Therefore, the estimated proportion of pyruvate transformed into lactic acid is 80%.

The gene regulation

The transcription and decomposition of lactose are tightly regulated at the genetic level: production of the enzymes can be decreased at the transcriptional level by regulatory protein LacR binding at appropriate DNA sites. We set M to stand for the transcript concentration of lacR gene. Since lacS and lacZ are located in the same operon, it is assumed that both genes are transcribed at same rate. We used N to describe the concentration of the two genes' mRNA. When describing the equations, the mRNAs are assumed to be transcribed at two distinct rates: basal one when there is no regulation and a higher/lower rate when being regulated [22].

The transcription of lacR gene

For M, it is described as function (1). The basal rate is given by Vand increased at a constant rate Vwhich is regulated by the pH value. GR stands for the repressor gene concentration. The degradation of the mRNA occurs at the rate of d. It is also diluted due to cell growth rate at μ. We assumed the decrease of mRNA is first-order proportion to the concentration of M.

As the mechanisms of pH-induced promoters have not been fully elucidated, we included a generic pH-dependent switch F(pH) which turns on when pH value is below the threshold. The switch takes the form of the smoothed step function [23], where 5.5 represents the threshold pH level around which the switch occurs [15]. n indicates the steepness of the smooth switch function.

Since lactic acid is the main product during fermentation, we assumed pH value is all controlled by the concentration of lactic acid. Besides, we needed to take the buffer capacity of milk medium into account. Furthermore, we assumed that all the lactic acid could spread to the medium even when the extracellular concentration was pretty high. We tested the pH value of milk medium when different amount of lactic acid was added. Three parallel experiments were taken and we fitted a function as the form of function (3). To the pH data, cand care constants. This function was then fed into the model via the switch function F(pH).

Transcription of lacS and lacZ gene

For N, it is described as function (4). The basal rate is given by Vand the transcription is inhibited by LacR. For the part to simulate the regulation rate of LacR, Tstands for the maximal rate of lacS and lacZ gene transcription, S represents the sensitivity of lacS and lacZ gene transcription to lactose permease and β-galactosidase [24]. G is the concentration of transcription gene. Again, we assumed first-order degradation rate at dand dilution rate at μ of the mRNA.

Cell growth

In milk medium, lactose is the main carbon source, so it is assumed that the cell growth rate is dependent on the extracellular lactose concentration. The product, lactic acid, which declines the pH value of milk medium performs inhibition on cell growth. The function to describe the production inhibition is the same as Concepcion presented [21]. In function (5), Kis the Monod constant for growth in extracellular lactose, and μis the maximum specific growth rate. Kis the maximum initial lactic acid concentration in which the microorganism growth is completely inhibited.

Combining the metabolism of lactose with gene regulation

All the definitions of the parameters are shown in Table 2. Combining with the enzyme reactions, gene regulations and cell growth discussed above, the process described in Figure 1 could be represented by:

Table 2

List of parameters used in the model

parameter	definition	value	reference/derivation
Km·LP	Michaelis constant for lactose permease	2.6 × 10-4 M	Patrick et al. (1997) [20]
Kcat·LP	rate constant for lactose transportation	1.29 × 105h-1	Patrick et al. (1997)
Km·GAL	Michaelis constant for β-galactosidase	9.8 × 10-4M	Moez et al. (2009) [25]
Kcat·GAL	hydrolysis rate for lactose	1.65 × 105h-1	Moez et al. (2009)
Km·GLYC	Michaelis constant for glycolysis	1 × 10-4M	L. Lactis Marcel et al. (2002) [26]
Vmax·GLYC	maximum reaction speed when enzymes was saturated by substrate	143.82M·h-1	L. Lactis Marcel et al. (2002)
Km·LDH	Michaelis constant for lactate dehydrogenase	0.1 M	L. Lactis Marcel et al. (2002)
Vmax·LDH	maximum reaction speed when lactate dehydrogenase was saturated by pyruvate	307.08M·h-1	L. Lactis Marcel et al. (2002)
KLS	LacS synthesis rate constant	564 h-1	Kennell et al. (1977) [27]
dLS	LacS decay rate constant	0.6 h-1	Kennell et al. (1977)
KLZ	LacZ synthesis rate constant	1.13 × 103h-1	Kennell et al. (1977)
dLZ	LacZ decay rate constant	0.6 h-1	Kennell et al. (1977)
KLR	LacR synthesis rate constant	564 h-1	Assumed
dLR	LacR decay rate constant	0.6 h-1	Assumed
dM , dN	mRNA decay rate constant	41.58 h-1	Varmus et al. (1970) [28] and Patrick et. al. (1997)
GR, G	lacR, lacSZ gene concentration	2.5 × 10-9 M	Cheng et al. (2001) [29]
n	steepness of the smooth switch function	1	Akyol et al. (2008) [15]
Kb	LacR/lactose binding rate constant	7.2 × 104 M-1·h-1	Assumed
V0	basal rate of lacR gene transcription for mRNA	1 h-1	Assumed
V1	increasing constant of lacR gene transcription for mRNA	20 h-1	Assumed
V2	basal rate of lacSZ gene transcription for mRNA	3.5 h-1	Obtained by 'fitting'
Tmax	maximum rate of lacSZ gene transcription	6 × 104 h-1	Assumed
S	sensitivity of lacSZ gene transcription to LacS and LacZ	1000	Guessed
KS	saturation constant for growth of L. bulgaricus on extracellular lactose	9.82 × 10-3 M	Concepcion et al. (2000) [21]
KLA	maximum lactic acid concentration to inhibit the growth of L. bulgaricus	0.448 M	Concepcion et al. (2000)
μ	maximum growth rate	0.8 h-1	Concepcion et al. (2000)

The final metabolic model is given by:

We set stoichiometric constant of Rwith two since two molecules of pyruvate are formed from one molecule of glucose. For the stoichiometric constant of R, it stands for the proportion of pyruvate which is converted to lactic acid.

Data estimated

Firstly, we did not take lacR gene into consideration. This is just the case of wild type L. bulgaricus which lacSZ operon is constitutive expression. Since most of the parameters in this reduced model have been reported, Vis the only one needs to be modified. We used SBToolbox in Matlab to set up the reduced model and simulated the metabolism process [30]. We compared the concentration changes of extracellular lactose and lactic acid within 8 hours with the data from Fatama, where the conditions used were 43 °C, 8% of skim milk concentration (without pH value control) and 4% of inoculum ratio [31]. After that, we added lacR gene. This time we did not have sufficient data for full parameterization of the model. However, we did have enough information to estimate their relative sizes and give the qualitative nature of our investigations. The parameter values are listed in Table 2.

Sensitivity and significance analysis of parameters

We still used SBToolbox to perform local sensitivity analysis and Sobol's method for global sensitivity analysis. The local sensitivity analysis investigates how small changes in a single parameter value could affect the model output. The method is based on the partial differentiation of the output with respect to the input parameters [32]. Herein, the partial differentiation is evaluated numerically by introducing a 1% increment from the specific parameter value. We chose Sobol sensitivity analysis method to calculate global sensitivity since the sensitivity index quantified the overall effects of a parameter, in combination with any other parameters, on the model output [33]. The number of simulation to carry out is 10000 times. We first set the number as 1000 according to Schmidt [30]. Since the results are varied among different simulations, we increased the number to 10000 and got stable outcomes. During the above sensitivity analysis, we chose lactic acid concentration as the model output.

After figuring out the sensitivity parameters to the model output, we carried out different simulations with each sensitivity parameter at two levels, a relative low value and a relative high value. The metabolism of lactic acid and changing of pH values was carefully considered. Those results indicate how pH values affect lactic acid production under different situations and help us to understand how to design LacR regulation in biology experiments in future.

Results and discussion

Reduced model without lacR gene

The results of modelled concentration changing of extracellular lactose and lactic acid are shown in Figure 2, along with the experimental data. The model fits the experimental data well. We have obtained the value of V, which indicates the basal rate of lacSZ gene transcription. The relative size of fitting value (V= 3.5 h-1) is determined according to Gustafsson's [24] report.

Figure 2

Simulation of batch fermentation conversion of lactose to lactic acid using wild type . The initial conditions were 43 °C 8% of skim milk concentration (without pH value control) and 4% of inoculum ratio. The dots and triangles indicated the measured concentration of lactose and lactic acid in the medium, respectively, every two hours. The curves were the simulated results of lactose and lactic acid concentration in 8 hours.

Sensitivity analysis of parameters

The results of local sensitivity analysis and Sobol sensitivity analysis are shown in Figure 3. The parameters with sensitivity from high to low in local sensitivity analysis are followed as Tand V. The other parameters perform little effect on the production of lactic acid. The reason is that when one of those parameters increased, it would only affect the production rate of the output but have no effect on the total amount of lactic acid. Thus, they are not sensitive parameters in this method. As for Sobol sensitivity analysis, the parameters with sensitivity from high to low is Tand d. The parameter Kshows the second most sensitivity in Sobol's method, just as Druesne has proven, the reduction of lactose permease's activity would greatly decline the production of lactic acid [11].

Figure 3

Local sensitivity analysis and Sobol's method of global sensitivity analysis. A Local sensitivity analysis B Sobol's method of global sensitivity analysis

Significant analysis of parameters

In sensitivity analysis, sensitivity does not mean importance, since sensitive parameters are always fixed inherently and unadjustable. We need to combine the sensitivity results with operability, especially for the design of biology experiments. Although Tand uperform high sensitivity in both methods, these parameters are usually inherently. Kvalue is adjustable, but to decrease its value would generate fatal disadvantage on the growth of the bacteria. So we do not include this parameter into consideration. Due to the same reason, we exclude K, which is the Michaelis constant for glycolysis. As for the parameter S, which means sensitivity of lacSZ gene transcription to LacS and LacZ proteins, it is difficult to measure and adjust in experiments. Therefore, the parameters which exhibit great influence on the production of lactic acid in both methods are K, Vand V, with Krepresenting RBS strength of lacR gene, and Vtogether with Vreflecting characters of the pH-induced promoter.

For the three parameters, the relative high levels of Kand Vare 1128 h-1, 1 h-1 and 21 h-1 respectively. The relative low levels of Kand Vare 564 h-1, 0.1 h-1, 10.5 h-1 respectively. The results are shown in Figure 4. From picture A, we could figure out that with the increase of promoter strength, the production rate and amount of lactic acid would be greatly decreased. However, for the relative high value of V, the finial pH value only declines to 4.83 which is higher than casein isoelectric point. Thus we need to find an appropriate promoter strength between the relative high and low values when it is turned on by pH value. From picture B, we could draw the conclusion that if the pH-induced promoter could be efficiently turned on when pH decreases to its threshold, leaking expression of the promoter would have little effect on the production of lactic acid. In picture C, the increase of Kleads to strong decline of both the production rate and total amount of lactic acid. By comparing the results of relative high level Kwith the relative high level Vin picture A, the two parameters almost perform the same influence on the production of lactic acid. This is due to the reason that the two parameters both affect the concentration of LacR. Therefore, we could adjust the strength of pH-induced promoter and/or RBS strength of lacR gene to make our switch perform well in controlling lactic acid production.

Figure 4

The effects of . The variations of lactic acid and pH in wild type are set as blank control. The curves indicated the simulated concentration of lactose and pH along with time under different values of Vand K. A The effects of V(10.5 h-1 or 21 h-1), where V= 1 h-1, K= 564 h-1. B The effects of V(0.1 h-1 or 1 h-1), where V= 10.5 h-1, K= 564 h-1. C The effects of K(564 h-1 or 1128 h-1), where V= 1 h-1, V= 10.5 h-1

Other two important features for the promoter are the threshold pH level around which it would start and whether the promoter would keep turning on when the pH value decreased to 3.5 or lower. Akyol [15] and Madsen [16] only described the status of pH-induced promoters along with the pH declined to 5.5. Further experiments are needed to test the characteristics of the promoters. This model neglect the effect of lower pH on the promoter, assuming that it would still turn on with the decrease of pH value.

Conclusions

In conclusion, we propose a new method to control the production of lactic acid by L. bulgaricus which is to build a pH-induced switch (promoter) to turn on repressor gene. The proposed method overcomes the disadvantage in bacterial growth by directly changing one of the enzyme related to lactic acid metabolism. Then, we build a reduced model of pH-induced promoter regulation on the production of lactic acid and adjust the model with data from wild type L. bulgaricus. After that we carry out sensitivity analysis of the parameters and figure out three significant ones in the model. To make our switch work well, we need to find an appropriate strength of promoter and/or an appropriate strength of RBS in lacR gene. The values of those two parameters should between the high level and low level we have set in the analysis, so that when the pH value declines to the promoter's starting point, the lacR gene could express moderately to inhibit the production of lactic acid.

Competing interests

The authors declare that they have no competing interests.

Authors' contributions

KL, LQ, XL, YY proposed the project and designed metabolic and gene regulation network. KL and XZ wrote the model and performed sensitivity analysis. CD, AH and DX guided the model network design and significance analysis of the parameters. KL wrote the manuscript. All authors have read and approved the final manuscript.