This tutorial demonstrates how to build and simulate a DFBA model corresponding to the growth of a single strain of Escherichia coli (based on the iJR904 genome-scale model) under aerobic and anaerobic conditions with glucose and xylose as limiting carbon substrates.
This notebook was adapted with minor changes from the dfba documentation.
from os.path import dirname, join, pardir
from cobra.io import read_sbml_model
from dfba import DfbaModel, ExchangeFlux, KineticVariable
A model will be build as the following ODE system: \begin{align} \frac{dX}{dt} &= v_{BM}X \\ \frac{dC_j}{dt} &= v_jX \end{align}
where $X$ is the Biomass (gDW); $v$ is the flux of a exchange reaction $j$ (mmol/gDW/h) or the exchange biomass reaction $BM$ (g/gDW/h, to measure growth).
First off, specify the path for loading file containing genome-scale metabolic model as cobra.Model object and set GLPK as LP solver of choice. After that, instantiate the object of class DfbaModel
with the cobrapy model.
path_to_model = join(pardir, "sbml-models", "iJR904.xml.gz")
fba_model = read_sbml_model(path_to_model)
fba_model.solver = "glpk"
dfba_model = DfbaModel(fba_model)
Note: A wall of warnings - which can be safely ignored - may be produced when loading this model with cobra.io.read_sbml_model
.
For this example, the metabolites $j$ in the medium will be Glucose, Xylose, Oxygen and Ethanol. Later, we will establish the kinetic rules for the uptakes of these metabolites. In addition, we need to set a variable to keep track of the biomass.
Instantiate kinetic variables to appear in model. The last command adds the kinetic variables to the model.
X = KineticVariable("Biomass")
Gluc = KineticVariable("Glucose")
Xyl = KineticVariable("Xylose")
Oxy = KineticVariable("Oxygen")
Eth = KineticVariable("Ethanol")
dfba_model.add_kinetic_variables([X, Gluc, Xyl, Oxy, Eth])
Instantiate exchange fluxes to appear in model, with ids corresponding to exchange reactions of the cobrapy model. The last command adds the exchange fluxes to the model.
mu = ExchangeFlux("BiomassEcoli")
v_G = ExchangeFlux("EX_glc(e)")
v_Z = ExchangeFlux("EX_xyl_D(e)")
v_O = ExchangeFlux("EX_o2(e)")
v_E = ExchangeFlux("EX_etoh(e)")
dfba_model.add_exchange_fluxes([mu, v_G, v_Z, v_O, v_E])
Provide symbolic expression for calculating the time derivative of each kinetic variable currently in the model. See how these correspond to our ODE system. \begin{align} \frac{dX}{dt} &= v_{BM}X \\ \frac{dC_j}{dt} &= v_jX \end{align}
dfba_model.add_rhs_expression("Biomass", mu * X)
dfba_model.add_rhs_expression("Glucose", v_G * X / 1000.0)
dfba_model.add_rhs_expression("Xylose", v_Z * X / 1000.0)
dfba_model.add_rhs_expression("Oxygen", v_O * X / 1000.0)
dfba_model.add_rhs_expression("Ethanol", v_E * X / 1000.0)
Add symbolic expressions for calculating lower/upper bounds of selected exchange fluxes currently in the model. Here, the convention is that both lower and upper bound expressions have positive signs, whereas lower bounds values are typically negative in cobrapy.
Here, the lower bounds will follow a Michaelis-Menten kinetics:
$$ lb = - \frac{C_jV_{max}}{ C_j + K_m} $$(The Kinetic parameters $V_{max}$ and $K_m$ can be found in the literature for the given exchange reactions)
dfba_model.add_exchange_flux_lb("EX_o2(e)", 15.0 * (Oxy / (0.024 + Oxy)), Oxy)
But we are not limited to this formula! Let's say Glucose is inhibited by Ethanol.
dfba_model.add_exchange_flux_lb(
"EX_glc(e)", 10.5 * (Gluc / (0.0027 + Gluc)) * (1 / (1 + Eth / 20.0)), Gluc
)
And Xylose in inhibited by both Ethanol and Glucose
dfba_model.add_exchange_flux_lb(
"EX_xyl_D(e)",
6.0
* (Xyl / (0.0165 + Xyl))
* (1 / (1 + Eth / 20.0))
* (1 / (1 + Gluc / 0.005)),
Xyl,
)
Initial values for each kinetic variable are provided in dictionary form.
The model is simulated using the simulate
method. This simulation covers the interval $[0.0, 25.0]$ hours, with results stored every $0.1$ hours. Results (trajectories of kinetic variables) will be returned as pandas.DataFrame. Optionally, the user can also provide a list of reaction ids whose flux trajectories will also be returned as a separate pandas.DataFrame, in this case three exchange fluxes in the model.
dfba_model.add_initial_conditions(
{
"Biomass": 0.03,
"Glucose": 4.0,
"Xylose": 5.0,
"Oxygen": 0.5,
"Ethanol": 0.0,
}
)
concentrations, trajectories = dfba_model.simulate(0.0, 16.0, 0.1, ["EX_glc(e)", "EX_xyl_D(e)", "EX_etoh(e)"])
creating library... adding model 140634011303832 to library... cpp file generated setting user data for model 140634011303832... compiling dynamic library... simulating...
Note: In order plot the results, one of plotly or matplotlib must be used. Check the original example in the documentation to see the matplotlib counterpart.
from dfba.plot.plotly import *
fig = plot_concentrations(concentrations)
fig.show()
fig = plot_trajectories(trajectories)
fig.show()
dfba_model = DfbaModel(fba_model)
# YOUR CODE...