Modeling a co-culture of Clostridium autoethanogenum and Clostridium kluyveri to increase syngas conversion to medium-chain fatty-acids

Graphical abstract


Introduction
One of the biggest challenges society faces nowadays is finding alternative processes for the sustainable production of fuels and chemicals. At present, the production of many commodities depends on fossil fuels, which is not sustainable or sugar crops, competing with human and animal food consumption [1]. To circumvent this, circular approaches are required, such as the conversion of lignocellulosic biomass or municipal waste as feedstocks to fuels and chemicals [2]. Although lignocellulosic biomass has been identified as a promising source for renewable energy and carbon [3], current technologies involving hydrolysis of this substrate result in a complex mixture of compounds that need further sepa-ration and individual processing [4]. However, gasification of these rigid materials allows for the conversion of the carbon in the original source to synthesis gas (syngas), consisting mainly of CO, H 2 and CO 2 . This energy-rich syngas can be further used as feedstock for chemocatalytic processes such as Fischer-Tropsh, but microbial fermentation of syngas is gaining more attention recently as a potential production platform [5,6]. Compared to chemical catalysts, microorganisms are more robust to variations of CO/H 2 ratio in syngas, and are also more resistant to the presence of certain impurities (e.g. sulfides), reducing the need for costly pretreatment of syngas [5].
Acetogenic clostridia are efficient microbial hosts for syngas fermentation as they can grow on CO and CO/H 2 via the Wood-Ljungdahl pathway [7]. However, the natural product range of most acetogens is limited to a mixture of acetate and ethanol [8]. Cocultivation of a syngas-fermenting organism with other organisms (that use the primary products of syngas fermentation) can be used https://doi.org/10.1016/j.csbj.2020. 10 to extend the range of possible products. Previously, a co-culture of Clostridium autoethanogenum and Clostridium kluyveri was described to produce medium-chain fatty acids (C4-C6) and their respective alcohols by assimilation of CO or syngas [9,10]. C. autoethanogenum is an acetogenic bacterium able to produce acetate and ethanol when growing on CO or syngas [11]. C. kluyveri grows on acetate and ethanol via reverse-b-oxidation, producing chainelongated acids like butyrate and hexanoate. When C. kluyveri is grown in co-culture with C. autoethanogenum on CO, it produces butyrate and hexanoate, which are further reduced by the acetogen to the corresponding alcohols, butanol and hexanol [9]. MCFA are used to produce pharmaceutical and personal care products, animal feed additives and lubricants, among other, and can be converted chemically or enzymatically into valuable biofuel molecules such as methyl esters, methyl ketones, alkenes and alkanes [12,13]. The theoretical maximum yield of hexanoate production in a co-culture of C. autoethanogenum and C. kluyveri is 0.056 mmol of hexanoate per mmol of CO, whereas the yield obtained in the most recent study [10] was 0.009 mmol hexanoate per mmol of CO, so there is substantial room for improvement and new strategies need to be developed.
Genome-scale, constraint-based metabolic models (GEMs) attempt to represent the complete set of reactions in a living organism, and have been used to gain better understanding of cellular metabolism, assessing theoretical capabilities or designing media and processes [14]. GEMs can be used to link the microbial consumption and production rates with cellular growth rates. Moreover, they enable linking these phenotypes with the genome content of the studied organisms and with internal phenotypes, such as metabolic fluxes that are usually difficult to measure experimentally. GEMs and their analysis with constraint-based techniques, such as flux balance analysis (FBA) for the calculation of steady-states, have been proven effective tools to devise strategies for increasing productivity of microbial fermentation processes [15][16][17]. Specifically, GEMs have been used to further understand the metabolism of clostridia. For instance, the GEM of Clostridium thermocellum allowed the design of metabolic strategies to increase ethanol production after identification of bottlenecks in central carbon metabolism that were inhibiting its production [18]. Stolyar and collaborators [19], generated a multi-species GEM by combining the GEMs of bacterium Desulfovibrio and archeon Methanococcus maripaludis S2 into a single model with a shared extracellular environment, bringing the use of GEMs to a next level. Since then, this type of community models have been used to describe metabolic interactions among community members and inter-species fluxes [20]. Li and Henson [21], recently used GEMs to compare mono-culture and co-culture systems to produce butyrate from carbon monoxide. They applied dynamic flux balance analysis (dFBA) [22] to analyze a community GEM to cover the changes in community composition over time and to assess the relative performance of these mixed cultures. The availability of GEMs for C. autoethanogenum [23] and C. kluyveri [24], enables the use of community modeling as a potential method to help optimizing the performance of this co-culture for syngas fermentation to elongated acids and alcohols.
In this study, we present a multi-species model built by combining the GEMs of C. autoethanogenum [23] and C. kluyveri [24]. The model accounts for experimental measurements informing on relative species abundances and steady-state production rates of syngas fermentation products obtained in chemostat runs under different conditions for mono-culture and co-culture [10]. In order to test the model, experimental values were introduced as environmental constraints by employing community flux balance analysis (cFBA) [20,25]. cFBA implicitly assumes equal abundances of the species when exchange fluxes are expressed on a per gDW basis. To circumvent this, and considering that in microbial communities different species can have distinct abundances, we have scaled fluxes by volumes in this study. Additionally, cFBA also assumes equal growth rates of the members of the community. In the current analysis of a co-culture in a chemostat, equal growth rates are achieved as the dilution rate ensures same growth rate for each organism [26]. Subsequently, the model was used to identify and assess strategies to optimize desired products, specifically butyrate and hexanoate.

GEMs of C.autoethanogenum and C. kluyveri
To represent the metabolism of C. autoethanogenum, the previously described GEM, iCLAU786, was retrieved in SBML (XML) format from the supplementary material provided by Valgepea et al. [23]. This model was amended with an exchange reaction to simulate acetate uptake when this is used as additional substrate (EX_AC_c). eQUILIBRATOR [27] was used to manually verify reaction directionality: Gibbs energy released (DG) at pH 7.0 and ionic strength (0.1 M) was computed. Reactions with DG 2 ½À30; 30 kJ/ mol were considered reversible.
The GEM of C.kluyveri, iCKL708, was downloaded in table format from its publication [24]. An additional reaction was added to excrete biomass, which was first included as new metabolite and as additional product in the biomass reaction BOF. Minor changes were applied affecting the reversibility of few reactions and addition of protons. Pyruvate synthase (Rckl119) was set to non-reversible in the direction of pyruvate production [28,29]. Pyruvate formate lyase (PFL) was set to non-reversible, allowing only the production of formate and acetyl-CoA. Protons were added in the exchange of heptanoate reaction (Rckl870). eQUILI-BRATOR [27] was used to manually verify reaction directionality as in previous model. The updated model was converted to SBML level 3 version 1 standardization [30].

