LSE Logo MBoC Logo

Single-cell variation leads to population invariance in NF-κB signaling dynamics

Published Online:


The activation dynamics of nuclear factor (NF)-κB have been shown to affect downstream gene expression. On activation, NF-κB shuttles back and forth across the nuclear envelope. Many dynamic features of this shuttling have been characterized, and most features vary significantly with respect to ligand type and concentration. Here, we report an invariant feature with regard to NF-κB dynamics in cellular populations: the distribution—the average, as well as the variance—of the time between two nuclear entries (the period). We find that this period is conserved, regardless of concentration and across several different ligands. Intriguingly, the distributions observed at the population level are not observed in individual cells over 20-h time courses. Instead, the average period of NF-κB nuclear translocation varies considerably among single cells, and the variance is much smaller within a cell than that of the population. Finally, analysis of daughter-cell pairs and isogenic populations indicates that the dynamics of the NF-κB response is heritable but diverges over multiple divisions, on the time scale of weeks to months. These observations are contrary to the existing theory of NF-κB dynamics and suggest an additional level of control that regulates the overall distribution of translocation timing at the population level.


The nuclear factor (NF)-κB signaling network plays a critical role in innate immune signaling (Hayden et al., 2006). The canonical NF-κB transcription factor, a heterodimer of p65 and p50, is normally held inactive in the cytoplasm by its association with inhibitor of κB (IκB) proteins. In response to any of a range of stimuli, IκB is phosphorylated and degraded, allowing NF-κB to translocate to the nucleus to induce expression of hundreds of genes (Pahl, 1999), a number of which inhibit NF-κB activity in negative feedback loops. The effects of these loops are evident in the beautiful oscillatory dynamics of NF-κB nuclear–cytoplasmic translocation that can be observed and quantified in single cells (Nelson et al., 2004). This dynamics drives NF-κB–dependent gene expression (Nelson et al., 2004; Tay, Hughey, et al., 2010) and thus has been intensively studied both experimentally and theoretically (Hoffmann, Levchenko, et al., 2002; Lipniacki et al., 2004; Ashall, Horton, Nelson, Paszek, et al., 2009; Lee et al., 2009; Cheong et al., 2011; Pękalski et al., 2013; Shinohara, Behar, et al., 2014; Sung, Li, Lao, et al., 2014).

Although the oscillations of NF-κB have been observed in multiple types of cultured mammalian cells, they have only been well characterized in response to a single stimulus—tumor necrosis factor (TNF). Signaling to NF-κB can be activated by many different stimuli, and it is unclear how different stimuli might change the oscillations. Recent theoretical work suggests that the period of the oscillations may be stimulus independent (Longo, Selimkhanov, et al., 2013). More fundamentally, how noisy are the oscillations in a single cell, and how long will the oscillations last in individual cells? Such questions have been addressed for p53 (Geva-Zatorsky, Rosenfeld, et al., 2006) and the circadian clock (Leise et al., 2012), two other mammalian oscillators, but only partially for NF-κB (Nelson et al., 2004; Sung et al., 2009).

In earlier work, we measured the nuclear localization of a fluorescent protein fused to NF-κB subunit p65 in response to TNF concentrations spanning five orders of magnitude and for thousands of individual cells (Tay, Hughey, et al., 2010). At that time, we examined many properties of p65 nuclear translocation dynamics and measured how key properties varied with TNF concentration (Figure 1A, left). Most features changed according to TNF concentration, but the interpeak time (i.e., the period), which averaged ∼80 min, did not. This observation of concentration-independent interpeak times was also made by others using a different cell type (Turner et al., 2010).


FIGURE 1: The time between p65-dsRed nuclear localization peaks is constant across stimuli and concentration. (A) The fraction of active cells (top row), as well as the time to (second row) and amplitude of (third row) the first peak, vary, depending on stimulus and concentration, whereas the time between the first and second peaks (bottom row) is constant. (B) Distributions of all interpeak times for five different environmental conditions: TNF and LPS (concentrations as shown in A), LPS (concentrations as shown in A), together with soluble TNF receptor II (sTNFR; 5 μg/ml), Pam3CSK4 (1 μg/ml), and Sendai virus (10 and 100 U/ml). See also Supplemental Figure S1.

The aim of this study was to quantitatively characterize the oscillations of NF-κB localization across stimuli and individual cells and over generations. We report that not only the mean, but also the entire distribution of the interpeak times is conserved—across concentrations and several stimuli. By making substantial improvements in our ability to observe and track single cells, we further show that the conservation of interpeak times holds only at the population level and not for individual cells. This finding is not explained by existing computational models of NF-κB, which focus on transcriptional variability within cells as opposed to variability across cellular populations. Finally, we show that over time and many generations, the clonal population derived from a single cell largely recapitulates the conserved distribution of interpeak time.


Distribution of interpeak times is conserved across concentrations and stimuli

