Model

The modeling team worked on developing full models of the circuit independently with the aim of combining the “best” parts of each model together and compile one overall model. Because SimBiology lacks the capability to model repression interactions, different team members implemented different workarounds to simulate lac repressor headpiece binding and the inhibited downstream effects. We also placed our emphasis on different portions of the model pathways; for example, one member delved into the influence of the presence of G6P and FBP on CcpA and HPrSerP complex binding affinity, while others focused on GlcT/GlcT-P kinetics or HPrK/P behavior.

During this preliminary phase of model development, we also consulted two modeling experts—Dr. Anders Nelson and Austin Rivera, both of whom were members of the 2016 Virginia iGEM team. With their feedback and suggestions, we developed new iterations of our individual models; for example, we incorporated scaled degradation following the triggering of a repression event, and we also removed some extraneous intermediates to narrow the scope of the models.

In coming together and forming one model, we considered the details of rate laws and decided upon consuming versus non-consuming interactions between different molecules (see Figure 6). We also discussed the assumptions and simplifications that formed the basis of our models, deciding upon which were most logical for the overall model (such as the grouping of intermediates in the glycolysis/metabolism pathway).

Figure 6. A preliminary diagram of the overall model following our joint brainstorming session. We outlined every single molecule we agreed to include, as well as the direction(s) that each reaction could proceed in. Solid blue circles indicate forward-only reactions, while open blue circles indicate reversible reactions. Dotted orange lines indicate situations in which reactions are dependent on the amount of the preceding molecule but do not consume them, whereas solid blue lines indicate reactions in which the reactants are consumed. In instances of the former, the preceding molecule is denoted in parentheses.

Ultimately, we decided that the effects of G6P and FBP on CcpA and HPrSerP binding affinity were significant enough to warrant incorporation into our overall model.

In defining the scope of our model, we established numerous assumptions and simplifications. It is important to clearly establish boundary conditions under which our system of ODEs can be solved so that we know under exactly what circumstances our modeling aims and questions are answered.

Assumptions

Firstly, we assumed that rate constants do not change with time. Additionally, if a promoter sequence is bound by an inhibitory transcription factor, we assumed that it experiences zero transcription. Conversely, if it is unbound, it is constantly transcribing. For the expression of T7 RNA polymerase, lipase, and procolipase, we assumed constant transcription and translation rates. Initial concentrations of most molecules within the full BLISS model were taken to be at 0 molecules.

Two types of reactions were incorporated into our model: Michaelis-Menten kinetics and mass action kinetics. Inherent in the former are the following assumptions:

1. There is no product present at the start of the reaction; thus, there is no reverse reaction of enzyme and product into the enzyme-substrate complex.
2. The rate of enzyme-substrate complex formation is equal to the rate of its breakdown; thus, the concentration of the complex is constant.
3. The concentration of the substrate is significantly larger than that of the enzyme; thus, the amount of substrate can be taken to be constant throughout the course of the reaction.
4. Only the initial velocity of the reaction is measured.
5. The total amount of enzyme remains constant; the enzyme is present either as free enzyme or bound in the enzyme-substrate complex.

In the latter, reaction kinetics are governed by the law of mass action. The rate of reactions that fall under this category are assumed to be proportional to the concentrations of the reactants. The probability of a reaction, then, is only dependent on the probability of the reactants colliding, and all other interactions are largely ignored.

In the literature, few concrete values for rate constants were found, and most measured values were determined in vitro. We had to assume that these values would translate over to an in vivo system. We also had to make educated guesses for certain rate constants based on literature for similar species (such as E. coli) in cases where numbers for B. subtilis had not been previously elucidated.

Simplifications

Due to the complexity of BLISS and the native pathways of B. subtilis it takes advantage of, particularly in the CCR pathway, we established simplifications to our model where we deemed it acceptable. For example, the metabolic pathway of glycolysis was largely ignored in our model—only the metabolites of G6P, FBP, and ATP were key molecules that we were concerned with. As such, most of the glycolysis intermediates have been simplified in our “FBP to ATP” umbrella reaction.

Research and Design

The modeling team performed a thorough literature review on the carbon catabolite repression (CCR) pathway, focusing on key metabolites. Based on this analysis, the team identified potential challenges and took the following actions to address them:

