"
Anatomy of an LBS model
To understand a full LBS model in action, lets walk through the code we wrote to model beta-carotene production in VitaYeast. A model for l-ascorbate (vitamin C) production has been produced in a similar fashion.
LBS source code:
In addition, we translated the LBS models into Matlab's SimBiology Toolbox to conduct more advanced analysis. Fully-interoperable SBML files generated by the Toolbox are provided here.
SBML:
File:Beta-carotene pathway.zip | File:L-ascorbate pathway.zip
Compartments and species
First we declare two compartments under the root "world" compartment. First is the cytosol, inside of which goes the nucleus.
We specify ribosomes and RNA polymerases globally because they are shared across all gene expression processes. We then declare the crtE gene which codes for the enzyme GGPP synthase, the crtYB gene which codes for the bifunctional phytoene synthase / lycopene cyclase, and the crtI gene which codes for carotene desaturase. Note that crtB and crtY are usually separate enzymes: phytoene synthase and lycopene cyclase respectively. crtYB is essentially crtB with an active domain from crtY, allowing it to catalyze both these two non-sequential reactions, in addition to a third reaction which converts neurosporene into dihydro-beta-carotene, resulting in a loss of carbon from this pathway.
Farnesyl diphosphate is the first metabolite that we consider in our synthetic pathway because it is the last metabolite found normally in yeast. We list our other metabolites in their order in the pathway until our final product beta-carotene. We conclude with a list of enzymatic kinetic parameters (Km and Kcat per Michaelis-Menten) and a list of expression kinetic parameters. These lengthy rate declarations are not shown below, but can be seen in the full code.
// Compartments
comp cytosol = new comp;
comp nucleus = new comp inside cytosol;
// Transcriptional machinery
spec Ribo = new{mrna:binding};
spec RNAP = new{dna:binding,mrna:binding};
// Gene-protein pairs
spec crtE = new{bind:binding}; spec GGPP_synth_exo = new{};
spec BTS1 = new{bind:binding}; spec GGPP_synth_endo = new{};
spec crtYB = new{bind:binding}; spec carotene_desat = new{};
spec crtI = new{bind:binding}; spec phytoene_synth_lycopene_cyc = new{};
// Metabolites
spec farnesyl_PP = new{};
spec GGPP = new{};
spec phytoene = new{};
spec neurosporene = new{};
spec lycopene = new{};
spec beta_carotene = new{};
spec hh_beta_carotene = new{};
Expression
The first three modules comprise our standard expression model.
module transcribe(spec DNA:{bind}, mRNA:{down}) {...};
module translate (spec mRNA:{up}, Prot) {...};
module express (comp nuc; spec DNA:{bind}, Prot) {...};
Simple reaction module
We create a module to handle simple reactions in which an enzyme converts a single substrate to a single product. All metabolic reactions that we are interested in follow the Michaelis-Menten pattern outlined below. Note that the rate here is not a rate constant as it would be in mass action kinetics, but is instead the Michaelis-Menten rate expression.
module mmRXN(spec sub, enz, prod; rate km, kcat) {
rate rxn_rate = kcat*enz/(sub+km);
sub + enz ->[rxn_rate] enz + prod
};
Expression and reaction module invocations
This is the heart of the code where we tell LBS to actually simulate something. We invoke the expression module on our two gene-enzyme pairs and the reaction modules for the metabolites that our enzymes process. We then explicitly degrade our enzymes at a uniform rate. The vertical pipes are parallel operators which tell LBS to evaluate these reactions in parallel. In the actual code, all the code below is encapsulated in cytosol[...].
It is important to note that our model greatly oversimplifies the complexity of the reactions. Each conversion between metabolites involves between 2 and 4 separate reactions. However, these reactions tend to be very fast. While the enzymes in our pathway are considered to catalyze all steps of these reactions, the overall rate is controlled by the enzyme's ability to catalyze a single rate-limiting reaction. In fact, the other minor reaction rates are often not even measurable. You will note one exception below in the case of carotene desaturase. Carotene desaturase converts phytoene to neurosporene at a rate of about 120 per second. It then converts neurosporene to lycopene at a rate of about 10 per second. Even though we model both reactions explicitly, the observation that one reaction is an order of magnitude slower than the other supports the single rate-limiting step simplification.
express(nucleus, crtE:{bind}, GGPP_synth_exo) |
express(nucleus, BTS1:{bind}, GGPP_synth_endo) |
express(nucleus, crtYB:{bind}, phytoene_synth_lycopene_cyc) |
express(nucleus, crtI:{bind}, carotene_desat) |
mmRXN(farnesyl_PP, GGPP_synth_endo, GGPP, GGPP_synth_endo_Km,
GGPP_synth_endo_Kcat) |
mmRXN(farnesyl_PP, GGPP_synth_exo, GGPP, GGPP_synth_exo_Km,
GGPP_synth_exo_Kcat) |
mmRXN(GGPP, phytoene_synth_lycopene_cyc, phytoene, phytoene_synth_Km,
phytoene_synth_Kcat) |
mmRXN(phytoene, carotene_desat, neurosporene, carotene_desat_neurosporene_Km,
carotene_desat_neurosporene_Kcat) |
mmRXN(neurosporene, carotene_desat, lycopene, carotene_desat_lycopene_Km,
carotene_desat_lycopene_Kcat) |
mmRXN(neurosporene, carotene_desat, hh_beta_carotene,
phytoene_synth_hh_beta_carotene_Km, phytoene_synth_hh_beta_carotene_Kcat) |
mmRXN(lycopene, phytoene_synth_lycopene_cyc, beta_carotene,
phytoene_synth_beta_carotene_Km, phytoene_synth_beta_carotene_Kcat) |
GGPP_synth_exo ->{prot_deg} |
GGPP_synth_endo ->{prot_deg} |
carotene_desat ->{prot_deg} |
phytoene_synth_lycopene_cyc ->{prot_deg} |
Initial conditions
We establish initial conditions for our system, namely 5 copies of each gene and a very limited number of ribosomes and RNA polymerases. We take the numbers here to mean "the number of ribosomes or RNA polymerases made available to the expression process of our genes of interest". There are of course many times more ribosomes and polymerases in the cell, but they are working on making other proteins. Only a small portion is ever free to deal with our synthetic construct.
In addition, the first line below is a black-box input: a certain number of our initial substrate is dumped into our system each turn by the cell. We are not interested in modeling the process of how it comes about, we only care about the result and what our system does to this exogenous input.
->{1} farnesyl_PP |
init Ribo 10 |
nucleus[
init crtE 5 |
init crtYB 5 |
init crtI 5 |
init RNAP 1
]
