Difference between revisions of "Team:KU Leuven/Modeling/Internal"

Transcription is the first step of gene expression. It involves the binding of RNA polymerase to the promoter region of DNA and the formation of mRNA. The transcription rate is dependent on a number of conditions like the promoter strength and the used strain. Our system has both constitutive promoters, which are always active, and inducable promoters, which can be activated or repressed.

Maximum Transcription rate

First we will try to find the maximum transcription rate. The prediction of the transcription rate has been an important hold-back in the past, even though Polymerases per Second was introduced as a unit. This is the amount of Polymerases that passes through a given position in the DNA per time unit and is essentially the transcription rate at a particular location on the DNA ^[1]. We can use this value as the transciption rate of the whole gene because it is the value of the slowest and thus rate-defining step, the binding of polymerase on the promoter. It remained difficult to measure in vivo but Kelly et al. have introduced a way to measure the activity of promoters using an in vivo standard, promoter J23101. This gave rise to a new unit: Relative activity of promoter (RPU)$=\frac{PoPS_{phi}}{PoPS_{J23101}}$ with $phi$ being the promoter of interest. Through the use of relative units, the variability due to equipment and conditions was drastically decreased. The PoPS of J23101 was found to be 0.03^[2].

As our constitutive promoter, we used J23114. On the iGEM website we found that the whole J constitutive promoter family had been characterized. From these measurements we conclude that J23114 has a PoPS of 0.00429. The strength of the pLuxPR promoter has been measured relatively to pLacIQ by iGEM Tokyo_Tech 2010. Using these measurements, we calculate the PoPS of pLuxPR to be 0.00169 PoPS. Next we sought values for pCI. iGEM NYMU-Taipei 2009 measured the pLux promoter relative to pCI, which allows us to calculate the PoPS for pCI: $\frac{1}{0.43} {\cdot} 0.00169=0.0394$ PoPS. They also measured pPen from which we are able to calculate the PoPS: $0.3106 {\cdot} 0.0394=0.01224$.

These values all depend on measuring activities and since the strains, media and antibiotic markers won’t be completely the same, the real values will diverge from the values found with these simple calculations. Nevertheless, these values can give us an idea about the relative strength of the different promoters. We could measure the real values of mRNA concentration in our cells by executing a qPCR, but for now we will use the found values for our model. Another important parameter for transcription is the copy number in which the genes are present in the cell. The plasmid that incorporates the genes, has an ORI with a copy number between 100 and 300. In our model we will use a mean copy number of 200.

Inducable promoters

For inducable promoters the maximum transcription rate has to be multiplied with a factor to integrate the effect of the transcription factors. Our system includes two types of DNA binding proteins: repressors and activators. These proteins bind certain regions in the promoter region of the DNA and can either repress or activate the promoter activity, affecting transcription rate. The repressor proteins used in our system are the cI protein of the lambda phage and the PenI protein. The activator protein is LuxR, which has to be activated with AHL. The proteins first need to form dimers before they are able to bind their target DNA.

DNA binding will be simulated using the Hill function : ${\theta} = \frac{{[L]}^n}{{K_d}^n + [L]^n} = \frac{1}{1+(\frac{K_{d}}{[L]})^{n}}$. $\theta$ is the amount of DNA bound by the protein, $L$ is the amount of protein, $K_d$ is the dissociation constant and $n$ is the Hill coefficient ^[3]. This Hill function has values between 0 and 1 and acts like an ON/OFF switch.
For the repressors, we are interested in how much of the DNA is still unbound and active: $\frac{1}{1 + ({\frac{[Repressor]}{K_{d}}})^{n}}$
For the activators, we are interested in how much of the DNA is bound and active: $\frac{1}{1+(\frac{K_{d}}{[Activator]})^{n}}$
The Hill coefficient gives an idea about the strength about the interaction between the DNA and the protein and is sometimes used to estimate the amount of proteins that bind the DNA, but caution is needed in using it. Indeed, the Hill function has known shortcomings, but it is still very useful because of its simplicity ^[3].

The constants needed for our model are $K_d$ and $n$. :

Table 1: The constants of the Hill function for cI, LuxR and PenI

