8. Advanced Topics
This is the Advanced Topics section of our Python documentation. Here, we explore how you can manipulate the simulation during it’s run and we explore the
KineticStruct
which can be used for Michaelis Menten kinetics. Next we discuss the Kinetics function
These tools enable you to add more biological information in the modelling.
Let’s dive in!
Manipulating the simulation
Each simulation for objects of the TimeStepDynamicFBABase
base class can be given an extra function enabling the user to manipulate the simulation on the go.
By manipulating you can think of adding or removing external metabolites from the system. But following some simple rules nearly everything can be changed during simulation.
Here we will give an example where we grow two strains of E. coli on Pyruvate and add L-Glutamate when the concentration of Pyruvate is half it’s initial concentration.
First we build a DynamicJointFBA
model:
import numpy
from cbmpy.CBModel import Model
from dcFBA.Models import CommunityModel
from dcFBA.DynamicModels import DynamicJointFBA
from dcFBA.DefaultModels import read_default_model
# Build the model and set all bounds
model1: Model = read_default_model("e_coli_core")
model1.getReaction("R_GLUt2r").setUpperBound(11)
model1.getReaction("R_PYRt2").setUpperBound(6)
model_2 = model1.clone()
model_2.getReaction("R_PYRt2").setUpperBound(10)
model_2.getReaction("R_GLUt2r").setUpperBound(8)
combined_model = CommunityModel(
[model1, model_2],
["R_BIOMASS_Ecoli_core_w_GAM", "R_BIOMASS_Ecoli_core_w_GAM"],
["ecoli_1", "ecoli_2"],
)
combined_model.getReaction("R_EX_pyr_e").setLowerBound(-numpy.inf)
# Set initial concentrations
dynamic_fba = DynamicJointFBA(
combined_model,
[0.1, 0.1],
{"M_pyr_e": 10, "M_glu__L_e": 0, "M_glc__D_e": 0},
)
Next we will define our deviation function. Remember that the function, no matter if you use it or not, needs to accept the DynamicJointFBA
model, the time steps
of the simulation
and a condition
which we will explain in a minute. Furthermore, you always need to return an int
.
def deviate_func(DFBA: DynamicJointFBA, used_time, condition) -> int:
if DFBA.metabolites["M_pyr_e"][-1] <= 5.0 and condition < 1:
combined_model.getReaction("R_EX_glu__L_e").setLowerBound(-numpy.inf)
DFBA.metabolites["M_glu__L_e"][-1] = 30
return 1
return condition
Here we define the deviate_func
as follows: run the function when the concentration of M_pyr_e runs under 5.0, run it only once. If the concentration is
not below or equal to 5.0 return the condition. The condition
parameter is by default set to 0 at the start of the simulation and is always passed to the deviate function.
This is done such that for example you can rerun the function for consecutive time steps. The value returned by the function will always be added to the global condition variable.
By doing so we can for example add Glucose to the system for N consecutive runs by modifying the conditional statement.
We can now pass the function to the simulation method and plot the results:
import matplotlib.pyplot as plt
dynamic_fba.simulate(0.1, deviate=deviate_func)
T = dynamic_fba.get_time_points()
metabolites = dynamic_fba.get_metabolites()
biomasses = dynamic_fba.get_biomasses()
plt.figure(1)
plt.plot(T, metabolites["M_glu__L_e"], color="blue", label="[glu__L]")
plt.plot(T, metabolites["M_pyr_e"], color="orange", label="[pyr]")
plt.xlabel("Time")
plt.ylabel("Concentration")
plt.legend()
plt.figure(2)
plt.plot(T, biomasses["ecoli_1"], color="blue", label="Biomass model 1")
plt.plot(T, biomasses["ecoli_2"], color="orange", label="Biomass model 2")
plt.xlabel("Time")
plt.ylabel("Concentration")
plt.legend()
plt.show()
The Kinetics Object for Michaelis Menten
All methods permit the passing of a KineticStruct
object during initialization. By doing so, the upper bounds of the reactions defined in this object
are not statically set, but are calculated using the Michaelis-Menten equation.
A KineticStruct
van be defined using the reaction id followed by the limiting substrate and the Km and Vmax
from cbmpy.CBModel import Model
from dcFBA.DynamicModels import DynamicSingleFBA
from dcFBA.DefaultModels import read_default_model
from dcFBA.Models import KineticsStruct
import matplotlib.pyplot as plt
model: Model = read_default_model("e_coli_core")
# Set bounds on the glucose import
model.getReaction("R_GLCpts").setUpperBound(10)
initial_biomass = 0.1
initial_concentrations = {"M_glc__D_e": 10}
#Define kinietcStruct
kin = KineticsStruct({"R_GLCpts": ("M_glc__D_e", 5, 10)})
ds = DynamicSingleFBA(
model,
"R_BIOMASS_Ecoli_core_w_GAM",
initial_biomass,
initial_concentrations,
kinetics=kin,
)
ds.simulate(0.15)
T = ds.get_time_points()
metabolites = ds.get_metabolites()
biomass = ds.get_biomass()
ax = plt.subplot(111)
ax.plot(T, biomass)
ax2 = plt.twinx(ax)
ax2.plot(T, metabolites["M_glc__D_e"], color="r")
ax.set_ylabel("Biomass", color="b")
ax2.set_ylabel("Glucose", color="r")
plt.show()
The Kinetics function
Under Construction
This section is currently under construction. Check back later for updates.
Reducing the memory usage of EndPointFBA
The bottleneck for EndPointFBA
is the amount of RAM required to store the joint stoichiometric matrix. We introduce two methods here to lower the total amount of RAM needed.
1. Scan unused reactions and species
By initializing the EndPointFBA
model, we define the medium. Based on the specified medium and its associated metabolites, certain reactions might never be activated, and some external metabolites may never be produced. To identify these inactive reactions, invoke the dcFBA.Helpers.ScanUnusedReactions
method.
This function runs consecutive FBA’s , each time setting the a reaction as an objective. Reactions that remain inactive get omitted from the initial CommunityModel
.
Additionally, species and reagents that are never involved are also pruned. This ensures that these reactions aren’t duplicated at each time point, effectively minimizing the total number of reactions and species.
2. Sparse Matrix
When invoking the EndPointFBA.simulate()
method, the stoichiometric matrix is constructed.
Setting the sparse option to True during this call will result in the creation of a scipy.sparse.csr_matrix
.
This structure is employed to depict a sparse matrix in the compressed sparse row (CSR) format, thereby conserving a significant amount of system RAM.