Time-varying Toggle Switches in Living Cells: Modeling, Analysis, and Validation

The genetic toggle switch is a well-known model in synthetic biology that represents the dynamic interactions between two genes that mutually repress each other. In fact, the switching behavior in biological networks interests many researchers across a variety of disciplines beyond synthetic biology, including epigenetics and dynamical systems. Complex biological networks can implement mixtures of switching behavior and nonlinear smooth system response, which makes them a key factor in understanding how network structure and system logic collectively map structure to function. Since the toggle switch provides one of the simplest biological realizations of multistable behavior and exists as a subnetwork motif in many natural biological networks, it serves as a good starting point to showcase developments in gene network architecture [2]. With these facts in mind, we aimed to design and implement a time-varying genetic toggle switch with dynamics that change over the life cycle of the cell.

The main idea behind the classical genetic toggle switch in bacteria cells is that there are two genes (\(x\) and \(y\)) such that if either one is turned on, then the other one toggles off — hence the name. In this way, we say that \(x\) represses \(y\) and vice versa. The equations of the classical genetic toggle switch are as follows:

\[\dot{x} = \frac{\alpha}{1+y^{n_y}}-\delta_1 x\]

\[\dot{y} = \frac{\beta}{1+x^{n_x}}-\delta_2 y.\]

Here, \(x,y \geq 0\) respectively denote the concentrations of repressors 1 and 2; \(\alpha, \beta \geq 0\) signify the effective synthesis rates of repressors 1 and 2; \(\delta_1, \delta_2 \geq 0\) represent the decay rates of repressors 1 and 2; and \(n_x, n_y \geq 0\) are cooperativity coefficients. For simplicity, we consider all parameters and state variables to be dimensionless [2]. Figure 1 provides a phase portrait of the classical genetic toggle switch with pre-chosen parameters to illustrate the bistable behavior.

To develop a time-varying toggle switch, we had to first create a promoter that ensures the feasibility of such a system. In this context, a promoter is simply a DNA sequence that governs the expression rate of its corresponding gene. Different types of promoters can encode for gene activation while the cell population increases (exponential growth phase) or becomes saturated (stationary phase). By stitching together a growth phase and stationary phase promoter, we crafted a hybrid promoter that exhibits different types of gene expression throughout the life cycle of the cell population. We then generated several mutated variants of this promoter to produce a library of hybrid promoters, each of which has a distinct gene expression profile during the two targeted phases of the cell population’s life cycle. Next, we selected a pre-screened promoter from the library to elicit distinctive switching behavior in cells. Finally, we replaced a promoter in an existing toggle switch with our new, hybrid promoter to obtain a hybrid toggle switch that is given by the following time-varying model: 

\[\begin{equation}
 \begin{split}
    \dot{x} = \frac{\alpha(t)}{1+y^{n_y}}-\delta (t) x \\
    \dot{y} = \frac{\beta(t)}{1+x^{n_x}}-\gamma (t) y 
  \end{split}
\quad \left| \quad
  \begin{split}
     \alpha (t) &=
       \begin{cases}
         \alpha_1 & \mathbf{if} \quad t \in [0,\tau] \\
         \alpha_2 & \mathbf{if} \quad t \in (\tau, \infty]
       \end{cases} \\
    \beta (t) &=
      \begin{cases}
        \beta_1 & \mathbf{if} \quad t \in [0,\tau] \\
        \beta_2 & \mathbf{if} \quad t \in (\tau, \infty]
      \end{cases} \\
    \delta (t) &=
      \begin{cases}
        \delta_1 & \mathbf{if} \quad t \in [0,\tau] \\
        \delta_2 & \mathbf{if} \quad t \in (\tau, \infty]
      \end{cases} \\
    \gamma (t) &=
      \begin{cases}
        \gamma_1 & \mathbf{if} \quad t \in [0,\tau] \\
        \gamma_2 & \mathbf{if} \quad t \in (\tau, \infty].
      \end{cases}
  \end{split}
\right.
\end{equation}\]

In addition to constructing the experimental hybrid toggle switch, we also wanted to build predictive time series models for the genetic circuit. Existing mathematical models of the classical genetic toggle switch help to describe the gene-circuit dynamics in rapidly dividing cells, assuming fixed or time-invariant kinetic rates. Since we can biologically engineer promoters to have varying rates of synthesis as cells transition from one growth phase to another, then logically, we can model these rates with time-varying parameters. Thus, we have provided a novel class of parameter-varying nonlinear models to capture the dynamics of genetic circuits—including the toggle switch—as they transition between different phases of growth (see Figure 2) [3]. Using Koopman operatory theory [4, 5], we then performed stability analysis on this class of systems to ensure the biophysical feasibility of our new models [1].

