# Consideration of maintenance in wine fermentation modeling

Alain Rapaport, Robert David, Denis Dochain, Jérôme Harmand, Thibault Nidelet

► To cite this version:

Alain Rapaport, Robert David, Denis Dochain, Jérôme Harmand, Thibault Nidelet. Consideration of maintenance in wine fermentation modeling. Foods, 2022, 11, 10.3390/foods11121682. hal-03671014

HAL Id: hal-03671014

https://hal.inrae.fr/hal-03671014

Submitted on 18 May 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Consideration of maintenance in wine
fermentation modeling

A. Rapaport<sup>a,*</sup>, R. David<sup>b</sup>, D. Dochain<sup>c</sup>, J. Harmand<sup>d</sup> and
T. Nidelet<sup>e</sup>

<sup>a</sup> MISTEA, Univ. Montpellier, INRAE, Institut Agro, France.

<sup>b</sup> TECHNORD, Tournai, Belgium.

<sup>c</sup> ICTEAM, Univ. Louvain-la-Neuve, Belgium.

<sup>d</sup> LBE, Univ. Montpellier, INRAE, Narbonne, France.

<sup>e</sup> SPO, Univ. Montpellier, INRAE, Institut Agro, France.

## Abstract

We show that a simple model with a maintenance term can satisfactorily reproduce the simulations of several existing models of wine fermentation from the literature, as well as experimental data. The maintenance describes a consumption of the nitrogen that is not entirely converted into biomass. We show also that considering a maintenance term in the model is equivalent to write a model with a variable yield that can be estimated from data.

*Key words:* Wine fermentation, nitrogen, mathematical modeling, population model, maintenance, variable yield.

## 1 Introduction

The overall principle of wine fermentation consists in the conversion of sugar into ethanol by yeast. It has been observed from a long time that nitrogen consumed during the yeast growth is also playing an important role. The fermentation can be indeed modeled by a two-steps process where the yeast first grows

<footnote>* Corresponding author.
Email addresses: alain.rapaport@inrae.fr (A. Rapaport), r.david@technord.com (R. David), denis.dochain@uclouvain.be (D. Dochain), jerome.harmand@inrae.fr (J. Harmand), thibault.nidelet@inrae.fr (T. Nidelet).</footnote>
on nitrogen as a limiting resource and then degrade the non-limiting sugar into
ethanol and carbon dioxide. However, experimental observations have shown
that the consumed nitrogen was not entirely converted into biomass. Several
mathematical models have been proposed to take into consideration this char-
acteristics. For instance, in [9, 15], the biomass growth follows a logistic law
whose carrying capacity depends on the initial quantity of nitrogen. In [6], a
model that distinguishes part of nitrogen used for yeast growth from another
part responsible of the synthesis of proteins (hexose transporters [8]) has been
developed. Both models have been calibrated with different sets of experimen-
tal data and provide satisfactory fitting. However, both models present some
drawbacks. The dependency of the dynamics on the initial condition of the
first model makes it sensitive to the precise knowledge of the initial quantity
of nitrogen (that needs to be "memorized" in the dynamical equations of the
model). Moreover it does not allow to consider non-batch operations or contin-
uous addition of nitrogen, such as in [3] for instance. The second model relies
on the knowledge of the time-varying concentration of transporters, which is
in general not easily accessible to experimental measurements, and several
assumptions have been necessary to estimate it from biomass measurements.

The objective of the present work is to propose a new model that reconciles
both approaches in a single one.

The observation of the ratio of produced biomass over nitrogen consumption
along the whole fermentation, determined on experimental database or nu-
merical simulations of models [6, 15], shows that this ratio is non-constant
and depends on the initial quantities. This highlights that the conversion of
nitrogen into biomass can be viewed as a variable yield process. The experi-
mental evidence that nitrogen is not entirely converted into biomass therefore
advocates for the consideration of a maintenance term in the modeling (see
for instance [10]), without necessarily requiring a detailed representation of
the internal mechanism or cells.

Indeed, different mechanisms in the internal functioning of the cells have been
investigated in the literature, in particular the role of carboxylate accumu-
lation [21, 23, 24] that could explain that the growth dynamics of yeast in
wine fermentation does not follow the classical mass-balanced models [12, 13].
However, the measurements of these biochemical compounds is experimentally
very difficult and almost impossible in an industrial framework.

The rationale of the results presented here is to test if the introduction of a
maintenance term (see [17, 18, 19] or [1, 11, 26]) can improve the wine fermen-
tation modeling. One of the originality of the proposed approach is to view
the nitrogen consumption as a global consumption for growth by considering
a variable yield. This allows to avoid to consider a specific structure to model
the maintenance. The purpose of the present work is thus to investigate the
ability of a simpler model with a maintenance term to reproduce and predict wine fermentation kinetics.

We hereafter propose a new modeling approach based on a maintenance term (which gives rise to a variable yield), a feature that has not been yet considered in the wine fermentation literature, to the best of our knowledge.

