Team:Aberdeen Scotland/internal/SimBiology

Gene Regulatory Network Modelling in SimBiology

Introduction

SimBiology is the MATLAB toolbox used by many systems biologists and iGEM teams alike to model gene regulatory networks. We decided that using this software would be a way of verifying what we have found from our own C models, but we do so with caution. In SimBiology everything on the user-interface is modelled in terms of a series of mass action chemical reactions, as opposed to a system of differential equations. The rates of the forward and backward reactions can be defined by a normal rate equation, or by any function such as the Hill-Function. The actual solvers and processes are out of the control of the modellers, and we are left to trust SimBiology, which is why we prefer to program our own things in C.

The advantage of Simbiology is that a slightly more complex model could be created more easily than in C. In the literature about the Lac operon there is a detailed analysis of the on and off rates of the DNA going between bound and unbound states. In SimBiology it is possible to model the change between unbound and bound DNA because unlike in the C differential equations it is necessary to have a level of DNA from which mRNA is produced.

SimBiology also comes with sensitivity analysis software, however as you will later see, the results of SimBiology compared to the data that comes out of our Monte-Carlo simulations is greatly inferior.

Description of Model

A schematic diagram of the model as it works in SimBiology is displayed below:

The dashed lines represent an enzymatic reaction, and the arrows represent the direction of chemical reaction. As stated before all mRNA productions must be modelled by an enzymatic production from DNA, and all protein productions must be modelled as enzymatic production from mRNA. This is useful as the repression and activation of gene networks can be modelled by mass action reactions between the unbound and bound form of the DNA or mRNA.

Our SimBiology model can be downloaded on the download page if you want to see how we modelled it exactly, and if you want to see the set of parameters we used. Below is a graph of the system turning on in response to input signals of quorum sensing and IPTG, and below that is one of the system reaching steady state without inputs:

Sensitivity Analysis

SimBiology comes with sensitivity analysis software, and so we proceeded to see what results it might come out with. However, unlike our Monte-Carlo simulations we cannot include active events such as the turning on and off of signal chemicals and quorum sensing or test for functionality such as the production of glue before lysis or ensuring the system does not automatically activate. Therefore the results that we expect this to give will only be a general indication of which parameters have a greater affect than others on system behaviour.

General Sensitivities:

By initially selecting all species and all parameters it became apparent that the most sensitive parameters were the on and off rate constants for the dimerisation of LacI. These had a sensitivity of many thousand as can be seen from the graph below. The reason for LacI being singled out as a sensitive species and that certain parameters would make it sensitive is that the level of LacI dictates whether the system is latched on or off. However these are insignificant results, as we know the dimer and tetramer of LacI to be very stable. The units in all the graphs are unit-less and have been normalised and are time integrals used to compare one sensitivity to another.

Degradation Rate Sensitivity:

Most sensitive, as can be seen in the graph, is the degradation rate of the mRNA which if too great or small can disrupt the behaviour of the cell, as all mRNA here has been assumed to have the same degradation rate.

The other sensitive degradation rates, more clearly illustrated in the above graph, are that of lambda CI (3), LacI (4) and the holin-antiholin complex (16). This is what should be expected, as lambda CI and LacI affect each other and both form parts of the latching mechanism. In particular the degradation rate of lambda CI is important, this would be because the system could never latch on if it is too high, and might always turn on if it is too low.

Production Rate Sensitivity:

As can be seen in the graph, the anti-lysis mRNA (10) is very sensitive to its production rate. This is something of insignificance, as of course it will be sensitive to this value, but as long as the production rate is not so high that lysis never occurs then the system should work.

Again LacI (2) is the species that is most sensitive to parameters. It is sensitive to, not surprisingly, the minimal production rate and the production rate of the latch mRNA which makes lambda CI and is capable of decreasing the levels of LacI and the translation rates of proteins, gamma protein.

K Value Sensitivities:

As we have seen with the other sensitivity analyses, LacI is sensitive. Here it is only sensitive to the hill coefficient controlling its production, K_ci.

Conclusions and Results

The only thing to conclude from using SimBiology is that it demonstrates that with our parameter set, and data from literature the system should work, and within certain levels of input signals it exhibits some AND gate behaviour and latches on successfully. The sensitivity analysis has produced very little important information as we already knew the system was very sensitive to the levels of LacI and lambda CI as these two chemicals form the basis of the AND gate and latching mechanism.

References

[1] A.B. Goryachev, D.J. Toh, T. Lee. “Systems analysis of a quorum sensing network: Design constraints imposed by the functional requirements, network topology and kinetic constants”. BioSystems 83 (2006) 178–187

[2] Michail Stamatakis and Nikos V. Manttaris. “Comparison of Deterministic and Stochastic Models of the lac Operon Genetic Network” Biophysical Journal Volume 96 February 2009 887-906 887