Design and Construction of a Versatile Synthetic Network for Bistable Gene Expression in Mammalian Systems

Abstract

We constructed and modeled a novel synthetic network which may be able to exhibit bistable expression of a reporter gene in mammalian cells. This network is based on an aptamer-fused short-hairpin RNA (shRNA) directed against a single mRNA encoding both a EGFP reporter gene and the repressor tTR-KRAB, which, in turn, represses transcription of the shRNA. The activity of the shRNA can be controlled by an inducer molecule (theophylline) which prevents the aptamer-fused shRNA to be properly processed. Repression of the tTR-KRAB can be relieved by treatment with doxycyline. This reciprocal negative feed-back loop can exhibit a bistable response, as shown through the mathematical analysis performed here. Specifically, the network can be controlled to induce sustained expression of a shRNA, or the reporter gene, with a transient input of two different inducer molecules.

1. Introduction

Biological “toggle-switches” are synthetic, or naturally occurring, gene regulatory networks, which enable sustained expression, or sustained repression, of a gene of interest. These two states, often referred to the ON and OFF state. To be a proper toggle-switch the ON, or OFF, state has to be maintained indefinitely and inherited by daughter cells. In a synthetic toggle-switch, the switch from the ON to the OFF state, or vice-versa, can be obtained by transient application of exogenous molecules. From a mathematical point of view, the properties that a network must have in order to work as a toggle-switch have been well characterized (Cherry and Adler, 2000).

The first synthetic “toggle-switch” was constructed in bacteria using two transcription factors reciprocally repressing each other's transcription, and whose activity could be induced via external inducer molecules (Gardner et al., 2000). More complex networks have been built in bacteria (Elowitz and Leibler, 2000) and in simple eukaryotes (Cantone et al., 2009), however, attempts to construct toggle switches in mammalian organisms have met with limited success. Only one synthetic toggle-switch in mammalian cells was shown to work by Kramer et al. (2004). This circuit is based on the original scheme proposed in bacteria, consisting of two TFs repressing each-other. This “toggle-switch” was constructed in a plasmid and stable cell-lines were obtained following transient transfection.

In this work, we constructed and partially characterized a new “toggle-switch” based on an original design consisting of a shRNA directed against a TF repressing the shRNA transcription. A previous work used intronically encoded siRNA towards the construction of a tighter and less prone to leakiness switch in mammalian cells (Greber et al., 2008). The major novelty of our design is the use of an aptamer fused shRNA, which is inducible with theophyllline, to silence transcription (An et al., 2006). The use of this molecule also allows for the whole circuit to be contained in a single lentiviral vector (Fig. 1), thus enabling its use both in vitro, such as in primary cells, and in vivo, in model organisms. The potential uses of such a toggle-switch are several, from understanding the function of disease-genes to innovative gene-therapy approaches.

FIG. 1.

Schematic diagram of the synthetic gene regulatory network for bistable gene expression.

2. Methods

2.1. Cell culture, transfection, and infection

CHO AA8 Tet-Off Cell Line (Clontech) stably expressing the tetracycline-controlled repressor (tTR-KRAB) were maintained at 37°C in a 5% CO²-humidified incubator. CHO cells were cultured in alpha-MEM (GIBCO BRL) supplemented with 10% heat-inactivated fetal bovine serum (FBS) (Invitrogen) and 1% antibiotic/antimycotic solution (GIBCO BRL). Cells were transfected 1 day after seeding using Lipofectamine 2000 (Invitrogen) according to manufacturer's instructions with 4μg of empty plasmid (pGEM-Z) or a plasmid containing the aptamer-fused shRNA. Samples were either treated with theophylline at a final concentration of 2nM, or without theophylline. One day post-transfection, the media and ligand were replaced and induction reiterated. Transfected cells were collected 60h post-transfection for RNA extraction and subsequent analysis. Experiments were performed in biological triplicates.

