Process Engineering of the Acetone-Ethanol-Butanol (ABE) Fermentation in a Linear and Feedback Loop Cascade of Con- tinuous Stirred Tank Reactors: Experiments, Modeling and Op- timization

The production of butanol, acetone and ethanol by Clostridium acetobutylicum is a biphasic fermentation process. In the first phase the carbohydrate substrate is metabolized to acetic and butyric acid, in the following second phase the product spectrum is shifted towards the economically interesting solvents. Here we present a cascade of six continuous stirred tank reactors (CCSTR), which allows performing the time dependent metabolic phases of an acetone-butanol-ethanol (ABE) batch fermentation in a spatial domain. Experimental data of steady states under four operating conditions with variations of the pH in the first bioreactor between 4.3 and 5.6 as well as the total dilution rate between 0.042 h-1 and 0.092 h-1 were used to optimize and validate a corresponding mathematical model. Beyond a residence time distribution representation and substrate, biomass and product kinetics this model also includes the differentiation of cells between the metabolic states. Model simulations predict a final product concentration of 8.2 g butanol L-1 and a productivity of 0.75 g butanol L-1 h-1 in the CCSTR operated at pHbr1 of 4.3 and D = 0.092 h-1, while 31 % of the cells are differentiated to the solventogenic state. Aiming at an enrichment of solvent-producing cells, a feedback loop was introduced into the cascade sending cells from a later state of the process (bioreactor 4) back to an early stage of the process (bioreactor 2). In agreement with the experimental observations, the model accurately predicted an increase of butanol formation rate in bioreactor stages 2 and 3, resulting in an overall butanol productivity of 0.76 g L-1 h-1 for the feedback loop cascade. The here presented CCSTR and the validated model will serve to investigate further ABE fermentation strategies for a controlled metabolic switch.


Introduction
The acetone-butanol-ethanol (ABE) process is the production of the eponymous solvents by fermentation with solvent-forming Clostridia. In 1862 Louis Pasteur discovered motile, rod-shaped bacteria to produce butanol under anaerobic conditions [1,2]. It was 1915, when Chaim Weizmann established an industrial process for the biotechnological production of acetone (and butanol), which was used all around the world in the following decades before being replaced by synthesis of these solvents from mineral oil [3,4]. Nowadays we await a comeback of the industrial ABE process as not only the demand for sustainable produced solvents and bulk chemicals is increasing, but butanol is further discussed as promising biofuel. Beside the exploration of new non-food and inexpensive waste materials as substrates for the fermentation process [5,6], an important setscrew to improve the economics of the process is the development of new fermentation technologies [7,8]. Those should allow effective raw material conversion with high yield, accelerate the onset of solvent formation in a biphasic ABE process and stabilize the bacteria in this metabolic phase. While in the first phase of the ABE fermentation, called acidogenesis, acetic and butyric acid are the main fermentation products, the formation of butanol, acetone and ethanol under partial re-assimilation of the organic acids occurs in the second phase, called solventogenesis. The research of the past decades has shown that low pH and/or the high concentration of the undissociated organic acids trigger the metabolic switch from acidogenesis to solventogenesis [9,10]. Investigations with external electron donors, such as methyl viologen or neutral red, have further shown that the redox state of the cells influences the bias for acidogenesis and solventogenesis, respectively [11][12][13]. Moreover, it has been proposed that small concentrations of butanol might promote the metabolic shift [14]. Clostridium acetobutylicum was used as model organism in the majority of the studies investigating the metabolic switch. However, under certain conditions C. acetobutylicum as well as other clostridial species fails to switch its metabolism to the solventogenic state, resulting in a hyper-acidification of the fermentation broth (pH ≤ 3.8) and a complete arrest of metabolic activity. The molecular mechanism behind this phenomenon often described as "acid crash" is still open. In some cases the spontaneous loss of the mega-plasmid pSOL1, encoding genes essential for solvent production, is a sound explanation. However, how this event can occur repeatedly under certain conditions, i.e. low buffer concentrations, while it rarely occurs under other conditions, i.e. high buffer concentrations, remains to be clarified. Explaining and controlling the metabolic shift has been the aim of numerous studies in the past [15]. In batch-experiments, controlling the pH at 4.5 to 5.0 has been described as successful strategy for a stable solventogenic phase achieving high butanol concentrations [10]. However, onset of solvent formation is observed not earlier than 15 h after start of fermentation under such conditions. The maximal product concentrations are only reached after about 36 h [10,16], which limits the productivity of the batch fermentation processes. Additions of acetate and butyrate have shown to increase the productivity of batch processes [17]. However, determining the optimal fermentation phase for acetate (or butyrate) co-feeding is challenging and laborious in a batch process as precise reproducibility of the time profile is limited. We here developed a model system that captures the metabolic phases of the ABE fermentation in different stages of a continuously operated bioreactor cascade at steady state. On the one hand this system allows analyzing the physiological characteristics of cells in the different metabolic phases, including the otherwise only transient occurring phases. On the other hand we aim at accelerating the solvent production and hence increasing the productivity of the fermentation process by varying the operating conditions of the cascade. In order to reduce the experimental work for this second task, we started developing a mathematical model, that describes the continuous, multi-stage ABE fermentation process taking into account the presence of different subpopulation types in each bioreactor tank. Using experimental data from four operating conditions of the linear cascade, the unknown kinetic parameters of the model were estimated and the predictive capacity of the steady state simulations was validated.
To further validate the model and to test a potential strategy for faster solvent formation, the simulator as well as the experimental set-up was modified by introducing a feedback loop of from bioreactor 4 to bioreactor 2 of the cascade. In that way fermentation broth with solventogenic cells and a higher butanol concentration is introduced into an early phase of the ABE process.