In previous work, we found that the time between two subsequent p65 nuclear translocations was conserved across concentrations of TNF (Tay, Hughey, et al., 2010). We found this invariance compelling and therefore tested whether the period of oscillation was also conserved in the response to other stimuli. We first exposed cells to a range of concentrations of lipopolysaccharide (LPS), a component of the outer membrane of Gram-negative bacteria. Different preparations of LPS can signal through different receptors and induce different dynamics of NF-κB activity (Hirschfeld et al., 1999; Gutschow, Hughey, et al., 2013). For this study, we used a highly purified preparation of Escherichia coli LPS that acts only through Toll-like receptor 4 (TLR4), as we verified previously (Lee et al., 2009). We then determined the values of key features of the dynamics, as before (Figure 1, right). As with TNF, most features of the dynamics of NF-κB varied with LPS concentration. In addition, cells stimulated with LPS showed nuclear–cytoplasmic oscillations of NF-κB whose average interpeak time remained essentially constant across concentrations. Of great interest, the distributions of interpeak times induced by TNF and LPS, pooled across all concentrations, were statistically indistinguishable (Figure 1B, first and second rows; p = 0.53 by two-sample t test).

We and others previously reported that cells stimulated with certain preparations of LPS may secrete TNF, which can activate NF-κB in a paracrine and autocrine manner (Covert, Leung, et al., 2005; Lee et al., 2009). We therefore wondered whether the similarity in interpeak time distribution between TNF and LPS might simply be due to LPS-induced secretion of TNF. To the contrary, the distribution of interpeak times of cells stimulated with LPS was not affected by the presence of soluble TNF receptor (Figure 1B, third row), indicating that conservation of the interpeak time is not caused by secreted TNF.

To determine whether the interpeak time distribution was the same for other stimuli, we measured the period of nuclear localization in response to Pam3CSK4 and Sendai virus, which activate NF-κB through different receptors from those that are bound by either TNF (TNFR) or LPS (TLR4). Pam3CSK4, a synthetic triacylated lipoprotein that mimics bacterial lipopeptides, is believed to activate NF-κB through a complex of TLR1 and TLR2 (Ozinsky, Underhill, et al., 2000). Unlike TLR4, which signals through adapter proteins myeloid differentiation primary response 88 (MyD88) and TIR (Toll/Interleukin-1 receptor)-domain-containing adapter-inducing interferon-beta (TRIF), TLR1/2 is believed to signal only through MyD88 (Kawai and Akira, 2006). Whereas TNF, LPS, and Pam3CSK4 bind to extracellular receptors, Sendai virus is a single-stranded RNA virus that activates innate immune signaling primarily through the intracellular receptor retinoic acid-inducible gene 1 (RIG-I; Kato et al., 2006) and possibly melanoma differentiation-associated protein 5 (MDA-5; Gitlin, Benoit, et al., 2010).

Although Pam3CSK4 and Sendai virus activate NF-κB via distinct molecular pathways, the distribution of interpeak times of p65 nuclear localization was similar to our previous measurements (Figure 1B, fourth and fifth rows, and Supplemental Figure S1A). Taken together, our results suggest that the distribution of oscillation period is constant across multiple concentrations and diverse stimuli.

The interpeak time distribution varies from cell to cell

The distribution of interpeak times was constant, but it was also broad, with a coefficient of variation of ∼40%. Part of that large variability comes from a cell's first interpeak time, which tended to be longer and more variable than later interpeak times (Supplemental Figure S1, B and C). However, even including only the later oscillations, the distribution in interpeak times showed marked variability. Of importance, our results so far have been based on data pooled from many cells, which means the variability we measured was a property of the population of cells, not a property of NF-κB signaling in single cells. Given the cell-to-cell variability that we have observed in cells’ responsiveness to low concentrations of TNF (Tay, Hughey, et al., 2010), we wondered whether cell-to-cell variability also existed in the oscillatory period of NF-κB localization.

The variability in the population must arise from a combination of two sources—one intracellular and the other intercellular (Figure 2A). The distinction between intracellular and intercellular variability is similar to the concept of intrinsic and extrinsic noise (Elowitz et al., 2002). Intracellular variability describes the variation from one oscillation to another within a given cell. If intracellular variability produced all of the variability in the interpeak time that we observed across the population, the distribution of interpeak times for each cell would be the same and would resemble the distribution for the entire population (Figure 2A, left). One cause of intracellular variability, for example, might be a reaction in the network whose duration is random.


FIGURE 2: The distribution of interpeak times within individual cells varies from one cell to another. (A) Schematic illustrating variation within cells and between cells. (B) Nuclear localization traces and (C) interpeak time distributions of three cells (selected from 199 total across two experiments) monitored for 20 h upon stimulation with 10 ng/ml TNF. (D) Mean interpeak time for each cell plotted against CV. (E) Mean and SDs of interpeak time for each individual cell plotted together with the distribution that would be expected without intercellular variability (black). (F) Ratio of variance between cells to variance within cells (F statistic), as determined using bootstrap analysis.