HEK 293 cells were maintained at 37°C in a 5% CO²-humidified incubator. HEK 293 cells were cultured in Dulbecco's modified Eagle's medium (DMEM, GIBCO BRL) supplemented with 10% heat-inactivated fetal bovine serum (FBS, Invitrogen) and 1% antibiotic/antimycotic solution (GIBCO BRL). A retro-virus expressing the circuit was produced in HEK293 FT cells using the ViraPower Promoterless Lentiviral Gateway Expression System (Invitrogen) according to manufacturer's instructions. Cells from 35-mm culture plates were transduced with the newly produced retrovirus, maintained in doxycycline (1 μg/mL) for 5 days and then trypsinized, filtered, and subjected to Fluorescence-Activated Cell Sorting (FACS) analysis in a Becton Dickinson FACSAria to select only positive clones expressing the circuit. The stable clone expressing the Tet-ON TTR-Krab based system was treated with varying concentrations of doxycycline ranging from an initial concentration of 1 ng/mL to 1 μg/mL; after five days, cells were trypsinized, filtered and subjected to FACS analysis to evaluate the total amount of fluorescence. Measurements were made on a sample of 20.000 events and every experiment was performed in biological triplicate.

2.2. Real-time polymerase chain reaction

Total RNA extraction from 35-mm culture plates was performed using the Qiagen RNeasy Kit (Qiagen) according to manufacturer's instructions. Retro-transcription of 1 μg of the total RNA extracted was performed using the QuantiTectReverse Transcription Kit (Qiagen), according to manufacturer's instructions.

Quantitative real-time PCR was performed using a LightCycler (Roche Molecular Biochemicals, Mannheim, Germany) to analyze the amplification status of tTA. Amplification of the gene was performed from the cDNA obtained from the total RNA and using the LightCycler DNA Master SYRB Green I kit (Roche Molecular Biochemicals). Primer sequences for Chinese Hampster GAPDH (used as reference gene) were designed by Primer 3.0 http://frodo.wi.mit.edu/ (forward primer: ACCCAGAAGACTGTGGATGG; reverse primer: GGATGCAGGGATGATGTTCT). Primer sequences for tTR-KRAB were also designed with Primer 3.0 (forward primer: ACAGCGCATTAGAGCTGCTT; reverse primer: ACCTAGCTTCTGGGCGAGTT). The relative amounts of genes were compared with the reference gene GAPDH and calculated using the Principle of Relative Quantification Analysis according to the standard formula 2^−DCt. To confirm the specificity of the amplification signal, we considered the primer dissociation curve in each case.

2.3. Construction of the aptamer-fused shRNA

A siRNA oligonucleotide targeting the TTR-Krab sequence was designed against the coding sequence of tetracycline-controlled repressor (tTR) gene using the Ambion technology platform. Custom siRNA HPLC purified with sequence GGUUUAACAACCCGUAAACtt (conc. 100 M) was synthesized by Ambion on the target sequence AAGGTTTAACAACCCGTAAAC starting at position 57 from the ATG in the tTA gene coding sequence. The aptamer sequence responsive to theophylline (An et al., 2006) was then artificially added using the facility of a gene synthesis company (GENEART) and the resulting aptamer-fused shRNA was cloned in an empty vector provided by the company, between the restriction sites ApaI (GGGCCC) and ClaI (ATCGAT). The fragment was cut and pasted between ApaI and ClaI within a subcloning vector containing a fragment of the lentiviral vector pLVUTHshGATA comprising the two restriction sites BSTXI (CCANNNNNNTGG) and FseI (GGCCGGCC) inserted after PCR amplification in a pGEM-Z (Promega) SmaI blunt site (CCCGGG). This construct, expressing an aptamer fused shRNA targeting the TTR-Krab sequence, was then used in transfection experiments of CHO cells stably expressing TTR-Krab to test the efficiency of silencing and theophylline induction. The circuit was then finalized via cut and paste of the fragment between the two restriction sites BSTXI and FseI in the subcloning vector into the same restriction sites of the lentiviral vector pLVUTHshGATA described in (Szulc et al., 2006).

3. Findings

3.1. Construction of the synthetic network

The network we constructed is schematically shown in Figures 1 and 2. It was obtained by modifying a lentiviral backbone for inducible shRNA expression reported in (Szulc et al., 2006). The network contains a Tet-repressor-KRAB protein (tTR-KRAB) under the control of the Ubiquitin promoter (UBC), which transcribes both EGFP and tTR-KRAB species as a single mRNA molecule. In addition a TRE-H1 promoter, containing the Tet Responsive Element (TRE), drives the expression of an aptamer-fused shRNA (An et al., 2006). We designed the shRNA to target the tTR-KRAB coding sequence, thus allowing, in the future, to replace the EGFP protein with any other gene without impairing the function of the switch.

