Team:Freiburg/Model

Modeling of the CARTELTM AND Gate

In synthetic biology, modeling can be applied to a wide range of topics, for example modeling of genetic circuits. In this subfield of mathematical and computational biology, ordinary differential equations (ODEs) are used to describe the transcriptional and translational processes over time, predict the behavior of the desired circuit and also to support their further development such as optimizing for a desired output, if enough data is available (Chen et al., 1999).

Schematic representation of CARTELTM AND gate inputs and output

Finding the CARTELTM AND Gate

The hypoxia response element (HRE), cAMP response element (CRE) and CTLA4 promoter can be utilized for a more specific CAR expression (Fig. 1). Combining these three enhancers in an AND gate to allow a specific CAR expression (Brophy et al., 2014) proved to be quite challenging. Using endogenous systems came with a drawback - a knockout or knockdown of either HIF1, VEGFR-2 or TDAG8 had to be generated.

Generating a single knockout or knockdown is challenging enough but generating three to test all possibilities seemed quite impossible in six months. Therefore using modeling to find the best combination of enhancers and knockout was an attractive option instead of choosing a random AND gate design.

To create such a model, the changing concentration of CAR can be described with a system of coupled ODEs. For each component of the genetic circuit - proteins and mRNAs - rate equations describing its synthesis, degradation, association or other processes are added (Chen et al., 1999). Below, it is shown how these ODEs can be set up in their simplest form.

f'(t)=a(z)-k1*f(t)

(1)

y'(t)=k2*f(t)-k3*y(t)

(2)

The variables are functions of time t, where f(t) describes the mRNA concentration, a(z) the transcription function, k1 the mRNA degradation rate, y(t) the protein concentration, k2 translation rate and k3 protein degradation rate. This kind of ODEs can be solved via numerical integration to describe changing protein and mRNA levels. To obtain the rate constants, experimental data have to be produced which is suitable for a nonlinear regression (Ingalls et al., 2012).

Equations (1) and (2) were further simplified for the AND gate comparison into equation (3) by neglecting the mRNA level to describe the delay of protein translation.

x'(t)=k4+k5*b(t)-k6*x(t)

(3)

The variables are functions of time t, where x(t) describes the protein concentration over time, k4 the basal expression rate, k5 the maximum expression rate, b(t) the activation function and k6 the degradation rate.

Regarding our intention to raise the specificity of cell killing, different possible designs were conceivable. These are accompanied by the introduction of certain genes or their elimination through a knockout. Focusing on the controlled CAR expression, an AND gate with a clear defined allocation of inputs and output had to be generated.

The first possible system was a knockout of HIF1A or HIF1B and its reintroduction via lentiviral transduction under the control of either CRE or CTLA4 promoter. Considering the law of mass action, a knockout of HIF1B would have a kinetic disadvantage over a knockout of HIF1A because HIF1B is accumulated permanently and independently from external conditions (Pescador et al., 2005) while HIF1A is only accumulated under hypoxia. Therefore only a knockout of HIF1A was considered a possible candidate for an AND gate.

The second or third possibility would be a knockout of either TDAG8 or VEGFR-2. Both can be described by the same model as both receptors are activated by one input and trigger activation of the associated promoter afterwards.

Unfortunately, experimental data only exist for HEK cells which are not representative for T cells (Ausländer et al., 2014; Nguyen et al., 2012). Therefore, the different AND gate designs were modeled with the same rate constants to compare their performances (equations in following drop down menus, Fig. 1, 2).

Figure 1: Schematic representation of an HIF1A knockout AND gate utilizing either low pH or high VEGF concentration as first input and low O2 concentration as second input to express the CAR

kbasal describes the basal expression rate, kmax the maximum expression rate, kdeg the degradation rate, Km gives the concentration of the activator were the promoters half activity is reached and the Hill coefficient n the activational slope of the given promoter. kdim and kdis always describes the dimerization and dissociation rate of HIF1AB.


Figure 2: Schematic representation of a VEGFR-2 or TDAG8 knockout AND gate utilizing low oxygen concentration as a first input and either low pH or high VEGF concentration as a second input to express the CAR

