Introduction

With up to an estimated 20% of the general population experiencing Exocrine Pancreatic Insufficiency (EPI) [1] and unable to produce enough digestive enzymes as a result, it is imperative that any barriers to treatment are removed. However, the only current state-of-the-art treatment for EPI involves the intake of porcine-derived pancreatic enzymes, posing a significant obstacle to many patients seeking to manage their condition. Consequently, we proposed an alternative solution involving the secretion of the needed digestive enzymes in vivo via a bacterial host—Bacillus subtilis.

Called the Bacillus Lipase Integrated Secretion System, or BLISS, our device involves the expression and secretion of two recombinant human pancreatic proteins: the human pancreatic lipase (PNLIP) and procolipase (CLPS) enzymes. Lipase serves to break down lipids (such as triglycerides) into fatty acids and glycerols [2]. Procolipase, which is cleaved to become colipase, plays a role as a necessary coenzyme for lipase, particularly in the presence of bile salts [3]. Two other major digestive enzymes, amylase and protease, are naturally produced by B. subtilis, and they serve to break down starches and proteins, respectively.

BLISS is designed such that lipase and procolipase production will only occur in B. subtilis following glucose intake (akin to eating a meal) through modifications involving the terminator/anti-terminator interaction of GlcT (described in Figures 1 and 2) and the carbon catabolite repression (CCR) pathway (see Figure 3)—both native features of the host. We have modeled these mechanisms in what is termed the “full BLISS model.” In addition, we have developed a smaller model that directly reflects our tested combination of parts 1, 2, and 3A to predict quantities of secreted enzymes. This second model, termed the “experimental model,” involves coding sequences for lipase and procolipase, as well as the T7 RNA polymerase that transcribes these two genes.

Figure 1. Glucose-dependent binding of GlcT to RAT sequence in the presence of glucose (adapted from Schmalisch et al., 2003) [4]. When glucose enters the cell, it becomes phosphorylated—thus preventing the phosphorylation of GlcT. In its non-phosphorylated state, GlcT is activated and binds in dimers to the RAT sequence of the mRNA transcript of a DNA segment, leading to the unfolding of the terminator hairpin. This allows for translation of the gene of interest. This figure was created with BioRender.com.

Figure 2. Glucose-dependent binding of GlcT to RAT sequence in the absence of glucose (adapted from Schmalisch et al., 2003) [4]. When glucose is not present, GlcT becomes inactivated through phosphorylation into GlcT-P. In its inactive form, GlcT-P is unable to bind to the RAT sequence of the mRNA transcript, allowing a hairpin to form and terminating any further transcription. This figure was created with BioRender.com.

Figure 3. Carbon catabolite repression pathway in B. subtilis (adapted from Görke and Stülke, 2008) [5]. Following the uptake of glucose through the transmembrane protein EIIGlc, the metabolites of the glycolysis pathway, i.e. glucose-6-phosphate (G6P) and fructose 1,6-bisphosphate (FBP), activate HPrK/P activity. This promotes the phosphorylation of histidine phosphocarrier protein (HPr) on its respective Ser46 residue. The HPrSerP binds to and forms an activated complex with CcpA. These complexes can then bind at catabolite repression elements (cre sites) located in promoter regions to inhibit transcription of a gene of interest (GOI). Inversely, when no glucose is present, transcription of the GOI will be activated; HPr does not become phosphorylated, leading to no formation of the HPrSerP-CcpA complex that would otherwise bind cre. This figure was created with BioRender.com.

ODE-Based Modeling

We modeled the pathway of our device with MATLAB SimBiology to computationally predict quantitative yields of lipase and procolipase corresponding to the level of glucose introduced to the system. Additionally, we verified that our repression mechanism involving both the CCR pathway and the GlcT terminator/anti-terminator mechanics would work as intended when low-to-no glucose was present.

We chose to develop our model in SimBiology for the program’s capabilities in building, simulating, and analyzing ODE-based models—as well as its intuitive interface. Even with little-to-no prior MATLAB experience, we were able to quickly grasp how SimBiology worked. Its powerful tools enable us to solve ODEs, visualize pathway dynamics graphically, and capture the intricacies of multi-step interactions and feedback loops inherent in BLISS. These functionalities are essential for accurately modeling such a complex biological network. Coupled with experimental data, we can adjust and refine our model such that we can realistically simulate the dynamics of BLISS.