FIG. 2.

Schematic diagram of the function of the gene regulatory network. (a) In the presence of doxycycline, the H1-TRE promoter and the UBC promoter are active and can produce both the shRNA, and the EGFP and tTR-KRAB proteins. When the shRNA production reaches high levels of concentration, doxycycline can be removed from the system because the shRNAs will keep the tTR-KRAB concentration low and the system will stay ON. (b) In the absence of doxycycline, the tTR-KRAB protein stops the transcription from both H1-TRE and UBC promoters. (c) The variables of the mathematical model of the synthetic network APU.

The network is reversibly controlled via the external treatment with doxycycline and/or theophylline. The switch transition is provoked either by theophylline, which inhibits the activity of the aptamer-fused shRNA on the degradation of the mRNA encoding for the EGFP and the tTR-KRAB, or by doxycycline, which relieves the repression by tTR-KRAB protein on the TRE-H1 promoter, thus allowing transcription of the shRNA. In the presence of doxycycline, tTR-KRAB is inhibited from binding to the promoter, hence the repressive epigenetic modification is relieved, and the shRNA is produced (Fig. 2a), along with the EGFP/tTR-KRAB mRNA. In the absence of doxycycline, the tTR-KRAB protein binds to the H1-TRE promoter and, by creating a local heterochromatic state extending over regions of a radius of 2–3 kilobases, it prevents transcription of both the shRNA and the EGFP/tTR-KRAB mRNA (Fig. 2b). We tested the ability of the tTR-KRAB to bind the H1-TRE promoter and stop transcription of both the shRNA and the EGFP/tTR-KRAB mRNA. To this end, we transduced HEK 293 cells with a virus produced from the lentiviral vector described in Szulc et al. (2006), which contains the tTR-KRAB controlling expression of a non-specific shRNA from the H1-TRE. Stable clones were obtained via cell sorting flow cytometry (FACS) and treated with different concentrations of doxycycline to induce EGFP expression as shown in Figure 3.

FIG. 3.

Characterization of the tTR-KRAB doxycycline inducible expression (Switch-ON, in the TetON-Krab system). Hek293 cells were infected with a cassette containing a constitutive promoter (UBC) expressing the tTR-KRAB. EGFP fluorescence was quantified by FACS analysis on 20,000 events in Hek293 cells stably expressing the Tet-ON plasmid (pLVUTHshGATA) following treatment for 5 days with the indicated concentration of the inducer molecule doxycycline.

These results confirm that the tTR-KRAB system is working as expected and allowed us to identify an appropriate concentration of doxycycline to achieve maximal activation. In Figure 6 below, we report the dynamic behavior of this system following removal of doxycycline or induction, respectively, which we obtained from Szulc et al. (2006). We used these data to fit the model's parameters as reported below.

The efficiency of theophylline on the newly synthesized aptamer-fused shRNA directed against the tTR-KRAB coding sequence was also verified. Figure 4 shows that we can achieve an excellent knock-down of tTR-KRAB expression in a stable clone expressing the protein (CHO cells) when we transfect with a plasmid expressing the shRNA, and restoration of expression in presence of the inducer theophylline.

FIG. 4.

Characterization of the novel aptamer-fused shRNA. Real-time data showing relative expression of tTR-KRAB in CHO cells stably expressing tTR-KRAB. (Left) Relative expression of tTR-Krab in cells transfected with an empty plasmid, respectively, nontreated (black) and treated with the inducer molecule theophylline (gray). (Right) Relative expression of tTR-Krab in cells transfected with the plasmid expressing the aptamer-fused shRNA, respectively, nontreated (black) and treated with the inducer molecule theophylline (gray).

We also built a lentiviral vector containing both constructs, i.e., the shRNA expressed by the H1-TRE promoter and the tTR-KRAB expressed by the UBC promoter, as shown in Figure 2 and produced a virus from it. We are currently testing the function of the complete circuit, to assess whether it exhibits a bistable behavior as expected.

3.2. Mathematical analysis

The mathematical model of the network is schematically illustrated in Figure 2c. The model presented here is a quantitative model based on ordinary differential equations (ODEs) (Polynikis et al., 2009). The model's parameters can be identified from the experimental data and, hence, reproduce the network behavior in silico. The parameters used for the simulations, unless otherwise specified, can be found in Figure 7a below.