Microorganism and Medium
The strain used in this study was Clostridium acetobutylicum DSM 792 obtained from the DSMZ Germany. All fermentations were started from cryo-preserved spore aliquots, which were exposed to a 10 min heat shock at 75 °C directly before germination. First pre-cultures were cultivated in Clostridia Growth Medium (CGM) [18] at a scale of 50 mL in a 500 mL bottle with gas outlet. After 20 to 40 hours incubation at 37 °C, 5 mL of this first pre-culture were used to set-up a second pre-culture in the scale of 50 mL containing 25 g L -1 glucose, 1 g L -1 KH2PO4, 5 g L -1 yeast extract, 0.2 g L -1 MgSO4x7H2O, 0.01 g L -1 FeSO4x7H2O and 0.01 g L -1 MnSO4xH2O The medium was adjusted to pH 6.9. After incubation for 13 hours at 37 °C, the second pre-culture was used to start a batch fermentation with a total volume of 0.45 L in the first bioreactor of the cascade (Figure 1). The medium for the batch fermentation was the same as in the second pre-culture except for a lower level of glucose of about 5 g L -1 and a pH of 5.6. The cultivation temperature in the cascade was controlled at 32 °C during the entire process. The phosphate-limited feed medium of the continuous process contained 60 g L -1 glucose, 0.1 g L -1 KH2PO4, 5 g L -1 yeast extract, 0.2 g L -1 MgSO4 x7H2O, 0.01 g L -1 FeSO4 x7H2O and 0.01 g L -1 MnSO4 xH2O. To prevent contamination, the medium in the feed bottle was adjusted to pH 2.5 with 32 % HCl prior to autoclaving. During the process the pH was adjusted with sterile 1.5 M KOH. The feed medium was further supplemented with 0.1 % (v/v) Polypropylene glycol P2000 (Merck KGaA, Darmstadt/Germany) to prevent foam formation in the cascade. All media were gassed with nitrogen for at least 1 hour to ensure anaerobic conditions. In the cascade the fermentation broth was pressurized with nitrogen during the entire process.

Continuous Fermentation in a Six-Stage Bioreactor Cascade
The continuous bioreactor cascade consisted of six linearly coupled stirred tanks( Figure 1A), allowing for additional modifications of the configuration such as feedback loops ( Figure 1B). For starting the process, a batch fermentation in the first bioreactor was performed until the biomass had reached an OD600 of 2, then the system was switched to the continuous mode with a feeding rate of 0.1 or 0.22 L h -1 (Table 1). Now being in continuous mode, the liquid volume in the first bioreactor was controlled at 0.4 L and bioreactor stages 2 to 6 were subsequently filled with fermentation broth (0.4 L each) within 9 h or 20 h, depending on the feeding rate. By daily sampling we monitored the establishment of a steady state throughout the cascade of continuously stirred tank reactors (CCSTR) (see supplement Figure SIa). When detecting constant metabolite and cell concentrations for at least three residence times, corresponding to 33 h and 72 h at the given feeding rates (compare Table 1), one of the operating parameters -the setpoint pH value in the first bioreactor pHbr1 or the feeding rate Fin -was modified. After an adaption period of three residences times, sampling of the new steady state was started. For the configuration with feedback loop ( Figure 1B), fermentation broth of bioreactor 4 was pumped to bioreactor 2 with a rate Fback of 0.055 L h -1 . In total we performed six fermentations lasting up to 43 days. Means and standard deviations of the steady state concentrations along the cascade were summarized in Table SIIa.

Analytical Procedures
Acids and solvents were analyzed with GC-FID (GC-14B, Shimadzu Inc., Kyoto Japan) equipped with a ZB-FFAP column (Phenomenex Inc., Torrance/USA). Glucose was measured with the Reflectoquant Glucose Test (Merck KGaA, Darmstadt/Germany). Optical density was obtained from a UV-VIS-spectrophotometer at a wave length of 600 nm. The cell dry mass was calculated from the OD600 using the experimentally determined factor of 0.4 g L -1 OD600 -1 [6]. For monitoring the morphology of cells and sterile control, a microscope (Axiostar plus, Zeiss, Göttingen Germany) with 100x magnification was used. Sporadically the phosphate concentrations were measured along the cascade using the phosphate analysis method DIN EN ISO6878 2004/09.