This approach helps us simulate and optimize key interactions driving the production of therapeutic enzymes of lipase and procolipase, crucial for clinical applications. The model bridges the gap between the biological complexity of metabolic pathways and practical therapeutic development.

Why ODE Modeling?

Metabolic pathways are governed by complex, time-dependent interactions between metabolites, enzymes, and regulators. ODE models effectively capture these non-linear dynamics, tracking how component concentrations change over time. This allows us to model the intricate feedback loops and regulatory mechanisms that static models can't, giving us a more accurate and dynamic understanding of the reactions that comprise BLISS. By applying ODE models, we can achieve the following:

  1. Quantitative Precision: ODEs allow us to quantify how concentrations of molecules like HPr, CcpA, and their complexes evolve over time. This helps us accurately model how quickly enzymes like lipase and procolipase are produced, providing critical insights into optimizing bacterial strains for therapeutic production.
  2. Predictive Power: The models enable us to predict outcomes under various conditions, such as changes in glucose levels or enzyme kinetics. We can estimate how many colony-forming units, or CFUs, are required to produce therapeutic doses of lipase and procolipase comparable to current treatments, guiding us in designing effective dosages.
  3. Parameter Variability: By exploring a range of kinetic parameters from the literature, we can understand how variations affect enzyme production. This allows us to refine our model and validate it against experimental data, ensuring its robustness for real-world applications.
  4. Time-Course Analysis: ODE models help predict how long it takes for B. subtilis to produce therapeutic enzymes. This time-course analysis is crucial for optimizing production timelines and ensuring consistent, efficient output.

The Bigger Picture

Our ODE modeling aims to provide a comprehensive understanding of the CCR and terminator/anti-terminator pathways in B. subtilis and optimize it for therapeutic enzyme production. By comparing simulation predictions with experimental data, we can validate and refine our model, improving enzyme production consistency and scalability. Ultimately, this approach brings us closer to developing B. subtilis as a reliable platform for producing therapeutic proteins.

By integrating precision and adaptability, ODE modeling ensures our solution is scientifically sound, advancing synthetic biology and therapeutic production.

Full BLISS Model

Conditional Glucose Diagrams

The full construct of BLISS involves both the GlcT switch mechanism and the CCR pathway such that with glucose, lipase and colipase are secreted (see Figure 4). Without glucose, lipase and colipase production completely cease (see Figure 5).

Glucose Present:

Hover over each step to read about it below the diagram!

Read the full description below the diagram!

Full Diagram

Figure 4. Complete dynamic diagram of BLISS when glucose is present. 1, Glucose enters the cell via the EIIGlc transmembrane protein and becomes phosphorylated into glucose-6-phosphate (G6P). The phosphate group is donated by histidine phosphocarrier protein (HPr), phosphorylated on its respective His15 residue (HPrHisP). 2, Through a simplification of two enzyme-catalyzed reactions of the glycolysis pathway, G6P becomes fructose 1,6-bisphosphate (FBP). The downstream processes of metabolism are not pictured. 3, HPr becomes phosphorylated at Ser46 (HPrSerP) by HPrK/P through the donation of a phosphate group from ATP (when ATP levels are high). 4, HPrSerP binds to carbon catabolite protein A (CcpA). 5–6, Two intermediate metabolites of glycolysis, G6P and FBP, increase the binding affinity of HPrSerP to CcpA. 7, The HPrSerP-CcpA complexes bind at catabolite repression elements (cre sites) located in the promoter regions of DNA to inhibit transcription. 8, In BLISS, the binding of cre sites prevents the lacI gene from being transcribed. 9, GlcT is activated and binds in dimers at the ptsG RAT sequence of the transcript of the T7RNAP gene, unfolding the terminator hairpin and allowing for translation of the mRNA. 10, The T7RNAP mRNA is translated into T7 RNA polymerase. 11, T7 RNA polymerase binds at the promoter regions of CLPS and PNLIP and transcribes the two genes. 12, The two transcripts are translated into colipase and lipase. 13, The signal recognition particles (SRP) attached to colipase and lipase target the two proteins for extracellular secretion by B. subtilis. These SRPs will be cleaved from the final exported protein. This figure was created with BioRender.com.

Glucose Absent:

Hover over each step to read about it below the diagram!

Read the full description below the diagram!

Full Diagram