The mRNA tTR-KRAB and the mRNA EGFP are represented by the same variable, X_mRNA, since they are both transcribed as a single molecule. The concentration of translated protein tTR-KRAB is denoted as P_TTR. A separate variable, P_ACT, represents the concentration of protein tTR-KRAB that is not bound by doxycycline and hence is active and can trigger the heterochromatin formation. Concentration of the shRNA is denoted by X_shRNA. \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$X_{s_free}$$\end{document} represents the free shRNAs that are not bound to theophylline. In addition, the time delay for forming the mature EGFP (which gives the green fluorescent color) by protein folding is modelled explicitly by two differential equations as in Tigges et al., (2009); one for the translation of mRNA to the unfolded EGFP protein, E_U, and one for the folded protein E_F. Finally, variables D and T stand for the concentrations of doxycycline and theophylline that are present in the cells, respectively. The ODE equations of the mathematical model of the circuit are presented in Figure 5. The parameters of the model are described in Table 1.

FIG. 5.

The ODE equations of the mathematical model of the circuit.

Table 1.

Parameters of the Mathematical Model of the Circuit

k_mRNA	Maximum transcription rate of UBC promoter
k_shRNA	Maximum transcription rate of H1-TRE promoter
k_T	mRNA translation rate
a_m	Basal activity of UBC promoter
a_s	Basal activity of H1-TRE promoter
d_mRNA	mRNA tTR-KRAB (or EGFP) degradation rate
d_TTR	Protein tTR-KRAB degradation rate
d_shRNA	shRNA EGFP degradation rate
d_EGFP	EGFP protein degradation rate
d_DOXY	Doxycycline degradation rate
d_THEO	Theophylline degradation rate
n_m, n_s	Hill coefficients for transcription functions
h	Hill coefficient for RNAi interaction
λ	Fold change in the mRNA degradation rate due to RNAi
k	Michaelis-Menten constant
θ_m	Affinity UBC promoter-inhibitor interaction
θ_s	Affinity HI-TRE promoter-inhibitor interaction
θ_D	Affinity doxycycline-tTR-KRAB interaction
θ_T	Affinity tetracycline-shRNA interaction
k_DIF	Diffusion rate of doxycycline and theophylline
k_FOLD	Fold rate of EGFP protein
D_EXT	External amount of doxycycline
T_EXT	External amount of theophylline

3.2.1. Fitting model's parameters

In order to test the efficiency of the model, the equations presented in Figure 5 were numerically fitted to the in vivo experimental results that were presented by Szulc et al. (2006) and reported in Figure 6.

FIG. 6.

Fitting of the model parameters to time-course experiments measuring EGFP fluorescence at different time points following addition (a) or removal (b) of doxycycline. (a) Numerical fitting of our mathematical model (equations 5) using the experimental results by Szulc et al. (2006). The cells were treated in doxycycline, which was removed after 5 days (day 0 in the x axis). (b) Numerical of fitting our mathematical model (equations 5) using the experimental results by Szulc et al. (2006). The cells were treated in the absence of doxycycline for 5 days, and then doxycycline was added to the cells (day 0 in the x-axis). The optimized values obtained are k_mRNA = 0.048, k_shRNA = 0.05, a_m = a_s = 0.001, k_T = 0.02, θ_m = θ_s = 0.001, n_m = n_s = 2.0, d_mRNA = 0.0173, d_shRNA = 0.054, d_TTR = 0.0231, d_EGFP = 4.4433 × 10⁻⁴, k_FOLD = 0.015, k_DIF = 0.002, θ_D = 0.02, d_DOXY = 6.418 × 10⁻⁴.

In these experiments, dynamic behavior due to the effects of doxycycline are reported. The difference between this circuit and ours, lies in the absence of the negative feedback due to RNAi. In fact in the citepd circuit, the shRNA produced acts against GATA1 (a zinc-finger transcription factor that plays an important role in gene regulation during development and differentiation of hematopoietic cells) and not against TTR as it is in our case. Nevertheless our model can be adapted to describe Szulc's circuit if one removes the negative feedback from RNAi. Namely in the mRNA equation, dXmRNA/dt (Fig. 5, first equation), we remove the last term. The numerical fitting of the model's parameters were performed using the Genetic Algorithm and Direct Search Toolbox of Matlab. The model was fitted to both data curves of experiments, with or without doxycycline experiments, and the optimization algorithm was minimizing the total cost error for the two curves.

