Mechanistic model of S. cerevisiae utilizing ligocellulosis in batch fermentation#
This model describes an S. cerevisiae fermentation in a batch reactor. It builds a stochiometric matrix to follow all conversions. The model accounts for the changes of Carbon sources, the development of biomass and ethonal over time.
Package import#
This portion of the code handles the import of all the relevant python packages.
from scipy.integrate import odeint
# Package for plotting
import math
# Package for the use of vectors and matrix
import numpy as np
import pandas as pd
import array as arr
from matplotlib.figure import Figure
import sys
import os
import matplotlib.pyplot as plt
from matplotlib.ticker import FormatStrFormatter
import glob
from random import sample
import random
import time
import plotly
import plotly.graph_objs as go
import json
from plotly.subplots import make_subplots
---------------------------------------------------------------------------
ModuleNotFoundError Traceback (most recent call last)
Cell In[1], line 1
----> 1 from scipy.integrate import odeint
2 # Package for plotting
3 import math
ModuleNotFoundError: No module named 'scipy'
Model definition#
The model is defined by using a class. This has several adavantages, one being that different classes can be defined independtly, e.g., one class represents the model definition while another will handle parameter optimization with experimental data. The same “optimization” class can be used with different model classes, so the code becomes more easily reusable.
The class Scerevisiae_Ligno includes several functions, each with a specific function.
init This function initialises the model class by defining all the relavant parameters and initial conditions.
rxn This function includes all the model equations. The model uses matrix notation to define the ODEs. In this part the all rates for different carbon sources are defined individually.
solve This function generates the timesteps for solving the ODEs. Depending on the initial conditions that were set, the development of the observables is simulated.
create_plot This function stores the results of the simulation on a dataframe which is used for plotting all relevant variables.
class SCerevisiae_Ligno:
#initialize model
def __init__(self, Control=False):
# define value of model parameters
self.nuMaxGlu = 2.348# h-1
self.nuMaxXyl = 1.622 # h-1
self.Ks_Glu = 0.565 # kg Glu m-3
self.Ks_Xyl = 18.1 # kg Xyl m-3
self.Ki_Glu = 283.7 # kg Glu m-3
self.Ki_Xyl = 18.1 # kg Xyl m-3
self.Ki_GluXyl = 10 # kg Glu m-3
self.Y_XGlu = 0.115 # kg X/kg Glu
self.Y_XXyl = 0.162 # kg X/kg Xyl
self.Ki_EtOHmaxGlu = 103 # kg Glu m-3
self.Ki_EtOHmaxXyl = 60.2 # kg Xyl m-3
self.Y_EtOHGlu = 0.47# kg EtOH/kg Glu
self.Y_EtOHXyl = 0.4 # kg EtOH/kg Xyl
self.gammaG = 1.42 # no unit
self.gammaX = 0.608 # no unit
# Acetate parameters
self.nuHAcMax = 0.04428 # h-1
self.Ks_HAc = 2.5 # kg HAc m-3
self.Ki_HAcGlu = 2.74 # kg HAc m-3
self.Ki_HAcXyl = 0.073 # kg HAc m-3
self.Y_HAcHMF = 0.234 # kg Ac/kg HMF
# Furfural parameters
self.nuFurMax = 0.16812 # h-1
self.Ks_Fur = 0.05 # kg Furfural m-3
self.Ki_FurGlu = 0.75 # kg Furfural m-3
self.Ki_FurXyl = 0.35 # kg Furfural m-3
self.Ki_FurHMF = 0.25 # kg Furfural m-3
self.Y_FurFA = 1.02 # kg FA/kg Fur
# Furfuryl alcohol parameters
self.Ki_FAGlu = 5# kg FA m-3
self.Ki_FAXyl = 6 # kg FA m-3
# HMF parameters
self.nuHMFMax = 0.31536 # h-1
self.Ks_HMF = 0.5 # kg HMF m-3
self.Ki_HMFGlu = 2 # kg HMF m-3
self.Ki_HMFXyl = 10 # kg HMF m-3
# Initial state variable
self.X0 = 0.5 # g/L
self.Glu0 = 40 # g/L
self.Xyl0 = 20 # g/L
self.EtOH0 = 0 # g/L
self.Fur0 = 1 # g/L
self.HAc0 = 1 # g/L
self.HMF0 = 0.5 # g/L
self.FA0 = 0 # g/L
#t and V conditions
self.t_end = 30
self.t_start = 0
self.V0 = 2
self.T0 = 30
#define the stoichiometric matrix
def rxn(self, C, t, u, fc):
# number of components
n = 8
m = 5
# initialize the stoichiometric matrix, s
s = np.zeros((m, n))
s[0, 0] = (self.Y_XGlu)
s[1, 0] = (self.Y_XXyl)
s[2, 0] = (0)
s[3, 0] = (0)
s[4, 0] = (0)
s[0, 1] = (-1)
s[1, 1] = (0)
s[2, 1] = (0)
s[3, 1] = (0)
s[4, 1] = (0)
s[0, 2] = (0)
s[1, 2] = (-1)
s[2, 2] = (0)
s[3, 2] = (0)
s[4, 2] = (0)
s[0, 3] = (self.Y_EtOHGlu)
s[1, 3] = (self.Y_EtOHXyl)
s[2, 3] = (0)
s[3, 3] = (0)
s[4, 3] = (0)
s[0, 4] = (0)
s[1, 4] = (0)
s[2, 4] = (-1)
s[3, 4] = (0)
s[4, 4] = (0)
s[0, 5] = (0)
s[1, 5] = (0)
s[2, 5] = (0)
s[3, 5] = (-1)
s[4, 5] = (self.Y_HAcHMF)
s[0, 6] = (0)
s[1, 6] = (0)
s[2, 6] = (0)
s[3, 6] = (0)
s[4, 6] = (-1)
s[0, 7] = (0)
s[1, 7] = (0)
s[2, 7] = (self.Y_FurFA)
s[3, 7] = (0)
s[4, 7] = (0)
# initialize the rate vector
rho = np.zeros((m))
# Glucose uptake process
rho[0] = self.nuMaxGlu * C[0] * (C[1] / (self.Ks_Glu + C[1] + ((C[1] ** 2) / self.Ki_Glu)) *
(1 - (C[3] / self.Ki_EtOHmaxGlu) ** self.gammaG) *
(1 / (1 + (C[4] / self.Ki_FurGlu))) *
(1 / (1 + (C[5] / self.Ki_HAcGlu))) *
(1 / (1 + (C[6] / self.Ki_HMFGlu))) *
(1 / (1 + (C[7] / self.Ki_FAGlu))))
# Xylose uptake process
rho[1] = self.nuMaxXyl * C[0] * (C[2] / (self.Ks_Xyl + C[2] + ((C[2] ** 2) / self.Ki_Xyl)) *
(1 - (C[3] / self.Ki_EtOHmaxXyl) ** self.gammaX) *
(1 / (1 + (C[4] / self.Ki_FurXyl))) *
(1 / (1 + (C[5] / self.Ki_HAcXyl))) *
(1 / (1 + (C[6] / self.Ki_HMFXyl))) *
(1 / (1 + (C[7] / self.Ki_FAXyl))) *
(1 / (1 + (C[1] / self.Ki_GluXyl))))
# Fur uptake process
rho[2] = self.nuFurMax * C[0] * (C[4] / (self.Ks_Fur + C[4]))
# HAc uptake process
rho[3] = self.nuHAcMax * C[0] * (C[5] / (self.Ks_HAc + C[5]))
# HMF uptake process
rho[4] = self.nuHMFMax * C[0] * (C[6] / (self.Ks_HMF + C[6])) * (1 / (1 + (C[4] / self.Ki_FurGlu)))
#Solving the mass balances
dXdt = s[0, 0] * rho[0] + s[1, 0] * rho[1] + s[2, 0] * rho[2] + s[3, 0] * rho[3] + s[4, 0] * rho[4]
dGludt = s[0, 1] * rho[0] + s[1, 1] * rho[1] + s[2, 1] * rho[2] + s[3, 1] * rho[3] + s[4, 1] * rho[4]
dXyldt = s[0, 2] * rho[0] + s[1, 2] * rho[1] + s[2, 2] * rho[2] + s[3, 2] * rho[3] + s[4, 2] * rho[4]
dEtOHdt = s[0, 3] * rho[0] + s[1, 3] * rho[1] + s[2, 3] * rho[2] + s[3, 3] * rho[3] + s[4, 3] * rho[4]
dFurdt = s[0, 4] * rho[0] + s[1, 4] * rho[1] + s[2, 4] * rho[2] + s[3, 4] * rho[3] + s[4, 4] * rho[4]
dHAcdt = s[0, 5] * rho[0] + s[1, 5] * rho[1] + s[2, 5] * rho[2] + s[3, 5] * rho[3] + s[4, 5] * rho[4]
dHMFdt = s[0, 6] * rho[0] + s[1, 6] * rho[1] + s[2, 6] * rho[2] + s[3, 6] * rho[3] + s[4, 6] * rho[4]
dFAdt = s[0, 7] * rho[0] + s[1, 7] * rho[1] + s[2, 7] * rho[2] + s[3, 7] * rho[3] + s[4, 7] * rho[4]
dVdt=0
dTdt = 0
return [dXdt, dGludt, dXyldt, dEtOHdt, dFurdt, dHAcdt, dHMFdt, dFAdt, dVdt, dTdt]
#solve the ODES
def solve(self):
t = np.linspace(0, 30) #generation of the time-points
u = 0
fc = 1
C0 = [self.X0, self.Glu0, self.Xyl0, self.EtOH0, self.Fur0, self.HAc0, self.HMF0, self.FA0, self.V0, self.T0] #initial conditions vector
C = odeint(self.rxn, C0, t, rtol=1e-7, mxstep=500000, args=(u, fc,)) #solve ODEs
return t, C
#generate the plot of model variables
def create_plot(self, t, C):
figure = make_subplots(rows=1, cols=2) #make figure with 2 subplots
#assign simulation results to variable for plotting
X = C[:, 0]
Glu = C[:, 1]
Xly = C[:, 2]
EtOH = C[:, 3]
Fur = C[:, 4]
HAc = C[:, 5]
HMF = C[:, 6]
FA = C[:, 7]
V = C[:, 8]
#collect all variables to plot in 1st subplot in a dataframe
df = pd.DataFrame({'t': t, 'Glu': Glu, 'X': X, 'Xly':Xly, 'EtOH': EtOH, 'Fur':Fur, 'HAc': HAc, 'HMF':HMF, 'FA':FA})
#add the different traces to 1st subplot
figure.add_trace(go.Scatter(x=df['t'], y=df['Glu'], name='Glucose'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['X'], name='Biomass'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['Xly'], name='Xylose'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['EtOH'], name='Ethanol'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['Fur'], name='Furfural'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['HAc'], name='Acetic acid'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['HMF'], name='HMF'), row=1, col=1)
figure.add_trace(go.Scatter(x=df['t'], y=df['FA'], name='Furfuryl alcohol'), row=1, col=1)
#add the title and axes labels
figure.update_layout(title=('Simulation of the model for the Scerevisiae in fedbatch using lignocellulosic'),
xaxis_title='time (h)',
yaxis_title='Concentration (g/L)')
#dataframe with varible to plot in 2nd subplot
df2 = pd.DataFrame({'t': t, 'V':V})
#add trace to the 2nd subplot
figure.append_trace(go.Scatter(x=df2['t'], y=df2['V'], name='Volume'), row=1, col=2)
return figure
Extract results and make plots#
This portion of the code extracts the results from the model class and creates the figures to display the simulation results.
model = SCerevisiae_Ligno() # Instantiate the SCerevisiae_Ligno class
t, C = model.solve() # Call the solve method to get the simulation results
fig = model.create_plot(t, C) # Call create_plot with the simulation results
fig.show()
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In[3], line 2
1 model = SCerevisiae_Ligno() # Instantiate the SCerevisiae_Ligno class
----> 2 t, C = model.solve() # Call the solve method to get the simulation results
3 fig = model.create_plot(t, C) # Call create_plot with the simulation results
5 fig.show()
Cell In[2], line 166, in SCerevisiae_Ligno.solve(self)
164 def solve(self):
--> 166 t = np.linspace(0, 30) #generation of the time-points
168 u = 0
169 fc = 1
NameError: name 'np' is not defined
Created on Thu Sep 6 13:34:32 2018
@author: simoca