On physiological multiplicity and population heterogeneity of biological systems

Inducer exclusion, by itself, cannot account for the glucose-mediated lac repression of Escherichia coli

The lac operon of Escherichia coli is repressed several 100-fold in the presence of glucose. This repression has been attributed to cAMP receptor protein-mediated inhibition of lac transcription and EIIA^Glc-mediated inhibition of lactose transport (inducer exclusion). The growing evidence against the first mechanism has led to the postulate that the repression is driven by inducer exclusion. Although inducer exclusion reduces the permease activity only 2-fold in fully induced cells, it could be more potent in partially induced cells. Here, we show that even in partially induced cells, inducer exclusion reduces the permease activity no more than 6-fold. Moreover, the repression is so small because these experiments are performed in the presence of chloramphenicol. Indeed, when glucose is added to a culture growing on glycerol and TMG, but no chloramphenicol, lac expression is repressed 900-fold. This repression is primarily due to reversal of the positive feedback loop, i.e., the decline of the intracellular TMG level leads to a lower permease level, which reduces the intracellular TMG level even further. The repression in the absence of chloramphenicol is therefore primarily due to positive feedback, which does not exist during measurements of inducer exclusion.

An entropy-like index of bifurcational robustness for metabolic systems

Natural and synthetic metabolic pathways need to retain stability when faced against random changes in gene expression levels and kinetic parameters. In the presence of large parameter changes, a robust system should specifically avoid moving to an unstable region, an event that would dramatically change system behavior. Here we present an entropy-like index, denoted as S, for quantifying the bifurcational robustness of metabolic systems against loss of stability. We show that S enables the optimization of a metabolic model with respect to both bifurcational robustness and experimental data. We then demonstrate how the coupling of ensemble modeling and S enables us to discriminate alternative designs of a synthetic pathway according to bifurcational robustness. Finally, we show that S enables the identification of a key enzyme contributing to the bifurcational robustness of yeast glycolysis. The different applications of S demonstrated illustrate the versatile role it can play in constructing better metabolic models and designing functional non-native pathways.

Kinetic approaches to lactose operon induction and bimodality

The quasi-equilibrium approximation is acceptable when molecular interactions are fast enough compared to circuit dynamics, but is no longer allowed when cellular activities are governed by rare events. A typical example is the lactose operon (lac), one of the most famous paradigms of transcription regulation, for which several theories still coexist to describe its behaviors. The lac system is generally analyzed by using equilibrium constants, contradicting single-event hypotheses long suggested by Novick and Weiner (1957). Enzyme induction as an all-or-none phenomenon. Proc. Natl. Acad. Sci. USA 43, 553-566) and recently refined in the study of (Choi et al., 2008. A stochastic single-molecule event triggers phenotype switching of a bacterial cell. Science 322, 442-446). In the present report, a lac repressor (LacI)-mediated DNA immunoprecipitation experiment reveals that the natural LacI-lac DNA complex built in vivo is extremely tight and long-lived compared to the time scale of lac expression dynamics, which could functionally disconnect the abortive expression bursts and forbid using the standard modes of lac bistability. As alternatives, purely kinetic mechanisms are examined for their capacity to restrict induction through: (i) widely scattered derepression related to the arrival time variance of a predominantly backward asymmetric random walk and (ii) an induction threshold arising in a single window of derepression without recourse to nonlinear multimeric binding and Hill functions. Considering the complete disengagement of the lac repressor from the lac promoter as the probabilistic consequence of a transient stepwise mechanism, is sufficient to explain the sigmoidal lac responses as functions of time and of inducer concentration. This sigmoidal shape can be misleadingly interpreted as a phenomenon of equilibrium cooperativity classically used to explain bistability, but which has been reported to be weak in this system.

Stochastic de-repression of Rhodopsins in single photoreceptors of the fly retina

The photoreceptors of the Drosophila compound eye are a classical model for studying cell fate specification. Photoreceptors (PRs) are organized in bundles of eight cells with two major types – inner PRs involved in color vision and outer PRs involved in motion detection. In wild type flies, most PRs express a single type of Rhodopsin (Rh): inner PRs express either Rh3, Rh4, Rh5 or Rh6 and outer PRs express Rh1. In outer PRs, the K₅₀ homeodomain protein Dve is a key repressor that acts to ensure exclusive Rh expression. Loss of Dve results in de-repression of Rhodopsins in outer PRs, and leads to a wide distribution of expression levels. To quantify these effects, we introduce an automated image analysis method to measure Rhodopsin levels at the single cell level in 3D confocal stacks. Our sensitive methodology reveals cell-specific differences in Rhodopsin distributions among the outer PRs, observed over a developmental time course. We show that Rhodopsin distributions are consistent with a two-state model of gene expression, in which cells can be in either high or basal states of Rhodopsin production. Our model identifies a significant role of post-transcriptional regulation in establishing the two distinct states. The timescale for interconversion between basal and high states is shown to be on the order of days. Our results indicate that even in the absence of Dve, the Rhodopsin regulatory network can maintain highly stable states. We propose that the role of Dve in outer PRs is to buffer against rare fluctuations in this network.

Complex networks of genetic interactions govern the development of multicellular organisms. One of the best-characterized networks governs the development of the fruit-fly retina, a highly organized, three-dimensional organ composed of a hexagonal grid of eight types of photoreceptor neurons. Each photoreceptor responds to a particular wavelength of light depending on the Rhodopsin protein it expresses. We present novel computational methods to quantify cell-specific Rhodopsin levels from confocal microscopy images. We apply these methods to study the effect of the loss of a key repressor that ensures each photoreceptor expresses only one Rhodopsin. We show that this perturbation has cell-specific effects. Our measurement of the cell-type specific Rhodopsin distributions reveals differences between photoreceptor cells, which could not otherwise be detected. Using mathematical models of gene expression, we attribute this variability to stochastic events that activate Rhodopsin production.