Figure 5. Complete dynamic diagram of BLISS when no glucose is present. 14, GlcT becomes phosphorylated into GlcT-P after accepting phosphate groups (donated by HPrHisP) in lieu of glucose. 15, The lack of GlcT leaves the RAT sequence unbound, allowing for the formation of the terminator hairpin in the leading 5’ region of the t7rnap mRNA (if transcription occurs). The hairpin prevents the subsequent translation of the transcript into T7 RNA polymerase. 3a, Due to the lack of glucose, G6P and FBP are not present, and ATP levels drop as inorganic phosphate levels rise. In these conditions, HPrK/P exhibits phosphorylase activity and dephosphorylates HPrSerP into HPr. 4a, CcpA only forms a complex with HPrSerP—not with HPr—which causes the cre sites to remain unbound. 16, The unbound cre sites allow for transcription of the downstream lacI gene, which produces the lac repressor headpiece (or LacI repressor protein). 17, The lac repressor headpiece binds to the lac operator sites. 18–20, The transcription of PNLIP (lipase), CLPS (colipase), and t7rnap will be inhibited due to the bound lac operator in their respective promoter regions. This figure was created with BioRender.com.

Our Process of Modeling BLISS

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.

Assumptions and Simplifications

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.

Predictive Modeling

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:

Problem Workaround
Rate Constant Variability: Many rate constants span a broad range (e.g., 10⁻⁴ to 10⁻³ s⁻¹). Two-Tier Modeling: We developed two models using the low-end and high-end values for the rate constants.
Plasmid Copy Number: The plasmid (pBS1C) has a variable copy number, ranging from 10-20 per cell. Plasmid Scaling: For reactions physically involving the plasmid, we scaled the rates by 10 and 20 for the low and high-end models, respectively.
Incomplete Data: Several necessary rate constants were not available in the literature. Informed Assumptions: Where data was unavailable, educated assumptions were made based on related studies.

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

Problem Workaround
Gene Suppression Implementation: SimBiology lacks a direct way to implement gene suppression based on introduced metabolites. Event-Based Suppression: Events were introduced to model gene suppression, triggered when key metabolites reached a specific threshold (e.g., when lacI = 0.003 mM, T7 RNA Polymerase = 0).
ODE Solver Limitations: Using molecules as the unit for rate constants caused the ODE solver (Stiff/NDF) to throw integration tolerance errors. Unit Adjustment: Rate constants were converted from molecules to millimolar (mM) and molar (M), which resolved the solver issues and improved the integration performance.

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.

The Core Framework

Metabolites and Variables

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

Ordinary Differential Equations

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

Reaction Rates for Low-Range Model

Reaction Rate
Glucose to G6P (Michaelis-Menten) Km = 0.1 mM = 2.77E5 molecules
Vmax = 12.3 µM/(min/mg protein) = 3.407E4 molecules/min
G6P to F6P (Michaelis-Menten) Km = 0.1 mM = 2.77E5 molecules
Vmax = 3 µM/min/mg protein* = 8.30E3 molecules/min
F6P to FBP (Michaelis-Menten) Km = 0.1 mM = 2.77E5 molecules
Vmax = 10 µM/min/mg protein* = 2.77E4 molecules/min
FBP to ADP (Michaelis-Menten) Km = 0.1 mM = 2.77E5 molecules
Vmax = 200 U/mg protein* = 200 µM/min/mg protein* = 5.54E5 molecules/min
G6P/F6P + CcpA to CcpA complex (Mass Action) Kf = 10E4 M-1 s-1 = 2.77E13 molecules-1 s-1
Kr = 10E-2 s-1 = 2.77E7 s-1
HPr + ATP to HPrSerP (Mass Action) Kon = 10E5 M-1 s-1 = 2.77E14 molecules-1 s-1
Koff = 10E-2 s-1 = 2.77E7 s-1
CcpA + HPrSerP to complex (Mass Action) Ka = 2E6 M-1 s-1 = 5.61E15 molecules-1 s-1
Kd = 1.6E-3 s-1 = 4.43E6 s-1
CcpA complex binding to cre site (Mass Action)*** Kon = 10E7 M-1 s-1 = 2.77E16 molecules-1 s-1
Koff = 10E-2 s-1 = 2.77E7 s-1
LacI gene (1083 bp) to mRNA (transcription and degradation rates)*** Ktranscription = 10 nt/s = 108.3 s-1
Kdegradation = 30 min = 5.6E-4 s-1
mRNA to lacI protein (translation rate) Kon = 0.1 s-1
Ktranslation = 30 nt/s = 36.1 s-1
Formation of lac complex (Mass Action) Kon = 10E7 M-1 s-1 = 2.77E16 molecules-1 s-1
Koff = 10E-4 s-1 = 2.77E5 s-1
GlcT-P to GlcT (Mass Action) Kon = 10E4 M-1 s-1 = 2.77E13 molecules-1 s-1
Koff = 10E-4 s-1 = 2.77E5 s-1
GlcT binding to RAT, transcribing t7rnap (Mass Action)*** Kon = 10E4 M-1 s-1 = 2.77E13 molecules-1 s-1
Koff = 10E-4 s-1 = 2.77E5 s-1
Ktranscription = 10 nt/s = 265.2 s-1
t7rnap gene (2652 bp) to T7 RNA polymerase (translation) Kon = 0.1 s-1
Ktranslation = 30 nt/s = 88.4 s-1
T7 RNA polymerase degradation (degradation) 0.077 min-1 = 1.28E-3 s-1
CLPS (578 bp) and PNLIP (1347 bp) transcription & translation Kon = 0.1 s-1 Ktranslation = 30 nt/s = 19.3 s-1 (CLPS), 44.9 s-1 (PNLIP)
Secretion Kcat = 0.01 s-1

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.

