Mechanistic model of S. cerevisiae in continous stirred tank reactor under aerobic conditions

Mechanistic model of S. cerevisiae in continous stirred tank reactor under aerobic conditions#

This model describes an aerobic S. cerevisiae fermentation with glucose as carbon source in a continous stirred tank reactor. It builds a stochiometric matrix to follow all conversions. The model accounts for the changes of glucose, ethanol, dissolved oxygen, biomass and volume.

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 includes several functions, each with a specific function.

  1. init This function initialises the model class by defining all the relavant parameters and initial conditions.

  2. rxn This function includes all the model equations. The model uses matrix notation to define the ODEs. In this part the settings for feed are specified as well, it includes terms for dilution, addtion and washing out effects.

  3. 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.

  4. create_plot This function stores the results of the simulation on a dataframe which is used for plotting all relevant variables.

class SCerevisiae_Cstr:
   #initialise the model
    def __init__(self):

      # define value of model parameters

        self.Yox_XG = 0.8
        self.Yred_XG = 0.05
        self.Yox_XE = 0.72
        self.Y_OG = 1.067
        self.Y_EG = 0.5
        self.Y_OE = 1.5
        self.q_g = 4.68
        self.q_o = 0.37
        self.q_e = 0.86
        self.t_lag = 4.66
        self.Kg = 0.17
        self.Ko = 0.0001
        self.Ke = 0.56
        self.Ki = 0.31
        self.O_sat = 0.00755
        self.kla = 1004

        #define initial conditions
        self.G0 = 18
        self.E0 = 0.34
        self.O0 = 0.00755
        self.X0 = 0.1
        self.T0 = 30
        self.V0 = 2

        #define feed parameters
        self.t_end = 100
        self.Cin = 100
        self.Fconst = 0.05
        self.t_feed_start = 15
        self.t_start=2
        self.steps = (self.t_end - self.t_start)*24


    #define the stoichiometric matrix
    def rxn(self, C,t):

        #matrix
        #number of components
        n = 4
        m = 4
        #initialize the stoichiometric matrix, s
        s = np.zeros((m,n))
        s[0,0] = -1
        s[0,1] = 0
        s[0,2] = -self.Y_OG
        s[0,3] = self.Yox_XG

        s[1,0] = -1
        s[1,1] = self.Y_EG
        s[1,2] = 0
        s[1,3] = self.Yred_XG

        s[2,0] = 0
        s[2,1] = -1
        s[2,2] = -self.Y_OE
        s[2,3] = self.Yox_XE

        s[3,0] = 0
        s[3,1] = 0
        s[3,2] = 1
        s[3,3] = 0
        #initialize the rate vector
        rho = np.zeros((4,1))
        ##initialize the overall conversion vector
        r=np.zeros((4,1))
        #Volume balance
        if (t>=self.t_feed_start):
            Fin = self.Fconst
            Fout = self.Fconst
        else:
            Fin = 0
            Fout = 0

        F = Fin - Fout

        rho[0,0] = ((1/self.Y_OG)*min(self.q_o*(C[2]/(C[2]+self.Ko)),self.Y_OG*(self.q_g*(C[0]/(C[0]+self.Kg)))))*C[3]
        rho[1,0] = ((1-math.exp(-t/self.t_lag))*((self.q_g*(C[0]/(C[0]+self.Kg)))-(1/self.Y_OG)*min(self.q_o*(C[2]/(C[2]+self.Ko)),self.Y_OG*(self.q_g*(C[0]/(C[0]+self.Kg))))))*C[3]
        rho[2,0] = ((1/self.Y_OE)*min(self.q_o*(C[2]/(C[2]+self.Ko))-(1/self.Y_OG)*min(self.q_o*(C[2]/(C[2]+self.Ko)),self.Y_OG*(self.q_g*(C[0]/(C[0]+self.Kg)))),self.Y_OE*(self.q_e*(C[1]/(C[1]+self.Ke))*(self.Ki/(C[0]+self.Ki)))))*C[3]
        rho[3,0] = self.kla*(self.O_sat - C[2])

        #Developing the matrix, the overall conversion rate is stoichiometric *rates
        r[0,0] = (s[0,0]*rho[0,0])+(s[1,0]*rho[1,0])+(s[2,0]*rho[2,0])+(s[3,0]*rho[3,0])
        r[1,0] = (s[0,1]*rho[0,0])+(s[1,1]*rho[1,0])+(s[2,1]*rho[2,0])+(s[3,1]*rho[3,0])
        r[2,0] = (s[0,2]*rho[0,0])+(s[1,2]*rho[1,0])+(s[2,2]*rho[2,0])+(s[3,2]*rho[3,0])
        r[3,0] = (s[0,3]*rho[0,0])+(s[1,3]*rho[1,0])+(s[2,3]*rho[2,0])+(s[3,3]*rho[3,0])


        #Solving the mass balances terms for dilution, addtion and washing out added
        dGdt = r[0,0] -F/C[4]*C[0] + Fin/C[4]*self.Cin - Fout/C[4]*C[0]
        dEdt = r[1,0] -F/C[4]*C[1] - Fout/C[4]*C[1]
        dOdt = r[2,0]
        dXdt = r[3,0] -F/C[4]*C[3] - Fout/C[4]*C[3]
        dVdt = F
        dTdt = 0

        return [dGdt,dEdt,dOdt,dXdt, dVdt, dTdt]

    #solve the ODES
    def solve(self):
        #time
      t = np.linspace(self.t_start, self.t_end, self.steps)  #generation of the time-points
      u=0
      C0 = [self.G0, self.E0, self.O0, self.X0, self.V0, self.T0]  #initial conditions vector
      C = odeint(self.rxn, C0, t, rtol = 1e-7, mxstep= 500000)  #solve ODEs

      return t, C


   #generate the plot of model variables


    def create_plot(self, t, C):
        fig = make_subplots(rows=1, cols=2) #make figure with 2 subplots
        #assign simulation results to variable for plotting
        G = C[:, 0]
        E = C[:, 1]
        O = C[:, 2]
        B = C[:, 3]
        V = C[:, 4]

        #collect all variables to plot in 1st subplot in a dataframe
        df = pd.DataFrame({'t': t, 'G': G, 'B': B, 'E': E, 'O':O, 'V':V})

        #add the different traces to 1st subplot
        fig.add_trace(go.Scatter(x=df['t'], y=df['G'], name='Glucose'), row=1, col=1)
        fig.add_trace(go.Scatter(x=df['t'], y=df['B'], name='Biomass'), row=1, col=1)
        fig.add_trace(go.Scatter(x=df['t'], y=df['E'], name='Ethanol'), row=1, col=1)
        fig.add_trace(go.Scatter(x=df['t'], y=df['O'], name='Oxygen'), row=1, col=1)

        #add the title and axes labels
        fig.update_layout(title=('Simulation of the model for the Scerevisiae in continuous'),
                          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
        fig.append_trace(go.Scatter(x=df2['t'], y=df2['V'], name='Volume'), row=1, col=2)

        return fig

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_Cstr() # Instantiate the SCerevisiae_Cstr 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_Cstr() # Instantiate the SCerevisiae_Cstr 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 108, in SCerevisiae_Cstr.solve(self)
    106 def solve(self):
    107     #time
--> 108   t = np.linspace(self.t_start, self.t_end, self.steps)  #generation of the time-points
    109   u=0
    110   C0 = [self.G0, self.E0, self.O0, self.X0, self.V0, self.T0]  #initial conditions vector

NameError: name 'np' is not defined

Created on Thu Sep 6 18:33:18 2018

By Simoneta Cano de las Heras and bjogut