Mathematical Modeling Approach
Our approach for simulations of the continuous, multi-stage ABE fermentation is based on an unstructured segregated model. Thus, while assuming empiric kinetic laws for the production and consumption of extracellular metabolites and cell growth without taking into consideration intracellular metabolites, we distinguished between three biomass subpopulations with individual metabolic properties, namely acidogenic (A), intermediate (I) and solventogenic (S) cells. A first version of the simulator, implemented by object-orientated programming in MATLAB R2011b, has been described [19]. The development of the metabolite concentrations in each bioreactor was simulated by iterative calculations of incremental changes of the metabolite concentrations in each bioreactor stage (eq. (1)), the development of the subpopulations (eq. (2)) and the specific metabolite conversion rates of the subpopulations and the total biomass depending on the present metabolite concentrations (eq. (6)-(12)) until a steady state is reached ( Figure 2). A steady state was assumed when the concentrations of biomass (x), glucose (glu), ethanol (eth), acetone (act), butanol (but), acetic acid (aa), butyric acid (ba) and potassium dihydrogen phosphate (PO4) changed less than 0.01 g L -1 in each bioreactor stage for the residence time of one bioreactor stage Dbr1 -1 . In comparison to the earlier model [19], the following modifications were introduced: (I) Addition of a butanol-and glucose-dependent lysis term, accounting for subpopulation-independent cell lysis in the presence of high butanol concentrations and glucose limitation (see eq. (7)); (II) Allowing growth of intermediate and solventogenic cells in addition to growth in the acidogenic phase described before; (III) Flexibilizing the carbon balance by introducing the parameter sk, which accounts for possible additional carbon sources introduced by the yeast extract or the flux of carbon into additional products (see eq. (12)). The quality of the simulation result was evaluated by calculating the distance of the simulated steady state concentrations from experimentally observed metabolite concentrations as weighted re-sidual square sum over all metabolites (except for PO4) and all bioreactor stages, hereafter referred as p (eq. in Figure 2). The simulated steady state depends on several inputs, which have to be assigned before starting the simulation ( Figure 2). The set-up parameters, i.e. the number of bioreactor, the culture volume in each bioreactor stage and the dilution rate in one bioreactor stage Dbr1, were the operating conditions in the experiment. The initial metabolite concentrations were taken from the experimental data set (Data in Table SIIa). However, these initial values only affected whether the simulation resulted in a stable steady state or not, in the steady state itself, the corresponding metabolite concentrations were independent from the applied initial metabolite concentrations. During the simulations, the pH values along the cascade were set to the experimentally observed pH values. Modifications of the pH values in the bioreactor stages had an influence on the simulated steady state concentrations. The feed concentration for glucose and potassium dihydrogen phosphate were 60 g L -1 and 0.1 g L -1 , respectively. The model equations for reaction kinetics are given below (eq. (3)-(5), (7)-(12)). A part of the kinetic model parameters were set to fixed values derived from literature, experimental investigations or estimations ( Table SIc). The remaining kinetic model parameters (Table SIb) were estimated in a parameter optimization process described in the results section. The optimization was aimed at minimizing the p-value (eq. in Figure 2). Table 2 gives an overview over all (fixed and optimized) kinetic parameters applied to obtain the simulation results presented in this work (Data in Table SIIb).