It focuses mainly on the new modeling of the growth of yeast on nitrogen.

This new model has been validated both on data generated by existing models (Section 4) and on experimental data (Section 5).

## 2 The proposed model

We denote by $N$, $S$, $E$, $CO_2$ and $X$ the concentrations of (total) nitrogen, sugar, ethanol, dioxide carbon and biomass, respectively. For simplicity, we derive here a model under isothermal conditions.

For the first step $N \to X$ (yeast growth on nitrogen), we propose the following equations

$$ \frac{dX}{dt} = \mu_N(N, X)X \qquad (1) $$

$$ \frac{dN}{dt} = - \frac{\mu_N(N, X)X}{Y} - m(N, X)X \qquad (2) $$

where $Y$ is the growth yield, $\mu$ the Contois growth function

$$ \mu_N(N, X) = \frac{\mu_N^{\max} N}{N + K_N X} $$

and $m$ a maintenance function, which is positive for $N > 0$ and $X > 0$. We choose here a ratio-dependent kinetics function $\mu_N$ to reproduce the observation that the growth is slowing down under an excess of yeast, with a Contois expression as in [6]. In the literature, the maintenance $m$ is often considered as constant [17, 18], which has been validated in continuous culture (chemostat). In general, continuous culture are intended to be operated at a stationary phase, very differently to batch operating mode. However, as already investigated in [26], maintenance terms have to depend on the level of available resources, say $R$ ($N$ here). In particular, a constant maintenance in a batch model would imply $\frac{dR}{dt} < 0$ when the resource is exhausted i.e. $R = 0$, and thus $R$ could take unrealistic negative values, as underlined in [19]. In [1, 11], the maintenance is directly related to the microbial activity which is stopped in absence of nutrients. This is why here we consider a maintenance function proportional to the growth activity, with a factor that might depend on the
nitrogen concentration (one may expect that it decreases when the substrate $N$ becomes rare)

$$m(N, X) = \alpha(N)\mu_N(N, X)$$

where $\alpha$ is a positive function equal to zero for $N = 0$. Then one can consider the function $y$ defined as follows

$$y(N) := \frac{Y}{1 + \alpha(N)Y}, \quad N \ge 0$$

Formally, model (1)-(2) can be rewritten equivalently as

$$\frac{dX}{dt} = \mu_N(N, X)X \qquad (3)$$

$$\frac{dN}{dt} = - \frac{\mu_N(N, X)X}{y(N)} \qquad (4)$$

where the function $y$ is playing the role of a *variable yield*. Identifying the function $m$ or the function $y$ is thus formally equivalent. However, we shall see in the next section that identifying the function $y$ instead of $m$ presents some practical advantages.

For the second step $S \to E + CO_2$, we follow the model proposed in the literature [6]

$$\frac{dE}{dt} = \frac{dCO_2}{dt} = [\mu_N(N, X) + \beta\nu_E(E)]\mu_S(S)X \qquad (5)$$

$$\frac{dS}{dt} = -k\frac{dE}{dt} \qquad (6)$$

where $\mu_S$ is a Monod function and $\nu_E$ a function inhibited by the ethanol

$$\mu_S(S) = \frac{\mu_S^{\max} S}{K_S + S}, \quad \nu_E(E) = \frac{1}{1 + K_E E} \qquad (7)$$

The inhibition by the consumption of sugar $S$ by ethanol $E$ has been reported many times in the literature [2, 5, 14, 20, 25]. The constant yield of production $k$ of $CO_2$ and consumption of $S$ follows a mass balance assumption, verified experimentally [7] and that can be determined with thermodynamics considerations [22].

Note that this model can be extended to anisothermal conditions, considering that the maximal specific rate parameters $\mu_N^{\max}$, $\mu_S^{\max}$ and affinity constants $K_S, K_E$ are temperature-dependent, as in [6].
### 3 Calibration of the model

From model equation (1), the parameters of the function $\mu_N$ can be identified independently of the yield and maintenance terms. To validate the hypothesis of ratio-dependency of the function $\mu_N$, one can first plot from experimental data the slope of the logarithm of $X$ versus the ratio $r = N/X$ and check if it follows qualitatively a function of the form

$$ \mu(r) = \frac{\mu_N^{\max} r}{K_N + r} $$

Then, a classical least-square method can be applied to fit parameters $\mu_N^{\max}$, $K_N$ on the data. Alternatively, one can plot the inverse of the slope of the logarithm of $X$ versus the inverse of the ratio $r$ to check if it follows qualitatively a linear dependency, as one gets from equation (1)

$$ \left( \frac{d \log X}{dt} \right)^{-1} = \frac{1}{\mu_N^{\max}} + \frac{K_N}{\mu_N^{\max}} \left( \frac{N}{X} \right)^{-1} \quad (8) $$

