In [1]:
from imp import reload
import re
import numpy as np
from scipy.integrate import ode
import NetworkComponents
Create a new network object and load the informations from the model Chassagnole2002. In the original model, the concentration of the coffactors depends explicitly on the time. To obtain steady state it is necessary to get rid of this explicit dependence, instead the concentrations of those coffactors are defined as constants. Also the Phosphotransferase system reactions has an unbalanced stoichiometry different from its actual stoichiometry. In the studied model, the stoichiometry is rectified but to maintain the rate to the author's choice, the $r_{max}^{PTS}$ is scaled by a factor 65. The file Chassagnole2002_info.csv contains the information about the number of carbons constituing each metabolite and it tells us wich metabolite exchange labelled carbons. A metabolite that do not echange labelled carbons behaves as sink, it is an exit for the system.
In [2]:
chassagnole = NetworkComponents.Network("chassagnole2002")
chassagnole.readSBML("./published_models/Chassagnole2002.xml")
chassagnole.readInformations("./published_models/Chassagnole2002_info.csv")
A Network object containts other objects stored in arrays:
To derive the tracer dynamics, one need to know the values of the forard and the backward values of the reactions. The function separateForwardBackwardFluxes perform this separation of a rate law from the original model into two new rate laws; one accounts for the forward rate and the second accounts for the backward rate.
The function updateNetwork compiles the network to assign an index and a formula to every reactions and species. After this step it is possible to create the derivative function for the concentration vector.
In [3]:
chassagnole.separateForwardBackwardFluxes()
chassagnole.updateNetwork()
chassagnole.generateDerivatives()
chassagnole.generateRates()
chassagnole.testCarbonBalance()
The following calls are required before generating the jacobians for the tracer and concentration perturbation.
In [4]:
Jtracer = chassagnole.generateTracerJacobian()
To find the jacobian that accounts for the tracers dynamics the algorithm first searches for the steady state of the model. At steady state the probability for a labelled molecule $A^t$ to be transformed through a reaction $v^+$ is proportional to the fraction of $A$ that is labelled. The tracer reaction releases labelled carbons that are shared between the substrate of the reaction proportionally to their stoichiometry and to the number of carbons they contain.
In [6]:
Jperturbation = chassagnole.generatePerturbationJacobian()
tauc,Tc = chassagnole.computeCharacteristicTimes("perturbation",method="integration")
taut,Tt = chassagnole.computeCharacteristicTimes("tracer",method="inverseJacobian")
print("tau_c = %f s"%(tauc))
print("tau_t = %f s"%(taut))
print("T_c = %f s"%(Tc))
print("T_t = %f s"%(Tt))
In [7]:
teusink = NetworkComponents.Network("Teusink2000")
teusink.readSBML("./published_models/Teusink2000.xml")
teusink.readInformations("./published_models/Teusink2000_info.csv")
teusink.separateForwardBackwardFluxes()
teusink.updateNetwork()
teusink.generateDerivatives()
teusink.generateRates()
teusink.testCarbonBalance()
Jtracer = teusink.generateTracerJacobian()
Jperturbation = teusink.generatePerturbationJacobian()
tauc,Tc = teusink.computeCharacteristicTimes("perturbation",method="integration")
taut,Tt = teusink.computeCharacteristicTimes("tracer",method="integration")
print("tau_c = %f s"%(tauc*60))
print("tau_t = %f s"%(taut*60))
print("T_c = %f s"%(Tc*60))
print("T_t = %f s"%(Tt*60))
In [8]:
mosca = NetworkComponents.Network("Mosca2012")
mosca.readSBML("./published_models/Mosca2012.xml")
mosca.readInformations("./published_models/Mosca2012_info.csv")
mosca.separateForwardBackwardFluxes()
mosca.updateNetwork()
mosca.generateDerivatives()
mosca.generateRates()
mosca.testCarbonBalance()
Jtracer = mosca.generateTracerJacobian()
Jperturbation = mosca.generatePerturbationJacobian()
tauc,Tc = mosca.computeCharacteristicTimes("perturbation",method="integration")
taut,Tt = mosca.computeCharacteristicTimes("tracer",method="inverseJacobian")
print("tau_c = %f s"%(tauc*60))
print("tau_t = %f s"%(taut*60))
print("T_c = %f s"%(Tc*60))
print("T_t = %f s"%(Tt*60))
In [9]:
curto = NetworkComponents.Network("Curto1998")
curto.readSBML("./published_models/Curto1998.xml")
curto.readInformations("./published_models/Curto1998_info.csv")
curto.separateForwardBackwardFluxes()
curto.updateNetwork()
curto.generateDerivatives()
curto.generateRates()
curto.testCarbonBalance()
Jtracer = curto.generateTracerJacobian()
Jperturbation = curto.generatePerturbationJacobian()
tauc,Tc = curto.computeCharacteristicTimes("perturbation",method="inverseJacobian")
taut,Tt = curto.computeCharacteristicTimes("tracer",method="inverseJacobian")
print("tau_c = %f s"%(tauc*60))
print("tau_t = %f s"%(taut*60))
print("T_c = %f s"%(Tc*60))
print("T_t = %f s"%(Tt*60))
In [0]: