Team:WHU-China/Model
1.Overview
The model part provides us an opportunity to attempt some interdisciplinary methods and one of the most important features of modeling is that we can design experiments reasonably and have a better knowledge of some quantitative features of the biological system after modeling. Interestingly, models can also help us to discover something novel and achieve goals that can’t be fulfilled with only traditional experiments.
In our project, the engineered bacteria (E. coli modified by design) limit the growth of the Propionibacterium acnes (P. acnes) in two ways: the engineered bacteria can secrete bacteriocins, which directly kill the Propionibacterium acnes; at the same time, the engineered bacteria compete with the Propionibacterium acnes for nutrient fatty acids.
We've developed two models to explain how our co-cultural system functions and how we accelerate fatty acids consumption of engineered bacteria. These are the growth simulator model and the flux accelerator model which consists of two parts.
Highlights
1. The growth simulator model predicts the trend of the change of our coculture system and gives
feedback to how to modify our product, considering it as a whole.
2. The first part of The flux
accelerator model reflects the metabolism of E. coli under different oxygen conditions to a certain
extent.
3. The second part of The flux accelerator simulates a more complex metabolic network and
successfully selects enzymes which have the most potential to promote the decomposition of higher
fatty acids.
In the growth simulator model, we not only explain how bacteria relieve acne on the skin surface, but also predict changes of the number of Propionibacterium acnes, which helps us to consider more influential factors that have not been considered in wet experiments. At the same time, this model simulates the final design of our project, reflecting it as a whole, which theoretically complements the experiments.
In the flux accelerator model, we first build an FBA model based on fatty acid β-Oxidation to analyze the enzymes that we need to overexpress to improve the fatty acids consumption ability of the engineered bacteria. Then, we use an existing molding tool to build a more complex FBA model for a more satisfying result. Deciding which enzymes had most potential to promote the decomposition of fatty acids saves us a considerable amount of work on the less effective enzymes.
2.The growth simulator model
2.1.Introduction
In the growth simulator model, we include Propionibacterium acnes and engineered bacteria into the same system, and construct an ordinary differential equation model simultaneously considering the processes of fatty acid consumption and bacteriocin secretion. The model reasonably explains the mechanism of the two pathways of growth inhibition of Propionibacterium acnes by engineered bacteria. Due to the limit of experimental conditions, the wet experimental group hasn't conducted coculture experiments of the strains. Therefore, the work of our mathematical model group is an extension of the wet experimental work, and at the same time illustrates the feasibility of the acne treatment idea to a great extent.
2.2.Assumptions
To simplify our model and emphasize the characters we want to investigate, we bring up following assumptions:
1. The promoter efficiency of bacteriocins is proportional to the concentration of fatty acids in the
environment[1], and the rate of the consumption of bacteriocins due to the exposure to the outside
environment is considered as a constant.
2. There only exists the competitive effect between the
engineered bacteria and the Propionibacterium acnes on fatty acids as well as the killing
effect through bacteriocins.
3. To emphasize the effect of fatty acids on engineered bacteria and
Propionibacterium acnes, we set the environmental capacity of both as a linear function of
fatty acids.
4. The killing effect of bacteriocins on P. acnes remains stable during the
whole process.
5. To eliminate the effect of units, we use the percentage of relative initial
time as data to fit the parameters.
2.3.Symbols and definitions
| Symbols | Definitions |
|---|---|
| S1 | Relative number of engineered bacteria |
| S2 | Relative number of Propionibacterium acnes |
| C1 | Concentration of bacteriocins |
| C2 | Concentration of fatty acid |
2.4.Parameters and description
| Symbols | Definitions |
|---|---|
| r1 | Growth coefficient of engineered bacteria |
| r2 | Growth coefficient of P. acnes |
| p1 | Bacteriocin production rate affected by engineered bacteria |
| p2 | Bacteriocin production rate affected by fatty acids |
| p3 | Decrease of bacteriocin per unit time |
| p4 | Increase of fatty acids per unit time |
| K | Ratio of reduction of Propionibacterium acnes to the number of bacteriocins |
| K1 | Rate of fatty acid consumption of engineered bacteria |
| K2 | Rate of fatty acid consumption of Propionibacterium acnes |
| a1 | Growth coefficient 1 of engineered bacteria affected by fatty acids |
| b1 | Growth coefficient 2 of engineered bacteria affected by fatty acids |
| a2 | Growth coefficient 1 of Propionibacterium acnes affected by fatty acids |
| b2 | Growth coefficient 2 of Propionibacterium acnes affected by fatty acids |
2.5.Model construction and explanation
We use the following symbol system[2] to represent the interrelationship of the four elements of the system (the bacteria, the Propionibacterium acnes, the bacteriocins and the fatty acids). Considering the indirect effect of fatty acids on bacteriocins, we use dashed lines to indicate indirect facilitation.
We use this symbol system to make a graphical representation of the relationship among engineered bacteria, P. acnes, fatty acids, and bacteriocins as follows.
This graph shows that both engineered bacteria and P. acnes consume fatty acid, engineered bacteria produce PctA( a kind of bacteriocins), and PctA inhibits P. acnes.
According to the graph, we transform the relationships into a system of ordinary differential equations as a simulation. After collecting data from online databases, we can work out the parameters, then predict the changes of the four parts.
Our ideas for designing each ordinary differential equation are as follows：
1. The first equation: In the modeling process, we have assumed that the growth of
the engineered bacteria is limited by Propionibacterium acnes and should roughly show an 'S-
shaped growth curve', so we refer to the famous logistic population size change theory to support
the growth model of the engineered bacteria.
2. The second equation: We assume
that bacteriocins are secreted only by engineered bacteria, so the p1S1 in the
equation represents this relationship, and because of the positive effect of fatty acid on promoter
efficiency of bacteriocins, we use p2C2 to represent this relationship.
Besides, the p3 indicates that the consumption of bacteriocins due to the environment is
stable through the whole process.
3. The third equation: The change in the
number of Propionibacterium acnes is the same as that of the engineered bacteria in
accordance with the Logistic growth model, and there is a corresponding decrease KC1 due
to the inhibition of bacteriocins, which for simplicity we assume to be proportional to the
concentration of bacteriocins.
4. The fourth equation: We consider that fatty
acids are only secreted by human skin at a constant rate p4, while fatty acid consumption
is caused by both engineered bacteria and Propionibacterium acnes, and the rate of
consumption is proportional to the number of both, and can be described by
k1S1 and k2S2 respectively.
2.6.Results
Based on the reasonable assumptions, we first gave a set of parameters, whose result can fit with our project well. Bringing them into the Growth Simulator, the final overall effect is shown as follows. And due to the lack of data from real experiments, the horizontal axis representing time doesn’t have real meanings but just indicates the relative time of the change.
By observing the trend of the four main elements in the system over time, we can initially derive the effect of each parameter on each of the four elements.
S1 (engineered bacteria) increases steadily with time, while we can also find a gradual decrease in the rate of growth (which is represented by the slope of the tangent line of the curve) with time. This indicates that in the co-culture system, the number of engineered bacteria increases over time, but this increase is limited by the size of its own population, which is also consistent with the Logistic model assumption.
The number of S2 (P. acnes) climbs with time at the very beginning but then drops suddenly and quickly, which shows the effect of bacteriocins and indicates that engineered bacteria limits the increase of P. acnes by competing for the nutrient substrate fatty acids.
C1 (the concentration of bacteriocins) gradually increases with time, with a strong correlation between its growth and the engineered bacteria. This indicates that the bacteriocins are secreted into the system by the engineered bacteria, and since the rate of secretion is assumed to be constant, a strong correlation ought to exist between the S1 curve and the C1 curve.
The decrease in C2 (concentration of the fatty acid) over time is easily explained by the fact that engineered bacteria and Propionibacterium acnes use fatty acids as nutrients.
2.7.Sensitivity Analysis
In this section, we change the parameters in our model separately by 5%, 10% and 20% with other parameters staying the same, to find the best character to improve our engineered bacteria.
Aimed at relieving acnes, we only choose to analyze the parameters relevant to the change of fatty acids and P. acnes. According to the description of all the parameters, we select r1, p1, p3, K and K1. All the parameters are increased by 5%, 10% and 20% in the following analysis except for p3, which is decreased by different percentages.
In the following figures, the vertical axis represents the value of y, which is calculated with the following definition:
y = y2 - y1
y1 represents the content of material when the content of P. acnes decreases to close to 0 after the change of a parameter. So y1 is a minimal positive value.
y2 represents the original content of material without the change of Parameters at the same time of y1 determined.
Fig. 5 above shows the difference in the change of percentage of fatty acids and P. acnes before and after the change of parameters shown in the horizontal axis by 5%. We can find that:
1. The changes of r1, p1and K have no significant effect on the relative content of
P. acnes. But when p3 or K1 is changed by 5%, the relative content of P.
acnes significantly changes, which means that if we manage to enhance the decrease of fatty
acids per unit time or the rate of fatty acids consumption by engineered bacteria, we can effectively
reduce the number of P. acnes.
2. Except for K1, changes of other parameters have
little effect on the relative content of fatty acids. The changes of r1, p1and K
even have negative effects because of the increase of the relative number of fatty acids.
When it comes to changes by 10%, all parameters have more significant effects on accelerating the reduction of P. acnes content and this kind of trend will continue when parameters are changed by 20%. The effects of changing p3 and K1 are still the most significant in both two cases.
The whole analysis process shows that all the five important parameters we choose will not cause obvious changes in fatty acids but they all have positive effects on the decrease of P. acnes in the whole system.
And surprisingly, for each degree of change, we can see that p3 and K1 both have quite apparent influences on the decrease of fatty acids and P. acnes, which means the system (including engineered bacteria, P. acnes, fatty acid and bacteriocins) represented by the growth simulator model is quite sensitive to these two parameters.
2.8.Conclusions
Overall, though we can only display the trend and relative value in the results of our model, the sensitivity analysis is really informative.
Thus, we can improve the effect of our ultimate product by reduce p3 (improving the stability of bacteriocins in the outside environment) and increase k1 (improving the ability of our engineered bacteria to decompose fatty acids). Especially if we can increase k1, the engineered bacteria will be more competitive for nutrient fatty acids, which is reasonable and is exactly what we want to do in the next model.
However, we have to notice that the positive effects resulting from the change of parameters are limited. According to the vertical axis from above figures, the degree of change from 5% to 10% will lead to almost doubling values, but this phenomenon cannot be seen in the change from 10% to 20%. Therefore, as long as we get a satisfying result, there is no need for us to chase the best effect.
3.The flux accelerator model
Part One
Introduction
Based on the analysis of the The growth simulator model, we propose the design idea of inhibiting the growth of P. acnes by accelerating the decomposition of fatty acids by engineered bacteria. Because of this, we designed the second model - The flux accelerator model. Targeted at improving engineered bacteria’s ability of consuming fatty acids, we build a two-part FBA (flux balance analysis) model[3] based on fatty acid β-Oxidation to see which enzymes should be overexpressed. Fatty acid β-Oxidation is the main process of engineered bacteria decomposing fatty acids, which involves the catalysis of a variety of enzymes. Overexpressing these enzymes will accelerate the decomposition rate of fatty acids, so as to improve the competitive effect of engineered bacteria on fatty acids. However, due to the fact that there is a limitation of metabolic pressure and the system should not be too complicated, it’s not a wise choice to overexpress all enzymes. We establish a flux balance analysis model and draw the conclusion that the enzymes most worthy of overexpression are FadD and FadE.
Assumptions
Apart from the basic assumptions of FBA method, we have identified some additional assumptions based on the characteristics of our project.
1. In any environment, engineered bacteria will finally reach a stable equilibrium under the constraint
of stoichiometry.
2. The purpose of this part is to consume fatty acids as much as possible, so
exogenic FA (fatty acids) is assumed to be excessive.
3. The rate of each reaction is relatively
stable in the Fatty Acid β-Oxidation process.
4. The reaction is complete.
Symbols and definitions
When we construct the metabolic network of Fatty Acid β-Oxidation process, the following symbols are used to represent the reaction process and substances:
| Symbols | Definitions |
|---|---|
| M1 | acyl-CoA |
| M2 | enoyl-CoA |
| M3 | 3-hydroxyacyl-CoA |
| M4 | 3-ketoacyl-CoA |
With these symbols, Fatty Acids β-oxidation can be expressed as:
Thus, the metabolic network can be described in the following metrological matrix form, with each row representing a substance and each column representing a one-step reaction:
| Reaction1 | Reaction2 | Reaction3 | Reaction4 | Reaction4* | Reaction5 | |
|---|---|---|---|---|---|---|
| M1 | 1 | -1 | 0 | 0 | 0 | 1 |
| M2 | 0 | 1 | -1 | 0 | 0 | 0 |
| M3 | 0 | 0 | 1 | -1 | -1 | 0 |
| M4 | 0 | 0 | 0 | 1 | 0 | -1 |
The metrology matrix is recorded as:
Set the reaction flux as 
, each component represents the reaction rate of the corresponding reaction.
Model construction and explanation
System Definition:
Reaction2,3,4,5 form a loop. Reaction1 is the cyclic input substance M1, and reaction4 outputs the substance acetyl coA. In Reaction4, each one-unit participates in the reaction and generates one-unit M4 and one-unit acetyl CoA at the same time. We determine the process of generating acetyl CoA in Reaction4 * as the boundary of reaction flow and Reaction1 as another boundary.
Goal:
Maximize the consumption speed of exogenic FA, namely, find the maximum value of Reaction1 rate under constraints.
Restraint Conditions:
1. In Reaction4, M4 and acetyl CoA have the same coefficient, which means that v4 = v4*
2. In the final state, the material reaction in the cycle reaches equilibrium, so the variation of M1, M2, M3 and M4 is 0, which means equation constraint：
3. The upper and lower bounds of all reaction Fluxes are shown as follow:
 component