Dynamic analysis of the KlGAL regulatory system in Kluyveromyces lactis: a comparative study with Saccharomyces cerevisiae

The GAL regulatory system is highly conserved in yeast species of Saccharomyces cerevisiae and Kluyveromyces lactis. While the GAL system is a well studied system in S. cerevisiae, the dynamic behavior of the KlGAL system in K. lactis has not been characterized. Here, we have characterized the GAL system in yeast K. lactis by developing a dynamic model and comparing its performance to its not-so-distant cousin S. cerevisiae. The present analysis demonstrates the significance of the autoregulatory feedbacks due to KlGal4p, KlGal80p, KlGal1p and Lac12p on the dynamic performance of the KlGAL switch. The model predicts the experimentally observed absence of bistability in the wild type strain of K. lactis, unlike the short term memory of preculturing conditions observed in S. cerevisiae. The performance of the GAL switch is distinct for the two yeast species although they share similarities in the molecular components. The analysis suggests that the whole genome duplication of S. cerevisiae, which resulted in a dedicated inducer protein, Gal3p, may be responsible for the high sensitivity of the system to galactose concentrations. On the other hand, K. lactis uses a bifunctional protein as an inducer in addition to its galactokinase activity, which restricts its regulatory role and hence higher galactose levels in the medium are needed to trigger the GAL system.

ELECTRONIC SUPPLEMENTARY MATERIAL

The online version of this article (doi:10.1007/s11693-011-9082-7) contains supplementary material, which is available to authorized users.

Design of the lac gene circuit revisited

The lactose (lac) operon of Escherichia coli serves as the paradigm for gene regulation, not only for bacteria, but also for all biological systems from simple phage to humans. The details of the systems may differ, but the key conceptual framework remains, and the original system continues to reveal deeper insights with continued experimental and theoretical study. Nearly as long lasting in impact as the pivotal work of Jacob and Monod is the classic experiment of Novick and Weiner in which they demonstrated all-or-none gene expression in response to an artificial inducer. These results are often cited in claims that normal gene expression is in fact a discontinuous bistable phenomenon. In this paper, I review several levels of analysis of the lac system and introduce another perspective based on the construction of the system design space. These represent variations on a theme, based on a simply stated design principle, that captures the key qualitative features of the system in a largely mechanism-independent fashion. Moreover, this principle can be readily interpreted in terms of specific mechanisms to make predictions regarding monostable vs. bistable behavior. The regions of design space representing bifurcations are compared with the corresponding regions identified through bifurcation analysis. I present evidence based on biological considerations as well as modeling and analysis to suggest that induction of the lac system in its natural setting is a monostable continuously graded phenomenon. Nevertheless, it must be acknowledged that the lac stability question remains unsettled, and it undoubtedly will remain so until there are definitive experimental results.

Noise contributions in an inducible genetic switch: a whole-cell simulation study

Stochastic expression of genes produces heterogeneity in clonal populations of bacteria under identical conditions. We analyze and compare the behavior of the inducible lac genetic switch using well-stirred and spatially resolved simulations for Escherichia coli cells modeled under fast and slow-growth conditions. Our new kinetic model describing the switching of the lac operon from one phenotype to the other incorporates parameters obtained from recently published in vivo single-molecule fluorescence experiments along with in vitro rate constants. For the well-stirred system, investigation of the intrinsic noise in the circuit as a function of the inducer concentration and in the presence/absence of the feedback mechanism reveals that the noise peaks near the switching threshold. Applying maximum likelihood estimation, we show that the analytic two-state model of gene expression can be used to extract stochastic rates from the simulation data. The simulations also provide mRNA–protein probability landscapes, which demonstrate that switching is the result of crossing both mRNA and protein thresholds. Using cryoelectron tomography of an E. coli cell and data from proteomics studies, we construct spatial in vivo models of cells and quantify the noise contributions and effects on repressor rebinding due to cell structure and crowding in the cytoplasm. Compared to systems without spatial heterogeneity, the model for the fast-growth cells predicts a slight decrease in the overall noise and an increase in the repressors rebinding rate due to anomalous subdiffusion. The tomograms for E. coli grown under slow-growth conditions identify the positions of the ribosomes and the condensed nucleoid. The smaller slow-growth cells have increased mRNA localization and a larger internal inducer concentration, leading to a significant decrease in the lifetime of the repressor–operator complex and an increase in the frequency of transcriptional bursts.

Expressing genes in a bacterial cell is noisy and random. A colony of bacteria grown from a single cell can show remarkable differences in the copy number per cell of a given protein after only a few generations. In this work we use computer simulations to study the variation in how individual cells in a population express a set of genes in response to an environmental signal. The modeled system is the lac genetic switch that Escherichia coli uses to find, collect, and process lactose sugar from the environment. The noise inherent in the genetic circuit controlling the cell's response determines how similar the cells are to each other and we study how the different components of the circuit affect this noise. Furthermore, an estimated 30–50% of the cell volume is taken up by a wide variety of large biomolecules. To study the response of the circuit caused by crowding, we simulate the circuit inside a three-dimensional model of an E. coli cell built using data from cryoelectron tomography reconstructions of a single cell and proteomics data. Correctly including random effects of molecular crowding will be critical to developing fully dynamic models of living cells.

Investigating heterogeneity of intracellular calcium dynamics in anterior pituitary lactotrophs using a combined modelling/experimental approach