kbasal describes the basal expression rate, kmax the maximum expression rate, kdeg the degradation rate, Km gives the concentration of the activator were the promoters half activity is reached and the Hill coefficient n the activational slope of the given promoter. kdim and kdisalways describes the dimerization and dissociation rate of HIFAB

Figure 3: Modeling of knockout possibilities.
In the following (1) stands for high and (0) for low induction level. The first number represents the induction level via VEGF or pH and the second via oxygen concentration. The colors symbolize the following conditions: Red (1/1), Blue (1/0), Yellow (0/1), Green (0/0). a) Depiction of relative CAR expression under all four conditions over the course of 24 h in a HIFA knockout. It was assumed that at 0 h relative CAR expression is at 0 %. b) Depiction of relative CAR expression under all four conditions over the course of 24 h in a HIF1A knockout. It was assumed that at 0 h relative CAR expression is at 100 %. c) Depiction of relative CAR expression under all four conditions over the curse of 24 h in a TDAG8 of VEGFR-2 knockout. It was assumed that at 0 h relative CAR expression is at 0 %. d) Depiction of relative CAR expression under all four conditions over the course of 24 h in a TDAG8 of VEGFR-2 knockout. It was assumed that at 0 h relative CAR expression is at 100 %.

The in silico induction of the different knockouts possibilities revealed that a TDAG8 or VEGFR-2 knockout would have a higher relative CAR expression under the conditions of high VEGF and low oxygen concentration than a HIF1A knockout (Fig. 3). Therefore the HIF1A knockout was chosen.

After deciding for an AND gate design, we got in contact with the group of Professor Jens Timmer. Cooperating with Svenja Kemmer from his group, we gained insights into the methodology of processing and evaluating experimental data by modeling. In order to model our system of interest, the previously found data of CRE and HRE characterizations in HEK cells were used (Ausländer et al. 2014; Nguyen et al. 2012). Subsequently, parameters as kinetic rates were estimated by the Maximum Likelihood method and used to further analyze the system. The software used for the modeling of ODEs is data2dynamics, which is applied in Matlab and has been developed by the group of Professor Jens Timmer.

The composition of the model is described below in Figure 4.

Cell Culture

Figure 4: Schematic representation of the steps describing the input processing in the AND gate

The respective equations are equivalent to those of the first model except for an additional equation describing the HIF1A mRNA since the CRE characterization data showed delay in the processing.

The mathematical model is able to describe dynamics as well as dose response data. However, experimental measurements were only available for a few components of the system. Thus, detailed predictions require further investigations. Performing promoter characterizations and measuring degradations rates will enable more reliable predictions.

The respective response curves were fitted and combined in the AND gate model, thus describing the functionality of our system.

Additionally, the dependency of the HIF1A concentration on the pH and the oxygen concentration is shown in a heat map (Fig. 6). All data refer to a time of 24 hours after exposure to the respective conditions. The HIF1A concentration depends on both a low pH and hypoxia as it is present in the tumor microenvironment. Since it regulates the CAR expression, it can be estimated that the promoter for the CAR expression will only be activated if both conditions are fulfilled.

Heat map

Figure 5: Heat map
Heat map illustrating the model predictions of the HIF1A concentration in dependence of the pH and the oxygen concentration in the tumor microenvironment. Different colors describe the relative amounts of HIF1A specified in the sidebar. The marked areas show the characteristic conditions in healthy tissues as well as in the tumor microenvironment (TM). It indicates that the CAR expression is largely restricted to the tumor microenvironment.

The mathematical description of our AND gate leads to more ideas of how to improve the design of our logic gate. For example an improvement of the AND gate could be achieved by varying the number of enhancer elements. Depending on whether there is any cooperation in the binding to the enhancer elements and the leakiness of the promoter, different constructions are conceivable. At this point, modeling becomes even more important to create an encompassing overview about the different possibilities and simulate the system’s behaviour.

If you are interested in trying out different functions of data2dynamics applied to our AND gate model, you can download the code for it here. To run the model, also Matlab is required.

Effective Range of Activated CARTELTM T cells

We often received the question 'How far can a T cell migrate after being activated in a tumor microenvironment?'. We consulted this year’s iGEM team from ETH Zürich due to their modeling expertise to find an answer to this issue.