Input
Initial parameters: Kinetic model and kinetic parameters:  Table SIc) kinetic model parameters for estimation (suppl. Table SIb) Simulation of metabolite concentrations cM and population fractions wP The incremental changes of the metabolite concentrations over time in each bioreactor stage were calculated considering the in-and out-flux of the respective metabolite and its conversion by the present biomass (eq.(1), an overview of all variables and parameters of the model is given in supplementary tables SIa, SIb and SIc): The changes in the size of a biomass-subpopulation xP in a bioreactor stage k were calculated on the basis of the in-and out-flux of respective cells, the growth of the respective subpopulation and the differentiation processes d (Eq.(2)): By Eqs.  Table 2.) The specific metabolite conversion rate of the total biomass in one bioreactor stage was calculated by Eq.(6) as sum of the specific metabolite conversion rates of the individually sized subpopulations delivered by the kinetic model (contributing expressions see Eqs. The following kinetic model equations (7)-(12) for biomass evolution, product formation and substrate uptake have been assumed: For the calculation of growth-dependent product formation and substrate consumption, only the growth-related part, but not the cell lysis-related part of the biomass evolution rate from Eq. (7) was taken into consideration in Eq. (8). Model parameters in Table 2, presented in italic were fixed to the given values according to literature information, experimental data or estimations (see also supplementary material Table SIc). All other parameters were retrieved from the parameter estimation process (see also Table SIb).
acid production and uptake r_aa_max_ (P)

Continuous fermentations in a linear CCSTR under four different operating conditions
ABE fermentations have been investigated for many decades now. However no ideal operating conditions enabling high space-time yields and long-term stability of the process have been found yet. In order to analyze in which way the switch from acid forming to solvent forming cells can be accelerated, we established a Cascade of Continuous Stirred Tank Reactors (CCSTR). This bioreactor conformation with six reactors operated in series allows separating the different metabolic stages in a spatial manner. In the initial set-up, two operating parameters were investigated: the feeding rate Fin and the setpoint pH value in the first bioreactor of the cascade pHbr1 (see Figure 1A). Long-term continuous fermentations in the bioreactor cascade were performed in order to characterize steady states with constant metabolite concentrations in all bioreactor tanks under four different operating conditions. To this end, the feeding rate was set to 0.1 L h -1 (D = 0.042 h -1 ) and to 0.22 L h -1 (D = 0.092 h -1 ) and the pHbr1 was set to pH 5.6 and pH 4. resulting four operating settings were repeatedly applied and the corresponding steady states were observed for periods corresponding to 3 to 18 residences times of the cascade (means of the reproduced steady states under the four operating conditions are shown in Figure 3; for an overview of the individually observed steady states see Table SIIa). Under all operating conditions we observed an accumulation of acetic and butyric acid in the first bioreactor stage ( Figure 3A, C, E and G). In the following bioreactors the organic acids were then metabolized together with glucose, shifting the product spectrum towards the organic solvents -butanol, acetone and ethanol. These dominate the metabolite spectrum at the end of the cascade, though their yields depend on the operating conditions ( Figure 3A, C, E and G). At a dilution rate of 0.042 h -1 , corresponding to a residence time of 4 h per bioreactor, and a pHbr1 of 5.6, the acetic and butyric acid concentration in the first bioreactor reached 10.3 ± 4.2 g L -1 and 8.8 ± 3.6 g L -1 , respectively ( Figure 3A). Despite these high acid concentrations the biomass grew to a concentration of 2.0 ± 0.8 g L -1 and stayed on this level along the cascade ( Figure 3B). The solvent concentrations are negligible in the first bioreactor, but increase from bioreactor 2 to 4. When the glucose concentration dropped below 5 g L -1 in bioreactor 5 and 6, the solvent formation stopped, yielding final concentrations of 6.7 ± 3.3 g L -1 butanol, 1.8 ± 0.9 g L -1 acetone and 0.4 ±0.3 g L -1 ethanol ( Figure 3A). Changing the pHbr1 to 4.3 while keeping the total dilution rate at 0.042 h -1 , reduced the accumulation of acetic and butyric acid in the first bioreactor to 7.5 ± 1.8 g L -1 and 3.1 ± 1.2 g L -1 , respectively ( Figure 3C). At the same time, the formation of the solvents as well as consumption of the glucose over the cascade was increased. Thus maximal concentrations of butanol (9.7 ± 0.5 g L -1 ), acetone (2.9 ± 0.4 g L -1 ) and ethanol (1.0 ± 0.2 g L -1 ) were already reached in bioreactor 3, while glucose was nearly depleted at this point ( Figure 3C). The glucose depletion was accompanied by a severe decrease of the biomass concentration from 3.0 ± 0.4 g L -1 in bioreactor 2 to 0.5 ± 0.2 g L -1 in bioreactor 6 ( Figure 3D).
To get a better resolution of the metabolic processes in bioreactor 1 to 3 we increased the dilution date to 0.092 h -1 which corresponds to a residence time of 1 h 49 min per bioreactor. Again the pHbr1 was set to 5.6 ( Figure 3E and F) and 4.3 ( Figure 3G and H). Under both conditions the acetic acid and butyric acid concentrations rose in bioreactor 2 (and 3) and a net uptake of the organic acids only appeared in the later bioreactors of the cascade. Thus maximal organic acids concentrations were observed in bioreactor 3 and 4 at pHbr1 = 5.6 ( Figure 3E , Table SIIa) and in bioreactor 2 at pHbr1 = 4.3 ( Figure 3G, Table SIIa), respectively. The minima of the pH profiles coincide with the maxima of the acid concentrations (see Figure 3F and H). Despite of this discontinuous development of organic acids along the cascade there was a formation of solvents from bioreactor 1 to bioreactor 6 as well as a nearly linear uptake of glucose under both pHbr1 conditions ( Figure 3E and G). At pHbr1 = 5.6 solvent formation and glucose uptake were slower than at a pHbr1 of 4.3 and resulted in concentrations of 4.8 ± 1.9 g L -1 butanol, 1.6 ± 0.6 g L -1 acetone, 0.3 ±0.1 g L -1 ethanol and 25.9 ±9.9 g L -1 glucose in bioreactor 6 ( Figure 3E). This roughly corresponds to the composition of bioreactor 3 in the cascade operated at a dilution rate of 0.042 h -1 and a pHbr1 of 5.6 ( Figure 3A, Table SIIa). In the cascade operated at pHbr1 = 4.3 (and D = 0.092 h -1 ) butanol, acetone, ethanol and glucose reached final concentrations of 8.2 ± 2.6 g L -1 , 2.8 ± 0.9 g L -1 , 0.8 ±0.4 g L -1 and 13.9 ±8.3 g L -1 , respectively ( Figure 3G). These solvent concentrations were significantly lower than those observed in bioreactor 3 of the cascade operated at D = 0.042 h -1 and pHbr1 = 4.3 amounting to.7 ± 0.5 g L -1 butanol, 2.9 ± 0.4 g L -1 acetone and 1.0 ±0.2 g L -1 ethanol ( Figure 3C, Table SIIa). However, the remaining glucose concentration of 13.9 ±8.3 g L -1 in the final bioreactor of the process at D = 0.092 h -1 which is observed after a total mean residence time of 10.9 h is also much higher than the glucose concentration of 2.5 ±3.4 g L -1 in bioreactor 3 of the process at D = 0.042 h -1 that is observed after a mean residence time of 12 h (3x 4 h).
These results indicate that process operating conditions were in a range, in which the fermentation time (residence time) is crucial for a complete transformation of the substrate to solvents. The fourth data set (D = 0.092 h -1 / pHbr1 = 5.6) was reserved for model validation. In total 30 parameters were estimated within rationally predefined boundaries ( Table 2, Table SIb). At the same time, the model included 25 fixed parameters, which had been derived from literature, experimental studies or calculations ( Table 2, Table SIc). As starting point of the optimization process we used the parameter set that had been found by educated guessing to deliver reasonable simulation results for a first version of the mathematical model [19]. The parameter for the maximal butyric acid uptake rate of intermediate cells k_BA_up_ (I) was initially set to 0.25 corresponding to the center of the solution space ( Table SIb). The parameter optimization was realized using the MATLAB function "fmincon" (with "active set algorithm"), which allows defining dependencies between the parameters. Here we used the restriction "mu_d_ (A) (differentiation rate of acidogenic to intermediate cells)" < "mu_d_ (I) (differentiation rate of intermediate to solventogenic cells)" to prevent accumulation of intermediate cells observed in first attempts of the parameter optimization without constraints. The estimation algorithm minimized the distance between the simulated and the experimentally observed metabolite (glucose, ethanol, acetone, butanol, acetic acid and butyric acid) and biomass concentrations in all six bioreactor stages of all three operating conditions. The differences between simulated and experimental values were weighted by the experimentally observed variances, which we defined as maximum of the standard deviations in the individual data sets, the standard deviation of the averaged data set and 0.1 (compare Table SIIa -note  § ). Thus, the predictions of the model for the experimentally accurately determinated solvent concentrations counted more for the quality of the respective parameter set than a good agreement of the model with data from the fluctuating measurements of organic acid concentrations. Furthermore, the discrepancy of a mismatch of 1 g L -1 in e.g. the modeled and experimentally measured glucose concentrations (very small difference) and of 1 g L -1 in e.g. the ethanol concentrations (very big difference) is balanced by the weighting with the experimental variance of the respective metabolite concentration.
The improvements in the p-value became negligible after 29 cycles of the fmincon function, and the optimization algorithm was stopped. The resulting parameter set OPT10e (see Table 2 and also Table SIb) was used to simulate the metabolite concentrations, shown on the y-axes of Figure 4 (see also Table SIIb). These were compared with the experimentally observed metabolite concentrations, shown at the x-axes of Figure 4 (compare Figure 3, Table SIIa). The closer a data point is to the diagonal (y=x) in this representation the better is the simulation of this metabolite concentration. The error bars correspond to the experimental variance, i.e. the larger the error bars the less important is the respective metabolite concentration for the quality of the parameter set (p-value). Or simply stated, as long as the error bars cross the y=x diagonal, it is still a good simulation result.
The simulations of the solvent concentrations are in very good agreement with the experimental data for all operating conditions, including those reserved for validation ( Figure 4A, D, G and J). Also the moderate variations of the organic acids concentrations, observed in the cascade at D = 0.092 h -1 are reproduced by the simulations (Figure 4H and K). However, the important changes in the acids concentrations along the cascade at D = 0.042 h -1 are not covered by the simulations (Figure 4B and E). Also, the glucose concentrations are not simulated accurately with respect to the experimental variance in glucose measurements ( Figure 4C, F, I and L). Still, the overall tendencies of the glucose concentrations along the cascade are qualitatively represented in the simulations. As shown in Figure 5, there is a high level of agreement between simulated and experimentally observed biomass concentrations for all operating conditions. The discrepancies in the glucose concentrations are the main contributor to the weighted residual square sums, depicted in Table 3. Globally seen, the simulations had an equal quality for all bioreactor stages and all operating conditions (Table 3). Particularly good results were retrieved for the simulation of the cascade at D = 0.092 h -1 and pHbr1 = 4.3. It should be mentioned that also the experimental data set for D = 0.092 h -1 / pHbr1 = 5.6), which has not been used in parameter estimation, was simulated with comparable quality, validating the model with its parameter values. Table 3 -Quality of simulation for individual bioreactor stages and data sets given by the weighted residual square sum over all metabolites p (see eq. (13)). Data set of D=0.092 h -1 /pHbr1=5.6 was not used for parameter optimization and serves for validation of the estimated parameter set.