Multi-species GEM reconstruction
The multi-species GEM of C. autoethanogenum and C. kluyveri was generated by combination of single species models: iCLAU786 [23] and the updated version of iCKL708 [24], respectively, following a compartmentalized approach [19] were each species is considered a single compartment. Therefore, we consider two internal compartments: 'cytosol auto' and 'cytosol kluy', with 'c' and 'ck' as their respective identifier (id). Intracellular metabolites were assigned to their corresponding compartment and the flag ' c' was added to the id of metabolites belonging to 'cytosol auto' and ' ck' to those belonging to 'cytosol kluy'. In addition to these two internal compartments, the model has an extracellular compartment that is unique for both species. To achieve this, all metabolites that were defined as extracellular (' e') in its original models, will be defined in the common extracellular compartment of the community model, id: ' e'. As some metabolites will appear in both species, names need to be unified and corrected to have the same naming system (namespace). Metabolites that are shared between species, will be exchanged through this extracellular compartment, being first transported from the corresponding intracellular compartment to the extracellular compartment, or vice versa. In principle, all metabolites that are present in both internal compartments and are defined in the extracellular compartment, can be exchanged, being the directionality of the associated reactions favorable to produce the exchange. However, some dependencies have been assumed in the model based on experimental evidences.
Each species has its own biomass synthesis reaction. An extra biomass metabolite was created and defined in the extracellular compartment for each species: 'BIOMASS c e' and 'BIOMASS ck e'. In addition, two extra reactions were added for each species, one to transport biomass from the intracellular to the extracellular compartment, and a second one to excrete biomass (exchange reaction). A reaction was included to distinguish the amount of H 2 excreted by C. kluyveri, from the amount of H 2 metabolized by C. autoethanogenum. The same was done for acetate. A reaction was included to distinguish the amount of acetate metabolized by C. autoethanogenum, from the amount of acetate produced by C. autoethanogenum. The model also contemplates the possible production of butanol and hexanol via butyrate and hexanoate uptake by C. autoethanogenum. The added reactions are a transport reaction from the external compartment to the internal compartment of C. autoethanogenum, reactions for production of butyraldehyde and caproaldehyde from the corresponding fatty-acids and reactions for production of alcohols from their corresponding aldehydes. Finally, the multi-species model was transformed into SBML level 3 version 1 (see supplementary material).

Multi-species modeling framework
In order to model the community, we have followed an approach similar to the one proposed by SteadyCom [31] and that is based on community FBA (cFBA) [25]. Environmental fluxes (mmol l À1 h À1 ) are integrated as model constraints instead of specific fluxes (mmol gDW À1 h À1 ), where gDW indicates grams of dry weight. The biomass reaction of each species incorporates as new term, the biomass of the relative species together with the growth rate term. In this way, we can account for species abundance in the community. The biomass of each species is calculated based on the community biomass and the species ratio. In addition to this, steady-state and equal growth rate of species are assumed.

Calculation of species abundances
The ratio between C. autoethanogenum and C. kluyveri in coculture, was estimated from partitioning RNAseq reads and confirmed via cell counting in microscopy observations as two independent methods. Transcriptomic data was obtained from steady-state co-cultures grown in chemostats [10]. The Genomes of C. autoethanogenum: DSM 10061 (GCA 000484505.1) [32] and C. kluyveri: DSM 555 (GCA 000016505.1) [33] were retrieved from the European Nucleotide Archive. The genomes have similar size with sequence length 4.352.205 and 4.023.800, respectively [34]. Reads were mapped to each genome using BWA-SW (Burrows Wheeler Aligner) [35] and the ratio was calculated based on the amount of reads associated to each species.
The second method consisted of direct cell counting under microscopy observations. This led to a proportion between cell numbers of 10 C.autoethanogenum by 1 C. kluyveri. This proportion was considered to calculate the accumulated dry weight. To calculate the respective dry weights, the cellular volume of each species was calculated based on their average size. C. autoethanogenum is a rod-shaped bacterium with an average size of 0.5 Â 3.2 lm [11] and C. kluyveri cells are curved rods, with average size of 12.5 lm in length and 1.5 lm in width [36]. Cell volume was calculated following a previously proposed formula for rod-shaped cells [37]: V ¼ ½ðw 2 Á p=4Þ Á ðl À wÞ þ ðp Á w 3 =6Þ, with l and w indicating length and width, respectively. The associated dry weight (DW) was then derived using: DW ¼ 435 Á V 0:86 [37]. Then, the dry weight of C. autoethanogenum was multiplied by 10 and the dry weight of C. kluyveri was multiplied by 1 (as the observed proportion). Finally, the biomass-species ratio was calculated based on the ratio of the accumulated dry weights. The proportion was observed to be constant among experimental conditions with CO and CO/H 2 , so we assume that the relative abundances are constant for the rest of conditions too.

Use of experimental values to constrain the model
Experimental measurements were converted to mmol l À1 h À1 . Product concentrations, measured in mM, were used to compute product secretion rates (mmol l À1 h À1 ) by using the same hydraulic retention time (HRT) that in laboratory settings. Dilution rate is the inverse of the HRT. In steady-state conditions, growth rate is considered equal to the dilution rate and therefore, it was calculated as the inverse of the HRT and expressed in h À1 . In co-culture, the growth rate of both species was assumed to be the same and equal to that of the community, since HRT was kept constant both, in mono-culture and co-culture experiments [26]. We have followed the usual convention in constraint-based modeling, so that uptake is represented by negative fluxes whereas production corresponds to positive fluxes. To model experimental conditions, we fix substrate uptake rates to the desired ones by setting the lower bounds of the corresponding exchange reactions to the measured values multiplied by À1 (as it corresponds to consumption). Biomass reactions were constrained with the growth rate multiplying the total biomass by the ratio of each species (gDW l À1 h À1 ). Similarly, product formation was set to be at least 80% of the calculated product formation by modifying the lower bound of the corresponding exchange reaction. ATP maintenance reactions of each species, ATPM auto and ATPM, were transformed to mmol l À1 h À1 from the pre-set values in mmol gDW À1 h À1 multiplying by the total biomass and species ratio. In cases where metabolites behave as products that are further metabolized by the other species, the transport reactions of these metabolites are forced to operate in the direction from the external compartment to the other species compartment.