${K_d}$ is the value which corresponds with the amount of repressor or activator necessary to repress or activate half of the promoters. Since our copynumber is 200, we multiply the Hill constants with 200 because this higher value seems more logical to us.

Leakiness

A last important value in transcription is the leakiness of the promoter. A shutdown promoter still has a very small transcription which is called the leakiness of the promoter. This value should be as low as possible, because there should be a big difference between ON and OFF promoters. The effect of leakiness of promoters is something that should be checked.

iGEM 2012 Upssala measured the leakiness of the pLuxPR promoter. Without AHL and LuxR, the pLuxPR was partially active and had a RPU of 10. This means that the leaking promoter has a PoPS of $\frac{10}{390} {\cdot} 0.040 {\cdot} 0.03 = 3.077 {\cdot} 10^{-5} PoPS$. We can use this value to estimate the leakiness of the other promoters.

The next step in gene expression is the translation of the synthesized mRNA to a protein. This involves the binding of a ribosome complex to the RBS site of the mRNA, followed by elongation and termination. The initiation step is the rate-limiting step of the three. The translation rate is thus tied to the thermodynamics of RBS-ribosome interaction, but other sequences (5’-UTR, Shine-Delgarno,repeats, internal startcodons) and possible secondary structures (hairpins) also play a role ^[8]. Salis et al have developed an algorithm to predict the translation initiation rate based upon thermodynamic calculations. This algorithm is used in the Ribosome Binding Site Calculator⁹. We used this calculator to predict the translation rates of our mRNAs.

The results of the calculator are shown in table 2. The program gives us results in au (arbitrary units). Since we know the translation rate of LuxI and LuxR, we can use these as a standard to calculate the other translation rates since the used scale is proportional^[10]. LuxI and LuxR have almost the same output from the calculator which corresponds with the paper where they also have the same translation rate. Our values are also normal values since values between 1 and 100 000 are possible^[7]. We did get warnings about the accuracy of the prediction (NEQ: not at equilibrium) which happens when mRNA may not fold quickly to its equilibrium state.

Table 2: The translation rates for the different proteins as found by the RBS calculator using LuxR and LuxI as standard

For Transaminase B we also need to add a production term because it is a protein which is also naturally occurring in E. coli. We estimate a value of 0.0033 for this constitutive production rate.

Before the proteins can bind the promoter region, they first have to form complexes. cI and PenI form homodimers, while LuxR first forms a heterodimer with AHL and afterwards reacts to form a homodimer with another LuxR/AHL molecule. This next part will describe the kinetics of such complexation. The kinetics will be modelled by simple mass action rates. We only need the parameters of LuxR and cI, because the Hill function found for PenI was already adapted for the protein in monomer form.

Table 3: Dissociation and association rates for the oligomerization of the transcription factors used in our circuit

Our models comprises two product forming enzymes: LuxI and Transaminase B. LuxI is the enzyme responsible for AHL production and Transaminase B is responsible for leucine formation. Therefore, we need to describe the following reactions:
$$ \text{acylACP} + \text{SAM} \rightarrow \text{ACP} + \text{MTA} + \text{AHL} $$ $$ \text{glutamate} + \text{aKIC} \rightarrow \text{aKG} + \text{Leucine} $$
For LuxI kinetics Schaefer et al. (1996) found a value of 1 molecule/minute as Vmax. Since we assume that LuxI is fully saturated with substrate, we take this value as the synthesis rate of AHL. ^[12]

Transaminase B synthesizes leucine from aKIC and glutamate. In this reaction glutamate is converted to aKG. The reaction kinetics is of the Ping-Pong Bi-Bi form. This model describes a mechanism in which the binding of substrates and release of products is ordered. The enzyme shuttles between a free and a substrate-modified intermediate state. ^[13]

Yang et al. (2005) studied these reactions and found a method to estimate the constants: the Lambda approximation ^[14]. Their constants of this reversible Ping-Pong Bi-Bi reactions are put in the next table :

Table 4: Constants of the Ping-Pong Bi-Bi reactions using Transaminase B as found by Yang et al. (2005) ^[13]