However, for the accurate identification of the parameters $\mu_N^{\max}$, $K_N$, a linear regression on equation (8) is expected to be less reliable than a nonlinear least square optimization on the solution $X(\cdot)$ of (3), because $\left(\left(\frac{d}{dt} \log X\right)^{-1}, \left(\frac{N}{X}\right)^{-1}\right)$ data might be far to be uniformly distributed.

Note from equations (1)-(2) that one has

$$ \lim_{t \to +\infty} N(t) = 0 $$

(because the derivative of $N$ cannot vanish when $N$ is not exhausted). In absence of the maintenance term $m$, one gets $\frac{dX}{dt} + Y \frac{dN}{dt} = 0$ which implies that one should have

$$ Y = \frac{X(+\infty) - X(0)}{N(0) - N(+\infty)} = \frac{X(+\infty) - X(0)}{N(0)} $$

To test the validity of the model with maintenance, one can plot from experimental data the ratio $\frac{X(+\infty)-X(0)}{N(0)}$ for different values of $N(0)$ to check that it is not constant. If it is the case, one can then look for identifying a non-constant function $y$. For this purpose, we write from equations (3)-(4)

$$ X(+\infty) - X(0) = - \int_{0}^{+\infty} y(N(t)) \frac{dN}{dt}(t) dt $$

and as $t \mapsto N(t)$ is a monotone decreasing function, one can make the change of variable $n = N(t)$ in this last integral to obtain

$$ X(+\infty) - X(0) = \int_{0}^{N(0)} y(n) dn $$
Therefore, if one fits a differential function f such that $f(0) = 0$ that satisfies

$$X(+\infty) - X(0) = f(N(0))$$

for experimental data with different values of $N(0)$, then one simply get $y = f'$.

Let us underline that identifying the function $y$ in this way can be done independently to the knowledge of the kinetics $\mu_N$, differently to the function $m$, what clearly presents some robustness advantages. Once the function $\mu_N$ is identified, the maintenance function can then be determined as

$$m(N, X) = \left( \frac{1}{y(N)} - \frac{1}{Y} \right) \mu_N(N, X)$$

where $Y = y(0)$ (to fulfill $\alpha(0) = 0$).

For model equations (5)-(6), the coefficient $k$ is kept from the literature, and the parameters $\beta, \mu_S^{max}, K_S, K_E$ are identified (with a least-square method) from experimental data of $CO_2$ production rate.

# 4 Validation of the model on synthetic data

We have used synthetic data generated by models of the literature that have been previously validated on experimental data [6, 15] for a range of initial conditions and operating conditions.

Fitting comparisons of the proposed model with the different data sets are reported in Section 6.

## 4.1 Validation on simulations of a model with transporter

We have considered the model with transporters developed in [6], which is more complex with two additional state variables: the concentrations of hexose transporters and the nitrogen dedicated to these transporters. Data have been generated by simulating this model with the parameters given in [6] and operating conditions given in Table 1.

This model distinguishes explicitly two forms of nitrogen, one available for the yeast $N_X$ and the other one $N_{tr}$ for the transporters. To compare with the variable $N$ of our model, we have considered the total nitrogen $N = N_X + N_{tr}$.
<table><caption>Table 1</caption><tr><td>X(0)</td><td>0.02 g.l<sup>-1</sup></td></tr><tr><td>N(0)</td><td>0.071 – 0.57 g.l<sup>-1</sup></td></tr><tr><td>S(0)</td><td>200 g.l<sup>-1</sup></td></tr><tr><td>time horizon</td><td>350 hours</td></tr><tr><td>temperature</td><td>constant equal to 24°</td></tr><tr><td>others</td><td>no initial transporter<br/>no nitrogen addition</td></tr></table>

Operating conditions for the simulation of the model with transporters

### 4.1.1 Estimation of the Contois function

We have used a nonlinear least square method based on a Newton algorithm with a finite difference approximation of the Jacobian matrix (function `leastsq` of `scilab`). Figure 1 shows a good fitting of the Contois function $\mu_N$ on data $(N/X, \frac{dX}{dt}/X)$ of the transporter model, with parameters given in Table 2.

<figure><img src="image_2.png" /><figcaption>Fig. 1. Result of the fitting of the Contois function on data from the model with transporters</figcaption></figure>

<table><caption>Table 2</caption><tr><td>$\mu_N^{max}$</td><td>0.103 h<sup>-1</sup></td></tr><tr><td>$K_N$</td><td>0.0381 g.l<sup>-1</sup></td></tr></table>

Parameters of the Contois function $\mu_N$
### 4.1.2 Estimation of the variable yield function

On Figure 2, data $X(T) - X(0)$ versus $N(0)$ from the model with transporters have been plotted for $T = 350$ hours (we have checked that $N$ is quasi-null at $T$ and that $X$ does no longer increase after $T$). One can see that the points are aligned but the line that passes through these points does not touch 0, which is not possible for a constant yield (for a constant yield, the points have to be aligned on a line that passes through 0 because when $N(0) = 0$ there is no biomass production).