Model-based predictions of biomass subpopulations regarding acidogenic, intermediate and solventogenic states
During the experiments, the total biomass concentration in each bioreactor stage was monitored. However, contrary to the long spread hypothesis that metabolically distinct states strictly correlate with different morphologies of Clostridium cells, more recent studies have clarified that solvent-forming cells can also appear with rod-shaped morphologies [16,20]. Thus, a straightforward microscopical characterization of the metabolic states of cells along the bioreactor cascade is not possible. Therefore, our model, postulating biomass subpopulations of acidogenic, intermediate and solventogenic states in each bioreactor ( Figure 5), may provide insights into experimentally so far inaccessible process characteristics. For fermentations with a pHbr1 of 5.6, the model predicts almost exclusively acidogenic cells in the first bioreactor of the cascade ( Figure 5 A and C). While at a dilution rate of 0.042 h -1 , the proportion of solventogenic cells rises to 42 % in bioreactor 6, the maximal proportion of solventogenic cells at a dilution rate of 0.092 h -1 (after a total mean residence time of 10.9 h) is only 26 %. However, this exceeds the ratio of 21 % solventogenic cells after 12 h mean residence time, predicted for the 3 rd bioreactor of the cascade at 0.042 h -1 ( Figure 5A and C, Table SIIb).
In fermentations with a pHbr1 of 4.3, the model predicts significant cellular differentiation occurring already in the 1 st bioreactor of the cascade ( Figure 5B and D). At a dilution rate of 0.042 h -1 the predicted proportion of solventogenic cells is 19 % in bioreactor 1 and increases to 36 % along the cascade ( Figure 5B, Table SIIb). It is striking that under the respective operating conditions, the part of the intermediate cells is very low in all bioreactor stages. Its maximum of 18 % is predicted for bioreactor 1 ( Figure 5B, Table SIIb). At the higher dilution rate, the proportion of solventogenic cells is significantly lower, with 4 % and 31 % in the 1 st and 6 th bioreactor, respectively ( Figure 5D, Table SIIb). The ratio of 31 % solventogenic cells in the 6 th bioreactor (after a total mean residence time of 10.9 h) corresponds roughly to the 29 % of solventogenic cells, calculated for the 3 rd bioreactor of the cascade at the low dilution rate. However, under this condition (mean residence time of 12 h at D = 0.042 h -1 ), the ratio of acidogenic (57 %) to intermediate (14 %) cells was still higher than in the last bioreactor stage at D = 0.092 h -1 . For the latter, ratios of 39 % acidogenic and 30 % intermediate cells were predicted ( Figure 5B and D, Table SIIb).
It should be noticed that under all conditions the fraction of acidogenic cells was predicted as equal or larger than the fraction of solventogenic cells ( Figure 5B, Table SIIb).