|
v1 | v2 | v3 | v4 | v4* | v5 |
|---|---|---|---|---|---|---|
| Upper boundary 1 | 7.5611 | 0.8610 | 0.8610 | 1.1456 | 0.7809 | 0.8610 |
| Upper boundary 2 | 8.6254 | 0.7651 | 0.7651 | 0.4257 | 0.6852 | 0.7651 |
| Upper boundary 3 | 8.2003 | 0.7467 | 0.7467 | 0.6312 | 0.0888 | 0.7467 |
| Lower boundary | 0 | 0 | 0 | 0 | 0 | 0 |
Remark: we use the data in the database[4] to obtain the relative data of a unified unit to set our upper and lower bound values, which can save the troublesome steps of measuring accurate values and enhance the operability. The three values in the upper bound represent the data after E. coli is exposed to air for 5 min, 10 min and 15 min. The data is gained from a survey of E. coli transcriptome dynamics during the transition from anaerobic to aerobic conditions.
Result
We use the software MATLAB to solve the optimization problem and compare the values of the target flow vector with the upper bound. It is found that the rate of Reaction2 will improve the restricted target when exposed to the air for a short time, but the target reaction rate is obviously limited by Reaction1 when exposed to the air for 15min. As a result, we should at least increase the expression of FadD, the coenzyme of Reaction1. However, according to the results obtained by previous teams, we found that our results were improvable, so we used the existing flow balance analysis tools to consider β-Oxidation pathway more comprehensively and carefully (Part 2).
Part Two
Introduction
For the fatty acid β-Oxidation, if the process is only roughly considered as one cycle, many metabolites and metabolic processes will inevitably be ignored. Therefore, in order to simulate the process more accurately, we need to take the whole reaction and metabolite network in E. coli into account and establish a more comprehensive FBA model.
We have obtained a genome-scale metabolic network (GEM) from the Biochemical, Genetic, and Genomic (BiGG) database. iJO1366[5] is a genome-scale reconstruction of the entire metabolic network in Escherichia coli str. K-12 substr. MG1655, a well-researched and well-documented strain that is assumed to be metabolically equivalent to E. coli K-12 TG1. And it accounts for 1366 genes, 2251 metabolic reactions, and 1136 unique metabolites.
Assumptions
Because this model is an improvement of the model before, the basic assumptions are consistent. In addition, there are the following assumptions:
1. Only the fatty acid β-Oxidation is considered.
2. The reaction flux of the decomposition of
palmitic acid is used as an index to measure the ability to decompose fatty acids.
Model and explanation
The fatty acid β-Oxidation model is constructed by using the BiGG Model iJO1366 of Escherichia coli str. K-12 substr. MG1655, and the following metabolic network diagram has been obtained:[6][7]
Aiming at the maximum growth rate of Escherichia coli, we uploaded metabolic data and the following initial metabolic network diagram can be constructed:
Note:
1. In this case, the growth rate of Escherichia coli is 0.982h−1.
2. Color of each
reaction:
Min: gray; Med, purple; Max, red, as shown in the Fig. below.
That is, the color represents the relative reaction flux, and the absolute reaction flux is the number on the right side of this step.
3. The direction of the arrow is the defined direction of the forward reaction.
4. Some reaction rates
seem to be zero, but in fact, they are due to the number of significant figures, which are actually not
zero.
Result
The optimization goal is to maximize the decomposition rate of fatty acids. The 4 results are as follows:
1. Under basic conditions, the maximum consumption of fatty acids is 88.324 mmol/ (gDW*hr). (The amount of the metabolite in millimoles in a gram of bacterial cells (dry weight) over the time of an hour.) That is, it is feasible to increase the consumption of fatty acids under limited conditions.
2. In addition, we set the reaction flux of the specific step as 1 to simulate the situation when FadE is overexpressed. With the maximum decomposition rate of fatty acids as the optimization objective, the results are as follows:
According to the result, the reaction rate of the decomposition of fatty acid under the circumstances is 86.610 mmol/ (gDW*hr).
3. When we set the reaction fluxes of the reaction that FadB and
FadJ are involved in to simulate the overexpression of these two enzymes (set the
reaction rate of these two steps as 1), taking the maximum decomposition rate of fatty acid as the
optimization objective as well, there is no solution, which indicates that the reaction flux of
these two steps cannot be larger if the decomposition rate of fatty acid aims to be the
fastest.
The maximum reaction fluxes of these steps are both 1.95*e-14 mmol/ (gDW*hr) if the
maximum reaction fluxes of these steps are taken as the optimization objective.
4.In terms of FadA and FadL (setting the reaction flux of these two steps as 1), taking the maximum decomposition rate of fatty acid as the optimization objective, the results are consistent with the FadB and FadJ in result3.
Conclusions of two parts
Based on the above results, the following conclusions can be drawn:
a)From a computational point of view, it is very feasible to increase the consumption of fatty acids.
b)Over-expression of FadD and FadE can improve the ability of our engineered bacteria to decompose fatty acids, while the over-expression in the other two steps of fatty acid β-Oxidation is of little significance and may increase metabolic pressure.
According to the results, only the effect of the most promising enzymes FadD and FadE are overexpressed in our engineered bacteria and tested experimentally.
4.Future
Due to the limit of time and experimental conditions, our models are not complete to some extent. If time permits, we plan to polish our models and investigate deeper into them in these aspects:
1. Adjust parameters of our growth simulator or add more details to it in accordance with experiment
data and try to make it appliable to more conditions that haven’t been considered yet.
2. Learn
more about the FBA method and bring up more reasonable and detailed constraints in order to get a
better solution.
Reference
[1]https://2019.igem.org/Team:NTHU_Taiwan/Parts
[2]Momeni, B. , X. Li , and W. Shou
. "Lotka-Volterra pairwise modeling fails to capture diverse pairwise microbial interactions."
eLife,6,(2017-03-18) 6(2017).
[3]Kenneth J Kauffman; Purusharth Prakash; Jeremy S Edwards (2003).
Advances in flux balance analysis. , 14(5), 491–496. doi:10.1016/j.copbio.2003.08.001
[4]Partridge
JD, Scott C, Tang Y, Poole RK, Green J. Escherichia coli transcriptome dynamics during the transition
from anaerobic to aerobic conditions. J Biol Chem. 2006 Sep 22;281(38):27806-15. doi:
10.1074/jbc.M603450200. Epub 2006 Jul 20. PMID: 16857675.
[5]BiGG
Model iJO1366 ( ucsd.edu )
[6]Rowe E, Palsson BO, King ZA. Escher-FBA: a web
application for interactive flux balance analysis. BMC Syst Biol. 2018 Sep 26; 12(1):84. doi:
10.1186/s12918-018-0607-5. PMID: 30257674; PMCID: PMC6158907.
[7]Heirendt L, Arreckx S, Pfau T,
Mendoza SN, Richelle A, Heinken A, et al. Creation and analysis of biochemical constraint-based models
using the COBRA Toolbox v.3.0. Nature Protocols. 2019 Mar 1; 14(3):639–702.
[8]https://2019.igem.org/Team:BNU-China/Model