Chemostat experimental data
Experimental data was collected from reactor run 3 and 4 of the recent study [10] on C. autoethanogenum in mono-culture and cocultivation of C. autoethanogenum and C. kluyveri grown on CO/H 2 and CO/acetate. In co-culture experiments, C. kluyveri was inoculated in the reactor on top of C. autoethanogenum in a 1:20 volume ratio. The organisms were cultivated in chemostat to control environmental conditions such as pH (6.2), temperature (37°), HRT (between 1.5 to 2 days) and medium composition during the entire reactor run. A reactor run starts with inoculation of C. autoethanogenum in mono-culture followed by co-cultivation with C. kluyveri after reaching stationary phase. Total reactor volume is 1.5 l. Working volume was set between 0.75 l to 1 l. Experiments were run on different conditions of CO/H 2 and CO/acetate as initial substrates. Concentrations of organic acids and alcohols in the reactor and gas composition in the outflow were tracked during the runs.

Model simulations
Model simulations were done using COBRApy, version 0.17.0 [38], IBM ILOG CPLEX 128 and Python 3.6. Simulations based on changes/addition of parameters were done by constraining the associated reactions with the mentioned values. When simulations required removal of a substrate or product, flux through the associated reaction was set to 0. Constraints on the profile of fermentation products were kept unchanged when simulations were based on substrate uptake ratios in C. kluyveri, unless stated otherwise. For each explored condition, the solution space and the set of fluxes compatible with the measured constraints were sampled using the sample function in the flux_analysis submodule COBR-Apy. Flux sampling is a method to get a distribution of fluxes [39] under specific conditions. Presented results are the average and standard deviation based on 15000 iterations generated at each condition. All additional assumptions taken into account during model simulations are listed in the supplementary material.

Genetic intervention strategy
OptKnock and RobustKnock [40,41] were applied as algorithms that suggest reactions to be knocked out that can potentially increase the yield of a target reaction. The algorithms were applied to increase ethanol production in the GEM of C. autoethanogenum [23]. Both algorithms were integrated in a python script adapted for COBRApy and CPLEX as solver. OptKnock identifies a set of reaction knockouts that allows high production of a target product under the constraint of optimal growth in wild type. RobustKnock guarantees a minimal production rate by considering alternative optimal solutions that produce less of the target product. This is achieved by employing a bi-level max-min optimization. The possible reactions to be modified were adjusted in order to avoid essential reactions, reactions associated to essential genes, extracellular reactions and reactions with no associated genes. The identified mutants were further implemented in the model of C. autoethanogenum deleting the corresponding reactions. Each mutant was assessed at each experimental condition and compare to the wild type. The media as well as the biomass reaction were constrained using the experimental data of mono-culture experiments for those conditions [10]. Fluxes are expressed following the modeling framework. For each explored condition, the solution space and the set of fluxes compatible with the measured constraints were sampled. Presented results are the average and standard deviation based on 15000 iterations generated at each condition.

Results
The objective of this study is to find optimization strategies for the production of medium-chain fatty-acids from syngas using the co-culture of C. autoethanogenum and C. kluyveri. The generated multi-species GEM, together with the GEM of C. autoethanogenum, were used to assess these strategies.

