Difference between revisions of "Team:SDU-Denmark/testvyff"

Approach

The Gillespie algorithm is a way to get calculate the evolution of stochastic functions; in this case cell concentrations. To use the algorithm, two things are needed:

1. reaction rates of the system (a1, a2, a3) at a given configuration.
2. a random number generator.

For each time step two things are calculated, using the random number generator: The time before next reaction and which reaction occurs. The time before next step is given by
( Figure (well... equation): reaction time insert here (tex doc))
Where S is the sum of the reaction rates and r1 is a random number between 0 and 1. This gives the time as if the system was one reaction with reaction rate S, using the random number to give a statistic distribution. The reaction is chosen proportionally to each individual reaction rate, using another random number. This way, reaction with high rates compared to other reactions will happen the most.
The reaction is carried out, the time is moved the calculated amount, and new reaction rates can now be calculated, to repeat the whole process. This continues until the time reaches the wanted limit or a specific number of reactions have occurred. It is necessary to have a limit on the number of reactions as it is possible for the time steps to grow smaller and smaller, in which case the calculation time quickly become either immense or impossible.

Approach

We are modelling the RelE RelB toxin-antitoxin system. RelE is a toxin restricting growth, by inducing a dormant state. This is inhibited by RelB, which forms complexes with RelE. Two different complexes are made: RelB₂RelE and RelB₂RelE₂, containing 1 and 2 RelE molecules respectively.
Both RelE and RelB are expressed from the same promoter, RelBE. When only small amounts of RelE is present, RelB and RelB₂RelE represses transcription of RelBE, by binding to the operator.
At higher concentrations of RelE, the toxin mitigates this repression, by reacting with bound complexes.Cataudella I., Trusina A., Sneppen K., Gerdes K., Mitarai N. Conditional cooperativity in toxin-antitoxin regulation prevents random toxin activation and promotes fast translational recovery. Nucleic Acids Res. 2012;40:6424–6434. doi: 10.1093/nar/gks297

For natural purposes the half life of RelB decreases significantly under starvation due to lon-protease, with shifts the equilibrium of RelB and RelE to a high state of RelE. The interactions with the promoter, keeps the amount of free RelE at a very low value outside starvation and stabilises the system. In our simulation the shift in equilibrium is made by introducing additional translation of RelE.
We used two model in two ways. First we saw how a given configuration of relB and RelE production increased the relE concentration and if it could cause dormancy within 2 hours. Second we investigated for how long each configuration

Rates and reactions

In the units off all reaction rates we use the approximation that in an E.coli. with a size of 1-2 μm, 1 molecule in the cell = 1nm. Thus we convert all units to be measured in molecules, as this fits the premises of the gillespie algorithm.
The RelB forms dimers at a high rate, so we assume all present RelB to be in dimers, capable of forming complexes with RelE.
RelB has a relatively low half-life at about 3-5 minutes, while RelE is stable and it’s half life is an effect from dilution due to growing bacteria (we use 43 min). During dormancy, growth is restricted and we increase RelE half life to 2000 min (around a day) as the dilution.
The transcription rates of RelE and RelB is based on the concentration of RelE and RelB under stable conditions. Here RelB is 10 times more prevalent than RelE (citation 2), so to make up for the higher half life of RelE, RelB has a much higher transcription rate than RelE (100 times)
The complexes are close to stable and given the same half life as RelE. However, to get free RelE to work RelB in complexes need to decay as well. The rate is set to a fourth of free RelB.

For the promoter bindings we let the operator be inhibited by binding with either RelB or 1-2 RelB₂RelE, given that the operator has to binding sites. We consider the cell to have four chromosomes with one promoter on each (less chromosomes would let the system work, but with more noise).
The initial values in the model are given by:

insert table with numbers and constants

Notice that all values are integers as the gillespie algorithm works with discrete numbers of molecules. The values were chosen based on a stable equilibrium found by letting the model run a simulation of the inherent system over 450 minutes, with different starting values.
The implemented total production rates shown in the model might seem too high as they range from 1-350 molecules pr. min, while the rates in the inherit system is effectively around 80-100 for relB and 2-5 for relE. The possibility of placing the system on high-copy plasmids, however, makes the high total production values reasonable, as the individual relE promoter only need a production of 0.01-2 (assuming a few hundred copies).
Considering the inherent toxin-antoxin system activating under starvation, we see that the the magnitude of relE copies around 50-80 molecules pr. cell. This makes it reasonable to believe that the cell enters dormancy when a few tens of free relE copies.

Running the model

The model is using the gillespie algorithm to give a stochastic view of the system and is run on matlab. The code uses an implementation by matlab user Nezar (https://se.mathworks.com/matlabcentral/fileexchange/34707-gillespie-stochastic-simulation-algorithm).
For the dormancy runs we input deterministic initial values, and then let the system run 30 min without activating the inserted toxin promoter. This results in a stochastic distribution of initial values mimicking variations between cells. Analysis show that 30 minutes is enough for the model to find a stable distribution, which is realistic considering the growth cycle of an e.coli.
For wakeup runs, we use the data generated at the end of a sleep run as initial value and remove the production of relE. We consider a cell to be woken up when the concentration of free relE decreases below 15 copies. This is perhaps a low value, but tests show marginal difference between 15 and 45 copies, where the low bound is chosen to decrease uncertainty of the cells state.
All runs simulate 1000 cells, which should be sufficient to get stable averages. The model assumes well-mixed conditions in each cell, but considers each cell independently. The model has no cut off for maximum values of relE, as we don’t know the exact relation between relE concentration and hibernation state, yet as a functional cutoff is found through wake-up times, it is not necessary.