Our time-varying models contain several parameters that can accurately fit the gene expression data from the hybrid toggle switch, unlike time-invariant models. To validate these models, we developed an algorithmic framework that both accurately fits time-varying parameters to our model and allows us to estimate the varying parameters of other models as well — including other genetic circuits that incorporate our hybrid promoter library. In the context of the time-varying toggle switch model, algorithms in our framework can estimate all of the varying \(\alpha\), \(\beta\), \(\delta\), and \(\gamma\) parameters and also identify and hold static parameters \(n_x\) and \(n_y\) as constant. 

To broadly overview our framework, let us define \(X_{data} \in \mathbb{R}^{m \times N}\) as the input data from our system of interest and \(X_{model} \in \mathbb{R}^{m \times N}\) as the predictive model data that approximates \(X_{data}\) upon the completion of an algorithm within our framework. These quantities are both matrices: \(m\) is the number of states and \(N\) is the number of data points (or snapshots) such that a row of either \(X_{data}\) or \(X_{model}\) is a state trajectory. 

As Figure 3 illustrates, the framework begins with the initialized continuous-time model, detects parameter switches in the input data, and then segments the data into \(n\) intervals. From there, it identifies model parameters \(p(t)\) and \(p_{static}\) to fit the model to the input data over each segmented time interval; specifically, it numerically integrates the segmented dynamic model to obtain \(X_{model}\) and then uses numerical optimization to minimize the distance (norm) between \(X_{data}\) and \(X_{model}\) where \(i \in \{1,...,n\}\). If the parameters are presupposed to be piecewise constant, then our algorithm outputs them as such. If the parameters are presupposed to be continuously varying functions, however, we then presuppose that a high number of parameter switches occur in the data and thus create highly segmented data intervals and parameters. The next step involves the use of dictionary-based sparse regression [6] to express the parameters \(p(t)\) as a linear combination of dictionary functions \(\mathfrak{D}\). Finally, we use the identified model parameters to generate a model trajectory \((X_{model})\) that accurately approximates \(X_{data}.\) Since this model trajectory is an approximation of the \(N\) sequential snapshots in \(X_{data}\), it may also be known as an \(N\)-step prediction. By providing the parametrized models for the hybrid toggle switches that accurately reconstruct their gene expression profiles, we have thus validated these models and their underlying constructions (see Figure 4).

In conclusion, our approach involves the strategic design of a hybrid promoter that causes shifts in gene expression dynamics throughout the cell population's life cycle. By amalgamating growth and stationary phase promoters, we introduced a library of hybrid promoters that are each characterized by a distinctive gene expression profile. We integrated one promoter from our vast library into the toggle switch, engendering an assemblage that we term the hybrid toggle switch; we then applied Koopman operator theory for stability analysis and the assurance of our model’s biophysical plausibility. Furthermore, we showcased the use of stability analysis in the time-varying Koopman framework. Finally, we implemented an algorithm for parameter estimation that can fit a hybrid toggle switch model to experimental data, and used this algorithm to parameterize the hybrid promoters in our library. Our algorithm is also generalizable for the estimation of varying parameters across a diverse set of models. In the intricate landscape of synthetic biology, our contributions serve as a testament to the perpetual dialogue between experimentation and theory, underscoring the need to unravel the inherent complexities of biological systems.

Jamiree Harrison delivered a minisymposium presentation on this research at the 2023 SIAM Conference on Computational Science and Engineering (CSE23), which took place in Amsterdam, the Netherlands, last year. He received funding to attend CSE23 through a SIAM Student Travel Award. To learn more about Student Travel Awards and submit an application, visit the online page. 

SIAM Student Travel Awards are made possible in part by the generous support of our community. To make a gift to the Student Travel Fund, visit the SIAM website.

References
[1] Budišić, M., Mohr, R., & Mezić, I. (2012). Applied Koopmanism. Chaos, 22(4), 047510.
[2] Gardner, T.S., Cantor, C.R., & Collins, J.J. (2000). Construction of a genetic toggle switch in Escherichia coli. Nature, 403(6767), 339-342.
[3] Harrison, J., & Yeung, E. (2021). Stability analysis of parameter varying genetic toggle switches using Koopman operators. Mathematics, 9(23), 3133.
[4] Mauroy, A., & Susuki, Y. (2018). Introduction to the Koopman operator in systems and control. In Proceedings of the 2018 57th annual conference of the Society of Instrument and Control Engineers of Japan (pp. 59-63). Nara, Japan: Institute of Electrical and Electronics Engineers.
[5] Mezic, I., & Surana, A. (2016). Koopman mode decomposition for periodic/quasi-periodic time dependence. IFAC-PapersOnLine, 49(18), 690-697. 
[6] Zhang, C., & Van Der Baan, M. (2022). Signal processing using dictionaries, atoms, and deep learning: A common analysis-synthesis framework. Proc. IEEE, 110(4), 454-475.