Reaction Rates for High-Range Model

Reaction Rate
Glucose to G6P (Michaelis-Menten) Km = 1.0 mM
Vmax = 26.0 µM/min/mg protein
G6P to F6P (Michaelis-Menten) Km = 0.2 mM
Vmax = 5 µM/min/mg protein
F6P to FBP (Michaelis-Menten) Km = 0.15 mM
Vmax = 10 µM/min/mg protein*
FBP to ATP (Michaelis-Menten) Km = 0.5 mM
Vmax = 600 U/mg protein* = 600µM/min/mg protein*
G6P/F6P + CcpA to CcpA complex (Mass Action) Kf = 10E6 M-1 s-1
Kr = 10E-1 s-1
HPr + ATP to HPrSerP (Mass Action) Kon = 10E7 M-1 s-1
Koff = 10E-1 s-1
CcpA + HPrSerP to complex (Mass Action) Ka = 4E6 M-1 s-1
Kd = 2.4E-3 s-1
CcpA complex binding to cre site (Mass Action) Kon = 10E8 M-1 s-1
Koff = 10E-1 s-1
LacI gene (1083 bp) to mRNA (transcription and degradation rates) Ktranscription = 40 nt/s = 108.3 s-1
Kdegradation = 30 min = 5.6E-4 s-1
mRNA to lacI protein (translation rate) Kon = 0.1 s-1
Ktranslation = 30 nt/s = 36.1 s-1
Formation of lac complex (Mass Action) Kon = 10E7 M-1 s-1 = 6.022E30 molecules-1 s-1
Koff = 10E-4 s-1 = 6.022E19 s-1
GlcT-P to GlcT (Mass Action) Kon = 10E4 M-1 s-1 = 2.77E13 molecules-1 s-1
Koff = 10E-4 s-1 = 2.77E5 s-1
GlcT binding to RAT, transcribing t7rnap (Mass Action) Kon = 10E4 M-1 s-1 = 6.022E27 molecules-1 s-1
Koff = 10E-4 s-1 = 6.022E19 s-1
t7rnap gene (2652 bp) to T7 RNA polymerase (translation) Kon = 0.1 s-1
Ktranslation = 30 nt/s = 88.4 s-1
T7 RNA polymerase degradation (degradation) 0.077 min-1 = 1.28E-3 s-1
CLPS (578 bp) and PNLIP (1347 bp) transcription & translation Kon = 0.1 s-1
Ktranslation = 30 nt/s = 19.3 s-1 (CLPS), 44.9 s-1 (PNLIP)
Secretion Kcat = 0.01 s-1

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.

MATLAB SimBiology View

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.

Results


"Low-End" Model: All Metabolites

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).

"Low-End" Model: Key Metabolites

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

"High-End" Model: All Metabolites

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).

"High-End" Model: Key Metabolites

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

Conclusions and Future Applications

"Low-End" Model

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.

"High-End" Model

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.

Extensions of the Model