Introducing a feedback loop for recirculation of fermentation broth from bioreactor 4 to bioreactor 2
Our experimental investigations have shown that there is a high solvent productivity for the cascade at D = 0.092 h -1 and pHbr1 = 4.3. However, complete substrate conversion and maximal solvent concentrations were not reached within the 10.9 h of total mean residence time under these operating conditions. Thus, we introduced a feedback loop from bioreactor 4 to bioreactor 2 into the cascade ( Figure 1B), aiming at an increase of the butanol productivity within the cascade. Butanol pulse experiments of Junne had raised the hypothesis that butanol can induce the solventogenic state of the cells [14], thus by recirculation of butanol-enriched fermentation broth, the differentiation process should be triggered according to this hypothesis. On the other hand, our here established model proposes that fermentation broth from bioreactor 4 is enriched with solventogenic cells, which can directly increase the butanol productivity in bioreactor 2.
The recirculation flow Fback from bioreactor 4 to 2 was set to 0.055 L h -1 , corresponding to a quarter of the feeding flow Fin. The outflow from the modified cascade, being equal to Fin, remained unchanged in comparison to the linear cascade ( Figure 1B). The steady states observed in the cascade with feedback loop showed a 50 % increase of butanol concentration in bioreactor 2 compared to the linear cascade ( Figure 5A, Table SIIa). Even though 94 % of this increase in the butanol concentration has to be attributed to transport of butanol from stage 4 to 2, 6 % of the increase resulted from an augmented metabolic activity in bioreactor 2. However, the beneficial effect of the feedback decreases in the following bioreactor stages. In bioreactor 4, the earliest termination point of the modified cascade, the butanol concentration was only 6 % higher than in the linear cascade. In bioreactor 6 there were no more than 1.6 % of butanol increase in the feedback loop system. Further, there was an important experimental variance, questioning the significance of the butanol increase ( Figure 5A, Table SIIa). However, simulations with the model, adapted to the feedback loop cascade, also showed an increase in butanol concentration for bioreactor 2, 4 and 6 of 87.8 %, 10.1 % and 2.9 %, respectively ( Figure  5B, Table SIIb). Thus, model simulation and experimental investigation indicate that a small increase of the butanol concentration can be achieved by introducing a recirculation into the cascade. Results from computer simulations quantifying the influence of a recirculation flow show that if Fback is raised to 0.5 or 1 times Fin, the butanol concentration in the 6 th bioreactor increases by 3.9 % or 4.8 % (in comparison to the linear cascade), respectively ( Figure 6A). At the same time the model predicts a reduction of the organic acid uptake. As shown in Figure 6B, the acetic acid concentrations in bioreactor 2 and 3 decrease at higher recirculation rates. However, as result of reduced organic acid uptake in bioreactor 5 and 6 the acetic acid concentrations at the outflow of cascades, run at recirculation ratios of 0.25, 0.5 and 1, increase by 0.7 %, 2.7 % and 5.6 % in comparison to the linear cascade, respectively ( Figure 6B, Table SIIb). Similar increases were modeled for the butyric acid concentrations (Table SIIb). The proportion of solventogenic cells on the total biomass is predicted to notably increase in bioreactor 2, when fermentation broth is recirculated ( Figure 6C). However, up to bioreactor 4 of the feedback-loop cascade, the subpopulation proportions adapted to the levels of the linear cascade. At bioreactor 6 the part of solventogenic cells is predicted to be even slightly lower in the cascade with recirculation ratio above a factor of 0.25 compared to the linear cascade ( Figure 6C, data in Table SIIb). The experimentally observed butanol concentrations are in good agreement with the predictions of our model, based on a recirculation of already differentiated cells. However, the declining beneficial effect of the feedback loop along the cascade does not substantiate a significant effect of a possible butanol-triggered differentiation of cells to the solventogenic state proposed by Junne [14].
3.5 Product concentrations, productivities and yields of ABE processes run in the continuously operated, multi-stage bioreactor cascades.
High solvent concentrations as discussed above are essential for an efficient downstream process. However, productivity and yield, i.e. the time needed to produce the solvents and the efficiency in using the carbon and the energy stored in the substrate, are equally important criteria for an economic process. Looking at the productivities and yields at the outflow of the cascade (Table 4), one can note that although butanol and total solvent concentrations are lower at high dilution rates (D = 0.092 h -1 ), the productivities for butanol and solvents as well as the substrate yields are higher compared to D = 0.042 h -1 Table 4). However, maximal product concentrations at D = 0.042 h -1 / pHbr1 of 4.3 and 5.6 were already reached in bioreactor 3 and 5, respectively (Table 4, Figure 3, Table SIIa). Thus, the cascade should be truncated to operate efficiently at this dilution rate. Regarding the maximal productivities instead of the maximal product concentrations, the optima are already reached in bioreactor 2 and 3 at D = 0.042 h -1 / pHbr1 of 4.3 and 5.6, respectively (Table 4). Also the yields were maximal in this range of the cascade with D = 0.042 h -1 , i.e. at bioreactor 2 for pHbr1 of 4.3 and at bioreactor 4 for pHbr1 of 5.6 ( Table 4). The maximal butanol and solvent productivities observed under the operating conditions applied in this study were 0.93 ±0.14 g butanol L -1 h -1 and 1.29 ±0.21 g total solvents L -1 h -1 at the outflow of bioreactor 2 in the cascade run at D = 0.042 h -1 and pHbr1 of 4.3 (Table 4). Yet the butanol and solvent concentrations at this state (7.4 ±1.1 g butanol L -1 and 10.3 ±1.7 g total solvents L -1 ) were rather low ( Figure 3, Table SIIa).
Also at a dilution rate of 0.092 h -1 and a pHbr1 of 4.3 the productivity was found to be maximal in a cascade shortened to 5 stages. However, at these operating conditions, the productivities at the outflow of bioreactor 6 (0.75 ±0.23 g butanol L -1 h -1 and 1.08 ±0.36 g total solvents L -1 h -1 ) were only slightly lower than the maximal values at the outflow of bioreactor 5 (0.76 ±0.23 g butanol L -1 h -1 and 1.08 ±0.36 g total solvents L -1 h -1 ) ( Table 4).
It is interesting to see, that the maximum for the solvent production rate could be shifted to bioreactor 4 by introducing a feedback loop into the cascade (Table 4). But the product concentrations in bioreactor 4 of the feedback loop cascade (5.7 ±1.7 g butanol L -1 and 8.0 ±2.0 g total solvents L -1 ) were significantly lower than those in bioreactor 6 of the same cascade (8.3 ±0.9 g butanol L -1 and 11.9 ±1.5 g total solvents L -1 ; Figure 3, Table SIIa). Furthermore, the productivities in bioreactor 6 of the feedback loop cascade (0.76 ±0.08 g butanol L -1 h -1 and 1.09 ±0.14 g total solvents L -1 h -1 ) were slightly above those observed in the linear cascade ( Table 4). The conversion yields of 0.18 g butanol g -1 glucose and 0.26 g total solvents g -1 glucose did not change between the feedback loop and the linear cascade operated at D = 0.092 h -1 , pHbr1 of 4.3 (Table 4).
When operating the linear cascade at D = 0.092 h -1 and pHbr1 of 5.6, product concentrations, productivities and yields increased until the outflow of the cascade ( Table 4). The low butanol and solvent concentrations of 4.8 ±1.9 g butanol L -1 and 6.7 ±2.7 g total solvents L -1 together with the comparatively low productivities and yields indicated that for this combination of dilution rate and pHbr1 the number of bioreactors in the cascade was not sufficient for effective conversion.