<figure><img src="image_3.png" /><figcaption>Fig. 2. Result of the fitting of the function $f$ on data from the model with transporters</figcaption></figure>

Then we have fitted a $C^2$ function $f$ such that $f(0) = 0$ with the following expression

$$f(N) = \begin{cases} aN + b \left(1 - \left(\frac{N_{\text{f}} - N}{N_{\text{f}}}\right)^3\right), & N < N_{\text{f}} \\ aN + b, & N \ge N_{\text{f}} \end{cases}$$

whose parameters are given in Table 3.

The calibration of the parameters $a, b$ of the function $f$ has been performed with a linear regression (function reglin of scilab).

<table><caption>Table 3</caption><tr><td>a</td><td>7.55</td></tr><tr><td>b</td><td>0.808 g.l<sup>-1</sup></td></tr><tr><td>N<sub>f</sub></td><td>0.176 g.l<sup>-1</sup></td></tr></table>

Parameters of the variable yield function $y$
Then, we obtain the variable yield function y as the C¹ function

$$y(N) = f'(N) = \begin{cases} a + b \frac{3(N_{\downarrow} - N)^2}{N_{\downarrow}^3}, & N < N_{\downarrow} \\ a, & N \ge N_{\downarrow} \end{cases}$$

and the function α which describes the maintenance as

$$\alpha(N) = \frac{1}{y(N)} - \frac{1}{y(0)} = \begin{cases} \frac{N_{\downarrow}^3}{aN_{\downarrow}^3 + 3b(N_{\downarrow} - N)^3} - \frac{N_{\downarrow}}{3b + aN_{\downarrow}}, & N < N_{\downarrow} \\ \frac{3b}{a(3b + aN_{\downarrow})}, & N \ge N_{\downarrow} \end{cases}$$

that are both depicted on Figure 3.

<figure><img src="image_4.png" /><figcaption>Fig. 3. Graphs of the obtained variable yield function y and of the function α</figcaption></figure>

Note that the model with transporters has been validated only for N(0) in the interval [0.071, 0.57] g.l⁻¹, and that we have no a priori information about the behavior of the yield for values of N(0) smaller than 0.071 g.l⁻¹. The threshold parameter N₄ has been simply chosen so that the simulations of the variables X and N of the model (3)-(4) were the closest from the ones of the transporter model.

### 4.1.3 Estimation of the other parameters and comparison of the models

For the model of the second step S → E + CO₂, the stoichiometric parameter k has been taken from the literature, while the other parameters β, μₘₐₓ, Kₕ, Kₕ have been estimated with a least-square optimization on the CO₂ chronicles only (the CO₂ production rate being a variable that is usually measured in experiments), starting from values in [6]. Values are given in Table 4.

Here also, we have used a nonlinear least square method based on a Newton algorithm with a finite difference approximation of the Jacobian matrix (function leastsq of scilab). All data have been re-normalized to 1 (i.e. for each variable, the figures have been divided by the largest one).

Finally, we present on Figures 4, 5, 6 simulations of the new model for three largely different initial values of nitrogen from 0.170 g.l⁻¹ to 0.567 g.l⁻¹. The
<table><caption>Table 4</caption><tr><td>k</td><td>2.17</td></tr><tr><td>β</td><td>2.41</td></tr><tr><td>μ<sub>S</sub><sup>max</sup></td><td>0.197 h<sup>-1</sup></td></tr><tr><td>K<sub>S</sub></td><td>21.1 g.l<sup>-1</sup></td></tr><tr><td>K<sub>E</sub></td><td>72.7 g.l<sup>-1</sup></td></tr></table>

Parameters for the second step $S \to E + CO_2$ model

evolution of the ethanol concentration $E$ has not been reproduced as it is proportional to the $CO_2$ concentration.

These simulations shows the ability of the new model to reproduce, with a single set of parameters, close simulations to the model with transporters, in terms of production of biomass and dioxide carbon, estimation of the peak of the $CO_2$ production rate and depletion of (total) nitrogen and sugar.

<figure><img src="image_5.png" /><figcaption>Fig. 4. Comparison with the model with transporters (in dashed) for $N(0) = 0.170 g.l^{-1}$ (constant temperature of 24°C)</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Fig. 5. Comparison with the model with transporters (in dashed) for $N(0) = 0.283 \text{ g.l}^{-1}$ (constant temperature of 24°C)</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 6. Comparison with the model with transporters (in dashed) for $N(0) = 0.567 \text{ g.l}^{-1}$ (constant temperature of 24°C)</figcaption></figure>
## 4.2 Validation on the SOFA model

The model proposed in [15] does not consider explicitly transporters with an additional state variable as the previous model, and present instead a more sophisticated expression of the dynamics that depend on the initial condition, with an additional latency term at the beginning of the simulations.

