Light Operon

1. Introduction

A deterministic model of differential equations over a continuous time interval is generally easy to implement using known numerical methods. However, the equations alone do not account for the variability in expression that we often observe at low molecular count. In other words, deterministic simulations fail to capture the actual physical basis of the reaction (Gillespie 1977; Wilkinson 2012). This is due to the intrinsic and extrinsic noise that occurs at low intracellular molecular counts; interactions do not follow traditionally defined constant parameters which limits the accurate characterization of a genetic system (Elowitz et al., 2002).

It then becomes essential when characterizing a genetic system to have a deterministic approach (to capture the large average overtime) in addition to a stochastic approach. A stochastic approach allows one to model the probabilistic trajectory based on the initial conditions. One commonly used stochastic algorithm is the Gillespie Algorithm— a simple but powerful approach to simulation that takes into account the initial state of the reactants, the reaction rate and the number of molecules present for a reaction with each timed step of the reaction drawn from a probability distribution (Gillespie 1977).

As a working example, we will be working with the pDawn and pDusk operon (BBa K161609 and BBa K1075044). Many teams have made either the pDawn/pDusk system as their inducer (IONIS Paris 2015, Wageningen 2016, NUS Singapore 2017, Cornell 2017, Kingsborough 2017). In 2017, NUS Singapore 2017 and Kingsborough 2017 used the pDawn/pDusk inducer to induce expression of MazF to induce cell death. To our knowledge, there is no literature or past iGEM project exploring a stochastic model of the pDawn/pDusk system at low molecule count.

2. Methods

In order to improve reproducibility and to build on past iGEM teams’ efforts, we used the differential equations and parameters described in the 2016 Wageningen iGEM Team’s kill switch design. Using the “smfsb” package in R for stochastic simulations and the NDSolve function in Mathematica, we modeled the 3 light-sensitive states of the Yf1 homodimer over a period of 20 hours.

Equations, Constants, & Parameters

\(\frac{dy_{DD}}{dt}=k_1+ 2\cdot k_2 \cdot y_{DL,LD} - 2 \cdot (N\cdot k_3) \cdot y_{DD} - \beta_1 \cdot y_{DD}\)
\(\frac{dy_{DL,LD}}{dt}=2 \cdot (N \cdot k_3) \cdot y_{DD} + 2 \cdot k_2 \cdot y_{LL} - 2 \cdot k_2 \cdot y_{DL,LD} - 2 \cdot (N \cdot k_3) \cdot y_{DL,LD} - \beta_2 \cdot y_{DL,LD}\)
\(\frac{dy_{LL}}{dt} = 2 \cdot (N \cdot k_3) \cdot y_{DL,LD} - 2\cdot k_2 \cdot y_{LL} - \beta_3 \cdot y_{LL}\)

Constant/Parameter	Value	Description
\(N\)	Variable of \(\frac{\mu \cdot mol}{m^2 \cdot h}\)	Concentration of light
\(k_1\)	\(2.6921 \frac{\mu \cdot mol}{hr}\)	Production rate of \(y_{DD}\)
\(k_2\)	\(0.0008 \frac{1}{hr}\)	Relaxtion rate of \(y_{DL,LD}\) and \(y_{LL}\)
\(k_3\)	\(0.4219 \frac{m^2}{\mu \cdot mol}\)	Conversion cross-section of light intensity activated production rate of \(y_{DL,LD}\) and \(y_{LL}\).
\(\beta_1\)	\(0.3049 \frac{1}{hr}\)	Degradation rate of \(y_{DD}\)
\(\beta_2\)	\(0.8406 \frac{1}{hr}\)	Degradation rate of \(y_{DL,LD}\)
\(\beta_3\)	\(0.1477 \frac{1}{hr}\)	Degradation rate of \(y_{DD}\)

Post-Reaction Pre-Reaction Matrix for Gillespie

3. Results

The different states of Yf1 are being hit with different intensities of light, producing alternating results unique to the state. Click on any of the images to view full size.

In the stochastic results, we observe a much greater variation in expression at 95 \(\mu mol\) compared to 948 \(\mu mol\). Contrast this to the deterministic model in which the low/high concentration plots are identical with respect to a scaling factor. Gillespie causes step-behavior in the concentration. Although the plot of 948 \(\mu mol\) appears to be much smoother, we would see the same “cliffs” that we observe at 95 \(\mu mol\) by stretching the y-axis enough. The impact of a fixed gain or loss in concentration is much more “felt” in a small system than it would be in a large system (here size is relative to the initial concentration of Yf1). If we only cared about modeling a large, stable concentration however, the deterministic model clearly suffices to capture the general behavior and is less computationally expensive.

1. Starting concentration of 95 Yf1 molecules (Gillespie)

yDD (Dark-Dark state); yDL/LD (Dark-Light state); yLL (Light-Light state)

2. Starting concentration of 948 Yf1 molecules (Gillespie)