Team:Baltimore BioCrew/Model

Model



Overview

We started by creating a predictive model that shows the population of B. theta over time when grown with/without the presence of MAFF. However, in order to create this model, we needed certain parameters. As B. theta has not been studied well in-depth, these parameters have not been well documented. Considering this, we created a model to find the growth rate, half-saturation rate, and cell yield of B.theta, with the help of published data (Nakajima et al 2018).

schematic



Part 1: Introduction to Monod Model

The Monod Model simulates a situation where microbial growth is limited by a nutrient. Typical cell population growth follows the following phases: the lag, exponential growth, stationary, and death phases. However, the Monod Model only incorporates and considers the exponential growth and stationary phases. The exponential phase is where growth is proportional to the cell number and is most rapid in the culture. The stationary phase occurs next, and is indicated by no net growth of the cells in the culture (death rate is balanced out by the growth rate).

Outside of the scope of the Monod Model , there are the lag and death phases. During the lag phase, the growth rate is basically 0 as cells adapt to culture conditions. This phase becomes the exponential phase after the initial population doubles. In the final death phase, the death rate overtakes the growth rate, resulting in a net loss of cells.

part 1 graph

A typical growth phase for bacterial growth in culture. Source: Maier, 2009.

The Monod Model can be described using a set of ordinary differential equations. Ordinary differential equations (ODEs) describe the rate of change of a dependent variable in relation to itself and the other variables. We use the Monod Model to predict how the amount of cells and substrate change over time based on the cell population and the substrate concentration.

part 1 graph
Parameters and Variables Representation
t Time, the main independent variable
St Substrate level over time
Nt Cell density over time
mumax Maximum specific growth rate for a culture, with units of 1/time
Ks Half-saturation constant with units of mass/volume
Z Cell yield, representing the unit amount of substrate over unit amount of cells produced



Part 2: Naive Optimization to Establish Initial Estimate of Parameters

What is optimization?

Optimization is the process by which we can determine the values of certain parameters that best satisfy defined certain output conditions. We optimized for mumax, Ks, and Z, such that the simulation closely approximates experimental data (below). To do so, we implemented used the lsqnonlin MATLAB function.

Our experimental data was taken from the following graph. Since the paper didn’t give us exact numbers for each point, we had to find the values of those points ourselves. To do this we placed the graph within WebPlotDigitizer,which allowed us to find the exact values of each point based on the axes.

Our experimental data comes from the following graph:

part 2 nakajima

Growth curves of wild-type B.theta and B.ThetaΔMAFF strains. Source: Nakajima, et al 2018.

We modeled all four of the lines. The red line represents the growth curve of Wild-Type B.Theta (that can produce MAFF) growing in a competitive media. The blue line represents the growth curve of B.Theta that cannot produce MAFF (through gene deletion) growing in a competitive media. The black line is the same as the red but in nutrient rich media. The green line is the same as the blue but in nutrient rich media.We are assuming that the competitive media will best reflect the gut.

Our Results:

Optimization results for MAFF growth

Parameters:

Parameters (units) WT-rich ΔMAFF-rich WT-poor ΔMAFF-poor
mumax (1/hr) 3.00 1.40 0.70 0.48
Ks (unitless) 4.58 2.12 0.69 0.0043
Z (unitless) 3.09 3.17 3.63 3.55


Conclusions:
    Overall, model simulations fit the data pretty well, but there were some deviations
    Decided to model lag phase because we knew the Monod Model did not model lag phase but it was present in data
    We expected that including the lag phase could help improve the accuracy



Part 3: Limitations of the Monod Model

As discussed earlier, the Monod Model is specifically designed to model the exponential growth and stationary phases. This creates an issue because the graph containing our experimental data includes a lag phase (since it’s based on an actual experiment) which is not included in the Monod Model. The lag phase is an important part of cell growth and has a fundamentally different graph shape from the exponential and stationary growth phases modeled by the Monod Model. Thus, we needed to modify our model to account for the lag phase.




Part 4: Improved Optimization

In order to develop an improved optimization procedure, we aimed to optimize for the lag phase and understand when to best take measurements during our own experiments.

To include the lag phase in our optimization, we added a new parameter to the code: t-shift (time shift). We used t-shift to move the starting point for the data used to fit a Monod Model. Our code will re-optimize the model parameters for progressively smaller and smaller data sets and then display the residuals, which is the sum of the difference between each experimental point and its nearest simulation point squared for each one.

