Skip to main content

Clostridium butyricum maximizes growth while minimizing enzyme usage and ATP production: metabolic flux distribution of a strain cultured in glycerol



The increase in glycerol obtained as a byproduct of biodiesel has encouraged the production of new industrial products, such as 1,3-propanediol (PDO), using biotechnological transformation via bacteria like Clostridium butyricum. However, despite the increasing role of Clostridium butyricum as a bio-production platform, its metabolism remains poorly modeled.


We reconstructed iCbu641, the first genome-scale metabolic (GSM) model of a PDO producer Clostridium strain, which included 641 genes, 365 enzymes, 891 reactions, and 701 metabolites. We found an enzyme expression prediction of nearly 84% after comparison of proteomic data with flux distribution estimation using flux balance analysis (FBA). The remaining 16% corresponded to enzymes directionally coupled to growth, according to flux coupling findings (FCF). The fermentation data validation also revealed different phenotype states that depended on culture media conditions; for example, Clostridium maximizes its biomass yield per enzyme usage under glycerol limitation. By contrast, under glycerol excess conditions, Clostridium grows sub-optimally, maximizing biomass yield while minimizing both enzyme usage and ATP production. We further evaluated perturbations in the GSM model through enzyme deletions and variations in biomass composition. The GSM predictions showed no significant increase in PDO production, suggesting a robustness to perturbations in the GSM model. We used the experimental results to predict that co-fermentation was a better alternative than iCbu641 perturbations for improving PDO yields.


The agreement between the predicted and experimental values allows the use of the GSM model constructed for the PDO-producing Clostridium strain to propose new scenarios for PDO production, such as dynamic simulations, thereby reducing the time and costs associated with experimentation.


The rising biodiesel industry has resulted in a major overproduction of glycerol as a byproduct, which now threatens the economic viability of this industry [1, 2]. This situation has spurred research into glycerol utilization as a carbon source [3,4,5] and for the generation of products such as 1,3-propanediol (PDO), a precursor of important commercial polymers, such as polyester and polyurethane [6, 7]. PDO can be biosynthesized from glycerol by bacteria such as Clostridium butyricum or Klebsiella spp. [3, 7]. Clostridium species are the more attractive alternative because they are safer and achieve higher yields than Klebsiella [8]. However, industrial PDO production using bacteria is still limited by insufficient yields, which presents a serious obstacle to the competitiveness of this process [9,10,11]. Therefore, strategies such as fed-batch cultures and random mutagenesis have been developed, resulting in improvements in PDO production of up to 137% and 78%, respectively [11,12,13]. A more detailed understanding of the metabolic pathways in species such as Clostridium butyricum could therefore shed light on a better ways to promote glycerol transformation to PDO in this organism.

Metabolism studies of glycerol by the anaerobic bacterium Clostridium butyricum have generally focused on its central metabolism, which is composed of oxidative and reductive branches [14]. The oxidative branch is mainly related to the production of ATP and reducing equivalents (NADH), with the formation of acetic and butyric acids as byproducts. By contrast, the reductive branch produces PDO while simultaneously regenerating reducing equivalents by conversion of NADH to NAD [7, 9, 15]. Bizukojc et al. [16] reported the most detailed metabolic model for a PDO producer Clostridium strain, indicating the functioning of 77 reactions and 69 metabolites. The model, in addition to the oxidative and reductive branches, also included simplified synthesis reactions for amino acids, macromolecules, and biomass. However, at present, metabolic models based on genome annotation information, also known as genome-scale metabolic (GSM) models [17, 18], are lacking for Clostridium butyricum.

A proteomics study of the native Colombian strain Clostridium sp. IBUN 158B cultured in glycerol [19] has provided experimental validation of the enzyme expression involved in PDO metabolic networks in this specie. The proteome contained 21 enzymes classified as follows: one from the reductive branch (PDO dehydrogenase), three from the oxidative branch, eleven from carbohydrate synthesis, four from amino acid synthesis, and two from nucleotide synthesis. Gungormusler et al. [20] also used proteomics for the experimental detection of 262 different enzymes expressed by Clostridium butyricum 5521 cultured in glycerol. Nevertheless, despite this experimental information and the computational tools available, the prediction of PDO production by Clostridium based on its metabolic behavior is still limited.

One computational tool commonly employed for metabolic modeling is flux balance analysis (FBA). FBA allows the use of a steady state assumption of defined culture conditions to predict the phenotype of one microorganism based on its GSM model [21,22,23,24,25]. However, a GSM model expressed as stoichiometric matrix is an undetermined system, that is, it has more reactions than metabolites. This creates a situation with infinite solutions, so an objective function is required to predict the microorganism phenotype. FBA then becomes an optimization process in which the constraints are the culture conditions, mass balances, and thermodynamic feasibilities [22, 25,26,27,28].

In general, predictions using GSM models assume biomass yield maximization as the objective function, based on the assumption that cells have evolved to select the most efficient pathways that achieve the best yields [29]. Nevertheless, predictions with biomass maximization do not always capture the cellular physiology, and alternative objective functions have been developed [28, 30,31,32,33]. Studies have included error minimization by bi-level optimization [30, 34, 35], objective function selection by Bayesian inference [31] or by Euclidian distance minimization [32], and linear combination of objective functions [28, 36]. The results, overall, highlight that a cell does not maximize biomass yield under scenarios like substrate excess, so that one single function is unable to predict all the evaluated scenarios [28, 32, 33, 37,38,39].

For these reasons, the initial purpose of the present research was to construct the first GSM model of a PDO producer Clostridium strain. The biological model selected was the Colombian strain Clostridium sp. IBUN 13A, a strain isolated by our Bioprocesses and Bioprospecting Group. This strain is a natural PDO producer and has been employed over the last 20 years in several studies aimed at understanding PDO production, including the annotation of its genome [40,41,42]. Additionally, as second objective, our intent was to predict the phenotypic states of this bacterium during culture in glycerol and in other substrates using the GSM model and FBA with the adequate objective function. Our overall aim was to evaluate the effect of perturbations in the constructed GSM model on PDO yield improvements.

Results and discussion

Genome-scale metabolic model iCbu641 reconstruction and curation

A draft metabolic model for the Clostridium sp. IBUN 13A strain was constructed based on RAST annotation [41]. The draft was composed of 641 genes, 365 enzymes, 671 reactions, and 606 metabolites. GapFind [43] analysis of the draft model identified 303 blocked metabolites, which were reduced to 63 by adding 59 reactions based on experimental fermentation evidence from Clostridium butyricum cultured in glycerol [15, 19, 44] and on curated GSM models from other solventogenic clostridia [16, 45,46,47,48,49,50,51,52]. The biomass reaction was adapted from C. beijerinckii GSM [45], which does not account for the proton formation associated with ATP hydrolysis during the growth-associated maintenance (GAM), as is also observed in C. acetobutylicum [49]. This excluded proton would accumulate in the biomass, thereby preventing stabilization of the biomass charge. By contrast, the GSM models of C. thermocellum [50], C. ljungdahlii [52], and C. cellulolyticum [51] include this proton production in their biomass reactions. Therefore, the proton formation was included in the present biomass reaction to resolve the inconsistency in the elemental composition of the biomass, as well as the charge balance. The elemental composition per C atom, calculated based only on stoichiometric consumption of precursors, was therefore CH1.624O0.456N0.216P0.033S0.0047.