The other possible source of variability is intercellular (also referred to as cell-to-cell), that is, between cells in the population. If intercellular variability were at least partly responsible for the variability in the interpeak time that we observed across the population, the distribution of interpeak times would be different from one cell to another (Figure 2A, right). In addition, the distribution for an individual cell would tend to be narrower than the population distribution. Intercellular variability could involve each cell having a different concentration of a given signaling molecule, thereby affecting an important rate constant.

Which type of variability is responsible for the interpeak time distributions we observed? To quantify accurately the intracellular and intercellular variability of the dynamics of p65 translocation, we needed to determine the distribution of interpeak time for many individual cells. This required measuring many more oscillations for each cell than had previously been reported. A number of technical improvements (see Materials and Methods) enabled us to observe p65 dynamics continuously in the same cells for 20 h (Figure 2B). Of the cells stimulated with 10 ng/ml TNF, almost all of them showed sustained nuclear–cytoplasmic oscillations of p65 throughout the time course. In total, we observed ∼200 cells that exhibited at least nine oscillations. To aid exploration of this rich data set, we developed an interactive figure that includes all of the individual cell data (Interactive Figure;

Our results showed that the individual cells varied considerably in their distributions of interpeak time (Figure 2, C–F). The mean interpeak time of the cells ranged from 50 to 100 min, with an interquartile range from 62 to 77 min (Figure 2D and Supplemental Figure S2). Consistent with later interpeak times tending to be shorter than the first (Supplemental Figure S1, B and C), the distribution of mean interpeak times is somewhat narrower and shifted to the left compared with the distributions of interpeak times from the various stimuli (Figure 1B). With regard to intracellular variability, the interquartile range of the coefficient of variation (CV) of the interpeak time for each cell was between 8 and 21%. The mean and the CV of the period in a given cell were uncorrelated (r2 = 2 × 10−5).

Finally, we tested for intercellular variability. We pooled all the interpeak times from all the cells and then randomly and repeatedly drew 10 interpeak times from that pool (with replacement) to create a null distribution corresponding to no intercellular variability. This randomization of the interpeak times produced a significantly different distribution than was actually observed in single cells (Figure 2E, p = 1.6 × 10−9 by two-sample Kolmogorov–Smirnov [KS] test), and the variance between cells was found to be about sixfold higher than that within cells (Figure 2F). We therefore concluded that the population variability of the period of p65 oscillations is driven largely by cell-to-cell variability.

Computational modeling suggests that stochastic transcription alone cannot reproduce intercellular variability

Given the intercellular variability in p65 oscillations that we observed, we next sought to examine the sources of variability in a computational model of the NF-κB signaling network. To represent the heterogeneity in single cells, a model must contain a stochastic or variable element. The model we used represents the binding of NF-κB to the promoters of its inhibitors (A20, IκBα, and IκBε) as a stochastic process, which leads to stochastic transcription of the corresponding mRNAs (Paszek et al., 2010). This has been a standard assumption for the source of variability in NF-κB signaling. The other reactions in the model are approximated as deterministic. The model also contains stimulus-independent A20 activity, which has been shown to correctly reproduce the invariance of the interpeak time across simulated TNF concentrations (Turner et al., 2010). The model was originally fitted to data from murine embryonic fibroblasts, which are similar to the 3T3s used in this study.

Using this model, we ran 500 simulations of a 20-h stimulation with 10 ng/ml TNF. Each simulation corresponds to the dynamics of the signaling network in an individual cell. These simulated cells show sustained oscillations of p65 localization (Figure 3A). We quantified the interpeak times for each simulation (Figure 3B) and found an overall interpeak time distribution of 69 ± 13 min.


FIGURE 3: Computational modeling captures variability within, but not between, cells. (A) Three model simulations of p65 nuclear localization. (B) Interpeak time distributions of three model simulations. (C) CV vs. mean interpeak time for each simulation (orange) and for each experimentally observed cell (blue). Ratio of variability between cells to variability within all simulations, as determined using bootstrap analysis. (D) Mean and SDs of interpeak time for each individual simulation, plotted together with the distribution that would be expected without any intercellular variability (black). (E) Ratio of variance between simulations to variance within all simulations (F statistic), as determined using bootstrap analysis, on the same scale as Figure 2F.

We then calculated the intracellular and intercellular variability by comparing interpeak times within and between individual simulations. A typical simulated cell had a mean interpeak time of between 66.5 and 71.1 min and a CV between 15 and 20% (interquartile ranges; Figure 3C). The ranges of both the mean and CV of the simulated data were about two times smaller than in the experimental data, which indicates a larger component of intercellular variability in the experimental data. Consistent with the experimental observations, the mean interpeak time and the CV of the interpeak time for a given simulated cell were uncorrelated (r2 = 0.08).

Finally, to test statistically for intercellular (i.e., intersimulation) variability, we performed the same randomization procedure that we did for the experimentally observed cells. We found that the distribution of mean interpeak times from the simulations differed only slightly from the null distribution corresponding to no intercellular variability (Figure 3D; p = 0.05 by two-sample KS test). Furthermore, the ratio of intracellular to intercellular variability for the model simulations and the randomized case showed only a small (26%), albeit significant (F statistic, p < 0.05), difference (Figure 3E). We concluded that the stochastic promoter binding in the model produces variability primarily only from one oscillation to another in every cell but does not create meaningful differences between cells.

Besides stochastic processes, another way to produce heterogeneity in a model is to vary the parameters from one simulation to another. With that in mind, we explored how varying the model parameters affected the oscillations of NF-κB. We independently varied each model parameter up and down twofold. For each set of parameters, we ran 50 simulations of a 12-h stimulation with TNF, then calculated the mean and CV of the interpeak time for those simulations. Finally, we calculated the sensitivity of the mean and CV of the interpeak time to each parameter (Supplemental Figure S3). Intriguingly, the results suggest that parameters can have widely varying effects on the oscillations of NF-κB. For example, an increase in A20's translation rate (parameter c2) increases the oscillatory period and strongly reduces the variability of the oscillations, whereas an increase in A20 mRNA's degradation rate (parameter c3) strongly reduces the oscillatory period but has no effect on the variability.

Our simulation results therefore suggest that many parameters could produce the intercellular variability that we observe if the parameters vary from one cell to another. Such parameters might be considered epigenetic factors, which could drift over time and generations.

From single-cell to population interpeak time distributions: the emergence of variability

The finding that the individual cell interpeak time distributions were so different from each other, especially given that the overall population distribution is so constant across stimuli and concentrations, led to questions about time scales. In particular, how long does it take for the progeny of an individual cell to generate the variability seen in the population? To answer this question, we considered cells at generations 0 (the founding cell) and 1 (the daughter cells), as well as after development of a clonal population (∼4–8 wk).

First, we considered the founding cell. Because our cells were actively growing during the experiment, it seemed possible that the cell cycle might affect the oscillations of NF-κB. However, we found that the average period of the population did not change significantly over time (Figure 4A). We also calculated the correlation of interpeak time with time after stimulation for each cell (Figure 4B). The distribution of these correlations was centered at zero, and no individual cell had a significant nonzero correlation (after Bonferroni correction). Together these data suggest that the period of p65 oscillations does not change appreciably on the time scale of our experiments (20 h). These data also imply that the cell cycle is likely not responsible for the intercellular variability that we observe.


FIGURE 4: Analysis of the time scale on which the period changes. (A) For every 50-min interval from 100 to 1100 min after TNF stimulation, the mean (black) and SD (gray) of all observed interpeak times whose left peak occurs in that interval. (B) Distribution, along with mean and SD, of the Spearman rank correlation of interpeak time with time for every cell. (C) Two examples of cell division that occurred during TNF stimulation: mother cell (black) and two daughter cells (light and dark blue). (D) Distribution of absolute difference in mean interpeak time between two daughter cells with the same mother and between any two daughter cells. (E) Distribution of interpeak times, in response to 10 ng/ml TNF stimulation, for each clonal cell line. (F) Same calculation as in D. Distribution of all pairwise comparisons of mean interpeak time within cells of a clone compared with the distribution of all pairwise comparisons across cells of any clone. See also Supplemental Figure S2.

Next we investigated whether the distribution of interpeak time was affected by cell division. To do this, we analyzed the oscillations of p65 in cells that divided during our experiments (Figure 4C). For each pair of daughter cells, we calculated the absolute difference between the mean interpeak times. We also calculated the absolute differences between mean interpeak times for all cells in the sample. We found that the typical absolute difference in mean interpeak time between two daughters was about half what would be expected if the mean interpeak times of the two daughters were independent (Figure 4D). Furthermore, after Bonferroni correction, none of the daughter cell pairs showed a significant difference in interpeak time. Thus a single cell division does not appear to create variability between cells with respect to the period of p65 oscillations.

Finally, to study how the period of NF-κB oscillations change on a longer time scale, we created multiple clonal lines of cells derived from our original cells. At the time of our experiments, we estimate that the cells of each clonal line were separated from their respective founding cell by at least 30 cell divisions. For each clonal cell line, we analyzed the dynamics of p65 nuclear localization in response to 10 ng/ml TNF and quantified interpeak times as before. We found that the interpeak time distributions for each clone resembled each other (Figure 4E). Moreover, the clone-based interpeak time distributions were similar to but slightly shifted from the distributions for long-term TNF stimulated cells (∼10 min difference in means; Figure 4E), as well from as distributions for all other stimuli and concentrations (Figure 1B).

We then calculated the absolute differences between mean interpeak times within clones and for all cells in the clone data sets and compared the distributions for each (Figure 4F). The clonal cell data were much more consistent with the pairwise average between cells than with the data we derived from daughters. Furthermore, the majority of the clones exhibited significant intercellular variability, and three of the six clones showed levels of intercellular variability approaching that of the original polyclonal cells (Supplemental Figure S4). These results suggest that after a period of weeks to months, at least some individual cells can regenerate most of the intercellular variability of the original population.


In summary, we have identified an invariant property in the dynamic response of NF-κB to multiple stimuli and concentrations: the period of oscillation of p65 between the cytoplasm and the nucleus. In addition, the distribution of oscillation periods is conserved across concentrations and for a variety of ligands, namely a cytokine, a virus, and bacterial wall components. These experimental results are consistent with recent theoretical work (Longo, Selimkhanov, et al., 2013). It has been unclear whether the oscillations in NF-κB signaling are specific to cytokine stimulation or are a general phenomenon of the network. Our data suggest that the oscillatory pattern of p65 localization is a general tendency of the NF-κB network that can be triggered by multiple receptor complexes. Just as highly conserved genetic sequences often have critical importance, our results thus provide compelling, if circumstantial, evidence that the stimulus independence of the oscillatory period of p65 and the time scale on which it oscillates are important for the physiological function of NF-κB.

Whether the TNF-induced dynamics of p65 localization is oscillatory on time scales longer than a few hours has been controversial. Extending the results of others (Sung et al., 2009), we now show that, in our cells, TNF can induce undamped oscillations in p65 localization for at least 20 h. In fact, the oscillations may be capable of continuing indefinitely, implying that in the presence of TNF, the network enters a limit cycle. Because our cells continue to grow and divide in the presence of TNF, the cells eventually become confluent, after which it is difficult to measure the fluorescence of individual cells. Although the precise mechanism for these sustained oscillations is unknown, the simulations of the mathematical model suggest that very low levels of active IKK are sufficient to maintain the oscillations. Consistent with this suggestion, phospho-IKK in our cells becomes undetectable by enzyme-linked immunosorbent assay after >30 min of TNF stimulation (unpublished data).

Much effort has gone into uncovering the forms and mechanisms of stimulus specificity in immune signaling, of which NF-κB is just one of several transcription factors involved (Cheng et al., 2011; Alves et al., 2014). Our findings suggest that specificity in NF-κB signaling arises from upstream signaling events, which can affect the first cycle of p65 translocation or how long the network stays active, or from the combinatorial effects of multiple transcription factors. We also cannot exclude the possibility that the posttranslational status of p65 varies with stimulus, but if it does, it does not perturb the oscillations.

Surprisingly, we also found that the invariance we observed in the population depends on variability created as single cells diverge over time. Although we previously explored the sources of variability affecting whether a cell responds to intermediate concentrations of TNF (Tay, Hughey, et al., 2010), the origins of variability in the oscillations (once a cell does respond) have remained unclear. We see no evidence that the period within a given cell changes on the time scale of our experiments, indicating that the cell-to-cell variability in the timing of the oscillations is not related to the cell cycle. It seems theoretically possible, albeit unlikely, that a stochastic event at the onset of activation determines the period of oscillation.

Most intriguingly, although the distribution of the period is conserved within the cell population, this conservation does not extend to individual cells. Our long time courses enabled us to measure accurately the distributions of oscillation period in single cells, and we found the intercellular variability to be significantly greater than the intracellular variability. These findings indicate that, although the population distribution of the period is not influenced by the multiple stimuli we tested, there are intracellular parameters that can vary among cells to tune the period. Previous work suggested that these parameters might relate to the expression of IkBa (Nelson et al., 2004) or RelA (Barken et al., 2005). This finding is also consistent with the fact that others using different cell lines reported different interpeak times for p65 oscillations, which nevertheless remain invariant across TNF concentrations (Turner et al., 2010).

Moreover, these data raise the possibility that the interpeak time is at least partially determined by epigenetic factors (e.g., protein levels or chromatin modifications), in that it is conserved for at least one division cycle but shifts over longer time scales. This phenomenon resembles what has been observed for cells undergoing apoptosis (Spencer, Gaudet, et al., 2009; Flusberg et al., 2013). At present it is difficult to determine whether, given more time, the clones would eventually show similar distributions of interpeak time to our original measurements, because another key difference is the polyclonal nature of the original cell line. Nevertheless, the combined findings that 1) the overall interpeak time distribution is conserved across populations, and 2) single cells with divergent interpeak times eventually give rise to a population with considerable intercellular variability suggest a control or feedback mechanism that warrants further investigation, characterization, and eventual inclusion in future mathematical models of NF-κB signaling.