Differently to the previous model which is built as a "mass-balanced" model, this one relies on an empirical dynamics of logistic shape for the biomass growth, with some parameters that depend on the initial concentration of nitrogen $N(0)$, instead of the two-dimensional model (3)-(4).

Therefore, this is not a Markovian model. It has been validated on different operating conditions, and has been encoded into the SOFA software exploited for decision making [9]. We have launched simulations of this model for the same operating conditions than for the previous model (Table 1). Although simulations look qualitatively similar, they do not overlap, especially for the biomass chronicle. This could be explained by the fact that this model is intended to predict a number of cells and not precisely a biomass (an average number of 4.15 $10^9$ cells for one $g$ of biomass has been used to have $X$ expressed in $g.l^{-1}$ as for the previous model). We have proceeded to a new validation of our model on these data.

### 4.2.1 Estimation of the Contois function

Figure 7 shows that the data $\left(\frac{N}{X}, \frac{dX}{dt}, \frac{dX}{dt}\right)$ do not follow precisely the graph of a function (this is most probably due to the fact that the model is not Markovian). Indeed, this happens mainly for large value $N_0$ of the initial nitrogen. We believe that this could be explained by the dynamics of the biomass $X$ of this model, which is a logistic law with a carrying capacity given by an heuristic expression that depends on $N_0$, and not a dynamics coupled with the dynamics of $N$ (indeed the interval of tested values of $N_0$ might be larger than the validity of this model). However, we have fitted the graph of a Contois function to these data with the parameters given in Table 5, which has been able to reproduce satisfactorily the trajectories of the model for a large amplitude of values of $N_0$, as we shall see later on.

As for the previous model, we have used a nonlinear least square method based on a Newton algorithm with a finite difference approximation of the Jacobian matrix (function `leastsq` of `scilab`). As one can see on Table 5, the values of $\mu_N^{max}$ and $K_N$ are significantly larger and smaller (respectively) than in Table 2, which is consistent with the observation that this model predict a faster convergence of the biomass to its maximal value, despite the latency term
(compare Figures 4, 5, 6 with Figures 10, 11, 12).

<figure><img src="image_8.png" /><figcaption>Fig. 7. Result of the fitting of the Contois function on data from the SOFA model</figcaption></figure>

<table><caption>Table 5
Parameters of the Contois function μ<sub>N</sub></caption><tr><td>μ<sub>N</sub><sup>max</sup></td><td>0.270 h<sup>-1</sup></td></tr><tr><td>K<sub>N</sub></td><td>0.00952 g.l<sup>-1</sup></td></tr></table>

### 4.2.2 Estimation of the variable yield function

Data X(T) - X(0) from the simulation of the SOFA model have been plotted on Figure 8 at T = 350 hours, for different values of N(0) in the interval [0.071, 0.57] g.l<sup>-1</sup> (here also we have checked that the fermentation was quasi-ended at T). On can see that the points follows an increasing concave curve and further increase very slowly, quite differently to the model with transporters (see Figure 2).

We have then fitted a C<sup>2</sup> function f with f(0) = 0 for the expression

$$f(N) = \begin{cases} bN - aN^2, & N < N_{\dagger} \\ bN - aN^2 + bN + \frac{A}{B}(e^{-BN_{\dagger}} - e^{-BN}) & N < N_{\dagger} \end{cases}$$

with

$$A = (b - 2aN_{\dagger})e^{BN_{\dagger}}, \quad B = \frac{2a}{b - 2aN_{\dagger}}$$

and parameters a, b, N<sub>†</sub> given in Table 6.

Parameters a and b have been determined with a linear regression (function reglin of scilab).
<figure><img src="image_9.png" /><figcaption>Fig. 8. Result of the fitting of the function f on data from the SOFA model</figcaption></figure>

<table><tr><td>a</td><td>15.1 l.g<sup>-1</sup></td></tr><tr><td>b</td><td>15.2</td></tr><tr><td>N<sub>†</sub></td><td>0.465 g.l<sup>-1</sup></td></tr></table>

Table 6

Parameters of the variable yield function y

Then, we obtain the expression of the variable yield function

$$y(N) = f'(N) = \begin{cases} b - 2aN, & N < N_{\dagger} \\ Ae^{-BN}, & N \ge N_{\dagger} \end{cases}$$

as well as the function α

$$\alpha(N) = \frac{1}{y(N)} - \frac{1}{y(0)} = \begin{cases} \frac{1}{\frac{b-2aN}{e^{b(N-N_{\dagger})}}} - \frac{1}{b}, & N < N_{\dagger} \\ \frac{e^{b(N-N_{\dagger})}}{b-2aN_{\dagger}} - \frac{1}{b}, & N \ge N_{\dagger} \end{cases}$$

whose graphs are drawn on Figure 9.

### 4.2.3 Estimation of the other parameters and comparison of the models