The optimization results are illustrated in Figure 6. It can be seen that the model can capture well the experimental data.

3.2.2. Nullclines

In order for the circuit to be bistable, it must possess two (or more) stable steady states when the inducer drugs are absent. The complete model of the circuit is composed of seven differential equations, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$d \bar{X} / dt$$\end{document} , where \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$\bar{X} = ( X_{mRNA} , X_{shRNA} , P_{TTR} , E_{U} , E_F , D , T ) $$\end{document} . The equilibria of the system are defined by \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$d \bar{X} / dt = 0$$\end{document} . Since, we are looking for a bistability in the absence of inducer drugs, we consider dD/dt = dT/dt = 0. Therefore, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$X_{s_free} = X_{shRNA}$$\end{document} and P_ACT = P_TTR. Additionally, quantities E_U, E_F do not influence the bistability nature of the circuit since they are uncoupled from X_mRNA, X_shRNA, P_TTR.

Therefore, the bistability nature of the circuit will be determined by the three differential equations for X_mRNA, X_shRNA and P_TTR. Let \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}\begin{align*} \frac { dX_ { mRNA } } { dt } = \frac { dX_ { shRNA } } { dt } = \frac { dP_ { TTR } } { dt } = 0.\end{align*}\end{document}

From dP_TTR/dt = 0 one can easily derive that the equilibrium of the protein tTR-KRAB, \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$\tilde{P}_{TTR}$$\end{document} , is linearly proportional to the mRNA equilibrium \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$\tilde{X}_{mRNA}$$\end{document} , namely \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$\tilde{P}_{TTR} = k_T / d_{TTR} \tilde{X}_{mRNA}$$\end{document} . Therefore, one can now study the two-dimensional dynamical system defined by variables X_mRNA, X_shRNA.

Figure 7a illustrates the nullclines defined by dX_mRNA/dt = 0 and dX_shRNA/dt = 0 for a given set of parameters. The intersection of nullclines dX_mRNA/dt = 0, dX_shRNA/dt = 0 define the equilibria of the system. It can be seen that the nullclines intersect at three different points which determine two stable and one unstable equilibria as indicated in the plot. This confirms that there are parameter regions at which the circuit can indeed exhibit bistability. The parameter values used here, were obtained after parameter optimization of experimental data reported in Szulc et al. (2006).

FIG. 7.

Bifurcation analysis of the model. (a) The nullclines of mRNA and shRNA equilibriums of the model. The nullcline \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$\dot{X}_{mRNA} = 0$$\end{document} is shown in black color when its solution is real, and when the solution is complex we plot in dotted line the real part. The nullcline \documentclass{aastex}\usepackage{amsbsy}\usepackage{amsfonts}\usepackage{amssymb}\usepackage{bm}\usepackage{mathrsfs}\usepackage{pifont}\usepackage{stmaryrd}\usepackage{textcomp}\usepackage{portland, xspace}\usepackage{amsmath, amsxtra}\pagestyle{empty}\DeclareMathSizes{10}{9}{7}{6}\begin{document}$$\dot{X}_{shRNA} = 0$$\end{document} is always real for the specific parameter region and is shown in red color. Parameter values are k_mRNA = 0.048, k_shRNA = 0.05, a_m = a_s = 0.0001, k_T = 0.02, θ_m = θ_s = 0.0316, n_m = n_s = 2, d_mRNA = 0.0173, d_shRNA = 0.054, d_TTR = 0.0089, k = 0.4, λ = 4.5, h = 4. (b) Numerical continuation of the circuit's equilibria when the strength of the H1-TRE promoter, k_shRNA, is varied. λ = 6, h = 4, k_FOLD = 0.15, k_DIF = 0.003 (other parameters of the model are common with the ones in plot a). (c) Numerical continuation in two parameter space defined by k_mRNA and k_shRNA (the parameters of the model are common with the ones in plot b).

3.2.3. Bifurcation and continuation analysis