It is expected that as data corresponding to the lag phase is removed, the residuals will decrease since the Monod Model is only intended to model the remaining exponential and stationary phases. However, in using this approach, we suspect that there would be a point where the optimization would improve, not from the removal of the lag phase, but due to it being easier to optimize a line to a smaller data set (since it would have fewer points increasing the residuals). To explore if this were the case, we generated data using known parameters for E.Coli growing in culture with the addition of 2 hours of minimal change from the initial population. Using this data allows us to check how well our approach to optimization is performing.







The graph above shows the optimization with the lowest residuals is the one which includes a lag phase of two hours. The residual at this time is significantly lower than any of the points around it, showing that this is where the Monod Model does the best at fitting a line to the data set.

After we could identify the lag phase and exclude it from our experimental data, we proceeded to analyze the effect of unintended errors at several time points in the experiment. This process is referred to as sensitivity analysis. There are two approaches to sensitivity analysis: systematic and unsystematic. Systematic errors often come from the measuring instruments, and errors are consistent throughout the entire experiment. Unsystematic errors, on the other hand, come from unpredictable and unintentional changes at some points in the experiment, and only some data points are faulty. Since we aimed to analyze which phase of bacterial growth was the most sensitive to error, we decided to go with the unsystematic error approach.

We started by introducing a 5% error into each data point in the Nakajima paper and analyzed the effect of the error by calculating the percent change in the residual values, which we referred to as residual cost. We had predicted that the residual cost would peak during the exponential phase because the error of data in the exponential phase would be exponentially different. However, the result (as shown in the graph below) indicated that changes at points later in the experiment would yield more errors.



To double-check our findings, we returned to the E. Coli hypothetical experiment that we generated previously to analyze t-shift. However, instead of 9 data points, we simulated 200 and reran the sensitivity analysis on them. The result of this simulation correlated with what we found from the Nakajima paper.



The graph above shows that the residual cost of the experiment will increase exponentially if the error is introduced later in the experiment and will plateau when we reach the stationary phase at around t =10.





In addition to the residual cost, we also analyzed the effect of errors on the optimized parameters. The graphs above suggest that the parameters generally do not experience significant changes in their values (peaks of mumax and Z barely reach 0.1%) except for Ks at the initial time point and during the 6-10 hour time frame, which is when the culture is transitioning from the exponential phase to the stationary phase. Additionally, the graphs show that all three parameters experience the same trend as the residual, where the error exponentially increases until we reach the stationary phase. This trend suggests that our model is more susceptible to errors at later time points in the experiment, meaning that experimental data should be collected more carefully and frequently at these times so that more accurate parameters could be optimized.



Conclusions and Future Application

At the outset, our goal was to model the effect that MAFF has on the growth of bacteria within the gut’s microbiome. Using data from the Nakajima paper (currently the only published paper to study MAFF’s effect on bacterial growth) and the Monod Model, we developed a rudimentary model of the system. From there, we worked to improve the quality of our model. Optimized parameters from this model gave us equations that accurately described the effect of MAFF on bacteria cell growth. We then rigorously tested the sensitivity of our model. The results from this process showed us the sections of our model that are most sensitive to a slight change in measurements of cell density.

We have done a few things that have practical applications. The first is that we have found the parameters and equations that model the growth of B. theta when it is in a system with and without MAFF. This allows us to accurately model a system that has only recently been discovered. It also allows us to understand the strength of the effect of MAFF on the growth of cells. The sensitivity analysis we did has practical applications as well. We can use it to figure out which sections of cell development errors in measurements would affect the parameters of the model the most. With a more accurate set of parameters, we can predict the necessary density in the probiotic in order to achieve a desirable level in the gut.

Citations:

Alfred B. Cunningham, John E. Lennox, and Rockford J. Ross, Eds. (2001-2010) Chapter 2 Biofilm Growth and Development (pp.1)

Maier, R. M. (2009). Bacterial Growth. In Environmental Microbiology (pp. 37-54). Elsevier Inc.. https://doi.org/10.1016/B978-0-12-370519-8.00003-1

Nakajima, A., Vogelzang, A., Maruya, M., Miyajima, M., Murata, M., Son, A., … Suzuki, K. (2018). IgA regulates the composition and metabolic function of gut microbiota by promoting symbiosis between bacteria. The Journal of experimental medicine, 215(8), 2019–2034. doi:10.1084/jem.20180427


</html>