After curation, elemental balancing, and loop deletion, the constructed iCbu641 GSM model included 641 genes, 891 reactions, and 701 metabolites. Table 1 summarizes the main features of the curated metabolic network, and Fig. 1 shows the pathway distribution of the cytosolic reactions. Comparison with other GSM models of solventogenic Clostridium strains showed 11 unique enzymes, including as PDO dehydrogenase (EC. and glycerol dehydratase (EC. These two enzymes function in lipid metabolism but are associated with the reductive branch of glycerol metabolism [14].

Table 1 Main features of the iCbu641 metabolic network
Fig. 1
figure 1

Distribution of cytosolic reactions in the iCbu641 GSM model by functional pathway. Notation (■) Gene-Associated reactions (■) Non Gene-Associated reactions

The iCbu641 metabolic network is the first GSM model curated for a PDO producing Clostridium strain [17, 18, 49] (See Additional files 1 and 2 for complete metabolic model at excel format and SBML format, respectively). The iCbu641 model also includes all the enzymes associated with glycolysis and the pentose phosphate pathway and most of the enzymes involved in the TCA cycle. The enzymes from the TCA cycle that were not included are malate dehydrogenase (EC., succinate-CoA ligase (EC. - EC. and fumarate reductase (EC., which were not detected in the genome. Therefore, additional experimentation is required to verify the presence or absence of genes encoding these three enzymes in the genome of Clostridium sp. IBUN 13A. The model is able to synthesize de novo all the precursors involved in the biomass reaction (e.g., amino acids, nucleotides, fatty acids, teichoic acid, and cofactors).

Flux distribution prediction using flux balance analysis with glycerol as substrate

The constructed iCbu641 GSM model and FBA were employed to predict the flux distribution of Clostridium butyricum cultured in glycerol. FBA was solved using linear programming (LP) with biomass maximization as an objective function. When compared with the experimental data, FBA predicted a biomass overestimation and no PDO production (Fig. 2 – Blue lines), giving biomass yield (YX/S) and PDO yield (YPDO/S) errors of 300% and 100%, respectively [44]. Therefore, taking into account only biomass yield, we could infer that, experimentally, Clostridium butyricum does not grow optimally in glycerol. This can be explained mathematically, because the objective function and all constraints (mass balances and thermodynamic feasibilities) used were linear. This means that the optimum found is indeed a vertex of the feasible solution space [38], where no PDO could be produced.

Fig. 2
figure 2

Robustness analyses in function of glycerol uptake flux. a Specific growth rate μ and b PDO secretion flux v PDO. Notation: Experimental data from Solomon et al. [44] at () glycerol limitation and () glycerol excess. FBA predictions using: biomass maximization (Blue Line), biomass maximization per enzyme usage (Red Line) and biomass maximization while minimizing both enzyme usage and ATP production (Green Line)

Contrary to FBA predictions, PDO production in Clostridium sp. IBUN 13A is related to dha operon expression and the enzyme activity of both PDO dehydrogenase (EC. and glycerol dehydratase (EC., which have been experimentally detected in the presence of glycerol as main carbon source [40, 53, 54]. The lack in predictions can be interpreted biologically based on redox balance together with the capability of C. butyricum to produce formic acid and hydrogen (H2) [44]. Clostridium butyricum is an anaerobic bacterium, so according to the redox balance, the substrate must act simultaneously as an acceptor and donor of electrons [55]. However, the aforementioned capability allows to LP optimization predicts that the substrate could be used mostly as an electron donor, thereby generating more ATP and biomass. This prediction is achieved due to the regeneration of reducing equivalents through formic acid and H2 formation, which would require no substrate as an electron acceptor to produce reduced compounds such as PDO. In other words, in the vertex predicted by LP optimization, the substrate is used mainly as an electron donor due to the formation of formic acid and H2, and no PDO is produced.

Similar results were observed using glucose as a substrate, where formic acid and H2 were overproduced instead of reduced products like butanol or ethanol, which have been detected experimentally (data not shown) [56]. Consequently, LP optimization could be considered as unsuitable for predicting the experimental yields of Clostridium butyricum, which is capable of producing formic acid or H2. Nevertheless, LP optimization is commonly employed in solventogenic Clostridium strains, although all of them are able to produce at least hydrogen [45, 47, 48, 50, 51]. This lack of prediction could be solved using experimental constraints of formic acid and H2 secretion in LP optimization; however, no linear trend was observed, especially in H2 secretion [44].

A new objective function of maximizing biomass yield per flux unit (Equation 1) was therefore employed to improve the FBA predictions. This objective function is based on the hypothesis that cells operate to maximize biomass yield while minimizing enzyme usage [32]. This non-linear programming (NLP) optimization is non-convex, but Schuetz et al. suggested that the predicted local optimum is indeed the global optimum [32]. In addition to the new objective function, a non-linear constraint was employed, which corresponds to the maximum acetic acid secretion flux and has an allosteric trend in the function of glycerol uptake flux [44]. Therefore the constraint was not used during simulations with other carbon sources, as glucose. This allosteric trend of acetate production has been previously reported for Clostridium butyricum cultured in glycerol as a mechanism to control acetyl-CoA/CoA and ATP/ADP ratios [57]. The new simulations (Fig. 2 – Red lines) predicted YX/S and YPDO/S errors of nearly 4.5%, and 1.5%, respectively, when compared with the experimental values obtained under limiting conditions of glycerol (<15 g/L) [44]. Therefore, the error reduction using NLP predictions suggests the resolution of the redox balance problems observed in LP optimization caused by the capability of producing formic acid and H2. It also suggests that Clostridium butyricum indeed minimizes enzyme usage and prefers short pathways (PDO production) to maximize its growth under limiting nutrient conditions.

$$ \mathit{\operatorname{Max}}\ \frac{\mu}{\sum_{i=1}^n{v}^2} $$

Identical phenotypic predictions were also observed using both LP and NLP optimizations when the production of hydrogen and formic acid were blocked as additional constraints (data not shown). This validates, on the one hand, the effect of these products in the lack of prediction using LP optimization. On the other hand, NLP optimization indeed predicts global optima, as suggested Schuetz et al. [32].

However, the new objective function overestimated YX/S under glycerol-excess conditions (without reaching inhibition conditions) [44], suggesting that Clostridium butyricum grows under sub-optimal conditions when the substrate is present in excess, as is also observed in E. coli [32]. This behavior is understandable if thermodynamics is considered. Growth, as with any reaction, is thermodynamically feasible if its driving force, expressed as Gibbs free energy, is negative (ΔGgrowth < 0). Growth Gibbs free energy can be expressed as a function of biomass yield and the Gibbs free energies for catabolism and anabolism, as shown Eq. 2 [55]. Therefore, an improvement in thermodynamic feasibility will be coupled to biomass yield reduction, as occurs in scenarios such as substrate competition with other microorganisms [55]. This suggests that a substrate excess condition could induce the above scenario in order for the organism to prevail in this environment.

$$ {\Delta G}_{growth}={\frac{1}{Y_{X/ S}}\Delta G}_{catabolism}+{\Delta G}_{anabolism} $$

The sub-optimal behavior was accounted for in the objective function of Eq. 1 by the use of a tunable weighting factor (w) that minimizes both the enzyme usage and ATP production of the network (Equation 3). The ATP incorporation in the objective function is based on experimental results [15, 44] and because a reduction in ATP production is related to the respective biomass yield reduction. Therefore, a weight w equal to 1 corresponds to the optimal conditions of Eq. 1. The use of weight factors has been described by Torres et al., who added ATP and NADH/NADPH minimization or maximization to the objective function, which improved S. cerevisiae growth predictions up to 98% [36, 39]. On this basis, simulations under excess conditions were made using a weight factor equal to 0.04 (w = 0.04), where the average experimental errors were minimal [44]. The average YX/S and YPDO/S errors were 5.3% and 2.5%, respectively, as shown Fig. 2 – Green Lines, confirming the ability of FBA to provide accurate predictions of these results through NLP optimization. However, this weight factor only applies to Clostridium butyricum cultured anaerobically in glycerol; additional experimental data would validate its usage.

$$ \mathit{\operatorname{Max}}\kern0.75em \frac{\mu}{w\left[\sum_{i=1}^n{v}^2\right]+\left(1- w\right)\left[{v_{ATP\ prod}}^2\right]}\kern.4em w\in \left(0,1\right) $$

We further validated the necessity of employing the weight factor at sub-optimal conditions of glycerol by comparing both predicted phenotypes with the results reported by Zeng [15]. Zeng found that the directionality of the reaction catalyzed by ferredoxin reductase (EC. or EC. changes if glycerol is in limitation or in excess. Under the limiting condition, ferredoxin reductase consumes reducing equivalents and produces more H2 and less PDO than is observed under the excess condition, where ferredoxin reductase produces reducing equivalents. Consequently, Zeng quantified ferredoxin reductase directionality by calculating the ratio between hydrogen and reduced ferredoxin formed (\( {\alpha}_{H_2/ Fd} \)). The experimental ratios are 1.1 and 0.4 under limiting and excess glycerol conditions, respectively. Similarly, the ratios predicted by FBA are 1.46 and 0.40 at optimal and sub-optimal conditions, respectively. Therefore, the predictions for both conditions agree with the reported values and validate the utility of the weight factor in predicting glycerol cultures.

We also compared the flux prediction for the 365 enzymes present in iCbu641 with the experimental expression of the 286 enzymes detected by Gungormusler et al. in the proteome of the Clostridium butyricum DSM 10702 strain cultured under limiting glycerol conditions [20]. An analysis of the 174 common enzymes using flux couple finding (FCF) revealed 80 enzymes that were partially or fully coupled to growth and 67 enzymes that were directionally coupled to growth [58]. The remaining 27 enzymes were blocked and were therefore excluded from the comparison. These 27 enzymes could be blocked due to a lack of gene annotation or exclusion from biomass synthesis (e.g., terpenoid synthesis). Terpenoids have been detected in the cell walls of Clostridium strains and may arise as a stress response to the acids formed during culture [59]. The terpenoid pathway could be unblocked if these metabolites are added to biomass synthesis, but their concentration levels first need accurate quantification [59].

Assuming a qualitative correlation between the expression and flux for the 147 enzymes included in the comparison, FBA was able to predict the expression of 123 of them (83.7%). The remaining 24 non-predicted enzymes correspond only to directionally coupled growth and most of them (21 enzymes) are involved in carbohydrate metabolism and the synthesis of nitrogenous compounds (amino acids and nucleotides). The absence of a prediction for these enzymes could be due the presence of alternate pathways and isozymes, as is the case for asparagine synthase (EC., isocitrate dehydrogenase (NADP) (EC., and glycerol-3-phosphate dehydrogenase (EC., where the alternative enzymes are asparagine synthetase (EC., isocitrate dehydrogenase (NAD) (EC., and glycerol kinase (EC., respectively. Another possible cause is the inability of FBA to predict regulation mechanisms [60], as is the case for the enzyme pyruvate phosphate dikinase (EC., which is involved in the gluconeogenesis pathway but appears to act in the place of pyruvate kinase, consuming AMP instead ADP [61, 62]. Other enzymes, such as nicotinate phosphoribosyltransferase (EC. or pyrimidine nucleoside phosphorylase (EC., show a lack of prediction because they are involved in RNA or DNA fragment recycling [63]. Finally, 1,4-alpha-glucan branching enzyme (EC. or starch synthase (EC. were not predicted as they are part of granulose synthesis, a process that is not included in the biomass reaction. Granulose is a polysaccharide employed as carbon source during sporulation; therefore, it is produced during exponential growth [64]. (See Additional file 3 for complete proteomic comparison results).

Flux distribution prediction using FBA and other carbon sources

The robustness of iCbu641 was further tested by comparing the experimental and simulated data using other substrates at optimal conditions, Eq. 1. We first compared the FBA predictions and experimental yields of Clostridium butyricum W5 cultured in glucose [56], as shown in Table 2. We observed an accuracy of nearly 97% for predicting the biomass yield using this substrate, confirming the ability to predict not only glycerol cultures but also glucose cultures. All the reported experimental yields also agreed with their respective predicted feasible ranges calculated using flux variability analysis (FVA).

Table 2 Comparison of the experimental and simulated yields (mol/mol) of Clostridium butyricum W5 cultured in glucose [56]

Junghare et al. [65] also evaluated biomass and hydrogen production of Clostridium butyricum TM-9A using different carbohydrates. Our comparison of the experimental yields with the predicted yields obtained through FBA is shown in Table 3. The trends in the predicted YX/S agree with the experimentally obtained values, showing smaller yields for pentoses and the highest yield for the trisaccharide raffinose, while the yields using monosaccharides were lower than those obtained using disaccharides. However, some differences are observed between the experimental and predicted values. First, ribose and xylose had considerably higher experimental than predicted yields. Second, the experimental yield for cellobiose was much lower than the predicted yield. Finally, although the simulations predicted the same yields for arabinose and ribose, their experimental yields differed. The first case could be a result of an incomplete curation of iCbu641 related to pentose consumption; therefore, more experimental information is needed. The second and third cases could be due to miscalculation of the experimental yields of arabinose and cellobiose, since these substrates were not consumed completely, as can be inferred from their pH reports [65]. The last can point to the possibility that Clostridium was not well adapted to these substrates and may have needed to undergo more generations to reach its optimal growth [33]. The hydrogen yields (YH2/S) showed experimental values that were mostly within their respective predicted feasible ranges. The only unpredicted YH2/S corresponded to arabinose, which supports the necessity of complementing iCbu641 for consumption of pentoses. Consequently, iCbu641 has the capacity to be employed to predict Clostridium butyricum growth using different carbohydrates as substrates.

Table 3 Comparison experimental and simulated yields (mol/mol) of Clostridium butyricum TM-9A cultured in different carbohydrates [65]

We also assumed a qualitative correlation between the enzyme flux and mRNA expression and compared FBA predictions and experimental transcriptomics results for the C. butyricum strain CWBI 1009 cultured in glucose [66]. Of the 288 enzymes shared by both systems, 51 were blocked according FCF analysis and excluded from comparison. Among the remaining 237 enzymes, FBA predicted the activity of the 123 enzymes as partially and fully coupled to growth and 56 enzymes directionally coupled to growth, for a total of 179 predicted enzymes (75.5%). Similar to the proteomics comparison, the last 58 non-predicted enzymes were also directionally coupled to growth, and their lack of prediction is consistent with the reasons mentioned in the proteomics comparison. However, the lack of prediction of enzymes involved in thiamine production is highlighted, which is because this cofactor is not included as a biomass precursor. A similar situation happens with holo-ACP synthetase (EC., an enzyme involved in the CoA hydrolysis required to synthetize acyl carrier proteins (ACPs) [67]. (See Additional file 3 for complete transcriptomic comparison results).

Qualitative comparisons between proteomic and transcriptomic data require the assumption that enzymes are active, but this depends on different factors, like post-translational modifications, allosteric control, etc. [68]. Such factors are also a problem even when quantitative omics data are used in mathematical approaches that are employed to predict phenotype states [69]. However, qualitative proteomic and transcriptomic comparisons with FBA predictions have been previously reported using E. coli K12 by Lewis et al. [70], who found predictions for up to 82% of the evaluated enzymes, and most of the unpredicted ones were isozymes, in agreement with our results.

Finally, the knockout mutant of Clostridium butyricum W5 obtained by ClosTron technology [71] was employed to evaluate the ability of iCbu641 to predict yields after perturbations in GSM. This mutant has butyric acid production blocked; its experimental yields are shown in Table 4. The wild type strain yields were predicted using FBA, while the mutant strain yields were obtained using regulatory on/off minimization (ROOM) [72]. ROOM was employed because it is better at predicting mutant phenotypic states when compared to other approaches, like minimization of metabolic adjustments (MOMA) [39]. This is because the ROOM approach seeks to maintain both the metabolic network and gene expression stabilities, as determined experimentally [72]. Simulations predicted an increased yield of ethanol using the mutant strain, and this was experimentally detected (Table 4). A biomass yield (YX/S) reduction was also predicted for the mutant strain. By contrast, the experimental reduction of hydrogen yield in the mutant was not predicted; however, this is not conclusive since the FVA ranges of the wild type and mutant strains did agree with their experimental values.

Table 4 Comparison experimental and simulated yields of wild type and mutant strains of Clostridium butyricum W5 cultured in glucose [71]

Modeling scenarios with PDO yield increment

Three scenarios were evaluated to predict an increase in YPDO/S using iCbu641. The first strategy was to use ROOM to predict single and double mutants through in silico enzyme deletion. The 145 enzymes associated with reactions directionally coupled to growth at culture conditions, as reported by Comba et al. [19], were considered in the analysis. Enzymes fully or partially coupled to growth were not considered, since their deletions would affect the culture time and therefore increase the fermentation costs. A total of 18 enzymes catalyzed at least two reactions that could block growth if simultaneously deleted, including dihydrodipicolinate reductase (EC. or proline oxidase (EC. Similar results were obtained for the 22 double enzyme deletions, such prephenate dehydrogenase (EC. and prephenate dehydrogenase (NADP) (EC. or shikimate dehydrogenase (EC. and quinate/shikimate dehydrogenase (EC. Most of the single and double deletions enhanced the YPDO/S by up to 1% in relation to the wild type strain. One of the best mutant predicted had simultaneous deletion of lactate dehydrogenases (EC. and EC., which increased YPDO/S only to nearly 1.2% (See Additional file 4 for complete mutant prediction results).

Single knockout mutants obtained from Colombian strain Clostridium sp. IBUN 158B [73] to improve the PDO production were used for validation. This is a different PDO producer strain, but research shows that the native strains currently sequenced share at least 99% of the genome (Article in preparation). Therefore, no significant differences were expected between the metabolic models of 13A and 158B, supporting the use of these mutants in validation. The inactivated enzymes were hydrogenase (ΔhydA-420 s), lactate dehydrogenase (ΔldhA-508 s), and 3-hydroxybutyryl-CoA dehydrogenase (Δhbd-414 s), corresponding the lack of production of hydrogen, lactic acid, and butyric acid, respectively. Montoya [73] reported that these three mutants were viable, as shown Table 5; however, he cultured only two of them due lack of time during his doctoral research. The biomass yield predictions for the mutants were overestimated, although a trend was predicted. This overestimation could be due to smaller glycerol uptake fluxes for the mutants, which would result in inadequate biomass yield normalization; however, the lack of measurements limits the elaboration of better comparisons.

Table 5 Comparison experimental and simulated yields of wild type and mutant strains of Clostridium sp. IBUN 158B cultured in glycerol [73]

The experimental values for PDO yields agreed with the predicted range, validating the ROOM simulations. Moreover, the FVA of the wild type strain indicated a PDO maximum flux that was 11.6% higher than the flux predicted by FBA; therefore, none of the mutants would show an increase in the PDO yield greater than this value without affecting biomass yield. This validates the single and double mutant predictions and suggests that mutant elaboration by knockout is an inadequate strategy for improving PDO yields. Similar results were obtained using Optknock, with up to 3 deletions [35, 74], where the maximum PDO yield predicted was 0.712 deleting hydrogen and butanol production. This agreed with the ROOM predictions shown in Table 5 and indicated a biomass yield reduction of nearly 28%. The Optknock maximum PDO yield value was 17.4% higher than the predicted value under glycerol limitation, but it was only 0.4% higher than the predicted value under glycerol excess, which reinforces that blocking reactions are useless. This can be understood by considering that PDO is a primary metabolite and its production is associated with growth [10, 75].

The second strategy evaluated was perturbation in the biomass composition by simultaneous random variation of stoichiometric coefficients of 44 precursors and 8 macromolecules. This strategy was evaluated since modifications in the culture media can affect biomass composition, such as accumulation of lipids during nitrogen starvation or protein accumulation with excess of nitrogen in the culture medium [76]. Normal distributions with relative standard deviations of 30% were employed for all stoichiometric coefficients considered in the perturbation. Figure 3 shows the YX/S and YPDO/S correlation values, indicating low correlation with the precursors (fatty acids, amino acids, nucleotides, polar lipids and cofactors) but a higher correlation with macromolecules, especially proteins (the main biomass component, accounting for 86.7% on a molar basis). This could suggest that a low protein content in the cell could improve YX/S and reduce the YPDO/S, due their negative and positive correlations, respectively. However, the relative standard deviations obtained for YX/S and YPDO/S were 3.2% and 0.45%, respectively, meaning that the model is sufficiently stable to perturbations in biomass composition. These predictions suggest that changes in culture media aimed at modifying biomass composition and enhancing PDO yields would be unnecessary. However, modifications in the culture media could reduce PDO yield, as Dabrock et al. found using iron in excess [77].

Fig. 3
figure 3

Pearson correlation coefficients of biomass and PDO yields to biomass precursors. Notation: (■) YX/S correlation. (■) YPDO/S correlation. Pearson values were calculated using a normalized random distribution of biomass precursors; a relative standard distribution of 30% was employed in all precursors. The covariance analyses were made with 9035 random combinations of precursor compositions

The third strategy evaluated was to use two substrates simultaneously: glucose and glycerol. Glucose is used as carbon source, while glycerol is used to maintain redox balance, therefore generating higher YPDO/S values than obtained using glycerol as single substrate [78]. Figure 4a shows the experimental values and FVA ranges, which suggest that the cell operates at optimal conditions to produce PDO when glucose is present in the medium. This validation allowed an evaluation of the YPDO/S at different glucose and glycerol uptake fluxes, as shown Fig. 4b, which shows no PDO production (YPDO/S = 0) in the absence of a glycerol uptake flux. The results also show the complete transformation of glycerol to PDO (YPDO/S = 1), using ratios of at least 0.375 between uptake fluxes of glucose and glycerol. Therefore, these predictions permit the proposal that co-fermentation is the best alternative for improving biomass and PDO yields. However, a priori prediction of the molar ratio between these substrates that could allow these flux ratios is difficult.

Fig. 4
figure 4

PDO yields using glucose and glycerol co-fermentation. a Comparison of experimental (scatter dots with standard deviation as error bars [78]) and FVA range prediction (Vertical boxes) of YPDO/S in the function of the glucose/glycerol uptake flux ratio. b Maximum YPDO/S predicted using FVA at different glucose and glycerol uptake fluxes. A ratio of glucose/glycerol uptake fluxes greater than or equal to 0.375 allows a complete conversion of glycerol to PDO (YPDO/S = 1). By contrast, PDO production is 0 without glycerol uptake flux (YPDO/S = 0)


We generated iCbu641, the first curated genome-scale metabolic model for a PDO producer Clostridium strain. During iCbu641 validation, we solved flux balance analysis using LP optimization; however, according to the experimental data, the model predicted errors of nearly 300% for biomass yield and failed to predict PDO production. Therefore, NLP optimization was employed in FBA simulations, and the new objective function maximized biomass yield per flux unit [32]. The validation allowed prediction of appropriate growth and PDO production of cultures under glycerol limitation, but it still overestimated the experimental yields of cultures under glycerol excess. Thus, sub-optimal growth predictions under glycerol excess were achieved through a second NLP optimization, where ATP minimization was added to objective function. Therefore, both objective functions were able to predict Clostridium butyricum growth and PDO production under limiting and excess glycerol conditions. Additional validations were developed using proteomics and transcriptomics data, as well data from knockout mutants, which allowed verification of the accuracy of predicting perturbations of iCbu641. All validations were completed using experimental data from different Clostridium butyricum strains and suggested that iCbu641 is an agnostic GSM model at the state steady, but the differences may be observed during dynamic predictions.

Subsequently, perturbations in the metabolic network and biomass composition were proposed to increase the PDO yield predictions. However, these perturbations predicted no significant increments. We also evaluated glucose-glycerol co-fermentation as a strategy to improve PDO yields. We found that a ratio of glucose and glycerol uptake fluxes greater than or equal to 0.375 would allow the complete conversion of glycerol to PDO; however, experimental analysis is needed to find the molar ratio that allows the achievement of this flux ratio. Finally, predictions of PDO production in state steady cultures using iCbu641 allows the proposal of this GSM model for predicting dynamic cultures (i.e. batch and fed-batch fermentations) capable of increasing PDO production, thereby minimizing the need for direct experimental efforts.


Genomic scale metabolic model iCbu641 reconstruction

The draft genome was obtained for the Colombian-native strain Clostridium sp. IBUN 13A, isolated and stored by the Bioprocesses and Bioprospecting Research Group from the Institute of Biotechnology of the Universidad Nacional de Colombia. This draft genome was previously sequenced and annotated and is available in GenBank with the accession no NZ_JZWG00000000.1 [41]. The strategy used for initial manual curation was reverse engineering proposed by Senger and Papoutsakis [47]; automated curation was also employed using GapFind and GapFill [43]. The GSM models of the solventogenic Clostridium strains C. acetobutylicum [46,47,48,49], C. thermocellum ATCC 27405 [50], C. beijerinckii NCIMB 8052 [45], C. ljungdahlii ATCC 55383 [52], and C. cellulolyticum H10 [51] were used as a database for the curation. The resulting network is based on KEGG nomenclature, whereas the SEED database [79] was used in mass and charge balances at pH 7. Finally, thermodynamically infeasible loops were eliminated according to the methodology of Schellenberger et al. [80].

The values for growth-associated maintenance (GAM) and non-growth–associated maintenance (NGAM) were reported for Clostridium acetobutylicum ATCC 824 by Lee et al. and are 40 mmol·ATP·g−1 and 5 mmol·ATP·g−1·h−1, respectively [46]. An allosteric model was also included as an upper bound constraint for the acetic acid secretion flux in the function of glycerol uptake flux. This trend was obtained using the experimental data reported by Solomon et al. [44] and Papanikolaou et al. [81]. The kinetic model was expressed as a logistic function, with 0.158 and 11.5·mmol·g−1 h−1 as initial and maximum values, respectively, and −0.0879·g·h·mmol−1 as the accumulation rate, as shown in Eq. 4.

$$ {v}_{acetic\ acid}\le \kern.75em \frac{11.5^{\ast }{0.158}^{\ast }{e}^{\left(-{0.0859}^{\ast }{v}_{glycerol}\right)}}{11.5+{0.158}^{\ast}\left({e}^{\left(-{0.0859}^{\ast }{v}_{glycerol}\right)}-1\right)} $$

Flux balance analysis

The dynamics of the mass balance of metabolite x i involved in N reactions is described in Eq. 5, where S ij is the stoichiometric coefficient of metabolite i in reaction j, and v j is the flux value in which this reaction occurs [25, 80]. Now, assuming a steady state, Eq. 5 can be expressed for M metabolites; however, since N > M, the prediction of fluxes v j can be achieved using FBA, which maximizes or minimizes an objective function Z (Equation 6). The constraints of this function are the mass balances for the M metabolites and the upper v j max, and lower v j min bounds of the N fluxes v j [28, 31]. Additionally, feasible ranges of fluxes predicted by FBA are calculated using FVA. Since the objective functions employed are non-linear, the objective function value Z calculated with FBA has to be relaxed by 5%, as suggested Mahadevan et al., Eq. 7 [82]. The mutant phenotypes were predicted using the ROOM approach, Eq. 8, where b j is the binary number of reaction j [72]. Also, v l j,wild and v u j,wild are the lower and upper confidence limits of wild type flux j. δ and ε are relative and absolute tolerance ranges, respectively.

$$ \frac{d{ x}_i}{ d t}=\sum_{j=1}^N{S}_{i j}{v}_i $$
$$ \begin{array}{c}\mathit{\operatorname{Max}}/\mathit{\operatorname{Min}}\kern1.25em Z= f\left({v}_j\right)\kern9.25em \\ {}\begin{array}{c} Subject\ to\kern13.75em \\ {}\sum_{j=1}^N{S}_{ij}{v}_j=0,\kern.75em \forall i\in 1,\dots, M\\ {}{v}_j^{min}\le {v}_j\le {v}_j^{max},\kern.75em \forall j\in 1,\dots, N\end{array}\end{array} $$
$$ \begin{array}{c}\begin{array}{c}\mathit{\operatorname{Max}}/\mathit{\operatorname{Min}}\kern1.25em {v}_j\kern15em \\ {} Subject\ to\kern15.75em \\ {}\sum_{j=1}^N{S}_{ij}{v}_j=0,\kern.3em \forall i\in 1,\dots, M\kern2.75em \end{array}\\ {}{v}_j^{min}\le {v}_j\le {v}_j^{max},\kern1.25em \forall j\in 1,\dots, N\kern.5em \\ {}{Z}_{relaxed}=\left\{\begin{array}{c}1.05 Z\kern0.5em if\ Z\ was\ minimized\\ {}0.95 Z\kern0.5em if\ Z\ was\ maximized\end{array}\right.\end{array} $$
$$ \begin{aligned} &\text{Min} \sum^{N}_{j=1}b_j \\ &\ \ Subject\ to \\ &\quad\sum^{N}_{j=1} S_{ij}v_{j,mutant}= 0, \quad \forall i\in 1, \ldots, M\ \\ & \quad v_{j,mutant}^{max}\Lambda v_{j,mutant}^{min}\left\{ \begin{array}{llll} 0 & \forall & v_j & catalyzed\ only\ by\ enzyme\ ec\\ v_{j,wild}^{max} \Lambda v_{j,wild}^{min} & \forall & v_j & \left\{ \begin{array}{l} not\ cat\ by\ enzyme\ ec\\ cat\ by\ isozyme\ ofec \end{array}\right. \end{array}\right. \\ & \quad v_{j,mutant}-b_j\left( v_{j,mutant}^{max} -v_{j,wild}^{u}\right)\leq v_{j,wild}^{u} \quad \forall\quad j\in1,\ldots, N\\ & \quad v_{j,mutant}-b_j\left( v_{j,mutant}^{min} -v_{j,wild}^{l}\right)\geq v_{j,wild}^{l} \quad \forall\quad j\in1,\ldots, N\\ & \quad v_{j,wild}^{u} =v_{j,wild}+ \delta \left|v_{j,wild}\right|+\varepsilon \qquad \forall\quad j\in1,\ldots, N\\ & \quad v_{j,wild}^{l} =v_{j,wild}- \delta \left|v_{j,wild}\right|-\varepsilon \qquad \forall\quad j\in1,\ldots, N\\ & \quad b_{j} = [0,1] \quad \forall\quad j\in1,\ldots, N \\ & \quad \delta\approx0.05\qquad \varepsilon \approx 0.001\\ & \quad Subject\ to\\ & \qquad \ FBA\ of\ wild\ type\ strain\\ \end{aligned} $$

Experimental validation

The experimental validation used two robustness analyses, as reported by Price et al. [83]: the former for the growth rate (μ) and the latter for PDO secretion flux (v PDO ), both according to the glycerol uptake flux (v Glycerol ). The robustness analyses were made using different objective functions for both glycerol limited and glycerol excess conditions. The objective function employed under glycerol limitation was biomass maximization per enzyme usage. Under glycerol excess, the biomass was maximized, while both enzyme usage and ATP production were minimized, where ATP production corresponds only to thermodynamically feasible reactions able to produce ATP. According to KEGG nomenclature, these reactions are R00156, R00158, R00200, R00315, R00332, R00512, R00570, R00722, R01512, R01547, R01665, R01688, R01724, R02090, R02093, R02094, R02098, R02326, R02331, R03005, R03035, R03530, and R03920. The predicted values were compared with the experimental values reported for the Clostridium butyricum DSM 5431 strain cultured in glycerol limited and glycerol excess conditions [44].

The prediction capability of enzyme expression was also evaluated by comparing the enzymes present in both metabolic model iCbu641 and the experimental proteome of the strain Clostridium butyricum DSM 10702 cultured in glycerol [20]. The enzymes in the model were classified as blocked or directionally, partially or fully coupled to growth using FCF, according to the methodology reported by Burgard et al. [58]. Blocked enzymes were excluded from the comparison; this comparison assumed that all the expressed enzymes were active and catalyzed some reaction. Therefore, the enzyme expression could be: a) predicted when the flux of some of the reactions catalyzed by such enzyme is different to zero; b) not predicted when all the fluxes of reactions catalyzed by the enzyme expressed in the proteome are equal to zero.

Validation was also obtained using data from experimental cultures of Clostridium butyricum strains in substrates other than glycerol, such as glucose [56] and other carbohydrates [65]. The transcriptome reported by Calusinska et al. [66] for strain C. butyricum CWBI 1009, was also used for validation; this organism had been cultured in glucose using batch fermentation with uncontrolled pH. The mRNA detected in the transcriptome coded a total of 913 enzymes (532 unique and 381 redundant), where values at the exponential growth phase were used in a similar way to those from the proteomic data.

In silico perturbations of the metabolic model

Regulatory on/off minimization (ROOM) [72] was employed as a strategy for prediction of the phenotypic states of mutants by knockout of enzymes directionally coupled to growth. In total, 145 single mutants were evaluated. Double deletion was also studied by simultaneously blocking two enzymes, excluding enzymes previously detected as essential in the single deletion; this led to evaluation of almost 7900 double mutants.

We also made perturbations in the biomass composition expressed by the variation of stoichiometric coefficients of 8 macromolecules, 7 fatty acids, 20 amino acids, 8 nucleotides, 3 polar lipids, and 6 cofactors. FBA was performed using different compositions randomly generated using normal distributions with standard deviations equivalent to 30% the values used in GSM model iCbu641. A total number of 10,000 simulations were made, where 965 of them were excluded because at least one of the concentrations was negative. Coefficients of correlation were calculated for biomass and PDO yields using the remaining 9035 biomass precursor combinations.

Technical implementation

FBA, FVA, and ROOM were computer simulated using GAMS (General Algebraic Modeling System, GAMS Development Corp., Washington, DC) software V.24.2.2 r44857 for Linux. Linear and Nonlinear Programming (LP and NLP) were developed with solver CONOPT v3.15 N, and Mixed Integer Programming (MIP) was developed with solver CPLEX Data were analyzed using Microsoft Excel® 2010.



Acyl carrier proteins


Flux balance analysis


Flux couple finding


Flux variability analysis


Growth-associated maintenance


Genome-scale metabolic


Linear programming


Mixed integer programming


Minimization of metabolic adjustments


Non-growth-associated maintenance


Non-linear programing




Regulatory On/Off Minimization


  1. Ayoub M, Abdullah AZ. Critical review on the current scenario and significance of crude glycerol resulting from biodiesel industry towards more sustainable renewable energy industry. Renew Sust Energ Rev. 2012;16(5):2671–86.

    Article  CAS  Google Scholar 

  2. EIA. International Energy Statistics: Total Biofuels Production. 2014. []. Accessed 19 May 2016.

  3. Almeida JRM, Fávaro LCL, Quirino BF. Biodiesel biorefinery: opportunities and challenges for microbial production of fuels and chemicals from glycerol waste. Biotechnol Biofuels. 2012;5 art 48:1–16.

    Google Scholar 

  4. Yazdani SS, Gonzalez R. Anaerobic fermentation of glycerol: a path to economic viability for the biofuels industry. Curr Opin Biotechnol. 2007;18(3):213–9.

    Article  CAS  PubMed  Google Scholar 

  5. Amaral PFF, Ferreira TF, Fontes GC, Coelho MAZ. Glycerol valorization: new biotechnological routes. Food Bioprod Process. 2009;87(3):179–86.

    Article  CAS  Google Scholar 

  6. Saxena RK, Anand P, Saran S, Isar J. Microbial production of 1,3-propanediol: recent developments and emerging opportunities. Biotechnol Adv. 2009;27(6):895–913.

    Article  CAS  PubMed  Google Scholar 

  7. Kaur G, Srivastava AK, Chand S. Advances in biotechnological production of 1,3-propanediol. Biochem Eng J. 2012;64:106–18.

    Article  CAS  Google Scholar 

  8. Drozdzyńska A, Leja K, Czaczyk K. Biotechnological production of 1,3-propanediol from crude glycerol. Biotechnologia. 2011;92(1):92–100.

    Article  Google Scholar 

  9. Dobson R, Gray V, Rumbold K. Microbial utilization of crude glycerol for the production of value-added products. J Ind Microbiol Biotechnol. 2012;39(2):217–26.

    Article  CAS  PubMed  Google Scholar 

  10. Kubiak P, Leja K, Myszka K, Celińska E, Spychała M, Szymanowska PD. Czaczyk K, Grajek W: physiological predisposition of various clostridium species to synthetize 1,3-propanediol from glycerol. Process Biochem. 2012;47(9):1308–19.

    Article  CAS  Google Scholar 

  11. Wilkens E, Ringel AK, Hortig D, Willke T, Vorlop KD. High-level production of 1,3-propanediol from crude glycerol by clostridium butyricum AKR102a. Appl Microbiol Biotechnol. 2012;93(3):1057–63.

    Article  CAS  PubMed  Google Scholar 

  12. González-Pajuelo M, Meynial-Salles I, Mendes F, Andrade JC, Vasconcelos I, Soucaille P. Metabolic engineering of clostridium acetobutylicum for the industrial production of 1,3-propanediol from glycerol. Metab Eng. 2005;7(5–6):329–36.

    Article  PubMed  Google Scholar 

  13. Otte B, Grunwaldt E, Mahmoud O, Jennewein S. Genome shuffling in clostridium diolis DSM 15410 for improved 1,3-propanediol production. Appl Environ Microbiol. 2009;75(24):7610–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Celińska E. Debottlenecking the 1,3-propanediol pathway by metabolic engineering. Biotechnol Adv. 2010;28(4):519–30.

    Article  PubMed  Google Scholar 

  15. Zeng AP. Pathway and kinetic analysis of 1,3-propanediol production from glycerol fermentation by clostridium butyricum. Bioprocess Eng. 1996;14(4):169–75.

    Article  CAS  Google Scholar 

  16. Bizukojc M, Dietz D, Sun J, Zeng AP. Metabolic modelling of syntrophic-like growth of a 1,3-propanediol producer, clostridium butyricum, and a methanogenic archeon, Methanosarcina mazei, under anaerobic conditions. Bioprocess Biosyst Eng. 2010;33(4):507–23.

    Article  CAS  PubMed  Google Scholar 

  17. Senger RS, Yen JY, Fong SS. A review of genome-scale metabolic flux modeling of anaerobiosis in biotechnology. Curr Opin Chem Eng. 2014;6:33–42.

    Article  Google Scholar 

  18. Dash S, Ng CY, Maranas CD. Metabolic modeling of clostridia: current developments and applications. FEMS Microbiol Lett. 2016;363(4):fnw004.

  19. Comba González N, Vallejo AF, Sánchez-Gómez M, Montoya D. Protein identification in two phases of 1,3-propanediol production by proteomic analysis. J Proteome. 2013;89:255–64.

    Article  Google Scholar 

  20. Gungormusler-Yilmaz M, Shamshurin D, Grigoryan M, Taillefer M, Spicer V, Krokhin OV, et al. Reduced catabolic protein expression in clostridium butyricum DSM 10702 correlate with reduced 1,3-propanediol synthesis at high glycerol loading. AMB Express. 2014;4:63.

    Article  PubMed  PubMed Central  Google Scholar 

  21. Kauffman KJ, Prakash P, Edwards JS. Advances in flux balance analysis. Curr Opin Biotechnol. 2003;14(5):491–6.

    Article  CAS  PubMed  Google Scholar 

  22. Orth JD, Thiele I, Palsson BO. What is flux balance analysis? Nat Biotechnol. 2010;28(3):245–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Chen Q, Wang Z, Wei DQ. Progress in the applications of flux analysis of metabolic networks. Chin Sci Bull. 2010;55(22):2315–22.

    Article  CAS  Google Scholar 

  24. Mahadevan R, Edwards JS, Doyle Iii FJ. Dynamic flux balance analysis of diauxic growth in Escherichia coli. Biophys J. 2002;83(3):1331–40.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Raman K, Chandra N. Flux balance analysis of biological systems: applications and challenges. Brief Bioinform. 2009;10(4):435–49.

    Article  CAS  PubMed  Google Scholar 

  26. Min Lee J, Gianchandani EP, Eddy JA, Papin JA. Dynamic analysis of integrated signaling, metabolic, and regulatory networks. PLoS Comput Biol. 2008;4(5):e1000086.

  27. Haggart CR, Bartell JA, Saucerman JJ, Papin JA. Whole-genome metabolic network reconstruction and constraint-based modeling. Methods Enzymol. 2011;500:411–33.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. García Sánchez CE, Vargas García CA, Torres Sáez RG. Predictive potential of flux balance analysis of Saccharomyces cerevisiae using as optimization function combinations of cell compartmental objectives. PLoS One. 2012;7(8):e43006.

  29. Westerhoff HV, Winder C, Messiha H, Simeonidis E, Adamczyk M, Verma M, et al. Systems biology: the elements and principles of life. FEBS Lett. 2009;583(24):3882–90.

    Article  CAS  PubMed  Google Scholar 

  30. Gianchandani EP, Oberhardt MA, Burgard AP, Maranas CD, Papin JA. Predicting biological system objectives de novo from internal state measurements. BMC Bioinformatics. 2008;9:43.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Knorr AL, Jain R, Srivastava R. Bayesian-based selection of metabolic objective functions. Bioinformatics. 2007;23(3):351–7.

    Article  CAS  PubMed  Google Scholar 

  32. Schuetz R, Kuepfer L, Sauer U. Systematic evaluation of objective functions for predicting intracellular fluxes in Escherichia coli. Mol Syst Biol. 2007;3:119.

    Article  PubMed  PubMed Central  Google Scholar 

  33. Ibarra RU, Edwards JS, Palsson BO. Escherichia coli K-12 undergoes adaptive evolution to achieve in silico predicted optimal growth. Nature. 2002;420(6912):186–9.

    Article  CAS  PubMed  Google Scholar 

  34. Burgard AP, Maranas CD. Optimization-based framework for inferring and testing hypothesized metabolic objective functions. Biotechnol Bioeng. 2003;82(6):670–7.

    Article  CAS  PubMed  Google Scholar 

  35. Chowdhury A, Zomorrodi AR, Maranas CD. Bilevel optimization techniques in computational strain design. Comput Chem Eng. 2015;72:363–72.

    Article  CAS  Google Scholar 

  36. Vargas García CA, García Sánchez C, Arguello Fuentes H, Torres Sáez RG. Balance de Flujos Metabólicos en Saccharomyces cerevisiae basado en Compartimentalización Intracelular. Rev Colomb Biotecnol. 2013;15(2):18-28.

  37. Feist AM, Palsson BO. The biomass objective function. Curr Opin Microbiol. 2010;13(3):344–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Gianchandani EP, Chavali AK, Papin JA. The application of flux balance analysis in systems biology. WIREs Syst Biol Med. 2010;2(3):372–82.

    Article  CAS  Google Scholar 

  39. García Sánchez CE, Torres Sáez RG. Comparison and analysis of objective functions in flux balance analysis. Biotechnol Prog. 2014;30(5):985–91.

    Article  PubMed  Google Scholar 

  40. Barragán CE, Gutiérrez-Escobar AJ, Castaño DM. Computational analysis of 1,3-propanediol operon transcriptional regulators: insights into clostridium sp. glycerol metabolism regulation. Univ Sci. 2015;20(1):129–40.

    Article  Google Scholar 

  41. Rosas-Morales JP, Perez-Mancilla X, Lopez-Kleine L, Montoya-Castano D, Riano-Pachon DM. Draft genome sequences of clostridium strains native to Colombia with the potential to produce solvents. Genome Announc. 2015;3(3):e00486–15.

    Article  PubMed  PubMed Central  Google Scholar 

  42. Montoya D, Spitia S, Silva E, Schwarz WH. Isolation of mesophilic solvent-producing clostridia from Colombian sources: physiological characterization, solvent production and polysaccharide hydrolysis. J Biotechnol. 2000;79(2):117–26.

    Article  CAS  PubMed  Google Scholar 

  43. Satish Kumar V, Dasika MS, Maranas CD. Optimization based automated curation of metabolic reconstructions. BMC Bioinformatics. 2007;8:212.

    Article  PubMed  PubMed Central  Google Scholar 

  44. Solomon BO, Zeng AP, Biebl H, Schlieker H, Posten C, Deckwer WD. Comparison of the energetic efficiencies of hydrogen and oxychemicals formation in Klebsiella pneumoniae and clostridium butyricum during anaerobic growth on glycerol. J Biotechnol. 1995;39(2):107–17.

    Article  CAS  PubMed  Google Scholar 

  45. Milne CB, Eddy JA, Raju R, Ardekani S, Kim PJ, Senger RS, et al. Metabolic network reconstruction and genome-scale model of butanol-producing strain Clostridium Beijerinckii NCIMB 8052. BMC Syst Biol. 2011;5:130.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Lee J, Yun H, Feist AM, Palsson BØ, Lee SY. Genome-scale reconstruction and in silico analysis of the clostridium acetobutylicum ATCC 824 metabolic network. Appl Microbiol Biotechnol. 2008;80(5):849–62.

    Article  CAS  PubMed  Google Scholar 

  47. Senger RS, Papoutsakis ET. Genome-scale model for clostridium acetobutylicum: part I. Metabolic network resolution and analysis. Biotechnol Bioeng. 2008;101(5):1036–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. McAnulty MJ, Yen JY, Freedman BG, Senger RS. Genome-scale modeling using flux ratio constraints to enable metabolic engineering of clostridial metabolism in silico. BMC Syst Biol. 2012;6:42.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Dash S, Mueller TJ, Venkataramanan KP, Papoutsakis ET, Maranas CD. Capturing the response of clostridium acetobutylicum to chemical stressors using a regulated genome-scale metabolic model. Biotechnol Biofuels. 2014;7:144.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Roberts SB, Gowen CM, Brooks JP, Fong SS. Genome-scale metabolic analysis of clostridium thermocellum for bioethanol production. BMC Syst Biol. 2010;4:31.

    Article  PubMed  PubMed Central  Google Scholar 

  51. Salimi F, Zhuang K, Mahadevan R. Genome-scale metabolic modeling of a clostridial co-culture for consolidated bioprocessing. Biotechnol J. 2010;5(7):726–38.

    Article  CAS  PubMed  Google Scholar 

  52. Nagarajan H, Sahin M, Nogales J, Latif H, Lovley DR, Ebrahim A, Zengler K: Characterizing acetogenic metabolism using a genome-scale metabolic reconstruction of Clostridium ljungdahlii. Microbial Cell Factories. 2013;12(1).

  53. Quilaguy-Ayure DM, Montoya-Solano JD, Suárez-Moreno ZR, Bernal-Morales JM, Montoya-Castaño D. Analysing the dhaT gene in Colombian clostridium sp. (clostridia) 1,3-propanediol-producing strains. Univ Sci. 2010;15(1):17–26.

    Article  CAS  Google Scholar 

  54. Cárdenas DP, Pulido C, Aragón OL, Aristizábal FA, Suárez ZR, Montoya D. Evaluación de la producción de 1,3-propanodiol por cepas nativas de Clostridium sp. mediante fermentación a partir de glicerol USP y glicerol industrial subproducto de la producción de biodiésel. Revista Colombiana De Ciencias Químico Farmacéuticas. 2006;35(1):120–37.

    Google Scholar 

  55. Von Stockar U, Maskow T, Liu J, Marison IW, Patiño R. Thermodynamics of microbial growth and metabolism: an analysis of the current situation. J Biotechnol. 2006;121(4):517–33.

    Article  CAS  PubMed  Google Scholar 

  56. Cai G, Jin B, Saint C, Monis P. Metabolic flux analysis of hydrogen production network by clostridium butyricum W5: effect of pH and glucose concentrations. Int J Hydrog Energy. 2010;35(13):6681–90.

    Article  CAS  Google Scholar 

  57. Saint-Amans S, Girbal L, Andrade J, Ahrens K, Soucaille P. Regulation of carbon and electron flow in clostridium butyricum VPI 3266 grown on glucose-glycerol mixtures. J Bacteriol. 2001;183(5):1748–54.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Burgard AP, Nikolaev EV, Schilling CH, Maranas CD. Flux coupling analysis of genome-scale metabolic network reconstructions. Genome Res. 2004;14(2):301–12.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Ladygina N, Dedyukhina EG, Vainshtein MB. A review on microbial synthesis of hydrocarbons. Process Biochem. 2006;41(5):1001–14.

    Article  CAS  Google Scholar 

  60. Covert MW, Schilling CH, Palsson B. Regulation of Gene expression in flux balance models of metabolism. J Theor Biol. 2001;213(1):73–88.

    Article  CAS  PubMed  Google Scholar 

  61. Feng X-M, Cao L-J, Adam RD, Zhang X-C, Lu S-Q. The catalyzing role of PPDK in Giardia lamblia. Biochem Biophys Res Commun. 2008;367(2):394–8.

    Article  CAS  PubMed  Google Scholar 

  62. Wood HG, O'Brien WE, Michaels G. Properties of carboxytransphosphorylase; pyruvate, phosphate dikinase; pyrophosphate-phosphofructikinase and pyrophosphate-acetate kinase and their roles in the metabolism of inorganic pyrophosphate. Adv Enzymol Relat Areas Mol Biol. 1977;45:85–155.

    CAS  PubMed  Google Scholar 

  63. Vinitsky A, Grubmeyer C. A new paradigm for biochemical energy coupling. Salmonella typhimurium nicotinate phosphoribosyltransferase. J Biol Chem. 1993;268(34):26004–10.

    CAS  PubMed  Google Scholar 

  64. Al-Hinai MA, Jones SW, Papoutsakis ET. The clostridium sporulation programs: diversity and preservation of endospore differentiation. Microbiol Mol Biol Rev. 2015;79(1):19–37.

    Article  PubMed  PubMed Central  Google Scholar 

  65. Junghare M, Subudhi S, Lal B. Improvement of hydrogen production under decreased partial pressure by newly isolated alkaline tolerant anaerobe, clostridium butyricum TM-9A: optimization of process parameters. Int J Hydrog Energy. 2012;37(4):3160–8.

    Article  CAS  Google Scholar 

  66. Calusinska M, Hamilton C, Monsieurs P, Mathy G, Leys N, Franck F, et al. Genome-wide transcriptional analysis suggests hydrogenase- and nitrogenase-mediated hydrogen production in clostridium butyricum CWBI 1009. Biotechnol Biofuels. 2015;8:27.

    Article  PubMed  PubMed Central  Google Scholar 

  67. Mootz HD, Finking R, Marahiel MA. 4′-Phosphopantetheine transfer in primary and secondary metabolism of Bacillus Subtilis. J Biol Chem. 2001;276(40):37289–98.

    Article  CAS  PubMed  Google Scholar 

  68. Winter G, Krömer JO. Fluxomics - connecting 'omics analysis and phenotypes. Environ Microbiol. 2013;15(7):1901–16.

    Article  CAS  PubMed  Google Scholar 

  69. Machado D, Herrgård M. Systematic evaluation of methods for integration of Transcriptomic data into constraint-based models of metabolism. PLoS Comput Biol. 2014;10(4):e1003580.

  70. Lewis NE, Hixson KK, Conrad TM, Lerman JA, Charusanti P, Polpitiya AD, et al. Omic data from evolved E. coli are consistent with computed optimal growth from genome-scale models. Mol Syst Biol. 2010;6:390.

    Article  PubMed  PubMed Central  Google Scholar 

  71. Cai G, Jin B, Saint C, Monis P. Genetic manipulation of butyrate formation pathways in clostridium butyricum. J Biotechnol. 2011;155(3):269–74.

    Article  CAS  PubMed  Google Scholar 

  72. Shlomi T, Berkman O, Ruppin E. Regulatory on/off minimization of metabolic flux changes after genetic perturbations. Proc Natl Acad Sci U S A. 2005;102(21):7695–700.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Montoya Solano JD: Metabolic engineering of the Colombian strain Clostridium sp. IBUN 158B in order to improve the bioconversion of glycerol into 1,3-propanediol. Germany: University of Ulm; 2012.

  74. Burgard AP, Pharkya P, Maranas CD. OptKnock: a Bilevel programming framework for identifying Gene knockout strategies for microbial strain optimization. Biotechnol Bioeng. 2003;84(6):647–57.

    Article  CAS  PubMed  Google Scholar 

  75. Chatzifragkou A, Aggelis G, Gardeli C, Galiotou-Panayotou M, Komaitis M, Papanikolaou S. Adaptation dynamics of clostridium butyricum in high 1,3-propanediol content media. Appl Microbiol Biotechnol. 2012;95(6):1541–52.

    Article  CAS  PubMed  Google Scholar 

  76. Lari Z, Moradi-kheibari N, Ahmadzadeh H, Abrishamchi P, Moheimani NR, Murry MA. Bioprocess engineering of microalgae to optimize lipid production through nutrient management. J Appl Phycol. 2016;28(6):3235–50.

    Article  CAS  Google Scholar 

  77. Dabrock B, Bahl H, Gottschalk G. Parameters affecting solvent production by clostridium pasteurianum. Appl Environ Microbiol. 1992;58(4):1233–9.

    CAS  PubMed  PubMed Central  Google Scholar 

  78. Malaoui H, Marczak R. Influence of glucose on glycerol metabolism by wild-type and mutant strains of clostridium butyricum E5 grown in chemostat culture. Appl Microbiol Biotechnol. 2001;55(2):226–33.

    Article  CAS  PubMed  Google Scholar 

  79. Henry CS, Dejongh M, Best AA, Frybarger PM, Linsay B, Stevens RL. High-throughput generation, optimization and analysis of genome-scale metabolic models. Nat Biotechnol. 2010;28(9):977–82.

    Article  CAS  PubMed  Google Scholar 

  80. Schellenberger J, Lewis NE, Palsson BØ. Elimination of thermodynamically infeasible loops in steady-state metabolic models. Biophys J. 2011;100(3):544–53.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  81. Papanikolaou S, Ruiz-Sanchez P, Pariset B, Blanchard F, Fick M. High production of 1,3-propanediol from industrial glycerol by a newly isolated clostridium butyricum strain. J Biotechnol. 2000;77(2–3):191–208.

    Article  CAS  PubMed  Google Scholar 

  82. Mahadevan R, Schilling CH. The effects of alternate optimal solutions in constraint-based genome-scale metabolic models. Metab Eng. 2003;5(4):264–76.

    Article  CAS  PubMed  Google Scholar 

  83. Price ND, Reed JL, Palsson BØ. Genome-scale models of microbial cells: evaluating the consequences of constraints. Nat Rev Microbiol. 2004;2(11):886–97.

    Article  CAS  PubMed  Google Scholar 

  84. Chelliah V, Juty N, Ajmera I, Ali R, Dumousseau M, Glont M, et al. BioModels: ten-year anniversary. Nucleic Acids Res. 2015;43(D1):D542–8.

    Article  PubMed  Google Scholar 

Download references


The authors would also like to thank Satyakam Dash and Juan Rosas for helpful discussions and suggestions. Finally, iCbu641 GSM model was deposited in BioModels [84] and assigned the identifier MODEL1704210000.


This work was supported by COLCIENCIAS, Colombia.

Availability of data and materials

All data generated or analyzed during this study are included in this published article [and its Additional files].

Authors’ contributions

LMSB performed simulations, analyzed data, and wrote the paper. LMSB, AFGB, and DMC conceived the study and participated in its design and coordination. CDM helped with key analyses of the data. All authors read and approved the final manuscript.

Competing interests

The authors declare that they have no competing interest.

Consent for publication

Not applicable.

Ethics approval and consent to participate

Not applicable.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Author information

Authors and Affiliations


Corresponding author

Correspondence to Dolly Montoya.

Additional files

Additional file 1:

Excel format of iCbu641Genome Scale Metabolic Model. (XLSX 203 kb)

Additional file 2:

SBML format of iCbu641Genome Scale Metabolic Model. (XML 2250 kb)

Additional file 3:

Proteomic and transcriptomic comparison results. (XLSX 30 kb)

Additional file 4:

Results of modeling scenarios with PDO yield increment. (XLSX 5507 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Serrano-Bermúdez, L.M., González Barrios, A.F., Maranas, C.D. et al. Clostridium butyricum maximizes growth while minimizing enzyme usage and ATP production: metabolic flux distribution of a strain cultured in glycerol. BMC Syst Biol 11, 58 (2017).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: