Single-Cell Model

The single-cell model of our system consists of three connected but easily separable modules: the sensor, the band detector and the filter. The sensor detects toxic substances from cigarette smoke. The band detector produces green fluorescent protein upon a detection of a certain range (band) of toxic substance concentration. After the toxic substance concentration crosses a certain threshold, the filter turns the system red (red fluorescent protein is expressed). The single-cell model helped us understand in detail whether and how our synthetic circuit works, especially to explore the characteristics of the GFP band and get a feeling of how it is related to concentration of our input substance.

Overview

Figure 1: Molecular mechanism of SmoColi with acetaldehyde sensor. The model with the xylene sensor is similar.

The single cell model is a mathematical model of the species involved in our synthetic system of a single SmoColi cell. The change of the species concentrations in time is given by non-linear ordinary differential equations (ODEs), most of which follow Hill kinetics. The parameters we used in the model are derived from literature supported by experimental evidence. Some of the parameters are slightly adjusted to meet our own experimental approaches. The single cell model was simulated using numerical Matlab solvers.

Our model has several important strengths which are:

1. MODULARITY
Our system consists of three modules:

• Sensor Module
  
• Band Detector Module
  
• Filter (Alarm) Module

The modules are almost independent, related between each other just by a single species. Most of the parameter values used within the modules do not depend on the other modules as well. The independence of the three modules enables their reusage in other systems. So if we want for example to add another sensor module to the existing band detector, the only thing we should change in the band detector is few parameter values for the linker species. This important property allows us to use the band detector + alarm with different sensors, making our synthetic cells a biological quantifiers of a range of substances.

2. ROBUSTNESS
We did a robustness analysis for our system and found out that it is indeed robust. This means that even upon parameter variation, the GFP band does not dissapear, just gets shifted or varies in thickness. Each of our species eventually ends up in a single stable steady state, there are no bifurcations nor oscillations in our system. So, we can freely vary the parameters found in literature to meet our experimental approaches and we can be sure that our system will still show the same global behaviour.

Initial States

The first step was to obtain the inital concentrations for all the species in our circuit. It is crucial to obtain biologically realistic values since we use them in the subsequent ODE solvers.These values describe the state of our cells at the start of the experiment, i.e. before connecting the channel to the acetaldehyde filled reservoir. When the cells are grown, they don't have any acetaldehyde in their environment and this affects the concentrations of intracellular species whose production depends on it.

Figure 1: Initial internal concentrations of cellular species

The obtained values are in line with our expectations based on the architecture of our engineered system. By analyzing the circuit we notice that when acetaldehyde is absent, the only proteins that are produced are TetR, LacI and LuxI, this being precisely what we see in the graph.

Sensors

The sensor module of our system detects toxic substance found in cigarette smoke. Upon binding of the substance to the sensor molecule, a series of downstream gene regulation events take place. We are dealing with two toxic substances (acetaldehyde and xylene). These two are detected by different sensor molecules with different mechanisms of action. Therefore, we designed two models, one for acetaldehyde detection, another for xylene detection. The models differ only by their sensor module, the band detector and the filter modules are same for both.

Acetaldehyde Sensor

As the name says, this sensor detects acetaldehyde from the air. Acetaldehyde binds to the transcriptional repressor AlcR, inhibiting TetR, which further inhibits CI and LacI production. Thus, the acetaldehyde binding to AlcR leads to CI and LacI expression (double inhibition equals activation).

There has been experimental evidence about the AlcR repression activity when it binds to acetaldehyde [5]. However, the exact mechanism of acetaldehyde binding to AlcR is not known. In Weber et al. [5] the acetaldehyde-inducible gene expression was assessed by measuring the expression of acetaldehyde adjustable SEAP for different acetaldehyde concentrations. We fitted the experimental data points for SEAP expression to the Hill equation and used the parameter values that were obtained from the fit to describe the repression kinetics of TetR. The ODE for TetR is the following:

Equation system 1: ODE for TetR

TetR gene is repressed by the AlcR-acetaldehyde complex (AlcR_AcA). The input acetaldehyde concentration is implicitly considered in the equation, due to the fitting procedure. The meaning and the values of the parameters are summarized in the Parameters section.

Figure 2: Repression of TetR by the AlcR-Acetaldehyde complex as a function of the input acetaldehyde concentration

We can see in the graph that TetR is expressed only at very small (close to zero) acetaldehyde concentrations, and it gets almost completely repressed already when the acetaldehyde concentration is 200 uM.

Xylene Sensor

The sensor that detects xylene is XylR. It is a hexamer which gets activated upon binding of xylene. Active XylR triggers CI and LacI expression.

At the beginning no measurement data was available, such that we modeled the xylene receptor activation based on the model published in [6]. Later during the course of creating SmoColi, we obtained own measurement data, based on which we redefined both the model structure as well as the parameter values. We here present both, our initial as well as our final model of xylene activation.

Initial Model Based on [6]

Our initial model incooperated both active and inactive Xylene (Xyla and XylI respectively). Xylene is found in two forms in this model: bound to XylR and free xylene. We consider the free xylene as a state (Xyl) and the initial xylene concentration as an input parameter (Xyl_ini). Thus, the bound xylene concentration can be found by substracting the free xylene concentration from the initial one (Xyl_ini - Xyl).

The ODEs for active and inactive XylR, as well as for xylene are given below and the parameter values are summarized in the Parameters section.

Equation system 2a: ODEs for the states involved in the xylene sensor. XylRi stands for inactive XylR, XylRa for active XylR, Xyl for unbound xylene

Model Based on own Measurement Data

Based on our measurement data (see Team:ETH_Zurich/Biology/Validation), we redefined our parameter values. In this process we also simplified the model as well as adapted it to our synthetic system. Different to the Equation System 2a, our new model assumes the concentration of xylene to be much higher than the concentration of XylR. Consequently, the concentration of free xylene can be considered to be approximately equal to the externally provided concentration (Xyl=Xyl_ini=const), and the last equation in the previous model disappears.

Furthermore, different to the system under consideration in [6], we could find no indication in the measurement data of our synthetic system for a special dissociation kinetic of the activated XylR (the Hill coefficient is nearly 1, see [6]), such that we assumed standard mass action association and dissociation kinetics. Finally we modeled both, activated and inactivated XylR, as hexamers, further simplifying the model structure.

Equation system 2b:ODEs for the states involved in the redefined xylene sensor. XylRi stands for inactive XylR, XylRa for active XylR, Xyl_ini for the provided xylene concentration

To compare the model to the measurement data, we assumed the model being in quasi steady state. In this case it is true that the measured value for the dissociation constant is equal to the ratio of the dissociation and the association rate constant:

K_D=r_RXylR / r_XylR=144μM.

Since we could not measure any temporal data with a resolution high enough to separate the value of the dissociation rate constant from the association rate constant, we kept the value obtained from [6] for the association rate constant and could then redefine the value for the dissociation rate constant based on our measurement data to

r_RXylR=144μM * 4*10^-5 1/μM min = 5.7*10^-3 1/min.

The rest of the parameter values were kept from the previous model (see parameters section).

Band Detector

The band detector is an amplitude band-pass filter that gives a GFP output when the toxic input substance is in a certain concentration range. The filter is a feed-forward system with 2 branches, interconnected by a NOR gate.

After the toxic substance is detected by the sensor module, the system splits into a slow branch and a fast branch.

• The fast branch is the amplitude low-pass filter. This means that it will repress the production of GFP when the toxic substance concentration goes above the high threshold. In other words, GFP is produced only at low toxic substance concentrations. High concentrations of toxic substance trigger LacIm1 expression, GFP production being thus inhibited.
• The slow branch is the amplitude high-pass filter which represses GFP when the toxic substance concentration falls below the low threshold or in other words GFP is produced only at high levels of toxic substance.

By combining these 2 filters, GFP is produced when neither LacI_M1 or LacI represses it, which is the region between the two thresholds. This region is the range of the toxic substance concentration which can trigger the GFP pulse. The ODEs for the species involved in the band detector are given below.

Equation system 3: ODEs for the states involved in the band detector module (CI, LacI and GFP)