Cell responses are commonly heterogeneous, even within a subpopulation. In the present study, we investigate the source of heterogeneity in the Ca(2+) response of anterior pituitary lactotrophs to a Ca(2+) mobilisation agonist, thyrotrophin-releasing hormone. This response is characterised by a sharp increase of cytosolic Ca(2+) concentration as a result of mobilisation of Ca(2+) from intracellular stores, followed by a decrease to an elevated plateau level that results from Ca(2+) influx. We focus on heterogeneity of the evoked Ca(2+) spike under extracellular Ca(2+) free conditions. We introduce a method that uses the information provided by a mathematical model to characterise the source of heterogeneity. This method compares scatter plots of features of the Ca(2+) response obtained experimentally with those made from the mathematical model. The model scatter plots reflect random variation of parameters over different ranges, and matching the experimental and model scatter plots allows us to predict which parameters are most variable. We find that a large degree of variation in Ca(2+) efflux is a likely key contributor to the heterogeneity of Ca(2+) responses to thyrotrophin-releasing hormone in lactotrophs. This technique is applicable to any situation in which the heterogeneous biological response is described by a mathematical model.

Comparison of deterministic and stochastic models of the lac operon genetic network

The lac operon has been a paradigm for genetic regulation with positive feedback, and several modeling studies have described its dynamics at various levels of detail. However, it has not yet been analyzed how stochasticity can enrich the system's behavior, creating effects that are not observed in the deterministic case. To address this problem we use a comparative approach. We develop a reaction network for the dynamics of the lac operon genetic switch and derive corresponding deterministic and stochastic models that incorporate biological details. We then analyze the effects of key biomolecular mechanisms, such as promoter strength and binding affinities, on the behavior of the models. No assumptions or approximations are made when building the models other than those utilized in the reaction network. Thus, we are able to carry out a meaningful comparison between the predictions of the two models to demonstrate genuine effects of stochasticity. Such a comparison reveals that in the presence of stochasticity, certain biomolecular mechanisms can profoundly influence the region where the system exhibits bistability, a key characteristic of the lac operon dynamics. For these cases, the temporal asymptotic behavior of the deterministic model remains unchanged, indicating a role of stochasticity in modulating the behavior of the system.

Bistability of the lac operon during growth of Escherichia coli on lactose and lactose+glucose

The lac operon of Escherichia coli can exhibit bistability. Early studies showed that bistability occurs during growth on TMG/succinate and lactose+glucose, but not during growth on lactose. More recently, studies with lacGFP-transfected cells show bistability during growth on TMG/succinate, but not during growth on lactose and lactose+glucose. In the literature, these results are invariably attributed to variations in the destabilizing effect of the positive feedback generated by induction. Specifically, during growth on TMG/succinate, lac induction generates strong positive feedback because the permease stimulates the accumulation of intracellular TMG, which in turn, promotes the synthesis of even more permease. This positive feedback is attenuated during growth on lactose because hydrolysis of intracellular lactose by beta-galactosidase suppresses the stimulatory effect of the permease. It is attenuated even more during growth on lactose + glucose because glucose inhibits the uptake of lactose. But it is clear that the stabilizing effect of dilution also changes dramatically as a function of the medium composition. For instance, during growth on TMG/succinate, the dilution rate of lac permease is proportional to its activity, e, because the specific growth rate is independent of e (it is completely determined by the concentration of succinate). However, during growth on lactose, the dilution rate of the permease is proportional to e2 because the specific growth rate is proportional to the specific lactose uptake rate, which in turn, proportional to e. We show that: (a) This dependence on e2 creates such a strong stabilizing effect that bistability is virtually impossible during growth on lactose, even in the face of the intense positive feedback generated by induction. (b) This stabilizing effect is weakened during growth on lactose+glucose because the specific growth rate on glucose is independent of e, so that the dilution rate once again contains a term that is proportional to e. These results imply that the lac operon is much more prone to bistability if the medium contains carbon sources that cannot be metabolized by the lac enzymes, e.g., succinate during growth on TMG/succinate and glucose during growth on lactose+glucose. We discuss the experimental data in the light of these results.

Bistable behavior in a model of the lac operon in Escherichia coli with variable growth rate

This work is a continuation from another study previously published in this journal. Both the former and the present works are dedicated to investigating the bistable behavior of the lac operon in Escherichia coli from a mathematical modeling point of view. In the previous article, we developed a detailed mathematical model that accounts for all of the known regulatory mechanisms in this system, and studied the effect of inducing the operon with lactose instead of an artificial inducer. In this article, the model is improved to account, in a more detailed way, for the interaction of the repressor molecules with the three lac operators. A recently discovered cooperative interaction between the CAP molecule (an activator of the lactose operon) and Operator 3 (which influences DNA folding) is also included in this new version of the model. The growth rate dependence on the rate of energy entering the bacteria (in the form of transported glucose molecules and of metabolized lactose molecules) is also considered. A large number of numerical experiments is carried out with this improved model. The results are discussed in regard to the bistable behavior of the lactose operon. Special attention is paid to the effect that a variable growth rate has on the system dynamics.

Steady state approach to model gene regulatory networks—Simulation of microarray experiments

Genetic regulatory networks (GRN) represent complex interactions between genes brought about through proteins that they code for. Quantification of expression levels in GRN either through experiments or theoretical modeling is a challenging task. Recently, microarray experiments have gained importance in evaluating GRN at the genome level. Microarray experiments yield log fold change in mRNA abundance which is helpful in deciphering connectivity in GRN. Current approaches such as data mining, Boolean or Bayesian modeling and combined use of expression and location data are useful in analyzing microarray data. However, these methodologies lack underlying mechanistic details present in GRN. We present here a steady state gene expression simulator (SSGES) which sets up steady state equations and simulates the response for a given network structure of a GRN. SSGES includes mechanistic details such as stoichiometry, protein-DNA and protein-protein interactions, translocation of regulatory proteins and autoregulation. SSGES can be used to simulate the response of a GRN in terms of fractional transcription and protein expression. SSGES can also be used to generate log fold change in mRNA abundance and protein expression implying that it is useful to simulate microarray type experiments. We have demonstrated these capabilities of SSGES by modeling the steady state response of GAL regulatory system in Saccharomyces cerevisiae. We have demonstrated that the predicted data qualitatively matched the microarray data obtained experimentally by Ideker et al. [Ideker, T., Thorsson, V., Ranish, J.A., Christmas, R., Buhler, J., Eng, J.K., Bumgarner, R., Goodlett, D.R., Aebersold, R., Hood, L., 2001. Integrated genomic and proteomic analyses of a systematically perturbed metabolic network. Science 292, 929-934]. SSGES is available from authors upon request.