For the second step, the same stoichiometric parameter k has been taken for the literature, and the other parameters β, μ<sub>S</sub><sup>max</sup>, K<sub>S</sub>, K<sub>E</sub> have been estimated with a least-square optimization on the CO<sub>2</sub> chronicles only, as for data generated by the model with transporters (see Table 7).

Figures 10, 11, 12 show the comparison between the SOFA model and our
<figure><img src="image_10.png" /><figcaption>Fig. 9. Graphs of the obtained variable yield function *y* and of the function *α*</figcaption></figure>

<table><tr><td>k</td><td>2.17</td></tr><tr><td>β</td><td>3.22</td></tr><tr><td>μ<sub>S</sub><sup>max</sup></td><td>0.197 h<sup>-1</sup></td></tr><tr><td>K<sub>S</sub></td><td>17.6 g.l<sup>-1</sup></td></tr><tr><td>K<sub>E</sub></td><td>36.4 g.l<sup>-1</sup></td></tr></table>

Table 7

Parameters for the second step $S \to E + CO_2$ model

calibrated model for the same initial condition than for the former comparison with the model with transporters. Here also, we see that the proposed model reproduces quite faithfully the simulations of the SOFA model, with the advantage of being a simpler Markovian model. Indeed, the difference between the model with transporters and the SOFA model can be translated into different maintenance terms (see Figures 3 and 9): for large values of nitrogen, the model with transporters behaves like a model with a maintenance proportional to the growth, while the SOFA model amounts to have a strongly increasing maintenance. Recall that the simulations for the largest value of $N(0)$ showed the most differences between these two models (for $N(0) = 0.567 \text{ g.l}^{-1}$, the model with transporters predicts a biomass production of $5.11 \text{ g.l}^{-1}$, while the SOFA model predicts $3.88 \text{ g.l}^{-1}$; see Figures 6 and 12). While the model with transporters has been validated experimentally for $N(0)$ in the interval $[0.170, 0.567] \text{ g.l}^{-1}$, we believe the validation of the SOFA model for initial concentrations of nitrogen larger than $0.4 \text{ g.l}^{-1}$ might need to be revisited (although our model once calibrated is able to reproduce the SOFA simulations).
<figure><img src="image_11.png" /><figcaption>Fig. 10. Comparison with the SOFA model (in dashed) for $N(0) = 0.170 \text{ g.l}^{-1}$ (constant temperature of 24°C)</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 11. Comparison with the SOFA model (in dashed) for $N(0) = 0.283 \text{ g.l}^{-1}$ (constant temperature of 24°C)</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>Fig. 12. Comparison with the SOFA model (in dashed) for $N(0) = 0.567 \text{ g.l}^{-1}$ (constant temperature of 24°C)</figcaption></figure>
5 Calibration of the model on real data

We have considered data from experiments conducted at SPO Lab (INRAE, Montpellier, France) in 2004, that have been used to calibrate the model with transporters and the SOFA model (see [6, 15]). The data consist in a set of three experiments with the same operating conditions given in Table 1 and different initial concentrations N(0) of nitrogen, exactly as for the simulations of Sections 4.1 and 4.2. For each experiment, one has

- height measurement points for X,

- no measurement point for N, S or E,

- about 400 measurement points for CO₂ and dCO₂/dt,

We have first calibrated a function $f(\cdot)$ to the data $(N(0), X(T)-X(0))$, with the same expression than in Section 4.2, to determine a yield function $y(\cdot)$ (see Figure 13), using a linear regression to estimate parameters $a$ and $b$.

<figure><img src="image_14.png" /><figcaption>Fig. 13. Results of the fitting of the function $f$ on the experimental data (left) and of the corresponding variable yield function $y$ (right)</figcaption></figure>

As we do not have measurements of $N$ over the time, we cannot estimate the Contois parameters independently of the $CO_2$ measurements, as we did with the synthetic data. All the parameters of the model have been fitted simultaneously with a least square method (values are given in Table 8), excepted for the sugar conversion yield for which we have used the value of the literature $k = 2.17$, as before.

The nonlinear least square method uses a Newton algorithm with a finite difference approximation of the Jacobian matrix (function `leastsq` of `scilab`), and the data set has been re-normalized to the maximal value of 1. Figures 14, 15, 16, show the results of the fitting for the three experiments. One can appreciate the goodness of fit for a unique set of parameters. In particular, the production of biomass and $CO_2$, as well as the height and date of the peak of $dCO_2/dt$ are well predicted with this model.
<table><caption>Table 8
Parameters fitted on the experimental data</caption><tr><td>μ<sub>N</sub><sup>max</sup></td><td>0.175 h<sup>-1</sup></td></tr><tr><td>K<sub>N</sub></td><td>0.0133 g.l<sup>-1</sup></td></tr><tr><td>β</td><td>1.622</td></tr><tr><td>μ<sub>S</sub><sup>max</sup></td><td>0.393 h<sup>-1</sup></td></tr><tr><td>K<sub>S</sub></td><td>19.2 g.l<sup>-1</sup></td></tr><tr><td>K<sub>E</sub></td><td>71.9 g.l<sup>-1</sup></td></tr></table>