As dynamic, single-cell data continue to be collected for other signaling networks, we expect that further invariant or recurrent features will be identified, providing new insight into the networks’ essential biological mechanisms and functions.


Cell culture and stimulation

Cell culture was performed as previously described (Gutschow, Hughey, et al., 2013). The cells were polyclonal p65/ mouse 3T3s infected with lentivirus to express p65-dsRed and H2B-GFP. Where specified, cells were stimulated with LPS-EB Ultrapure (tlrl-pelps; Invivogen, San Diego, CA), recombinant mouse TNF (11271156001; Roche, Indianapolis, IN), Pam3CSK4 (tlrl-pms; Invivogen), or Sendai virus (600503; Charles River, San Diego, CA). Solutions were prepared in imaging media (DMEM prepared without riboflavin, folic acid, or phenol red, with 1% fetal bovine serum [FBS]) and kept on ice, and then warmed to 37°C just before stimulation. Soluble TNF Receptor II (426-R2; R & D Systems, Minneapolis, MN) was used at 5 μg/ml. Clonal lines were derived from the parental strain by limiting dilution into a 96-well plate at a concentration of 1 cell/5 wells. Cells grew in 19 of the 96 wells, and based on a Poisson distribution, ∼17 (91%) of these 19 outgrowths should be derived from a single cell, that is, be clonal. We measured p65 nuclear localization in the clonal cell lines that expressed both p65-dsRed and H2B-GFP.