As shown in the previous subsection, the circuit has the potential to be bistable. However, it is of utmost importance to find parameter values that are close to the real biological system in order to understand if bistability can indeed be observed in real experiments. Additionally, the bistability of the circuit must be robust with respect to biological noise and parameter variations. For this reason, continuation and bifurcation analysis was used to analyse the relationship between the bistability properties and the parameter values of the circuit. For the common parameters of the model of our circuit and the circuit synthesized by Szulc et al. (2006), we chose values close to the ones obtained by optimization of Szulc's experiments as discussed in previous section. All numerical continuations are carried out using the continuation software package “Matcont” (Matcont, 2010).

Figure 7b plots the continuation of one of the equilibria of the circuit, as the parameter k_shRNA is varied. Recall that k_shRNA represents the maximum transcription strength of the H1-TRE promoter. We can observe that for a range of values of k_shRNA (see interval B), the concentration of X_shRNA has three possible equilibria, two stable and one unstable, which are connected via neutral saddle bifurcation points. Therefore, inside parameter region B, the circuit is bistable. When the H1-TRE promoter is too strong (see interval C), then the circuit has only one stable equilibrium, which corresponds to high concentration of shRNA (i.e., low EGFP). However, if H1-TRE promoter is very weak (see interval A), then the only equilibrium of the circuit will correspond to low concentration of shRNA (i.e., high EGFP).

Choosing a value of k_shRNA within the bistable region, we then continue the neutral saddle bifurcation point in two parameter space. As shown in Figure 7c, we follow the neutral saddle bifurcation point in the two parameter space defined by k_mRNA and k_shRNA. Region (B), is the region of bistability; for any pair of values of k_mRNA and k_shRNA that lie in region (B), the circuit has two stable equilibria. If H1-TRE promoter is weak (k_shRNA low), namely if we are in region (A), then the circuit has a unique equilibrium, which corresponds to low concentration of shRNA. The reason is that tTR-KRAB protein is able to inhibit transcription from both promoters. On the other hand, if H1-TRE is strong enough, and UBC promoter is weak (region C), then the circuit's unique equilibrium will correspond to high concentration of shRNA, which can inhibit the tTR-KRAB protein and prevent the formation of heterochromatin. The red line indicated in Figure 7c, is a one-dimensional section of the two-dimension bifurcation diagram. It corresponds to the one-dimensional bifurcation diagram illustrated in Figure 7b.

3.2.4. Transition between the two bistable states

Another important feature that the circuit must have, is the ability to switch between the two states. Figure 8 illustrates how this can be done experimentally by transient treatment with the inducer drugs doxycycline or theophylline.

FIG. 8.

Characterization of model's bistability. (a) Numerical continuation of shRNA equilibrium as the concentration variable of the external input of doxycycline, D_EXT, is varied. Other parameter values are k_mRNA = 0.048, k_shRNA = 0.05, a_m = a_s = 0.0001, k_T = 0.02, θ_m = θ_s = 0.0316, n_m = n_s = 2, d_mRNA = 0.0173, d_shRNA = 0.054, d_TTR =0.0089, k = 0.4, d_EGFP = 4.4433 × 10⁻⁴, k_FOLD = 0.015, k_DIF = 0.003, d_DOXY = 6.418 × 10⁻⁴, d_THEO = 6.418 × 10⁻⁴, T_EXT = 0. Circles represent neutral saddle bifurcation points. (b) Numerical continuation of shRNA equilibrium as the concentration variable of the external input of theophylline, T_EXT, is varied. D_EXT = 0 (all the other parameters are the same as in plot a). Circles represent neutral saddle bifurcation points. (c) In silico time course of switching between the two states of the circuit. At day 0, the system is on the OFF state (shRNA low) and doxycycline is added in the system. At day 2, doxycycline is removed. At day 4, Theophylline is added to the system, and is removed at day 6.

For example, if the circuit is at the equilibrium state where shRNA production is OFF, then if a high enough concentration of doxycycline is added to the system, the circuit will switch to the state where shRNA production is ON (Fig. 8a). It can be seen that if doxycycline concentration is high enough, then the high concentration of shRNAs will be able to prevent the production of tTR-KRAB protein and hence alleviate the epigenetic repression.