<figure><img src="image_15.png" /><figcaption>Fig. 14. Simulation for $N(0) = 0.170 \text{ g.l}^{-1}$ (experimental data in blue)</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>Fig. 15. Simulation for $N(0) = 0.283 \, gl^{-1}$ (experimental data in blue)</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. 16. Simulation for $N(0) = 0.567 \, gl^{-1}$ (experimental data in blue)</figcaption></figure>
# 6 Fitting comparisons

For the calibration of the variable yield function on both synthetic and experimental data (Sections 4 and 5), we have used a linear regression (function reglin of scilab) for the determination of parameters *a*, *b* of the function *f* (for the model with transporter) and *f'* (for the SOFA model and experimental data). The residual error is given in Table 9. This shows that the model

<table><caption>Table 9
Residual standard error (RSE) for the determination of *a* and *b*</caption><thead><tr><th><i>data</i></th><th>Tr model</th><th>SOFA model</th><th>exp.</th></tr></thead><tbody><tr><td><i>RSE</i></td><td>2.21 10<sup>-10</sup></td><td>0.199</td><td>0.225</td></tr></tbody></table>

with transporters behave very closely to a variable yield model. The fitting performances for the SOFA model and experimental data are more difficult to interpret, because the validity of the SOFA model for the large range of initial concentrations of nitrogen we considered is questionable, and the quantity of experimental data is quite poor compared to the synthetic data.

For the synthetic data, the calibration of the growth characteristics (parameters $\mu_N^{max}$, $K_N$ of the Contois function) has been done first independently to the CO₂ data. Then, parameters for the second step (parameters *k*, *β*, $\mu_S^{max}$, $K_S$, $K_E$ for the CO₂ production) have been calibrated. In both cases, a nonlinear least square method based on a Newton algorithm with a finite difference approximation of the Jacobian matrix (function `leastsq` of scilab) has been used. Table 10 shows a good fitting quality.

<table><caption>Table 10
Root Mean Square Error (RMSE) for the calibration of the growth function $\mu$ and the CO₂ chronicles</caption><thead><tr><th><i>data</i></th><th>Tr model</th><th>SOFA model</th><th>exp.</th></tr></thead><tbody><tr><td>RMSE(μ)</td><td>0.0414</td><td>0.292</td><td>-</td></tr><tr><td>RMSE(CO₂)</td><td>0.0543</td><td>0.0895</td><td>0.0519</td></tr></tbody></table>

We recall that for experimental data, we do not have measurement of *N* over time, so that it was not possible to estimate the growth function independently of the CO₂ measurements. The estimation of all the parameters has been made on the CO₂ measurements only. We have used the same nonlinear least square method, with data re-normalized to 1 (i.e. the figures have been divided by the largest one), so that all points have equal weight in the criterion. The errors shows a good fitting of the CO₂ curves with the model with maintenance.
# 7 Conclusion

In this work, we have demonstrated that the consideration of a maintenance term, or equivalently a variable yield, in wine fermenting modeling can satisfactorily replace more sophisticated models with a simpler structure. Indeed, the effects of the underlying mechanisms of transporters or carbohydrate accumulation, which are difficult to capture experimentally, are somehow encoded into a maintenance term, and are translated into a variable yield between biomass and nitrogen. We have shown that this variable yield, as a function of the nitrogen concentration, can be estimated from experimental data of biomass growth and nitrogen depletion, without the need to measure internal compounds. This consideration brings a flexibility to suit to different kind of models or experimental data (once calibrated) with a single common structure, that could correspond to different operating conditions or hypotheses in wine fermentation. This new approach opens new perspectives of control of fermentation with nitrogen addition, based on a simple Markovian model, as well as model extensions with aromatic compounds [16] or multi-strains [4].

# References

[1] H. Beeftink, R. van der Heijden, and J. Heijnen. Maintenance requirements: energy supply from simultaneous endogenous respiration and substrate consumption. *FEMS Microbiology Letters*, 73(3):203-209, 1990.

[2] I. Caro, L. Pérez, and D. Cantero. Development of a kinetic model for the alcoholic fermentation of must. *Biotechnology and Bioengineering*, 38(7):742-748, 1991.

[3] T. Clement, M. Perez, J. Mouret, J. Sablayrolles, and C. Camarasa. Use of a continuous multistage bioreactor to mimic winemaking fermentation. *International Journal of Food Microbiology*, 150(1):42-49, 2011.

[4] C. Conacher, N. Luyt, R. Naidoo, D. Rossouw, M. Setati, and F. Bauer. The ecology of wine fermentation: a model for the study of complex microbial ecosystems. *Applied Microbiology and Biotechnology*, 105, 04 2021.

[5] A. C. Cramer, S. Vlassides, and D. E. Block. Kinetic model for nitrogen-limited wine fermentations. *Biotechnology and Bioengineering*, 77(1):49-60, 2002.