Figure 1A is based on the data as described in our earlier work (Tay, Hughey, et al., 2010). For the other experiments, microscopy was performed as described (Gutschow, Hughey, et al., 2013). Cells were imaged on wells of a glass-bottom, 96-well plate (164588; Nunc [Thermo Scientific], Waltham, MA) that had been precoated with 10 μg/ml human fibronectin (FC010; Millipore, Billerica, MA). The day before imaging, ∼7000 cells/well were seeded onto the plate in DMEM with 10% FBS. One hour before stimulation, the wells were switched to imaging media. Microscopy was performed on a Nikon Eclipse Ti fluorescent microscope, using a 20× air/0.75 numerical aperture objective. The camera was a Photometrics Cascade II:1024 electron-multiplying charge-coupled device. Image acquisition was controlled by Micro-Manager (Edelstein et al., 2010). Images were acquired using fluorescein isothiocyanate (FITC) and tetramethylrhodamine isothiocyanate (TRITC) filter sets (Semrock, Lake Forest, IL) every 5–6 min. Temperature (37°C), CO2 (5%), and humidity were held constant. For the 20-h time courses, Breathe-Easy film was used to minimize evaporation from the wells. In addition, for some wells in the plate, two-by-two grids of fields of view overlapping by ∼15% were imaged, to be stitched together later. Because our cells move around during the experiment, this was designed to decrease the perimeter:area ratio and thereby decrease the proportion of cells lost due to their moving out of the field of view.