The current model provides a reasonable range for probiotic dosage but holds potential for more specific and finetuned applications. For instance, modeling glucose production in relation to an average day’s eating patterns, with corresponding glucose doses aligned to meal timings, could estimate probiotic responses to fluctuating inputs at varying time scales. Incorporating individual factors such as age, weight, and height would enhance the model’s applicability across diverse patient profiles, enabling personalized probiotic dosing tailored to physical characteristics.

As more experimental data becomes available, the model's applications could expand significantly. Understanding the proportion of bacteria that traverse the gastrointestinal tract versus those that successfully colonize the gut would clarify dosing frequency on daily and weekly scales. This information could also evaluate the BLISS therapeutic’s competitiveness with existing treatments in terms of dosage efficiency. Additionally, simulating bacterial responses within the duodenum's specific environment could provide insights into metabolic strain on the host organism, further elucidating the therapeutic's functionality under its intended physiological conditions.

Lipase and Colipase Experimental Model

To directly reflect our more specific construct of Part 1, 2, and 3A (which encode for lipase, procolipase, and T7 RNA polymerase, respectively), we also simulated the combination of these parts with a system of ODEs in MATLAB (see Figure 12).

Figure 12. MATLAB SimBiology representation of 1-2-3A construct. In this model, T7RNAP is constitutively expressed (as shown in top-right entity). Both the transcription and translation of T7 RNA polymerase are considered. Subsequently, CLPS and PNLIP are expressed and secreted.

In this model, we are modeling the constitutive expression of T7 RNA polymerase (with the P43 promoter), as well as subsequent lipase/procolipase production and secretion. The rate of T7 RNA polymerase mRNA production is associated with the “number” of copies of the t7rnap gene, and it is translated into the protein product or degraded. Then, the T7 RNA polymerase transcribes CLPS/PNLIP, or it gets degraded. Finally, the respective CLPS and PNLIP products of procolipase and lipase are secreted extracellularly.

Assumptions

We have assumed that the rate of transcription of t7rnap into its mRNA counterpart is at the maximal rate of 80 nt/s. Similarly, we have assumed that the rate of translation is at a maximum of 20 aa/s [6].

We have taken the time for the degradation of mRNA to be 5 minutes [7]. The half-life of T7 RNA polymerase has been elucidated to be 50 minutes [8]; we assumed first order kinetics to determine an appropriate degradation rate. The transcription and secretion rates of CLPS and PNLIP are assumptions based on approximate T7RNAP transcription rates. Secretion is taken to be a similar rate value (due to lack of literature). We will verify and adjust these values as needed as we perform the appropriate experiments.

Every reaction within the experimental model adopts mass action kinetics; thus, rates for each reaction are based solely upon the concentrations of the reactants, and other interactions are negligible.

Simplifications

For this model, we have simplified it to not include promoter region kinetics (ex. time needed for RNA polymerase initiation). The transcription and translation of CLPS and PNLIP have been simplified into one reaction, and the degradation of their respective mRNA transcripts and proteins has been ignored. The mechanisms behind secretion are taken to be negligible in this instance.

Initial Values

Molecule Initial Amount (in molecules)
t7rnap 2
T7 mRNA 0
T7RNAP 0
SP-CLPS 0
SP-PNLIP 0
Procolipase 0
Lipase 0

Rate Constants

Rate Constant Value
kf_9 (t7rnap mRNA transcription) 0.0302 (1/sec)
kf_8 (T7RNAP degradation) 2.31E-4 (1/sec)
kf_7 (t7rnap mRNA degradation) 0.0033 (1/sec)
kf_6 (t7rnap mRNA translation) 0.0226 (1/sec)
kf_4 (SP-CLPS, SP-PNLIP transcript. & translat.) 0.0275 (1/sec)
kf_1 (lipase secretion) 0.0257 (1/sec)
kf (procolipase secretion) 0.0257 (1/sec)

Reactions and Reaction Rates

Reaction Reaction Rate
t7rnap transcription: null → T7 mRNA kf_9*[B. subtilis].t7rnap*75
t7rnap mRNA degradation: T7 mRNA → null kf_7*[B. subtilis].mRNA*50
T7 mRNA translation: T7 mRNA → T7RNAP kf_6*[B. subtilis].mRNA*75
T7RNAP degradation: T7RNAP → null kf_8*[B. subtilis].T7RNAP*50
SP-CLPS transcription & translation: null → SP-CLPS kf_4*[B. subtilis].T7RNAP*10
SP-PNLIP transcription & translation: null → SP-PNLIP kf_4*[B. subtilis].T7RNAP*10
Procolipase secretion: SP-CLPS → procolipase kf*[B. Subtilis].[SP-CLPS]
Lipase secretion: SP-PNLIP → lipase kf_1*[B. Subtilis].[SP-PNLIP]