Effect of DNA looping on the induction kinetics of the lac operon

The induction of the lac operon follows cooperative kinetics. The first mechanistic model of these kinetics is the de facto standard in the modeling literature [Yagil, G., Yagil, E., 1971. On the relation between effector concentration and the rate of induced enzyme synthesis. Biophys. J. 11, 11-17]. Yet, subsequent studies have shown that the model is based on incorrect assumptions. Specifically, the repressor is a tetramer with four (not two) inducer-binding sites, and the operon contains two auxiliary operators (in addition to the main operator). Furthermore, these structural features are crucial for the formation of DNA loops, the key determinants of lac repression and induction. Indeed, the repression is determined almost entirely (>95%) by the looped complexes [Oehler, S., Eismann, E.R., Krämer, H., Müller-Hill, B., 1990. The three operators of the lac operon cooperate in repression. EMBO J. 9(4), 973-979], and the pronounced cooperativity of the induction curve hinges upon the existence of the looped complexes [Oehler, S., Alberti, S., Müller-Hill, B., 2006. Induction of the lac promoter in the absence of DNA loops and the stoichiometry of induction. Nucleic Acids Res. 34(2), 606-612]. Here, we formulate a model of lac induction taking due account of the tetrameric structure of the repressor and the existence of looped complexes. We show that: (1) The kinetics are significantly more cooperative than those predicted by the Yagil and Yagil model. The cooperativity is higher because the formation of looped complexes is easily abolished by repressor-inducer binding. (2) The model provides good fits to the repression data for cells containing wild-type tetrameric or mutant dimeric repressor, as well as the induction curves for 6 different strains of Escherichia coli. It also implies that the ratios of certain looped and non-looped complexes are independent of inducer and repressor levels, a conclusion that can be rigorously tested by gel electrophoresis. (3) Repressor overexpression dramatically increases the cooperativity of the induction curve. This suggests that repressor overexpression can induce bistability in systems, such as growth of E. coli on lactose, that are otherwise monostable.

Origin of bistability in the lac Operon

Multistability is an emergent dynamic property that has been invoked to explain multiple coexisting biological states. In this work, we investigate the origin of bistability in the lac operon. To do this, we develop a mathematical model for the regulatory pathway in this system and compare the model predictions with other experimental results in which a nonmetabolizable inducer was employed. We investigate the effect of lactose metabolism using this model, and show that it greatly modifies the bistable region in the external lactose (Le) versus external glucose (Ge) parameter space. The model also predicts that lactose metabolism can cause bistability to disappear for very low Ge. We have also carried out stochastic numerical simulations of the model for several values of Ge and Le. Our results indicate that bistability can help guarantee that Escherichia coli consumes glucose and lactose in the most efficient possible way. Namely, the lac operon is induced only when there is almost no glucose in the growing medium, but if Le is high, the operon induction level increases abruptly when the levels of glucose in the environment decrease to very low values. We demonstrate that this behavior could not be obtained without bistability if the stability of the induced and uninduced states is to be preserved. Finally, we point out that the present methods and results may be useful to study the emergence of multistability in biological systems other than the lac operon.

Noise in timing and precision of gene activities in a genetic cascade

The timing of events along the induction cascade of bacteriophage lambda is independent of UV dose and displays increased relative temporal precision with cascade progression.

This behavior is reproduced by a model of a cascade consisting of independent steps that shows that higher temporal precision can be attained by a cascade consisting of a large number of fast steps.

The observed cell-cell variability in cascade timing is not due to differences in uniform dilation of intervals between events among cells, but rather to the independent distribution of interval durations within the cascade, consistently with the modular architecture of the lambda genome.

The single-cell time lapse study reveals a bistable regime at low UV doses in which some cells are induced while others are not, evidence for a commitment point beyond which lysis will occur, and an unexpected shutoff of the lambda pR promoter.