Description and validation of the GEM of individual strains
The GEM of C. autoethanogenum, iCLAU786 is composed of 1108 reactions and 1094 metabolites. The model is able to simulate growth on CO or syngas as the sole carbon and energy source, producing acetate and ethanol as the main fermentation products.
The GEM of C. kluyveri, iCKL708 has been previously shown to predict growth on ethanol and one other organic acid (acetate, propionate, butyrate, or succinate), propanol and acetate, crotonate, and vinyl acetate, in accordance to published experimental data [36,[42][43][44]. The updated GEM of C. kluyveri, has 993 reactions and 811 metabolites. This updated model also simulates growth on acetate and ethanol uptake producing butyrate and hexanoate as the main chain-elongated products and H 2 .

Multi-species GEM
The multi-species GEM contains 2064 reactions and 1823 metabolites, from which 139 reactions correspond to extracellular reactions and 208 metabolites belong to the shared extracellular compartment. Fig. 1, shows the dependencies included in the model to describe the syngas fermentation process by the co-culture based on the experimental data [10]. H 2 ck ! reaction represents the amount of H 2 excreted by C. kluyveri. Reaction ! H 2 e represents the uptake of H 2 in C. autoethanogenum. Reaction ! AC c, represents the uptake of acetate by C. autoethanogenum in simulations where acetate acts as additional substrate. This serves to distinguish the fluxes between acetate feed rate (! AC c), acetate production rate (AC e !) and acetate consumed by C. kluyveri (AC e ! AC ck). These are special cases since H 2 and acetate can be shared, metabolized and produced in co-culture conditions.
Previous studies have shown that C. autoethanogenum is able to grow on CO, CO/H 2 producing ethanol and acetate as the main fermentation products [11]. Acetate and ethanol can further be taken up by C. kluyveri producing H 2 , butyrate and hexanoate. H 2 produced by C. kluyveri appears to be further metabolized by C. autoethanogenum [9,10]. Furthermore, the presence of aldehyde ferredoxin oxidoreductase and Ethanol:NAD + oxidoreductase enzymes in C. autoethanogenum allows for a potential two step conversion of butyrate and hexanoate, via the respective aldehdye, to butanol and hexanol, respectively. C. kluyveri is not able to utilize CO and its metabolism can be inhibited by this compound [9]. However, providing co-cultivation with C. autoethanogenum, dissolved CO concentration can be kept low. Naturally, dissolved CO will be dependent on the gas-liquid mass transfer and the CO uptake rate of C. autoethanogenum. Because kLa values for COwater are relatively low in stirred tanks [45], these systems are often not kinetically limited and low dissolved CO concentrations are expected in the culture broth. This low dissolved CO concentration will leave C. kluyveri metabolism unaffected [9]. In the model this is indicated by preventing flux through the reaction CO e ! CO ck as it is shown in Fig. 1.
Microscopy observation of the co-culture led to the estimation of a ratio of 10 cells of C. autoethanogenum per 1 cell of C. kluyveri. Analyses of transcriptome samples obtained by RNAseq of the community [10], were done to identify the fraction of RNA arising from each community member. The estimated relative abundances yielded between 90-95% of C. autoethanogenum and 5-10% of C. kluyveri. Differences in cell size and volume were considered as C. kluyveri cells have approximately 36 more volume that C. autoethanogenum [11,36,37]. The estimated volumes were used to esti- mate dry weight of each cell species, resulting the cell dry mass of C. kluyveri, 22 times more than C. autoethanogenum. Finally, the cell ratio (10:1) was taken into account resulting in a biomass ratio of 68.5% C. kluyveri and 31.5% C. autoethanogenum. This cell ratio was observed in CO and CO/H 2 conditions and it was assumed to be constant in the model simulations.

Multi-species GEM accurately predicts experimental results
The initial mono-culture experiments only involved C. autoethanogenum [10]. Fig. 2, shows the steady-state production rates of the fermentation products expressed in mmol l À1 h À1 . Experiments were run on CO/H 2 and CO/AC (acetate) as initial substrates. Acetate production rate ('EX AC e') refers to the sum of the acetate feed rate not consumed by C. autoethanogenum and the one directly produced by C. autoethanogenum. Fig. 2 shows that the model predictions match relatively well the experimental results for C. autoethanogenum. Accordingly, the model correctly predicts that ethanol production increases at higher H 2 feed rates and gradually with the addition of acetate. However, the model predicts slightly higher production rates for acetate in conditions with higher amounts of acetate in the background.
The co-culture experiments were run under same conditions as the mono-culture experiments. Fig. 3, represents steady-state production rates of fermentation products expressed in mmol l À1 h À1 by the co-culture. It shows the comparison between experimental results collected in co-culture experiments [10] and the results obtained via the multi-species model. Acetate production rate ('EX AC e') refers to the sum of the acetate in the feed, the acetate directly produced by C. autoethanogenum, minus the acetate consumed by C. kluyveri (reaction id Rckl835). Ethanol production rate ('EX ETOH e') refers to the ethanol produced by C. autoethanogenum minus the ethanol consumed by C. kluyveri (reaction id Rck-l837). The model correctly predicts production of medium-chain fatty acids upon introduction of C. kluyveri. Similarly to the mono-culture simulations, there is a slight mismatch between predicted and observed acetate production. The model correctly predicts the increase of medium-chain fatty-acids when more H 2 or acetate is added. When H 2 feed rate is equal to 5.3 mmol l À1 h À1 , butanol is also produced (0.075 mmol l À1 h À1 ). Also, ethanol accumulation is low in most co-culture conditions, suggesting most of it is metabolized by C. kluyveri. The steady-state production rates of acetate increases with increasing acetate uptake (see Fig. 2 and 3), but the amount of acetate expressively produced by C. autoethanogenum decreases with increasing acetate uptake, since the addition of acetate leads to more acetate converted to ethanol [10]. However, there is still a relative high level of acetate accumulated versus desired fatty acid products, which represents a loss of carbon into acetate that could be minimized if C. kluyveri could consume more acetate.

Assessing the distribution of metabolic fluxes with the multispecies model
After having shown that the model describes accurately the metabolic interactions between the two microbes, we used it to explore intracellular flux distributions that would be otherwise challenging to access. To study the metabolic fluxes in the coculture, we used a sampling approach that produces, for each reaction in the combined model, a distribution of possible fluxes. Fig. 4 provides an overview of selected reactions in the system (indicated by R#). Fluxes for all reactions can be found in the supplementary material.
CO or CO and H 2 are converted via the Wood-Ljungdahl pathway in C. autoethanogenum. In this pathway, CO is converted to CO 2 via CO dehydrogenase, providing reducing equivalents to the cell. Released CO 2 is shuttled into the Wood-Ljungdahl pathway via the bifurcating formate dehydrogenase [46]. H 2 taken up by C. autoethanogenum is used for redox generation in NADPdependent electron bifurcating hydrogenase (Hyt) reaction (R3) and in formate hydrogen lyase reaction (R2 in Fig. 4)) to produce formate. The flux through these two reactions increases with increasing H 2 supply. Part of the formate is excreted and part is further metabolized to acetyl-CoA (ACCOA) following the Wood-Ljungdahl pathway. Pyruvate is partly produced from acetyl-CoA via the pyruvate synthase (R8 in Fig. 4) for assimilation. The majority of the acetyl-CoA is converted to acetate via acetyltransferase (R6) and acetate kinase, yielding ATP. Ethanol can be formed in two ways [47]: from the reduction of acetate to ethanol via aldehyde ferredoxin oxidoreductase (R7) and alcohol dehydrogenase, or via reduction of acetyl-CoA to acetaldehyde (ACAL) and ethanol. Acetate and ethanol are secreted to the medium where it is partly taken up by C. kluyveri. Ethanol production rate by C. autoethanogenum (reaction R10) and ethanol uptake by C. kluyveri (R15) indicate that most of ethanol is removed by C. kluyveri and proves the metabolic change to solventogenesis due to addition of C. kluyveri. This can be observed comparing to mono-culture results showed in Fig. 2. Comparison of experimental ( e) and model ( m) results of the steady-state production rates of fermentation products in C.autoethanogenum in mono-culture under CO/H 2 and CO/AC (acetate) conditions. X axis represents the H 2 and acetate feed rate, respectively. Y axis represents the steady-state production rates of acetate and the secondary y axis represents the steady-state production rate of ethanol. In CO/H 2 conditions, CO feed rate = 4.8 mmol l À1 h À1 ; growth rate = 0.021 h À1 and working volume = 1 l. In CO/AC conditions, CO feed rate = 6.4 mmol l À1 h À1 ; growth rate = 0.028 h À1 and working volume = 0.75 l. Substrates feed rates and production rates are expressed in mmol l À1 h À1 . Fig. 2, where the steady-state production of ethanol was lower than in co-cultivation with C. kluyveri (flux through R10). In C. kluyveri ethanol is oxidized to acetyl-CoA (ACCOA) and part of acetyl-CoA is converted to acetate via acetyltransferase and acetate kinase (R16). Acetyl-CoA initiates the reverse b-oxidation pathway (R16-R19) to produce butyrate via acetoacetyl-CoA (AACCOA), then 3-hydroxybutyryl-coa (3HBUTCOA), crotonyl-CoA (CROCOA) and butyryl-CoA (C40COA). C40COA, transfers the CoA group to acetate, producing butyrate and acetyl-CoA. Some of the butyrate can be elongated further to hexanoate by reaction of butyryl-CoA together with hexanoyl-CoA (C60COA)(R19). Acetyl-CoA is also assimilated by fixing CO 2 to pyruvate via pyruvate synthase (R20) in C. kluyveri.
Succinate is converted to Crotonyl-CoA (CROCOA), yielding an additional 2 acetate (see Fig. 4), involving the pathway via succinyl-CoA, succinate semialdehyde, 4-hydroxybutyrate, and 4-hydroxybutyryl-CoA [33]. The model indicates that part of the acetate pool in C. kluyveri comes from uptake of succinate produced by C. autoethanogenum (see R28), which could explain the Fig. 3. Comparison of the experimental ( e) and model ( m) results of the steady-state production rates of the fermentation products in co-culture under CO/H 2 and CO/AC (acetate) conditions. X axis represents the H 2 and acetate feed rate, respectively. Y axis represents the steady-state production rates of acetate and the secondary y axis represents the steady-state production rate of ethanol, butyrate and hexanoate by the co-culture. In CO/H 2 experiments, CO feed rate = 4.8 mmol l À1 h À1 ; growth rate = 0.021 h À1 and working volume = 1 l. In CO/AC conditions, CO feed rate = 6.4 mmol l À1 h À1 ; growth rate = 0.028 h À1 and working volume = 0.75 l. Feed rates of substrates and production rates are expressed in mmol l À1 h À1 . Fig. 4. Schematic representation of the simulated metabolism of C. autoethanogenum and C. kluyveri under chemostat cultivation conditions. In CO/H 2 experiments, CO feed rate = 4.8 mmol l À1 h À1 ; growth rate = 0.021 h À1 and working volume = 1 l. In CO/AC conditions, CO feed rate = 6.4 mmol l À1 h À1 ; growth rate = 0.028 h À1 and working volume = 0.75 l. Blue arrows indicate the fluxes direction. Metabolites colored green are extracellular metabolites, partly exchanged between species and partly excreted or assimilated to/from the media. The heatmap on the right shows the fluxes of the R# reactions selected on the map for all conditions. Flux values are log transformed (log (Flux + 1)) for a better visualization. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) slight mismatch observed in Fig. 3 between acetate predicted by the model compared to experimental results. The model predicts a low amount of ethanol being oxidized to acetate (R16), supporting activity of the succinate pathway.
The model predicts reduction of CO 2 production by C. autoethanogenum with increasing H 2 feed rate (Fig. 4B). This was also observed in the experimental measurements [10], where CO 2 /CO ratio decreased linearly with increasing H 2 uptake. According to model predictions, CO 2 production rate drops to 0.28 mmol l À1 h À1 when H 2 is supplied, as compared to 2.5 mmol l À1 h À1 that is produced when CO is the only carbon source (see supplementary material). It is observed that more H 2 is metabolized by C. autoethanogenum with increasing H 2 feed rate (as shown in Fig. 4B) similar to what was found in experimental results [10]. When more H 2 is fed to the reactor (see R25), more protons are released. In contrast, less protons are released when more acetate is fed. ATP synthase increases in both species (R12, R22) when more H 2 or acetate is fed to the reactor.

Effect of the biomass ratio of C. kluyveri-C. autoethanogenum on the metabolic profiles of the culture
The model enables detailed inspection of production and uptake profiles. Therefore we investigated the sharing of metabolites between both organisms. Analysis of intracellular fluxes in Fig. 4 suggests that succinate is produced by C. autoethanogenum and metabolized by C. kluyveri producing part of the acetate and crotonyl-CoA pool needed for chain elongation and production of fatty-acids. Fig. 5, shows the different profiles depending on the biomass species ratio when we simulate chemostat cultivation experiments of C. autoethanogenum and C. kluyveri growing on CO and H 2 as carbon and energy sources. Succinate, acetate and ethanol uptake by C. kluyveri increases with more H 2 supply. Succinate uptake decreases when C. kluyveri is less abundant in the co-culture. On the contrary, the amount of ethanol and acetate that is available to be metabolized by C. kluyveri decreases along with the relative abundance of C. autoethanogenum.
To further investigate the role of succinate, similar simulations were performed but this time preventing the uptake of succinate by C. kluyveri (Fig. 6). The model shows that, without succinate uptake, a biomass-species ratio of 70-30% and CO as the only carbon source results in an unfeasible situation. In addition, the fluxes through reactions related to non-growth associated maintenance (ATPM, ATPM auto) decreased substantially when biomass ratios of 60-40%, 50-50% and 40-60% were considered. Thus, the in silico analyses show that a possible way to meet the experimentally observed constraints is through succinate uptake. Changes in the biomass ratio also affect acetate and ethanol exchange between the microbes. Acetate and ethanol uptake by C. kluyveri increases when more H 2 is fed to the system. Acetate/ethanol production ratios become higher when C. kluyveri is more abundant and decrease when there is more C autoethanogenum in the coculture. Based on these results, the species biomass ratio in cocultivation is estimated to be between 60-40% and 70-30% (C. kluyveri-C. autoethanogenum).

Strategies to increase production of medium-chain fatty-acids
We used the experimentally-validated model to simulate alternative scenarios that have so far not been explored experimentally. Accordingly, we present here an analysis of possible strategies to increase medium-chain fatty acid production by the co-culture.

Addition of succinate
As a strategy to increase the production of desired products, we simulated how the addition of succinate as extra carbon source would affect the production of butyrate and hexanoate at different H 2 feed rates and fixed CO feed rate (4.8 mmol l À1 h À1 ). Fig. 7, shows the product profile for the species biomass ratio 70-30% (C. kluyveri-C. autoethanogenum), respectively. The increased carbon availability has been considered and used to normalize the results, so they are represented as mmol of product per total substrate (CO and succinate) per carbon. As already indicated in Fig. 3, and in contrast to experimental results, the model predicts more hexanoate than butyrate production even when succinate uptake is 0. The model predicts a yield of 0.026 mmol of hexanoate per mmol of CO. This means, that this co-culture has the capacity to produce three times more hexanoate than what it is currently being produced (0.009 mmol per mmol of CO) [10], reducing the gap respect to the maximum theoretical yield (0.056 mmol of hexanoate per mmol of CO) to 3%. Increased hexanoate production is observed next to an increase in CO 2 uptake by C. kluyveri and changes in H+ balance and ATP synthase.
Furthermore, succinate addition allows increased production of both fatty-acids. Hexanoate production can increase up to four times when succinate is added and butyrate has potential to increase around five times with respect to the results obtained Fig. 5. The effect of changing biomass species ratios on the uptake of acetate, ethanol and succinate by C. kluyveri at different H 2 feed rates and fixed CO feed rate. CO = 4.8 mmol l À1 h À1 and growth rate = 0.021 h À1 , respectively. C. kluyveri-C. autoethanogenum biomass ratio are indicated on the x axis and y axis represents the fluxes of transport reactions from extracellular to C. kluyveri compartment of acetate, ethanol and succinate.
with no succinate addition. On the other hand, simulations show no relevant differences upon variations of H 2 feed rates.

Genetic intervention strategies
As a second strategy for co-culture optimization we have predicted and evaluated genetic interventions that could lead to higher medium-chain fatty-acid production.
The strain design algorithms OptKnock and RobustKnock [40,41] were applied to identify candidate reactions to be knocked out. These were subsequently evaluated through dedicated simulations. Fig. 8 shows the effect of knocking out three candidate reactions on the production of ethanol in C. autoethanogenum. The three reactions are located in the metabolism of C. autoethanogenum. Mutant 1 refers to the knock out of formate transport reaction via proton symport (FORt2, model id rxn05559 c0). Mutant 2 refers to the deletion of acetaldehyde:NAD+ oxidoreductase (ACALDx, model id rxn00171 c0). Mutant 3 refers to the knock out of formate dehydrogenase (ferredoxin) (FDH fer, model id rxn00103 c0). The impact of removal of these reactions is compared to wild type C. autoethanogenum, corresponding to the GEM without any modification. The deletion of each reaction results in increased ethanol production. Ethanol increases up to 83% with respect to the wild type in simulations with CO as the only carbon source. In conditions where H 2 acts as second substrate, ethanol production increases up to 150%. Acetate decreases up to 11% in simulations with only CO and decreases up to 30% in CO/H 2 simulations. Mutant 1 seems to have a higher ethanol yield when CO is the only carbon source compare to CO/H 2 conditions while mutant 2 and mutant 3 provide a higher effect on ethanol/ acetate ratio in CO/H 2 simulations. As seen in Fig. 6, the production of ethanol by C. autoethanogenum and therefore, the uptake of ethanol by C. kluyveri, increases with increasing H 2 feed rate, resulting in an increased production of medium-chain fatty-acids. This suggests that these deletions in C. autoethanogenum can improve the production of medium-chain fatty-acids in this respective co-culture, since ethanol production can potentially be increased. Simulations of the effect of these mutations in monoculture and co-culture are presented in the supplementary material.

Discussion
We present here a constraint-based model of the co-culture of C. autoethanogenum and C. kluyveri in the context of CO/syngas fer- Fig. 6. Effect of changing biomass species ratios on the uptake of acetate and ethanol by C. kluyveri at different H 2 feed rates without succinate uptake and fixed CO feed rate. CO = 4.8 mmol l À1 h À1 and growth rate = 0.021 h À1 . C. kluyveri-C. autoethanogenum biomass ratio are indicated on the x axis and y axis represents the fluxes of transport reactions from extracellular to C. kluyveri compartment of acetate and ethanol. Fig. 7. Effect of succinate addition on the production of butyrate and hexanoate under different H 2 feed rates when a biomass ratio of 70-30% is considered. X-axis represent succinate feed rate and y-axis represent mmol of butyrate or hexanoate normalized per mmols of total substrate per carbon. mentation to produce medium-chain fatty-acids. A model with similar characteristics had already been used to simulate CO to butyrate conversion by bacterial co-culture systems [21]. Our model extends previous efforts, and is calibrated and tested with a battery of available experimental measurements. Modeling bacterial communities using flux balance analysis and GEM is complicated by the fact that special attention has to be paid to the biomass abundances of the microbial species in order to achieve balanced growth of the co-culture. Previous efforts used dynamic flux balance analysis to consider biomass growth [21]. Here, steady-state conditions were used of which additional data was available, allowing to overcome the challenge of estimating relative abundance of each species by combining microscopy observations and RNA seq in terms of cell numbers. These were subsequently converted to biomass ratios by considering the relationships between cellular dimensions, cell volume and biomass dry weight. This enabled the application of community flux balance analysis [25], which has been shown to accurately predict flux distributions and exchange fluxes between species and the community environment when analyzing stable communities in chemostat experiments. In this model it is assumed that the exchange of metabolites occurs indirectly using the culture medium as an intermediate. Recently, a direct exchange of electrons and metabolites due to direct cell-to-cell interactions have been observed in a co-culture of C. ljungdahlii and C. acetobutylicum [48]. Further research on the latter co-culture also demonstrates a high exchange of proteins, showing persistence of cells with exchanged cellular components [49]. An extension of our model could include more detailed description of the mechanism of metabolite exchange. This could shed light on additional interactions that might take place such as a possible exchange of amino acids.
The presented co-culture model is versatile and can simulate the CO/syngas fermentation process leading to medium-chain fatty-acid production. The presented model describes both the behaviour of a C.autoethanogenum mono-culture thriving on syngas as well as the behaviour of a co-culture of C. autoethanogenum and C.kluyveri. The model accurately reproduces the steady-state production rates of fermentation products obtained in chemostat experiments, and predicts the shift of the metabolism of C. autoethanogenum towards solventogenesis in co-cultivation with C. kluyveri [10]. In addition, model predictions on production/consumption rates (see Fig. 4) agree reasonably with previous literature on CO/syngas fermentation [9][10][11]33]. The product profile has shown the relative high level of carbon loss in acetate, compared to the desired elongated fatty-acids. This could be improved by in-line product removal, pH adjustments, an increase of CO pressure to obtain a higher conversion of acetate to ethanol or genetic engineering.
Analyses of metabolic fluxes in the model surprisingly suggested succinate production by C. autoethanogenum as an intermediate in the co-culture. Accumulation of succinate has not been experimentally observed in the calibration experiments [10], and is not described as major physiological end-product of C. autoethanogenum when grown on syngas [11]. However, its production by C. autoethanogenum in the presence of C. kluyveri could take place as it has been reported to be an overflow product of acetogenic metabolism [50]. Here, succinate could be produced to overcome the temporal overflow of C/electrons, potentially in conditions where too much reduction equivalents are provided. Succinate is described as a possible substrate for C. kluyveri [33]. Presence of succinate could slow down consumption of ethanol/acetate by C. kluyveri. This would also affect co-culture compositions by limiting C. kluyveri abundances. The suggested exchange of succinate, can explain the slight mismatch observed in the acetate production simulated compared to the laboratory results. Mono-culture results of C. autoethanogenum however, do not show succinate production. Thus, this difference can be derived from the biomass drop observed with increasing acetate feed rate affecting the ATP maintenance requirements [10] and redox balance.
Model simulations have shown that omitting succinate uptake resulted in unfeasible growth conditions when grown on only CO (see Fig. 6, ratio 70-30%). This dependency is additionally shown in the case with both CO and H 2 , where ATPM decreased substantially to sustain growth. An alternative explanation for this dependency could be related to C. autoethanogenum cell size assumptions, being potentially bigger than the average size, since it was reported to have considerable variations [11]. Therefore, the co-culture might operate in ratios where the biomass of C. autoethanogenum is more abundant, with relative values between 60-40% and 70-30%, leading to possible changes in the relative species abundance among different conditions.
As the model closely predicted obtained experimental values, it can be used to design potential strategies for improved production. Here, we explored a series of strategies to optimize the production of hexanoate. The first strategy relied on the role of succinate, as the model predicts it increases the pool of acetate and crotonyl-CoA, a precursor of the desired fatty-acids, in C. kluyveri (see   4). The flux through the reverse b-oxidation pathway increases up to four times when succinate is added (see supplementary material). The increase of crotonyl-CoA results in a higher butyryl-CoA pool. The presence of more butyryl-CoA initiates the chain-elongation process to produce hexanoyl-CoA in the same way as butyryl-CoA is formed (see 4). According to the model, most of butyrate formed from butyryl-CoA reacts with hexanoyl-CoA producing hexanoate. This results in an increase in hexanoate production up to three times (see Fig. 7). Moreover, an increase in ethanol uptake by C. kluyveri and a decrease in acetate production by C. autoethanogenum and subsequent uptake by C. kluyveri appears to cause an additional boost in hexanoate production. According to the model, addiction of succinate raises hexanoate production up to 0.067 mmol per carbon of fed substrate (CO and succinate) and would possibly lead to a further increased production of MCFA and alcohols in conditions with higher H 2 influx (>5 mmol l À1 h À1 ), as it has been previously observed. Furthermore, it has already been proven that succinate leads to an increase of MCFA in C. kluyveri [43], so its addition in co-culture experiments could potentially confirm model results.
The second strategy aims to increase hexanoate by increasing ethanol production by C. autoethanogenum. An increasing ratio of ethanol/acetate ratio has been shown to result in increased hexanoate production in C. kluyveri [51,36]. In cases where butyrate and hexanoate are not constrained (see Fig. 7), the predicted ethanol/acetate ratio (around 6:4) is higher than when butyrate is more prominent (see Fig. 6). The model thus confirms previous experimental results and highlights the potential of an increased ethanol/acetate ratio in stimulating the production of hexanoate. It should be bear in mind that model predictions are based on optimality principles and assumptions. The model suggests that increased ethanol production leads to increased medium-chain fatty-acids (see supplementary material).
C. autoethanogenum is known to increase production of ethanol under acidic or redox overloading conditions [52,53,10]. Additionally, partial inactivation of the adhE cluster or knock out of one of the AOR genes has been shown to result in increased ethanol production in C. autoethanogenum [47]. Using the GEM herein developed for C. autoethanogenum we predicted that the individual knock out of one of the following three reactions could increase the production of ethanol: acetaldehyde oxidoreductase (ACALDx), formate transport (FORt2), and the bifurcating formate dehydrogenase (ferredoxin) (FDH fer) (see Fig. 8). The acetaldehyde oxidoreductase reaction (ACALDx, id rxn00171 c0) is associated with several isoenzymes encoded by genes: CAETHG-RS16140, CAETHG-RS08865, CAETHG-RS08810, CAETHG-RS18400 and CAETHG-RS18395. The same isoforms are associated to aldehyde ferredoxin oxidoreductase reaction (CODH-ACS) (leq000004) and two of them (CAETHG-RS18400 and CAETHG-RS18395) are also involved in ethanol oxidoreductase (ALCDx) (rxn00543 c0) reaction. The affinity of each isoenzyme to each reaction has to be studied in order to fully eliminate acetaldehyde oxidoreductase activity. An acetaldehyde oxidoreductase (ACALDx) mutant has previously been shown to indeed have increased ethanol production up to 180% [47], making the deletion of this reaction seems a promising application to increase fatty-acids production in the co-culture system. The knock out of formate-related reactions in C. autoethanogenum is not described previously, but model simulations done here, suggest that they contribute to ethanol production. Formate transport in via proton symport (FORt2, model id rxn05559 c0) is catalyzed by an enzyme encoded by only one gene -CAETHG-1601, which allows relatively easy removal of this activity. The model predicts that inactivation of this reaction forces more flux through the Wood-Ljungdahl pathway, increasing the amount of acetyl-CoA. Due to the increase in acetyl-CoA pool, the fluxes through aldehyde ferredoxin oxidoreductase and acetalde-hyde oxidoreductase are increased, thus producing more ethanol. As third option the model suggests to knock out the formate dehydrogenase (ferredoxin) (FDH fer, model id rxn00103 c0) activity. This reaction has three associated isoenzymes encoded by genes: CAETHG-RS00400, CAETHG-RS13720 and CAETHG-RS14690. The inactivation of this reaction forces the production of formate mostly via formate hydrogen-lyase from H 2 (FHL, rxn08518 c0). In mono-culture conditions where there is no H 2 supply, H 2 is produced via an NADP-dependent electron-bifurcating hydrogenase reaction (Hyt) (model id leq000001). This functionality of Hyt seems to occur in situations where redox mediators get too reduced [46]. The model shows a higher conversion of CO converted to CO 2 via CO dehydrogenase (CODH4, model id rxn07189 c0) which forces more flux through Wood-Ljungdahl pathway, producing more acetyl-CoA. Also, more CO 2 is fixed producing pyruvate via pyruvate synthase (rxn05938 c0) which is shuttled back via pyruvate formate lyase (rxn00157 c0) to produce more acetyl-CoA and formate. The extra pool of acetyl-CoA forces more flux through aldehyde ferredoxin oxidoreductase and acetaldehyde oxidoreductase which leads to more ethanol.
When the reactions of fatty-acids production are not constrained, the model always predicts more hexanoate as compared to butyrate production measured in the actual chemostat experiments [10]. A potential reason for this is the pH of the co-culture and related toxicity effects of medium-chain fatty acids. Hexanoate production in C. kluyveri has been reported to be better at higher pH [54]. Thus, potentially the pH of 6.2 in the co-culture limits its production in the actual experiments. In addition to toxicity effects, the function of membrane proteins such as ATP synthase, electron transport chains or transporters can be affected by the change of proton motive force at different pH [55]. This is reflected by the observation that model predictions show differences in ATP synthase and proton balance under different butyrate/hexanoate production conditions (see supplementary material and Fig. 4). However, the acid stress response is difficult to simulate in GEMs, which possibly results in the differences observed between the model prediction and experimental results. In addition it has been observed oscillations in gas uptake rates and extracellular byproducts synchronized with biomass levels in C. autoethanogenum [56]. This could lead to thermodynamic changes affecting the reversibility of reactions and thus, the product range. An extensive thermodynamic and metabolic flux analysis study even extending the component contribution method [57] could have helped to better identified those changes.
We observe an increase of CO 2 uptake by C. kluyveri in cases where hexanoate is more abundant (see Fig. 4 and supplementary material) compared to simulations of experimental conditions, where butyrate is more abundant (see Fig. 4). CO 2 is essential for growth and C1 intermediate production in C. kluyveri [28,29,58,59]. In line with model predictions, CO 2 is converted to formate via a cyclic mechanism [29]. CO 2 is first fixed via pyruvate synthase (model id Rckl119) producing pyruvate that is further converted to formate via formate lyase (model id PFL), for assimilatory purposes. Formate is then assimilated via the tetrahydrofolate pathway to, subsequently, be transformed to various amino acids needed for growth [28,29,58,59]. Model predictions show that part of the CO 2 is also metabolized following phosphoenolpyruvate carboxylase reaction (PPC), where CO 2 is fixed together with phosphoenolpyruvate (pep) producing oxaloacetate (oaa). Oaa produces aspartate via aspartate aminotransferase (Rck-l310). Aspartate is a precursor in the synthesis of threonine involving aspartate kinase (Rckl334), aspartate semialdehyde oxidoreductase (Rckl323), homoserine kinase (Rckl335) and threonine synthase (Rckl336). Then, threonine produces acetaldehyde and glycine via threonine aldolase (Rckl341). Acetaldehyde is converted to acetyl-CoA following aldehyde alcohol dehydrogenase (ADH) reaction. So, the increase in CO 2 assimilation could lead to more acetyl-CoA and thus, more fatty-acids. Around 30-40% of CO 2 is also converted to carbonic acid. Carbonic acid produces oxaloacetate via pyruvate carboxylase (Rckl014), which follows the aforementioned route to acetyl-CoA. Simulations made with higher CO 2 uptake rates than the ones predicted when hexanoate is more abundant (>0.75 mmol l À1 h À1 ) however, did not lead to higher production rates of hexanoate or butyrate. This suggest that C. kluyveri metabolizes CO 2 up to a maximum value. In fact, this is supported by the observed correlation between growth and the maximum CO 2 fixed by C. kluyveri [59]. The use of lower hydraulic retention time and high pressure bioreactors, could possibly increase the uptake of CO 2 , close to its maximum capacity.
The maximum hexanoate predicted by the multi-species GEM is reached when succinate is added into the system in combination with CO and H 2 (see Fig. 7). Table 1, shows a comparison of electron yields for the hexanoate production predicted in this study compared to hexanoate production in other studies with similar culture systems [60,61,10].
Electron yield is calculated based on the amount of electrons going to hexanoate per the total amount of electrons going into the system as carbon and energy sources. Co-cultures of C. autoethanogenum and C. kluyveri yielded more hexanoate growing on CO/H 2 or CO/AC compared to a co-culture of C. ljundahlii and C. kluyveri, potentially as in the latter relatively more alcohols and C8 acids were produced as well. A mixed culture enriched in Acinetobacter, Alcaligenes, and Rhodobacteraceae growing solely on CO [61], increased the electron yield with respect to previously mentioned co-cultures up to 0.32. However, the addition of succinate in co-cultivation of C. autoethanogenum and C. kluyveri grown on CO and H 2 (this study), is here predicted to increase the yield of hexanoate up to 0.6, reflecting the potential of this approach to produce medium-chain fatty-acids.