Additionally, the modeling team encountered several challenges during the implementation of the model in MATLAB SimBiology and applied the following workarounds to address them:

With these issues addressed and the workarounds in place, the modeling team developed three predictive models based on low, mid-range, and high-end rate constant values, all built on the same core framework for consistency.

Metabolites were assigned variables arbitrarily, except for the time, kf, and kr variables.

Ordinary differential equations (ODEs) were generated using MATLAB and then transcribed for further analysis.

Sources for low ranges can be found in the table in the references section.
*Assumed 1 mg of protein was present.
**The volume of the system was taken to be 4.6 fL.
***Multiply reaction rate equation by 10 for plasmid copy number.

Sources for low ranges can be found in the table in the references section.
*Assumed 1 mg of protein was present.
***Multiply reaction rate equation by 20 for plasmid copy number.

Figure 7. This diagram shows the simplified carbon catabolite repression pathway in B. subtilis. Metabolites that interact with one another are highlighted in the same color. Red circles indicate non-reversible mass action reactions, blue circles represent reversible mass action reactions, and green circles illustrate Michaelis-Menten reactions.

Figure 8. All variable metabolites for the low-end model are graphed above from t=0 hours to t=5 hours with one dose of glucose mimicking the amount of glucose that reaches the duodenum of the small intestine in a standard meal (1.664M).

Figure 9. After five hours of consuming an average meal, the BLISS low-end predictive model estimates 139 mM lipase and 30 mM colipase in the gastrointestinal tract. After two hours, the model shows lower levels of 5.35 mM lipase and 1.28 mM colipase. These changes illustrate the rapid enzymatic response to dietary glucose.

Figure 10. All variable metabolites for the high-end model are graphed above from t=0 hours to t=5 hours with one dose of glucose mimicking the amount of glucose that reaches the duodenum of the small intestine in a standard meal (1.664M).

Figure 11. After five hours of consuming an average meal, the BLISS low-end predictive model estimates 0.34 mM lipase and 0.14 mM colipase in the gastrointestinal tract. After two hours, the model shows lower levels of 0.03 mM lipase and 0.01 mM colipase.

Current dietician specialists recommend 70,000 to 80,000 U of both lipase and colipase per dose, equivalent to 293,000 to 338,000 U/L per cup of probiotic yogurt. This highlights the critical need for substantial enzymatic activity to effectively break down triglycerides into usable components like fatty acids and glycerols. Our low-end model indicates that B. subtilis can produce 0.564 U/L of lipase and 0.547 U/L of colipase per bacterial cell after five hours, and 0.530 U/L of lipase and 0.452 U/L of colipase after two hours. These values are derived from assumed parameters of Vmax = 940 nmol/(s•L) and Km = 0.5 mM. (sources 4-7) To match the required enzymatic activity, at least 654,586 CFUs of Bacillus subtilis would need to be present in each serving of our therapeutic formulation. This is essential for ensuring that the probiotics can deliver adequate enzyme levels during the digestive process, which typically takes 2 to 5 hours for a meal to pass through the duodenum. In contrast, most probiotics typically contain CFU counts that exceed this number by several orders of magnitude, often reaching billions or trillions of CFUs per dose. The lower CFU requirement in our model may be attributed to a targeted enzyme production strategy. By focusing on enhancing the metabolic pathways specifically responsible for lipase and colipase synthesis in Bacillus subtilis, we can achieve greater enzymatic output without the need for excessively high cell counts. This approach allows us to optimize enzyme production, ensuring that the necessary enzymatic activity is reached with fewer CFUs. Targeting specific metabolic routes can also improve the overall efficiency of enzyme synthesis, making it possible to produce adequate amounts of lipase and colipase with a streamlined number of bacterial cells.

Our high-end model indicates that B. subtilis can produce 0.223 U/L of lipase and 0.121 U/L of colipase per bacterial cell after five hours, and 0.031 U/L of lipase and 0.011 U/L of colipase after two hours. These values are derived from assumed parameters of Vmax = 940 nmol/(s•L) and Km = 0.5 mM. (sources 4-7). To match the required enzymatic activity, at least 26,897,544 CFUs of Bacillus subtilis would need to be present in each serving of our therapeutic formulation. As also expected with the low-end model, the lower CFU requirement in our model may be attributed to a targeted enzyme production strategy.