Image analysis

Flat-fielding (correcting for uneven illumination of the field of view) and time-lapse registration (correcting for small imprecision of the stage movement) of images were performed with custom Matlab code. Image stitching was performed with custom Matlab code, using an iterative algorithm to stitch together adjacent fields of view. For each grid of fields of view at each time point, the stitching program first calculated the overlaps using the FITC (H2B-GFP) image and then applied those overlaps to the FITC and TRITC (p65-dsRed) images. Segmentation of nucleus and cytoplasm was performed in Matlab and in CellProfiler (Kamentsky et al., 2011). Nucleus tracking, peak finding, and curation of the data were performed using custom Matlab and Python code.

Data analysis

Except for the data in Figure 1, the first interpeak time of the p65 oscillations (whether for a cell or a simulation) was always excluded, as it tends to be longer (Supplemental Figure S1, B and C). For the analysis of the 20-h time courses, cell divisions were excluded, unless otherwise noted, and only cells that showed at least nine peaks of nuclear p65 were used for analysis. For both experiment and modeling, a bootstrapped (randomized) distribution of mean interpeak times was generated by pooling all interpeak times from the respective data set and then drawing n samples (with replacement) 10,000 times, where n is the mean number of interpeak times per cell or per simulation. A bootstrapped distribution for the F statistic (between-group variance divided by within-group variance) was calculated in the same way. Where noted, cell divisions (from the same 20-h experiments) were identified manually. Because the identity of the two daughter cells is interchangeable, a correlation metric, which is normally between two distinct variables, is not applicable. Thus we used the absolute difference between two cells. For the analysis of the daughter cells, for which we generally had fewer observed peaks, only cells with at least six peaks were used, and at most one outlier (typically an error in the peak finding) was removed according to the Thompson tau method.

Mathematical modeling

We used the stochastic model of the NF-κB pathway as described. We first ran the model without any TNF to find the median initial conditions and then ran the 500 simulations using those initial conditions. To calculate the sensitivity of the interpeak time to the model parameters, each parameter was independently increased or decreased twofold, and 50 simulations of 20-h stimulation with TNF were run. For each parameter, the sensitivity of the mean interpeak time was calculated as the difference between the mean interpeak time (across all simulations) for the twofold increase and mean interpeak time for the twofold decrease, divided by the mean interpeak time for the baseline value of the parameter. The sensitivity of the CV was calculated analogously.