Conclusions
The generation of the multi-species GEM of C. autoethanogenum and C. kluyveri has provided insights into the fermentation of CO/ syngas to medium-chain fatty acids by this co-culture. The prediction of intracellular flux distribution in this consortium enabled to uncover the potential importance of succinate uptake via C. kluyveri to produce butyrate, and suggested an effect of the biomass species ratio on the substrate profile of C. kluyveri. Simulations indicated that succinate addition might result in a substantial increase in hexanoate yield from syngas. In addition, the model of C. autoethanogenum shows that the deletion of reactions FORt2 or ACALDx or FDH fer in C. autoethanogenum potentially increase ethanol production, suggesting a potential increase in hexanoate production when these deletions were to be applied in co-culture experiments. Altogether, our model-driven approach has set a good basis for the systematic design of strategies to modulate and optimize the production of valuable chemicals from syngas.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Table 1 Comparison of electron yield obtained for hexanoate production between predicted results by the multi-species GEM and other co-culture/mixed culture. Electron yield is expressed by the amount of electrons going to hexanoate per total amount of electrons entering the system. a Volumetric consumption rate mmol l À1 h À1 b Volumetric production rate mmol l À1 h À1 c hexanaote e À /total e À in; CO = 2 e À per mol, H 2 = 2 e À per mol; AC(acetate)= 8 e À per mol; SUCC (succinate)= 14 e À per mol and hexanoate = 32 e À per mol d Assuming 90% gas consumption