Team:Edinburgh/Modeling/Whole cell model

Whole cell model

Description

We chose to base our work on a whole cell model originally developed by Tobias Bollenbach[1], tailored to capture the influence of different antibiotics on the growth rate and metabolism of E.coli. This model was suitable because it is a simple high-level model, and it was possible to understand and re-implement it within the time course of one summer. At the same time, it captures essential features of the cell metabolism which we might be interested in, such as the resources availability, the growth rate, and protein levels inside the cell.

The model as is it in the original paper assumes a constant influx of resources from the outer environment inside the cell. This resource influx is the only input to the model. As a function of this resource level, the model can predict population average values for the following parameters: • Growth rate • Resource concentration inside the cell • Protein number • Amount of DNA • Number of ribosomes • Speeds of synthesis of proteins, ribosomes and DNA • The cell volume • Number of replication origins and forks in the cell

The model doesn’t model RNA specifically, nor does it take account of the different proteins, their different expression levels and their different functions. Additionally, the model doesn’t model degradation of proteins, only of resources (which is supposed to capture resource usage by other factors which are not accounted for by the model).

The model consists of a system of 15 differential and algebraic equations which track the evolution of the parameters over time. The predictions given by the model are found by finding the steady-state solution of the differential-algebraic system. This solution gives parameter averages which are expected to be found in the cell in a state of steady exponential growth.

The model makes the assumption that the ribosome synthesis rate is optimized for highest growth rate. Hence, the model runs an optimization algorithm over this parameter to find what is the best synthesis rate for the given resource influx rate before running the actual simulation for this resource environment.

Solving the system

We tried solving the system with different solvers from Mathematica and Matlab but we couldn’t get any of them to work. The reason is that the system of equations is quite complex, includes algebraic as well as differential equations and a lot of self-dependences between the variables. Some of the solvers didn’t allow for algebraic equations in the differential system, while others failed to solve the system altogether.

We decided to program our own solver based on time step integration. Firstly, we implemented the simplest possible model – Forward Euler – but we got a mathematical problem. For some levels of resource influx, the system never reached steady state, but kept oscillating between different sets of values for the parameters. Since the output of the system is the steady-state solution, the model didn’t give any predictions in these cases.

With help from Nicky from the SDU-Denmark iGEM team, we learned about stiff differential systems and explicit and implicit methods for solving them. We started implementing different explicit and implicit methods. With all explicit methods we had the same oscillation problem, but implicit methods solved the system. So we ended up solving the system with Backward Euler.

Bollenbach, T., Quan, S., Chait, R., and Kishony, R. (2009). Non-optimal microbial response to antibiotics underlies suppressive drug interactions. Cell 139, 707-718.