This article was published online ahead of print in MBoC in Press ( on December 3, 2014.

*These authors contributed equally to this work.

Abbreviations used:

inhibitor of κB




melanoma differentiation-associated protein 5


myeloid differentiation primary response 88


nuclear factor κB


retinoic acid-inducible gene 1


Toll-like receptor


tumor necrosis factor (receptor)


TIR (Toll/Interleukin-1 receptor)-domain-containing adapter-inducing interferon-β.


We are grateful for support from the National Institutes of Health, including a Director's Pioneer Award (5DP1LM01150), a National Cancer Institute Pathway to Independence Award (R00CA125994), an R21 grant (1R21AI104305), and a Systems Biology Center grant (P50 GM107615), in addition to a Paul Allen Family Foundation Distinguished Investigator Award to M.W.C., Bio-X and ARCS Fellowships to J.J.H., and Weiland and Rensselaer Engineering Fellowships to M.V.G. We thank Derek Macklin and Harendra Guturu for developing the interactive figure and thank the Covert lab for helpful discussions and advice.


  • AlvesBN, TsuiR, AlmadenJ, ShokhirevMN, Davis-TurakJ, FujimotoJ, HoffmannA (2014). IκBε is a key regulator of B cell expansion by providing negative feedback on cRel and RelA in a stimulus-specific manner. J Immunol 192, 3121–3132. Crossref, MedlineGoogle Scholar
  • AshallL, HortonCA, NelsonDE, PaszekP, HarperCV, SillitoeK, RyanS, SpillerDG, UnittJF, BroomheadDS et al. (2009). Pulsatile stimulation determines timing and specificity of NF-κB-dependent transcription. Science 324, 242–246. Crossref, MedlineGoogle Scholar
  • BarkenD, WangCJ, KearnsJ, CheongR, HoffmannA, LevchenkoA (2005). Comment on “Oscillations in NF-kappaB signaling control the dynamics of gene expression.”. Science 308, 52. Crossref, MedlineGoogle Scholar
  • ChengCS, FeldmanKE, LeeJ, VermaS, HuangDB, HuynhK, ChangM, PonomarenkoJV, SunSC, BenedictCA et al. (2011). The specificity of innate immune responses is enforced by repression of interferon response elements by NF-κB p50. Sci Signal 4, ra11. Crossref, MedlineGoogle Scholar
  • CheongR, RheeA, WangCJ, NemenmanI, LevchenkoA (2011). Information transduction capacity of noisy biochemical signaling networks. Science 334, 354–358. Crossref, MedlineGoogle Scholar
  • CovertMW, LeungTH, GastonJE, BaltimoreD (2005). Achieving stability of lipopolysaccharide-induced NF-kappaB activation. Science 309, 1854–1857. Crossref, MedlineGoogle Scholar
  • EdelsteinA, AmodajN, HooverK, ValeR, StuurmanN (2010). Computer control of microscopes using μManager. Curr Protoc Mol Biol Chapter 14, Unit 14.20. Crossref, MedlineGoogle Scholar
  • ElowitzMB, LevineAJ, SiggiaED, SwainPS (2002). Stochastic gene expression in a single cell. Science 297, 1183–1186. Crossref, MedlineGoogle Scholar
  • FlusbergDA, RouxJ, SpencerSL, SorgerPK (2013). Cells surviving fractional killing by TRAIL exhibit transient but sustainable resistance and inflammatory phenotypes. Mol Biol Cell 24, 2186–2200. LinkGoogle Scholar
  • Geva-ZatorskyN, RosenfeldN, ItzkovitzS, MiloR, SigalA, DekelE, YarnitzkyT, LironY, PolakP, LahavG, AlonU (2006). Oscillations and variability in the p53 system. Mol Syst Biol 2, 2006.0033. Crossref, MedlineGoogle Scholar
  • GitlinL, BenoitL, SongC, CellaM, GilfillanS, HoltzmanMJ, ColonnaM (2010). Melanoma differentiation-associated gene 5 (MDA5) is involved in the innate immune response to Paramyxoviridae infection in vivo. PLoS Pathog 6, e1000734. Crossref, MedlineGoogle Scholar
  • GutschowMV, HugheyJJ, RuggeroNA, BajarBT, ValleSD, CovertMW (2013). Single-cell and population NF-κB dynamic responses depend on lipopolysaccharide preparation. PLoS One 8, e53222. Crossref, MedlineGoogle Scholar
  • HaydenMS, WestAP, GhoshS (2006). NF-κB and the immune response. Oncogene 25, 6758–6780. Crossref, MedlineGoogle Scholar
  • HirschfeldM, KirschningCJ, SchwandnerR, WescheH, WeisJH, WootenRM, WeisJJ (1999). Inflammatory signaling by Borrelia burgdorferi lipoproteins is mediated by toll-like receptor 2. J Immunol 163, 2382–2386. MedlineGoogle Scholar
  • HoffmannA, LevchenkoA, ScottML, BaltimoreD (2002). The IκB-NF-κB signaling module: temporal control and selective gene activation. Science 298, 1241–1245. Crossref, MedlineGoogle Scholar
  • KamentskyL, JonesTR, FraserA, BrayMA, LoganDJ, MaddenKL, CarpenterAE (2011). Improved structure, function and compatibility for CellProfiler: modular high-throughput image analysis software. Bioinformatics 27, 1179–1180. Crossref, MedlineGoogle Scholar
  • KatoH, TakeuchiO, SatoS, YoneyamaM, YamamotoM, MatsuiK, UematsuS, JungA, KawaiT, IshiiKJ et al. (2006). Differential roles of MDA5 and RIG-I helicases in the recognition of RNA viruses. Nature 441, 101–105. Crossref, MedlineGoogle Scholar
  • KawaiT, AkiraS (2006). TLR signaling. Cell Death Differ 13, 816–825. Crossref, MedlineGoogle Scholar
  • LeeTK, DennyEM, SanghviJC, GastonJE, MaynardND, HugheyJJ, CovertMW (2009). A noisy paracrine signal determines the cellular NF-kappaB response to lipopolysaccharide. Sci Signal 2, ra65. Crossref, MedlineGoogle Scholar
  • LeiseTL, WangCW, GitisPJ, WelshDK (2012). Persistent cell-autonomous circadian oscillations in fibroblasts revealed by six-week single-cell imaging of PER2::LUC bioluminescence. PLoS One 7, e33334. Crossref, MedlineGoogle Scholar
  • LipniackiT, PaszekP, BrasierAR, LuxonB, KimmelM (2004). Mathematical model of NF-κB regulatory module. J Theor Biol 228, 195–215. Crossref, MedlineGoogle Scholar
  • LongoDM, SelimkhanovJ, KearnsJD, HastyJ, HoffmannA, TsimringLS (2013). Dual delayed feedback provides sensitivity and robustness to the NF- κB signaling module. PLoS Comput Biol 9, e1003112. Crossref, MedlineGoogle Scholar
  • NelsonDE, IhekwabaAE, ElliottM, JohnsonJR, GibneyCA, ForemanBE, NelsonG, SeeV, HortonCA, SpillerDG et al. (2004). Oscillations in NF-kappaB signaling control the dynamics of gene expression. Science 306, 704–708. Crossref, MedlineGoogle Scholar
  • OzinskyA, UnderhillDM, FontenotJD, HajjarAM, SmithKD, WilsonCB, SchroederL, AderemA (2000). The repertoire for pattern recognition of pathogens by the innate immune system is defined by cooperation between toll-like receptors. Proc Natl Acad Sci USA 97, 13766–13771. Crossref, MedlineGoogle Scholar
  • PahlHL (1999). Activators and target genes of Rel/NF-kB transcription factors. Oncogene 18, 6853–6866. Crossref, MedlineGoogle Scholar
  • PaszekP, RyanS, AshallL, SillitoeK, HarperCV, SpillerDG, RandDA, WhiteMRH (2010). Population robustness arising from cellular heterogeneity. Proc Natl Acad Sci USA 107, 11644–9. Crossref, MedlineGoogle Scholar
  • Pe˛kalskiJ, ZukPJ, Kochan´czykM, JunkinM, KelloggR, TayS, LipniackiT (2013). Spontaneous NF-κB activation by autocrine TNFα signaling: a computational analysis. PLoS One 8, e78887. Crossref, MedlineGoogle Scholar
  • PurvisJE, KarhohsKW, MockC, BatchelorE, LoewerA, LahavG (2012). P53 dynamics control cell fate. Science 336, 1440–1444. Crossref, MedlineGoogle Scholar
  • ShinoharaH, BeharM, InoueK, HiroshimaM, YasudaT, NagashimaT, KimuraS, SanjoH, MaedaS, YumotoN et al. (2014). Positive feedback within a kinase signaling complex functions as a switch mechanism for NF-κB activation. Science 344, 760–764. Crossref, MedlineGoogle Scholar
  • SpencerSL, GaudetS, AlbeckJG, BurkeJM, SorgerPK (2009). Non-genetic origins of cell-to-cell variability in TRAIL-induced apoptosis. Nature 459, 428–432. Crossref, MedlineGoogle Scholar
  • SungMH, LiN, LaoQ, GottschalkRA, HagerGL, FraserID (2014). Switching of relative dominance between feedback mechanisms in lipopolysaccharide-induced NF-κB signaling. Sci Signal 7, ra6. Crossref, MedlineGoogle Scholar
  • SungMH, SalvatoreL, De LorenziR, IndrawanA, PasparakisM, HagerGL, BianchiME, AgrestiA (2009). Sustained oscillations of NF-κB produce distinct genome scanning and gene expression profiles. PLoS One 4, e7163. Crossref, MedlineGoogle Scholar
  • TayS, HugheyJJ, LeeTK, LipniackiT, QuakeSR, CovertMW (2010). Single-cell NF-kappaB dynamics reveal digital activation and analogue information processing. Nature 466, 267–271. Crossref, MedlineGoogle Scholar
  • TurnerDA, PaszekP, WoodcockDJ, NelsonDE, HortonCA, WangY, SpillerDG, RandDA, WhiteMRH, HarperCV (2010). Physiological levels of TNFalpha stimulation induce stochastic dynamics of NF-kappaB responses in single living cells. J Cell Sci 123, 2834–2843. Crossref, MedlineGoogle Scholar