[6] R. David, D. Dochain, J.-R. Mouret, A. Vande Wouwer, and J.-M. Sablayrolles. Nitrogen-backed modeling of wine-making in standard and nitrogen-added fermentations. *Bioprocess and Biosystems Engineering*, 37(1):5-16, Jan 2014.

[7] N. El Haloui, D. Picque, and G. Corrieu. Alcoholic fermentation in wine-making: On-line measurement of density and carbon dioxide evolution. *Journal of Food Engineering*, 8(1):17-30, 1988.
[8] K. Elbing, C. Larsson, R. M. Bill, E. Albers, J. L. Snoep, E. Boles, S. Hohmann, and L. Gustafsso. Role of hexose transport in control of glycolytic flux in saccharomyces cerevisiae. *Applied and Environmental Microbiology*, 70(9):5323–5330, 2004.

[9] A. Goelzer, B. Charnomordic, S. Colombié, V. Fromion, and J. Sablayrolles. Simulation and optimization software for alcoholic fermentation in winemaking conditions. *Food Control*, 20(7):635–642, 2009.

[10] G. Goma, R. Moletta, and M. Novak. Comments on the "maintenance coefficient" changes during alcohol fermentation. *Biotechnology Letters*, 1(10):415–420, 1979.

[11] J. Heijnen and J. Roels. A macroscopic model describing yield and maintenance relationships in aerobic fermentation processes. *Biotechnology and Bioengineering*, 23(14):739–763, 1981.

[12] D. Henriques and E. Balsa-Canto. The monod model is insufficient to explain biomass growth in nitrogen-limited yeast fermentation. *Applied and environmental microbiology*, 87(20):e0108421, 2021.

[13] D. Henriques, R. Minebois, S. Mendoza, L. Macías, R. Pérez-Torrado, E. Barrio, B. Teusink, A. Querol, and E. Balsa-Canto. A multiphase multiobjective dynamic genome-scale model shows different redox balancing among yeast species of the saccharomyces genus in fermentation. *mSystems*, 6(4):e00260–21, 2021.

[14] C. Leão and N. van Uden. Effects of ethanol and other alkanols on the glucose transport system of saccharomyces cerevisiae. *Biotechnology and Bioengineering*, 24(11):2601–2604, 1982.

[15] S. Malherbe, V. Fromion, N. Hilgert, and J.-M. Sablayrolles. Modeling the effects of assimilable nitrogen and temperature on fermentation kinetics in enological conditions. *Biotechnology and Bioengineering*, 86(3):261–272, 2004.

[16] J. Mouret, V. Farines, J. Sablayrolles, and I. Trelea. Prediction of the production kinetics of the main fermentative aromas in winemaking fermentations. *Biochemical Engineering Journal*, 103(15):211–218, 2015.

[17] S. Pirt. The maintenance energy of bacteria in growing cultures. *Proceedings of the Royal Society of London. Series B*, 163(991):224–231, 1965.

[18] S. J. Pirt. *Principles of microbe and cell cultivation*. Blackwell Science Ltd, New edition, 1985.

[19] A. Rapaport, T. Nidelet, S. El Aida, and J. Harmand. About biomass overyielding of mixed cultures in batch processes. *Mathematical Biosciences*, 322:108322, 2020.

[20] J. M. Salmon, O. Vincent, J. C. Mauricio, M. Bely, and P. Barre. Sugar transport inhibition and apparent loss of activity in saccharomyces cerevisiae as a major limiting factor of enological fermentations. *American Journal of Enology and Viticulture*, 44(1):56–64, 1993.

[21] U. Schulze, G. Lidèn, J. Nielsen, and J. Villadsen. Physiological effects of nitrogen starvation in an anaerobic batch culture of saccharomyces cerevisiae. *Microbiology*, 142(8):2299–2310, 1996.
[22] J. VanBriesen. Evaluation of methods to predict bacterial yield using thermodynamics. *Biodegradation*, 13(3):171–190, 2002.

[23] C. Varela, F. Pizarro, and E. Agosin. Biomass content governs fermentation rate in nitrogen-deficient wine musts. *Applied and environmental microbiology*, 70(6):3392–3400, 2011.

[24] F. Vargas, F. Pizarro, J. Pérez-Correa, and E. Agosin. Expanding a dynamic flux balance model of yeast fermentation to genome-scale. *BMC Systems Biology*, 5(775):1–12, 2011.

[25] C. Viegas, I. Sa-Correia, and J. Novais. Synergistic inhibition of the growth of saccharomyces bayanus by ethanol and octanoic or decanoic acids. *Biotechnology letters*, 7(8):611–614, 1985.

[26] G. Wang and W. M. Post. A theoretical reassessment of microbial maintenance and implications for microbial ecology modeling. *FEMS Microbiology Ecology*, 81(3):610–617, 09 2012.