We simulated the models for 5000 minutes in order to obtain the steady state GFP value. We repeated this simulation for a range of input acetaldehyde concentrations (from 0 to 5000 μM). In the following figures we can see how the steady state concentration of GFP changes for different input concentrations of acetaldehyde. For the acetaldehyde model the band is formed when the acetaldehyde concentration is approximately in the range of 500-1500 μM. The peak of the GFP pulse is observed when the acetaldehyde concentration is 1000 μM.

Figure 2: GFP vs. increasing acetaldehyde input (both in μM). The band appears when the acetaldehyde concentration is in the range of 500-1500 μM

Toxicity and parameter tuning

The acetaldehyde concentration range where the band is formed might be reasonable in our model, but not in reality. We tested the toxicity of acetaldehyde by adding 10 mM to a culture of E.coli cells and we observed that at this concentration the cells die. By taking into account the fact that a fraction of the acetaldehyde, due to its volatility, will diffuse out of the culture medium, we concluded that the maximum concentration the cells can tolerate is around 2 mM. This means that although in our model the band appears at a lower concentration: 500-1500 uM that does not affect the cells, the bacteria that are closer in the channel to the acetaldehyde source will still die.

Clearly, this tells us that we need to explore the parameter space for some of the constants to find suitable values. This is not surprising since for some of the parameters we made assumptions based on values for similar molecules. We carefully studied the equations involved in the band detector module and decided to make the following changes:

• TetR production rate (α_TetR): From the curve fitting [1], we got the TetR production rate to be 3.93 uM/min and we used this in our simulation. However, the measurements for which we did the fitting were done in mammalian cells, not in "E. coli". Since we express TetR on the same plasmid with other genes, for which we consider production rate to be 1 uM/min, we also changed α_TetR to 1 uM/min.
• TetR repression coefficient (β_TetR): We didn't have any supporting evidence for this, but we assumed it to be 0.01 uM. This was the value used in Basu et al. [1], but there however the repressor is LuxR, not TetR. We increased TetR repression coefficient to 0.1 uM since we assume that TetR is a less strong repressor than LuxR. With the modified value we get more realistic results.

With these 2 changes, the acetaldehyde concentration range where the band appears is much lower than in the previous case. We can see that the band forms between 25 μM and 100 μM of acetaldehyde input, with the peak at 46 μM (see Figure 4).

Figure 4:The band for acetaldehyde after parameter fine tuning is now shifted to the left. It appears when acetaldehyde concentration range is 30-100 uM

In both cases, the shape of the band looks similar, but it is approximately 10 times narrower in the case where we used the modified parameters. However, in both cases, its left side is steeper than its ride side, that is it appears more quickly than it dissapears.

Filter (Alarm)

This module acts like a high pass amplitude filter for the toxic input substance, since it produces RFP only when the toxic substance concentration goes above a certain threshold. At high acetaldehyde concentrations, TetR is repressed, so the repression of CI by TetR can not be established. At high xylene concentrations, XylR is activated and it triggers CI expression.

In both models, at high concentrations of the toxic substance CI is expressed and it inhibits the production of LuxI, which means that AHL is not produced. Once AHL is not there, it can not bind to LuxR, so LuxR can not repress RFP anymore. This results in RFP expression; all the cells turn red. As mentioned, this module of the model acts as a high pass amplitude filter for the toxic substance, or low pass amplitude filter for AHL, since AHL is produced only when the toxic substance levels are low.

The ODEs for the states involved in the second module of the system are given below:

Equation system 4: ODEs for the states involved in the Filter module. AHLi represents the intracellular AHL concentration, AHLe is the extracellular one. R represents the complex between LuxR and the intracellular AHL. Short explanation of the parameters is given in the Parameter section.

Figure 4: The intracellular AHL and the RFP concentration plotted against Acetaldehyde

The RFP steady state concentrations plotted against increasing toxic substance concentrations are presented in the above figures. For the acetaldehyde model, RFP appears after acetaldehyde crosses the threshold of 500 μM.

This plot shows the behavior of a singe cell. However, in our final system we have a population of cells immobilized in agarose that communicate whith eachother through the diffusing AHL. Thus, in the final setup the channel will turn red only when all the cells stop producing AHL and all the existing AHL is degraded in the channel.