Proteins, mRNA and other metabolites have a turnover rate. They are degraded over time. The degradation rate will be described with first order kinetics where d represents the degradation constant. Not every molecule has the same degradation rate since some molecules are more stable than others. We can influence the stability of the molecules. For example, in cell B it is important that there is a fast switch between conditions and a fast turnover of CheZ and RFP is necessary. This is why we add a LVA-tag to these proteins. This tag destabilizes the protein and makes them degradade faster. For Transaminase B we choose a very high degradation rate, because we did not include degradation terms for the Transaminase B bound to substrate. The degradation rates used in the model are put in the next table:

Table 5: Degradation rates for the different biomolecules used in our circuit

Our model has two types of diffusion. Diffusion over the cell membrane from the inside of the cell to the outside and diffusion in the external medium. The diffusion over the cell membrane is more complicated because proteins play a role in it and the membrane is not equally permeable for every molecule. For diffusion we need to change our unit because diffusion is usually used with concentration as unit. This is why we added a correction term for the volume of the cell and the external volume.

We assume that the volume of the cell has a shape of a cilinder. For this simplified volume Goryachev et al. found a value of $5.65 {\cdot} 10^{-16}l$ (radius=0.3µm, length=2µm)^[10]. We can take this volume as a constant, since cell growth is very small compared to the diffusion. The outside compartment will be modeled as half a sphere with the cell as center and a radius equal to $\sqrt{2{\cdot}D{\cdot}t} + r_0$. We only take half a sphere because there is no upward diffusion. For the initial value of the outside compartment we take a volume (and thus $r_0$ slightly bigger than the cell volume (radius). We choose a value of $5.7 {\cdot} 10^{-16}$l which accords to a radius of $5.128{\cdot}10^{-6} dm$. Because there is also already a volume taken by the cell, the total initial outside compartment has an initial value of $5{\cdot}10^{-18} l$. The volume will increase per time step with $\frac{2{\cdot}\pi{\cdot}D}{\sqrt{2{\cdot}D{\cdot}t}}{\cdot}{(\sqrt{2{\cdot}D{\cdot}t} + r_0)^2}$

With this approximation we get more logical results. Now we can calculate the diffusivity of AHL and leucine. Goryachev et al. used this formule $$D_{AHL,in}=\frac{P{\cdot}{S_{cell}}}{V_{cell}}$$ with P equal to $(1-5){\cdot}10^{-6} cm^2/s$ and $S_{cell}=4.34{\cdot}10^{-8} cm^2$^[10]. This formula gives us $D_{AHL,in}$=0.08-0.4 1/s. We take the diffusion rate for AHL equal to 0.23 1/s. For leucine we make a difference between inward and outward diffusion. Since leucine is an amino acid, a cell does not want to lose its leucine but it will want to take leucine in. Inward diffusion is thus facilitated by transporters, so we choose this value to be the largest. The inward diffusion rate is estimated at 0.05 1/s and the outward diffusion 0.0005 1/s.

Table 6: Diffusion rates of the molecules in our circuit

The designed circuit in cell A is under control of a temperature sensitive cI repressor. Upon raising the temperature, cI will dissociate from the promoter and the circuit is activated. This leads to the initiation of the production of LuxR and LuxI. LuxI will consecutively produce AHL, which binds with LuxR. The newly formed complex will then activate the production of leucine and Ag43. Leucine and AHL are also able to diffuse out of the cell into the medium. Ag43 is the adhesine which aids the aggregation of cells A, while leucine and AHL are necessary to repel cells B.

We can extract the following ODE's from this circuit:

We visualize these ODE's in the Simbiology Toolbox which results in the following diagram:

Figure 1: The internal network of cell A visualized by Simbiology

The system of cell B is also under control of the cI repressor and is activated similar as cell A. The activation by the temperature raise, leads to the production of LuxR. AHL of the medium can diffuse into the cell, binding LuxR and activating the next component of the circuit. This leads to the production of CheZ and PenI. CheZ is the protein responsible for cells to make a directed movement, governed by the repellent leucine. PenI is a repressor which will shut down the last part of the circuit which was responsible for the production of RFP.

We can extract the following ODEs for cell B from this sytem:

We visualize these ODE's in the Simbiology toolbox. This gives us the following diagrams:

Figure 2: The internal network of cell A visualized by Simbiology