Experimental Model Results

Figure 13. Graph of experimental model results on a 10-hour timescale. The constitutive expression of the T7RNAP gene leads to the production of ~24 molecules of lipase and procolipase after 10 hours.

According to this model, given 2 copies of the T7RNAP gene, around 24 molecules of lipase and procolipase each will be secreted by one B. subtilis bacterium over a 10-hour period. If the plasmid copy number of pBS1C (which is unknown) is approximately 10 copies, then 120 molecules of lipase and procolipase will be produced within the same time period.

Future Applications

Our current ODE-based model, which predicts lipase and colipase production from B. subtilis, provides a strong framework for expanding into experimental validation. Once experimental data are collected, the model can be refined to more accurately represent the enzyme production dynamics within a biological context. By iterating between simulation and experimental results, we can fine-tune key parameters such as rate constants and enzyme-substrate affinities to better reflect real-world conditions.

Key elements for model refinement could include:

Rate Constant Adjustments: If experimental results show significant deviations in enzyme output, particularly in high or low-end CFU counts, rate constants will need to be recalibrated. This includes modifying the enzyme-substrate interaction rates and degradation rates based on observed kinetic behavior.

Substrate Availability and Environmental Factors: Experimental conditions, such as substrate concentration fluctuations and varying gastrointestinal pH, can also influence lipase and colipase production. Incorporating these factors into the model would lead to more accurate predictions.

Incorporating Experimental Feedback: New experimental insights could reveal nonlinearities or feedback mechanisms that were not previously accounted for in the model. Adding such complexity to the model will further ensure its predictive power.

BioBlueprint AI

The Problem

The current iGEM Parts Registry poses significant challenges in usability and accessibility:

Lack of Specific Search Features:

  • The registry lacks functionality to search for parts with precise characteristics (e.g., a strong yeast promoter with positive galactose regulation).
Incomplete Catalog Pages:
  • Category-specific catalogs are not comprehensive, meaning they do not list all parts within a category. For example, the terminator catalog excludes many terminators in the registry and offers no method to identify all terminators available.
Inaccessible Database Formats:
  • The “database dump” files are user-unfriendly. The XML file exceeds 12 million lines and lacks a search mechanism, while the SQL file, though organized into tables, requires technical proficiency and a standardized search mechanism for effective querying.
Compatibility Challenges:
  • The registry provides no intuitive method to identify parts that are compatible for constructing functional BioBrick circuits.

Our Solution: BioBlueprint AI

BioBlueprint AI addresses these issues by creating a streamlined platform that simplifies navigation and usability. The software compiles compatible parts and enables users to construct comprehensive BioBrick circuits based on specific inputs.

Workflow

Data Extraction:

  • Scrape individual part XML files to extract key information.
Data Integration:
  • Organize extracted data into a standardized data warehouse for accessibility and searchability.
AI Integration:
  • Fine-tune the Gemma LLM API with our curated data, optimizing it to meet the software's specific requirements and improve user experience.

Deployment and Progress

BioBlueprint AI is currently in the developmental stage, with significant potential to revolutionize how iGEM teams navigate and utilize the Parts Registry. While initial progress has been made, further training data is essential to enhance the model's accuracy and functionality. Due to the evolving priorities of our team, focus was shifted toward wet lab efforts and experimental validation. This adjustment reflected the immediate needs of our project but underscores the development of BioBlueprint AI.

This software could serve as an invaluable resource for iGEM teams, particularly during the ideation and development stages of their projects. By simplifying the process of identifying compatible parts and building functional circuits, BioBlueprint AI has the potential to streamline project workflows and spark innovation across diverse iGEM initiatives.

For those interested in exploring the progress made so far, all code and documentation are publicly accessible on our GitLab repository. Teams are encouraged to view the repository to understand the current state of development and envision how this tool might evolve into a game-changing asset for the synthetic biology community.

References

[1] Lewis, D. An Updated Review of Exocrine Pancreatic Insufficiency Prevalence Finds EPI to Be More Common in General Population than Rates of Co-Conditions. J. Gastrointest. Liver Dis. JGLD 2024, 33 (1), 123–130. https://doi.org/10.15403/jgld-5005.