Additionally, as shown in Figure 8b, if the system is at equilibrium with high concentration of shRNA (shRNA ON), then if enough theophylline is added to the system, the circuit will switch onto the second equilibrium (shRNA OFF). For example, for the specific parameter region in Figure 8b, if theophylline concentration is greater than 0.031, then the system jumps to the second equilibrium. The reason is that, tTR-KRAB is now produced and can inhibit both promoters by triggering the heterochromatin formation which stops their transcription.

Another important issue, especially for experimental applications, is to know for how long cells must be treated with doxycycline or theophylline, in order to force the circuit to switch from one state to the other. Additionally, it is important to ensure that with the removal of the drug, the circuit will stay on the desired equilibrium. Figure 8c illustrates the in silico time course of switching between the two states. At day 0, the system is at the OFF state (shRNA low) and doxycycline is added into the system. We observed, that two days of treatment with doxycycline is enough to turn the production of shRNA on. Additionally, it is shown that the system remains on after the removal of doxycycline (days 2–4). To switch to the other the state, theophylline is added into the system at day 4. Treatment for two days is enough to switch to the OFF state where shRNAs are not produced. This stays OFF, even without the presence of theophylline, which is removed at day 6.

4. Conclusion

We designed, modeled, and constructed a new synthetic network in a lentiviral vector that in principle could behave as a bistable “toggle-switch” in mammalian systems. We tested the correct function of both the tTR-KRAB inducible repressor, and of the new aptamer-fused shRNA directed against the tTR-KRAB coding sequence. The major novelty of our work is the use of an aptamer-fused shRNA to obtain the reciprocal inhibition required to have a toggle-switch (Cherry and Adler, 2000). This allows the whole network to be embedded in a lentiviral vector which greatly expands its potential usefulness both in systems biology and for gene therapy applications. Through mathematical analysis, we were able to make in silico predictions for the behaviour of a this new synthetic network. We observed that there are realistic parameter regions in which the network can indeed exhibit bistability. Namely, the network can behave as a toggle-switch, which can be controlled by external inducer molecules, which can switch ON and OFF the expression of a reporter gene. We are currently testing whether this network exhibits a bistable behaviour in vitro and in vivo, as predicted by the mathematical model.

Footnotes

Acknowledgments

We would like to thank Dr. Yokobaiashi for sending the aptamer-fused shRNA targeting EGFP and for technical advice, and Dr. B-Cienz for great technical support. This work was supported by the ITALBIONET 2007–2010 grant from the Italian Ministry of Research to D.d.B.

Disclosure Statement

No competing financial interests exist.

References

, Trinh

V.B.

, Yokobayashi

2006. Artificial control of gene expression in mammalian cells by modulating RNA interference through aptamersmall molecule interaction. RNA, 12:710–716.

Cantone

, Marucci

, Iorio

et al. 2009. A yeast synthetic network for in vivo assessment of reverse-engineering and modeling approaches. Cell, 137:172–181.

Cherry

J.L.

, Adler

F.R.

2000. How to make a biological switch. J. Theor. Biol, 203:117–133.

Elowitz

M.B.

, Leibler

2000. A synthetic oscillatory network of transcriptional regulators. Nature, 403:335–338.

Gardner

T.S.

, Cantor

C.R.

, Collins

J.J.

2000. Construction of a genetic toggle switch in Escherichia coli. Nature, 403:339–342.

Greber

, El-Baba

M.D.

, Fussenegger

2008. Intronically encoded sirnas improve dynamic range of mammalian gene regulation systems and toggle switch. NAR, 36:1–16.

Kramer

B.P.

, Viretta

A.U.

, Daoud-El Baba

et al. 2004. An engineered epigenetic transgene switch in mammalian cells. Nat. Biotechnol., 22:867–870.

Matcont . 2010. Continuation software in Matlab. www.matcont.ugent.be/. 2010 November 15.

Polynikis

, Hogan

S.J.

, Di Bernardo

2009. Comparing different ODE modelling approaches for gene regulatory networks. J. Theor. Biol., 261:511–530.

10.

Szulc

, Wiznerowicz

, Sauvain

M.O.

et al. 2006. A versatile tool for conditional gene expression and knockdown. Nat. Methods, 3:119–116.

11.

Tigges

, Marquez-Lago

T.T.

, Stelling

et al. 2009. A tunable synthetic mammalian oscillator. Nature, 457:309–312.