Discussion
Here we analyzed continuous ABE fermentations with suspended cells in a six-stage bioreactor cascade. Under the investigated operating conditions, a good compromise of high product concentration and high volumetric productivity was found for the 3-stage (truncated) cascade, operated at pHbr1 of 4.3 and Fin of 0.1 L h -1 as well as for the 6-staged cascade, operated at pHbr1 of 4.3 and Fin of 0.22 L h -1 , i.e. for total mean residence times in the cascade of 12 h and 10.9 h, respectively. Under these conditions, volumetric productivities and final butanol concentrations of 0.81 ±0.04 g butanol L -1 h -1 / 9.7 ±0.5 g butanol L -1 and 0.75 ±0.23 g butanol L -1 h -1 / 8.2 ±2.6 g butanol L -1 were reached, respectively.
These productivities are more than twice times higher than those observed in a two-stage ABE process [21], being the prototype of a multi-stage, continuous ABE fermentation. In that former work, residence times of 8 h in the first bioreactor (pH 4.3) and 33-42 h in the second bioreactor (pH 4.3) yielded final butanol concentrations of 10.5-12.6 g butanol L -1 , but productivities of only 0.30-0.32 g butanol L -1 h -1 [21]. Since then a high number of continuous ABE fermentation strategies has been explored and published, being summarized in a matrix of final butanol concentrations and volumetric butanol productivities by Setlhaku et al. [22,23]. This matrix shows that product concentrations above 13 g L -1 can only be reached in systems with integrated downstream processing, i.e. with multiple stages separating biological solvent production and technical product concentration (e.g. by gas stripping or pervaporation) [23][24][25][26][27]. Furthermore, this matrix illustrates that butanol productivities of 1 g L -1 h -1 and above are exclusively reached in systems with immobilized cells or cell retention [28][29][30][31]. A four-stage bioreactor cascade with Clostridium cells immobilized in a biofilm has been reported to yield up to 10.8 g butanol L -1 with a productivity of 9.2 g butanol L -1 h -1 [22].
However, the focus of our work was the establishment of a system that captures the temporal development of a batch process in a spatial dimension. Such an approximation of a plug flow bioreactor has been realized in the here presented six-stage CCSTR. In contrast to biofilm or other fixed-bed bioreactors, cells from a certain metabolic phase can be retrieved from the CCSTR for further characterization or for a transfer into another environment. The characterization of metabolic phases in the ABE fermentation based on batch processes is linked to sophisticated prearrangements and poor reproducibility [16]. In the future, experimental characterization of individual cell populations in the different metabolic phases using flow cytometry [32][33][34][35], electrooptical measurements [36] or single cell RNA sequencing [37,38] can be applied to steady state samples of the CCSTR. Already the here presented model, showing good agreement with the experimental data, delivers prospects of the composition and the characteristics of the individual cell populations.
For instance, it is worth to mention that the option for growth of intermediate and solventogenic cells, which had been newly introduced as an extension to an earlier model [19] resulted in maximal growth rate parameters of 0.46 h -1 and 0.41 h -1 for intermediate and solventogenic cells, respectively, after parameter estimation. These are lower than the maximal growth rate parameter of acidogenic cells, being 0.73 h -1 . This result is in good agreement with the experimental analyses from other groups [20,39], showing that also under solventogenic conditions a stable steady state can be reached in single-stage continuous fermentation. However, this is characterized by a lower optical density than under acidogenic conditions [9,15,39], which might result from a lower maximal specific growth rate of the solventogenic cells.
In the model of Millat [9], simulating the pH-induced metabolic shift in a one-stage chemostat, two subpopulations -acidogenic and solventogenic cells -had been considered. However, in that model the dimensions of the subpopulations develop independently from each other, i.e. the model functions of two subpopulations are not linked to each other. In the here presented model the development of the subpopulations is a differentiation process from acidogenic to intermediate and further to solventogenic cells which depends on the concentration of undissociated acids. Implementing the mechanism of cell differentiation in the model, we followed an example of modeling the segregation of dihydroxyacetone-producing Gluconobacter axydans subpopulations [40]. One result of our simulation was that only a comparatively small proportion of cells, i.e. 26-40 %, differentiated to solventogenic cells in the last bioreactor stage, while 29-54 % remained in the acidogenic state. This contrasts the model of Millat [9] where 100 % of the acidogenic cells were transformed to 100 % of solventogenic cells during the metabolic switch. However, shifting of subpopulation proportions rather than homogeneous differentiation agrees well with the recent insights into metabolic heterogeneity of bacterial communities, gained by now available single-cell analysis methods, e.g. single-cell RNA sequencing [41][42][43]. For instance, in antibiotic stress situations many bacterial populations tend to segregate into subpopulations of metabolic-active, antibiotic-sensitive cells and of persistent, metabolic-inactive, but antibiotic-resident cells. The latter group ensures the survival of the species and may profit from nutrients released by cell lysis of the first group once the antibiotic pressure is withdrawn [44]. Heterocyst formation of Nostoc punctiforme is another example for beneficial differentiation processes in bacterial populations [45].
Our model does not consider back differentiation from solventogenic cells to acidogenic cells, though shifting experiments in the single-stage continuous ABE fermentation propose that those are possible [46]. However, under the experimental conditions applied here we did not expect a back differentiation to occur. Further, the remaining part of acidogenic cells even under externally solventogenic conditions, proposed by our model, offers the perspective that not a back differentiation, but rather overgrowing of the solventogenic by the acidogenic subpopulation leads to a backward shift from overall solventogenic to overall acidogenic conditions. Validation and improvement of the here presented model predictions of the subpopulation composition will be enabled by future experimental investigations, e.g. by flow cytometry. So far, the model kept its predictive capacity even when switching op-erating conditions from the linear cascade to a cascade with feedback loop, which is another validation of the presented model.
The hypothesis from pulse experiments by Junne [14] proposes that butanol has an auto-inducing effect on the solvent production by promoting the differentiation process. Using the potential of our model, we simulated a feedback loop introducing solventogenic fermentation broth from bioreactor 4 (including ~21 % of solventogenic cells and 5.3 g of butanol) into an earlier phase of the process in bioreactor 2. On one hand, one would expect a beneficial effect of the recirculated butanol on differentiation to solventogenic cells, which should even amplify itself after bioreactor 2. On the other hand there should be a beneficial effect, resulting from the recirculation of already differentiated solventogenic cells, which is expected to be rather constant in all bioreactors after entry from the feedback loop. However, here we observed significant higher solvent production in bioreactor 2 and 3 in the feedback loop cascade, but the beneficial effect declined in the later bioreactors, resulting in final butanol concentrations only 0.2 g L -1 higher than in the linear cascade. So there is no evidence for an amplification of solventogenesis triggered by butanol. But the model prediction for the feedback loop cascade, matching the experimental results, indicates that a positive effect of a higher percentage of solventogenic cells in bioreactor 2 and 3 is compensated mainly due to reduced differentiation of new solventogenic cells caused by a decrease of the acetic and butyric acid concentrations.
As a new hypothesis, adding a feed of acetic acid to bioreactor 2 or 3 of the feedback loop cascade is one promising configuration that will be tested in the future with the here introduced modelling tool and experimental setup for studies of the ABE fermentation. Using the mathematical model not only cascade set-ups with additional feeding points, but also variations of the cascade with variable number of bioreactor stages, differently sized stages or cell recycling loops, can be evaluated and selected for experimental investigations. Thus, this work opens the door to a model-supported optimization of a continuous, multi-stage ABE process.