[2] Pirahanchi, Y.; Sharma, S. Biochemistry, Lipase. InStatPearls; StatPearls Publishing: Treasure Island (FL), 2024.

[3] Ben Bacha, A.; Karray, A.; Daoud, L.; Bouchaala, E.; Bou Ali, M.; Gargouri, Y.; Ben Ali, Y. Biochemical Properties of Pancreatic Colipase from the Common Stingray Dasyatis Pastinaca. Lipids Health Dis. 2011, 10, 69. https://doi.org/10.1186/1476-511X-10-69.

[4] Schmalisch, M. H.; Bachem, S.; Stülke, J. Control of the Bacillus Subtilis Antiterminator Protein GlcT by Phosphorylation. J. Biol. Chem. 2003, 278 (51), 51108–51115. https://doi.org/10.1074/jbc.M309972200.

[5] Görke, B.; Stülke, J. Carbon Catabolite Repression in Bacteria: Many Ways to Make the Most out of Nutrients. Nat. Rev. Microbiol. 2008, 6 (8), 613–624. https://doi.org/10.1038/nrmicro1932.

[6] Philips, R. M. & R. » What is faster, transcription or translation? https://book.bionumbers.org/what-is-faster-transcription-or-translation/ (accessed 2024-09-30).

[7] Philips, R. M. & R. » How fast do RNAs and proteins degrade? https://book.bionumbers.org/how-fast-do-rnas-and-proteins-degrade/ (accessed 2024-09-30).

[8] Arnold, S.; Siemann, M.; Scharnweber, K.; Werner, M.; Baumann, S.; Reuss, M. Kinetic Modeling and Simulation of in Vitro Transcription by Phage T7 RNA Polymerase. Biotechnol. Bioeng. 2001, 72 (5), 548–561. https://doi.org/10.1002/1097-0290(20010305)72:5<548::AID-BIT1019>3.0.CO;2-2.

Sources of high and low reaction rates

Reaction Source
Glucose to G6P https://microbialcellfactories.biomedcentral.com/articles/10.1186/1475-2859-7-10
https://bmcmicrobiol.biomedcentral.com/articles/10.1186/1471-2180-4-6
G6P to F6P (Michaelis-Menten) https://support.megazyme.com/support/solutions/articles/8000025954-phosphoglucose-isomerase-bacillus-subtilis-e-pgibs-data-sheets
https://iubmb.qmul.ac.uk/enzyme/EC5/3/1/9.html
F6P to FBP (Michaelis-Menten) https://support.megazyme.com/support/solutions/articles/8000025954-phosphoglucose-isomerase-bacillus-subtilis-e-pgibs-data-sheets
https://iubmb.qmul.ac.uk/enzyme/EC5/3/1/9.html
FBP to ATP (Michaelis-Menten) https://en.wikipedia.org/wiki/Pyruvate_kinase
https://bio.libretexts.org/Bookshelves/Biochemistry/Fundamentals_of_Biochemistry_(Jakubowski_and_Flatt)/02%3A_Unit_II-_Bioenergetics_and_Metabolism/15%3A_Glucose_Glycogen_and_Their_Metabolic_Regulation/15.04%3A_Regulation_of_Glycolysis
G6P/F6P + CcpA to CcpA complex (Mass Action) Educated assumption
https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1280314/
HPr + ATP to HPrSerP (Mass Action) Educated assumption
CcpA + HPrSerP to complex (Mass Action) https://febs.onlinelibrary.wiley.com/doi/10.1111/j.1742-4658.2005.04682.x
CcpA complex binding to cre site (Mass Action)*** Educated assumption
LacI gene (1083 bp) to mRNA (transcription and degradation rates)*** https://bmcgenomics.biomedcentral.com/articles/10.1186/1471-2164-11-259
mRNA to lacI protein (translation rate) https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1007070
Formation of lac complex (Mass Action) https://en.wikipedia.org/wiki/Lac_repressor
GlcT-P to GlcT (Mass Action) Reasonable assumption
GlcT binding to RAT, transcribing t7rnap (Mass Action)*** Reasonable assumption
t7rnap gene (2652 bp) to T7 RNA polymerase (translation) Reasonable assumption
T7 RNA polymerase degradation (degradation) https://jbioleng.biomedcentral.com/articles/10.1186/1754-1611-4-9
CLPS (578 bp) and PNLIP (1347 bp) transcription & translation Reasonable assumption