Stochasticity or noise, an inherent property of all biological networks, is often manifested by different phenotypic behaviors in clonal populations of cells (Raser and O'Shea, 2005). Noise can arise, for instance, from sources such as cell–cell variations in small numbers of regulatory molecules or from the stochastic nature of molecular interactions (Paulsson, 2005). Besides affecting the number of molecules in a cell, noise may also lead to variability in timing of particular events along a given pathway. In this work, we studied temporal noise in the induction cascade of phage lambda.

Infection of a bacterial cell by bacteriophage lambda can lead to two different fates (Ptashne, 2004; Dodd et al, 2005; Oppenheim et al, 2005): the phage can either multiply inside the host leading to its eventual lysis and the generation of progeny virions (the lytic pathway) or, alternatively, it can integrate its genome into the host's genome (prophage state), replicating passively with the latter (the lysogenic pathway). The prophage state is highly stable, being maintained by a phage-encoded repressor, which shuts off phage genes leading to lytic growth. However, the lytic pathway can be induced in a lysogenic cell, through the activation of the bacterial SOS response to DNA damage (Little, 1996), for example by UV irradiation. Once activated, the SOS response results in cleavage of the lambda repressor, leading to expression of the phage early and late genes, and culminating in the lysis of the host cell.

The lambda induction cascade has been extensively characterized over the years. We built upon this knowledge to tap the cascade at different points and quantitatively analyze the progressive loss of temporal coherence between cells, as different stages along the cascade are executed, following synchronous induction. Using time-lapse microscopy, we monitored the time of activation of early and late genes in individual cells using lambda pR and pR′-tR′ promoter-GFP fusions, respectively, by means of reporter plasmids, and finally the time of lysis. Sample results are shown in Figure 2.

At low UV levels (5 J/m²), the network exhibits bistability: only approximately 40% of the bacteria lyse, whereas the others continue to divide, following a lag period. At high UV levels (20 J/m²), almost all bacteria lyse. We found that the timing of events in cells that lyse is independent of UV dose. This is in contrast to the known behavior of the SOS network (Friedman et al, 2005), indicating that these two networks proceed independently. Following induction, a surprising shutoff in the activity of the pR promoter is observed in all cells (see Figure 2). Furthermore, the data show that whereas early genes are expressed in all cells irrespective of cell fate, late genes are expressed only in the lysing cells, indicating that similar to infection, a specific commitment checkpoint is operating.

To characterize the temporal variability in a cell population, we used the coefficient of variation, defined as the non-dimensional ratio of the standard deviation and the mean time of occurrence of a particular event. We studied the changes in both standard deviation and coefficient of variation in timing of various events along the lambda induction cascade, from the expression of the early genes to the ultimate lysis of the cells. As shown in Figure 6, the absolute noise as measured by the standard deviation increases as the cascade progresses. In contrast, the coefficient of variation, which measures variability relative to the time of occurrence, decreases. Simple theoretical considerations described in the text yield a necessary and sufficient condition for a monotonic decrease in the coefficient of variation. Higher temporal precision can be achieved when the cascade is composed of a large number of fast steps.

Further support for the independence of network modules is furnished by a correlation analysis of the times of occurrence of different steps along the lytic cascade. This analysis also indicates that the variability in lysis time is not due to differences in the global rate of cascade progression, but probably to random fluctuations in the execution time of the various cascade stages. Indeed, phage lambda gene expression architecture is well known to have evolved from a number of independent regulatory modules (Hendrix, 2003).

Biological developmental pathways require proper timing of gene expression. We investigated timing variations of defined steps along the lytic cascade of bacteriophage λ. Gene expression was followed in individual lysogenic cells, after induction with a pulse of UV irradiation. At low UV doses, some cells undergo partial induction and eventually divide, whereas others follow the lytic pathway. The timing of events in cells committed to lysis is independent of the level of activation of the SOS response, suggesting that the lambda network proceeds autonomously after induction. An increased loss of temporal coherence of specific events from prophage induction to lysis is observed, even though the coefficient of variation of timing fluctuations decreases. The observed temporal variations are not due to cell factors uniformly dilating the timing of execution of the cascade. This behavior is reproduced by a simple model composed of independent stages, which for a given mean duration predicts higher temporal precision, when a cascade consists of a large number of steps. Evidence for the independence of regulatory modules in the network is presented.

Bacterial gene regulation in diauxic and non-diauxic growth

When bacteria are grown in a batch culture containing a mixture of two growth-limiting substrates, they exhibit a rich spectrum of substrate consumption patterns including diauxic growth, simultaneous consumption, and bistable growth. In previous work, we showed that a minimal model accounting only for enzyme induction and dilution captures all the substrate consumption patterns [Narang, A., 1998a. The dynamical analogy between microbial growth on mixtures of substrates and population growth of competing species. Biotechnol. Bioeng. 59, 116-121, Narang, A., 2006. Comparitive analysis of some models of gene regulation in mixed-substrate microbial growth, J. Theor. Biol. 242, 489-501]. In this work, we construct the bifurcation diagram of the minimal model, which shows the substrate consumption pattern at any given set of parameter values. The bifurcation diagram explains several general properties of mixed-substrate growth. (1) In almost all the cases of diauxic growth, the "preferred" substrate is the one that, by itself, supports a higher specific growth rate. In the literature, this property is often attributed to the optimality of regulatory mechanisms. Here, we show that the minimal model, which accounts for induction and growth only, displays the property under fairly general conditions. This suggests that the higher growth rate of the preferred substrate is an intrinsic property of the induction and dilution kinetics. It can be explained mechanistically without appealing to optimality principles. (2) The model explains the phenotypes of various mutants containing lesions in the regions encoding for the operator, repressor, and peripheral enzymes. A particularly striking phenotype is the "reversal of the diauxie" in which the wild-type and mutant strains consume the very same two substrates in opposite order. This phenotype is difficult to explain in terms of molecular mechanisms, such as inducer exclusion or CAP activation, but it turns out to be a natural consequence of the model. We show furthermore that the model is robust. The key property of the model, namely, the competitive dynamics of the enzymes, is preserved even if the model is modified to account for various regulatory mechanisms. Finally, the model has important implications for the problem of size regulation in development. It suggests that protein dilution may be the mechanism coupling patterning and growth.

In silico evolved lac operons exhibit bistability for artificial inducers, but not for lactose

Bistability in the lac operon of Escherichia coli has been widely studied, both experimentally and theoretically. Experimentally, bistability has been observed when E. coli is induced by an artificial, nonmetabolizable, inducer. However, if the lac operon is induced with lactose, the natural inducer, bistability has not been demonstrated. We derive an analytical expression that can predict the occurrence of bistability both for artificial inducers and lactose. We find very different conditions for bistability in the two cases. Indeed, for artificial inducers bistability is predicted, but for lactose the condition for bistability is much more difficult to satisfy. Moreover, we demonstrate that in silico evolution of the lac operon generates an operon that avoids bistability with respect to lactose, but does exhibit bistability with respect to artificial inducers. The activity of this evolved operon strikingly resembles the experimentally observed activity of the operon. Thus our computational experiments suggest that the wild-type lac operon, which regulates lactose metabolism, is not a bistable switch. Nevertheless, for engineering purposes, this operon can be used as a bistable switch with artificial inducers.

Predicting stochastic gene expression dynamics in single cells

Fluctuations in protein numbers (noise) due to inherent stochastic effects in single cells can have large effects on the dynamic behavior of gene regulatory networks. Although deterministic models can predict the average network behavior, they fail to incorporate the stochasticity characteristic of gene expression, thereby limiting their relevance when single cell behaviors deviate from the population average. Recently, stochastic models have been used to predict distributions of steady-state protein levels within a population but not to predict the dynamic, presteady-state distributions. In the present work, we experimentally examine a system whose dynamics are heavily influenced by stochastic effects. We measure population distributions of protein numbers as a function of time in the Escherichia coli lactose uptake network (lac operon). We then introduce a dynamic stochastic model and show that prediction of dynamic distributions requires only a few noise parameters in addition to the rates that characterize a deterministic model. Whereas the deterministic model cannot fully capture the observed behavior, our stochastic model correctly predicts the experimental dynamics without any fit parameters. Our results provide a proof of principle for the possibility of faithfully predicting dynamic population distributions from deterministic models supplemented by a stochastic component that captures the major noise sources.

Stochastic model of coliphage lambda regulatory network

The dynamic properties of the regulatory network governing the choice between lytic and lysogenic growths of coliphage lambda is studied using a Markov chain stochastic model. Our computer simulation confirms the finding by Li et al. [Proc. Natl. Acad. Sci. USA 101, 4781 (2004)] on the dynamics of budding yeast: that the biological stationary states are global attractors; the biological pathways of lytic and lysogenic growths are attracting trajectories; and the network functions are robustly designed against structural perturbations. In addition, our model shows the stochastic switch from lysogen to lytic growth, which has been observed in experiments. A definition of pseudoenergy is introduced in the network dynamics to reveal a transitionlike behavior in the system.

Biological significance of autoregulation through steady state analysis of genetic networks

Autoregulation of regulatory proteins is a recurring theme in genetic networks. Autoregulation is an important component of a genetic regulatory network besides protein-protein and protein-DNA interactions, stoichiometry, multiple binding sites and cooperativity. Although the biological significance of autoregulation has been studied before, its significance in presence of other mechanisms is not clearly enumerated. We have analyzed at steady state the significance of autoregulation in presence of other molecular mechanisms by considering hypothetical genetic networks. We demonstrate that autoregulation of a regulatory protein can impart amplification to the response. Further, autoregulation of an activator binding to the DNA as a dimer can introduce bistability, thus forcing the system to reside in two distinct steady states. In combination with autoregulation, cooperative binding can further increase the sensitivity and can yield a highly ultrasensitive response. We conclude that autoregulation with the help of other molecular mechanisms can impart distinct system level properties such as amplification, sensitivity and bistability. The results are further discussed in relation to various examples of genetic networks that exist in biological systems.

Comparative analysis of some models of gene regulation in mixed-substrate microbial growth

Mixed-substrate microbial growth is of fundamental interest in microbiology and bioengineering. Several mathematical models have been developed to account for the genetic regulation of such systems, especially those resulting in diauxic growth. In this work, we compare the dynamics of three such models (Narang, 1998a. The dynamical analogy between microbial growth on mixtures of substrates and population growth of competing species. Biotechnol. Bioeng. 59, 116-121; Thattai and Shraiman, 2003. Metabolic switching in the sugar phosphotransferase system of Escherichia coli. Biophys. J. 85(2), 744-754; Brandt et al., 2004. Modelling microbial adaptation to changing availability of substrates. Water Res. 38, 1004-1013). We show that these models are dynamically similar--the initial motion of the inducible enzymes in all the models is described by the Lotka-Volterra equations for competing species. In particular, the prediction of diauxic growth corresponds to "extinction" of one of the enzymes during the first few hours of growth. The dynamic similarity occurs because in all the models, the inducible enzymes possess properties characteristic of competing species: they are required for their own synthesis, and they inhibit each other. Despite this dynamic similarity, the models vary with respect to the range of dynamics captured. The Brandt et al. model always predicts the diauxic growth pattern, whereas the remaining two models exhibit both diauxic and non-diauxic growth patterns. The models also differ with respect to the mechanisms that generate the mutual inhibition between the enzymes. In the Narang model, mutual inhibition occurs because the enzymes for each substrate enhance the dilution of the enzymes for the other substrate. The Brandt et al. model superimposes upon this dilution effect an additional mechanism of mutual inhibition. In the Thattai and Shraiman model, the mutual inhibition is entirely due to competition for the phosphoryl groups. For quantitative agreement with the data, all models must be modified to account for specific mechanisms of mutual inhibition, such as inducer exclusion.

The dynamics of single-substrate continuous cultures: The role of ribosomes

When a chemostat is perturbed from its steady state, it displays complex dynamics. For instance, if the identity of the growth-limiting substrate is switched abruptly, the substrate concentration and cell density undergo a pronounced excursion from the steady state that can last several days. These dynamics occur because certain physiological variables respond slowly. In the literature, several physiological variables have been postulated as potential sources of the slow response. These include transport enzymes, biosynthetic enzymes, and ribosomes. We have been addressing this problem by systematically exploring the role of these variables. In previous work Shoemaker et al. (J. Theor. Biol., 222 (2003) 307-322), we studied the role of transport enzymes, and we showed that transients starting from low transport enzyme levels could be quantitatively captured by a model taking due account of transport enzyme synthesis. However, there is some experimental data indicating that slow responses occur even if the initial enzyme levels are high. Here, we analyse this data to show that in these cases, the sluggish response is most probably due to slow adjustment of the ribosome levels. To test this hypothesis, we extend our previous model by accounting for the evolution of both the transport enzyme and the ribosomes. Based on a kinetic analysis of the data in the literature, we assume that the specific protein synthesis rate is proportional to the ribosome level, and the specific ribosome synthesis rate is autocatalytic. Simulations of the model show remarkable agreement with experimentally observed steady states and the transients. Specifically, the model predictions are in good agreement with (1) the steady-state profiles of the cell density, substrate concentration, RNA, proteins, and transport enzymes, (2) the instantaneous specific substrate uptake, growth, and respiration rates in response to a continuous-to-batch shift, and (3) the transient profiles of the cell density, substrate concentration, and RNA in response to feed switches and dilution rate shifts. Time-scale analysis of the model reveals that every transient response is a combination of two fundamental (and simpler) dynamics, namely, substrate-sufficient batch dynamics and cell-sufficient fed-batch dynamics. We obtain further insight into the transient response by analysing the equations describing these fundamental dynamics. The analysis reveals that in feed switches or dilution rate shift-ups, the transport enzyme reaches a maximum before RNA achieves its maximum, and in dilution rate shift-downs the cell density reaches a maximum before RNA achieves a minimum.

Autoregulation of regulatory proteins is key for dynamic operation ofGALswitch inSaccharomyces cerevisiae

Autoregulation and nucleocytoplasmic shuttling play important roles in the operation of the GAL regulatory system. However, the significance of these mechanisms in the overall operation of the switch is unclear. In this work, we develop a dynamic model for the GAL system and further validate the same using steady-state and dynamic experimental expression data. Next, the model is used to delineate the relevance of shuttling and autoregulation in response to inducing, repressing, and non-inducing-non-repressing media. The analysis indicates that autoregulation of the repressor, Gal80p, is key in obtaining three distinct steady states in response to the three media. In particular, the analysis rationalizes the intuitively paradoxical observation that the concentration of repressor, Gal80p, actually increases in response to an increase in the inducer concentration. On the other hand, although nucleocytoplasmic shuttling does not affect the dynamics of the system, it plays a dominant role in obtaining a sensitive response to galactose. The dynamic model was also used to obtain insights on the preculturing effect on the system behavior.

Multistability in the lactose utilization network of Escherichia coli

Multistability, the capacity to achieve multiple internal states in response to a single set of external inputs, is the defining characteristic of a switch. Biological switches are essential for the determination of cell fate in multicellular organisms, the regulation of cell-cycle oscillations during mitosis and the maintenance of epigenetic traits in microbes. The multistability of several natural and synthetic systems has been attributed to positive feedback loops in their regulatory networks. However, feedback alone does not guarantee multistability. The phase diagram of a multistable system, a concise description of internal states as key parameters are varied, reveals the conditions required to produce a functional switch. Here we present the phase diagram of the bistable lactose utilization network of Escherichia coli. We use this phase diagram, coupled with a mathematical model of the network, to quantitatively investigate processes such as sugar uptake and transcriptional regulation in vivo. We then show how the hysteretic response of the wild-type system can be converted to an ultrasensitive graded response. The phase diagram thus serves as a sensitive probe of molecular interactions and as a powerful tool for rational network design.

Cell cycle- and age-dependent activation of Sod1p drives the formation of stress resistant cell subpopulations within clonal yeast cultures

Phenotypic heterogeneity describes non-genetic variation that exists between individual cells within isogenic populations. The basis for such heterogeneity is not well understood, but it is evident in a wide range of cellular functions and phenotypes and may be fundamental to the fitness of microorganisms. Here we use a suite of novel assays applied to yeast, to provide an explanation for the classic example of heterogeneous resistance to stress (copper). Cell cycle stage and replicative cell age, but not mitochondrial content, were found to be principal parameters underpinning differential Cu resistance: cell cycle-synchronized cells had relatively uniform Cu resistances, and replicative cell-age profiles differed markedly in sorted Cu-resistant and Cu-sensitive subpopulations. From a range of potential Cu-sensitive mutants, cup1Delta cells lacking Cu-metallothionein, and particularly sod1Delta cells lacking Cu, Zn-superoxide dismutase, exhibited diminished heterogeneity. Furthermore, age-dependent Cu resistance was largely abolished in cup1Delta and sod1Delta cells, whereas cell cycle-dependent Cu resistance was suppressed in sod1Delta cells. Sod1p activity oscillated approximately fivefold during the cell cycle, with peak activity coinciding with peak Cu-resistance. Thus, phenotypic heterogeneity in copper resistance is not stochastic but is driven by the progression of individual cells through the cell cycle and ageing, and is primarily dependent on only Sod1p, out of several gene products that can influence the averaged phenotype. We propose that such heterogeneity provides an important insurance mechanism for organisms; creating subpopulations that are pre-equipped for varied activities as needs may arise (e.g. when faced with stress), but without the permanent metabolic costs of constitutive expression.

Metabolic switching in the sugar phosphotransferase system of Escherichia coli

Bacteria grown in a mixture of multiple sugars will first metabolize a preferred sugar until it is nearly depleted, only then turning to other carbon sources in the medium. This sharp switching of metabolic preference is characteristic of systems that optimize fitness. Here we consider the mechanism by which switching can occur in the Escherichia coli phosphotransferase system (PTS), which regulates the uptake and metabolism of several sugars. Using a model combining the description of fast biochemical processes and slower genetic regulation, we derive metabolic phase diagrams for the uptake of two PTS sugars, indicating regions of distinct sugar preference as a function of external sugar concentrations. We then propose a classification of bacterial phenotypes based on the topology of the metabolic phase diagram, and enumerate the possible topologically distinct phenotypes that can be achieved through mutations of the PTS. This procedure reveals that there is only one nontrivial switching phenotype that is insensitive to large changes in biochemical parameters. This phenotype exhibits diauxic growth, a manifestation of the winner-take-all dynamics enforced by PTS architecture. Winner-take-all behavior is implemented by the induction of sugar-specific operons, combined with competition between sugars for limited phosphoryl flux. We propose that flux-limited competition could be a common mechanism for introducing repressive interactions in cellular networks, and we argue that switching behavior similar to that described here should occur generically in systems that implement such a mechanism.

The dynamics of single-substrate continuous cultures: the role of transport enzymes

A chemostat limited by a single growth-limiting substrate displays a rich spectrum of dynamics. Depending on the flow rate and feed concentration, the chemostat settles into a steady state or executes sustained oscillations. The transients in response to abrupt increases in the flow rate or the feed concentration are also quite complex. For example, if the increase in the flow rate is small, there is no perceptible change in the substrate concentration. If the increase in the flow rate is large, there is a large increase in the substrate concentration lasting several hours or days before the culture adjusts to a new steady state. In the latter case, the substrate concentration and cell density frequently undergo damped oscillations during their approach to the steady state. In this work, we formulate a simple structured model containing the inducible transport enzyme as the key intracellular variable. The model displays the foregoing dynamics under conditions similar to those employed in the experiments. The model suggests that long recovery times (on the order of several hours to several days) can occur because the initial transport enzyme level is too small to cope with the increased substrate supply. The substrate concentration, therefore, increases until the enzyme level is built up to a sufficiently high level by the slow process of enzyme induction. Damped and sustained oscillations can occur because transport enzyme synthesis is autocatalytic, and hence, destabilizing. At low dilution rates, the response of stabilizing processes, such as enzyme dilution and substrate consumption, becomes very slow, leading to damped and sustained oscillations.

Modeling network dynamics: the lac operon, a case study

We use the lac operon in Escherichia coli as a prototype system to illustrate the current state, applicability, and limitations of modeling the dynamics of cellular networks. We integrate three different levels of description (molecular, cellular, and that of cell population) into a single model, which seems to capture many experimental aspects of the system.

Mitochondrial membrane potential differentiates cells resistant to apoptosis in hybridoma cultures

Previous research has implicated mitochondrial physiology and, by extension, respiratory capacity in the initiation and progress of apoptosis of cells in culture and tissue environments. This hypothesis was tested by separating a hybridoma cell population into subpopulations of varying mitochondrial membrane potential (MMP) using Rhodamine 123 stain and fluorescence-activated cell sorter analysis and subjecting them to two apoptosis inducers, rotenone and staurosporin. Apoptotic death was characterized morphologically through the determination of apoptosis-related chromatin condensation and biochemically through the measurement of caspase-3 enzymatic activity. We found dramatic differences in the apoptotic death kinetics for the subpopulations, with the high MMP cells showing higher resistance to apoptotic death. After incubation with 30 microM rotenone, the low MMP cells exhibited one-third of the viability of the high MMP cells and a three-fold increase in the capsase-3 enzymatic activity. No changes were observed in the DNA content or the cell cycle distributions of the two cell subpopulations, which maintained their mean MMP difference after 20 generations. These results suggest that heterogeneity exists in mammalian cell populations with respect to mitochondrial physiology, which correlates with resistance to apoptotic death.

Simulation of prokaryotic genetic circuits

Biochemical and genetic approaches have identified the molecular mechanisms of many genetic reactions, particularly in bacteria. Now a comparably detailed understanding is needed of how groupings of genes and related protein reactions interact to orchestrate cellular functions over the cell cycle, to implement preprogrammed cellular development, or to dynamically change a cell's processes and structures in response to environmental signals. Simulations using realistic, molecular-level models of genetic mechanisms and of signal transduction networks are needed to analyze dynamic behavior of multigene systems, to predict behavior of mutant circuits, and to identify the design principles applicable to design of genetic regulatory circuits. When the underlying design rules for regulatory circuits are understood, it will be far easier to recognize common circuit motifs, to identify functions of individual proteins in regulation, and to redesign circuits for altered functions.

Non-genetic population heterogeneity studied by in situ polymerase chain reaction

Expression of a lac operon in Salmonella typhimurium single cells was monitored using lac mRNA targeting in situ reverse transcription-polymerase chain reaction (RT-PCR). It is demonstrated that suboptimal induction of the lac operon in a culture of S. typhimuriuml/F'lac+ cells generates a subpopulation in which transcription of the lac operon occurs and another subpopulation in which transcription of the lac operon is repressed, whereas suboptimal induction of the lac operon in a culture of S. typhimuriuml/F'lacY cells generates a population with uniform levels of lac mRNA. The outcome of the single-cell lac mRNA detection assay was compared with the outcome of a single-cell beta-galactosidase assay. In cultures grown under different suboptimal lac induction conditions, the fraction of cells in which transcription of the lac operon occurred was concurrent with the fraction of cells showing beta-galactosidase activity. Besides supporting the hypothesis that the lactose permease has a role in generating non-genetic heterogeneity in suboptimally induced cultures of Lac+ cells, these results demonstrate the usefulness of in situ RT-PCR for the study of non-genetic population heterogeneities.

Multiplicity and stability analysis of microorganisms in continuous culture: effects of metabolic overflow and growth inhibition

Metabolic overflow (enhanced uptake of substrate and secretion of intermediates) is a phenomenon often observed for cells grown under substrate excess. Growth inhibition by substrate and/or product is also normally found for this kind of culture. An effort is made in this work to analyze the dynamic behavior of a continuous culture subject to metabolic overflow and growth inhibition by substrate and/or product. Analysis of a model system shows that in a certain range of operating conditions three nonwashout steady state solutions are possible. Local stability analysis indicates that only two of them are stable thus leading to multiplicity and hysteresis. Further analysis of the intrinsic effects of different terms describing the metabolic overflow and growth inhibitions reveals that for the model system and the parameters considered, the combined effects of product inhibition and an enhanced formation rate of product under substrate excess cause the multiplicity and hysteresis. Growth inhibition by substrate and/or an enhanced substrate uptake appear not to be necessary conditions. The combined effects of enhanced product formation and product inhibition can also lead to unusual dynamic behavior such as a prolonged time period to reach a steady state, oscillatory transition from one steady state to another, and sustained oscillations. Using the occurrence of multiplicity and oscillation as criteria, the operating regime of a continuous culture can be divided into four domains: one with multiplicity and oscillation, one with unique steady state but possible oscillatory behavior, the other two with unique and stable steady state. The model predictions are in accordance with recent experimental results. The results presented in this work may be used as guidelines for choosing proper operating conditions of similar culture systems to avoid undesired instability and multiplicity. Copyright 1998 John Wiley & Sons, Inc.