Conclusions
During the ABE process Clostridia go through different metabolic phases. With the here presented cascade of continuous stirred tank reactors (CCSTR), it is possible to separate these phases in a spatial manner. Investigation and application of this reactor system is supported by a mathematical model, which combines two interacting levels of the process:  Reactor system model, depending on the configuration of the reactors, describing residence time distributions and their influence on microbial population fractions and metabolite concentrations.  Kinetic model for microbial metabolism, allowing for populations with different metabolic activities, depending on the bioreactor environment. The mathematical model, which fits the experimental steady state data of all operating conditions under investigation, indicates that even in the later bioreactor tanks of the cascade only a relative small portion of the cells (about one third) is in the economically interesting solvent-forming state. Our attempt to increase the proportion of the solventogenic cells by introducing a feedback loop into the cascade from bioreactor 4 to bioreactor 2 resulted in a final butanol concentration of 8.3 g L -1 and a productivity of 0.76 g butanol L -1 h -1 , which were only slightly higher than in the corresponding linear cascade.
In the future the here established multi-stage, ABE laboratory process with its corresponding mathematical description will serve as tool for predicting and testing further fermentation strategies such as co-feeding of organic acids to different phases of the process. A collection of possible configuration modifications is shown in Figure 8.  Adapting individual mean residence time in bioreactors by using different working volumes.  Moving cells and metabolites between bioreactor environments by feed-back and feed-forward loops.  Feeding additional substrates along the cascade, e.g. organic acids  Biomass retention by separation and recycle.
Modeling metabolism in the cascade led to a more profound understanding of the shift from acidogenic to solventogenic state of the cells and the metabolic activity of physiologically different subpopulations. Together with the reactor model, allowing for simulation of any combination of configuration modifications, the performance of a cascade set-up with its operating conditions can be predicted.  Funding: This research was funded by the German ministry "Bundesministerium für Bildung und Forschung" under the ERA-Net EuroTransBio-7 initiative, reference number 031A231C.