A variable-gain stochastic pooling motif mediates information transfer from receptor assemblies into NF-κB

img ]

A myriad of inflammatory cytokines regulate signaling pathways to maintain cellular homeostasis. The IκB kinase (IKK) complex is an integration hub for cytokines that govern nuclear factor κB (NF-κB) signaling. In response to inflammation, IKK is activated through recruitment to receptor-associated protein assemblies. How and what information IKK complexes transmit about the milieu are open questions. Here, we track dynamics of IKK complexes and nuclear NF-κB to identify upstream signaling features that determine same-cell responses. Experiments and modeling of single complexes reveal their size, number, and timing relays cytokine-specific control over shared signaling mechanisms with feedback regulation that is independent of transcription. Our results provide evidence for variable-gain stochastic pooling, a noise-reducing motif that enables cytokine-specific regulation and parsimonious information transfer. We propose that emergent properties of stochastic pooling are general principles of receptor signaling that have evolved for constructive information transmission in noisy molecular environments.

Here, we develop genetically modified cells that endogenously express fluorescent protein (FP) fusions of NEMO and RelA, allowing same-cell measurements of CI-like structures and canonical NF-kB signaling from live-cell images. We establish differences between TNF and IL-1 responses in biophysical properties of NEMO complexes and demonstrate a continuum relating CI-like structures and downstream NF-kB responses in the same cell. By tracking single complexes, we demonstrate that (i) cytokine dosage and time-varying presentation modulate the timing and numbers of CI-like structures, (ii) single complexes have switch-like activation profiles where the aggregate of NEMO recruitment and time-varying properties of each complex are cytokine specific, and (iii) dynamics of formation and dissolution for single complexes during the primary cytokine response are independent of transcriptional feedback. Last, we characterize a signaling motif called a variable-gain stochastic pooling network (SPN) that encompasses our experimental observations. The variable-gain SPN (VG-SPN) motif has beneficial noise mitigation properties and provides a trade-off between information fidelity, ligand specificity, and resource allocation for intracellular signaling molecules. We propose that the VG-SPN architecture and its associated benefits to signal transmission are common mechanisms for receptor-mediated signal transduction.

When observed in single cells exposed to inflammatory stimuli, the RelA subunit of NF-κB encodes a dynamic transcriptional signal by translocating from the cytoplasm into the nucleus ( 2 , 3 , 22 – 24 ). Models calibrated to single-cell RelA data ( 23 – 26 ) have revealed numerous transcriptional mechanisms and emergent properties that place the NF-κB pathway among exemplars of dynamical biological systems ( 27 , 28 ). Key to these findings are two mediators of negative feedback, IκBα and A20, which are transcriptionally regulated by NF-κB. IκBα restores NF-κB to its baseline cytoplasmic localization through nuclear export and sequestration, whereas A20 diminishes kinase activation upstream of NF-κB through disassembly of CI-like structures in addition to noncatalytic mechanisms ( 10 , 23 , 25 , 26 , 29 ). Dynamical regulation of transcription and feedback via NF-κB is strongly recapitulated between models and experiments; however, there is a dearth of quantitative single-cell data at the level of cytokine detection and dynamical properties of CI-like complexes to substantiate our understanding of upstream signal transmission.

Ligation of TNF to the TNF Receptor 1 (TNFR1) recruits adaptor proteins and enzymes to form a large multiprotein complex near the plasma membrane ( 6 – 9 ). Ubiquitin-modifying enzymes are critical components that assemble linear, branched, and mixed polyubiquitin scaffolds around the multiprotein complex ( 10 – 13 ). The NF-κB essential modulator (NEMO) subunit of the cytoplasmic IκB kinase (IKK) complex is rapidly recruited via direct interaction with the polyubiquitin scaffold and accessory proteins, where IKK is activated through induced proximity with regulatory kinases ( 4 , 14 – 16 ). The fully assembled TNFR1 complex, referred to as “complex I” [CI; ( 6 )], is a master regulator of inflammation-dependent NF-κB signaling. Although other inflammatory molecules such as IL-1 signal through CI-like complexes using different receptors, adaptor proteins ( 17 , 18 ), and varying compositions of ubiquitin chain scaffolds ( 19 , 20 ), all regulate NF-κB through IKK activation mediated by induced proximity with other signaling mediators that reside on CI ( 21 ).

A limited number of transmembrane receptors expressed on the cell surface mediate crucial transmission of information between extracellular and intracellular signaling molecules. Key questions are understanding the mechanisms and limitations that underlie signal transmission, in particular for cytokine receptor signaling that is often deregulated in disease. The nuclear factor κB (NF-κB) signaling pathway is an archetypal molecular communication channel that transmits information about extracellular cytokines to regulate cellular adaptation through activation of the RelA transcription factor ( 1 – 3 ). When ligated with inflammatory molecules, such as tumor necrosis factor (TNF) and interleukin-1β (IL-1), many activated receptors converge on NF-κB signaling ( 4 , 5 ).

To test predictions of the VG-SPN, we first simulated R max values expected for TNF and IL-1 responses estimated from experimental G L and CI max values ( Figs. 2 and 7D ). We found that these values were consistent with AUC i values for NEMO in cells exposed to 1000 ng/ml of either cytokine ( Fig. 7D ). Next, we similarly estimated channel capacity values expected for TNF and IL-1 responses. Even though the response magnitude of IL-1 is significantly greater than for TNF ( Fig. 7D ), the VG-SPN model predicted that the system would be robust to variation in G L and CI max values, and the channel capacity for TNF will be only marginally smaller ( Fig. 7E ). Because throughput for NEMO imaging experiments is insufficient for channel capacity calculations, we instead measured the channel capacity downstream using dynamics of nuclear RelA (fig. S13) as a system output ( 2 , 3 ). Model predictions were again similar with the experimental channel capacity ( Fig. 7E ). Together, the VG-SPN motif generates trade-offs, and although the motif results in hard limits for the information content of a signaling pathway, it is also parsimonious and robust to stochastic noise and cell-to-cell variability in protein abundances.

Next, we explored the information transmission properties for different VG-SPN configurations by calculating their channel capacities ( 2 ). We assessed the VG-SPN model using parameters obtained from our single-cell IKK data and simplifying assumptions about distributions for CI properties (see Materials and Methods for details, and hyperparameter tuning in fig. S12). For these simulations, CI max per cell and amplifier gain G L were allowed to vary independently. Consistent with our previous analysis for noise propagation in VG-SPNs, when G L is greater than 10 molecules per complex, the channel capacity increases with the number of CI max per cell and rapidly approaches saturation ( Fig. 7C , bottom). The maximum channel capacity was found to be 2.5 bits with the model calibrated to live cell data (fig. S12), suggesting a theoretical maximum for information transmission at the level of CI. However, for G L values of 10 or lower, the system requires orders of magnitude more CI max per cell to reach comparable channel capacity values. We also calculated the R max that shows an inverse linear relationship, where different configurations of CI max per cell and G L can achieve activation of the same number of response molecules up to arbitrarily high numbers ( Fig. 7C , top).

( A ) Different configurations of the VG-SPN. For a mathematically controlled comparison, CI and “G L ” are adjusted so that each configuration can achieve the same R max . ( B ) Detection introduces shot noise (top), where a switch activates erroneously with a probability given by the SNL. Because each CI is independent, shot noise is mitigated through parallelization via multiple CI until a “noise floor” is achieved. The noise floor is determined by the shot noise value. Variability in the gain from each CI switch (bottom) introduces noise during transmission. Gain noise is effectively averaged and approaches zero via parallelization of independent switches. ( C ) R max (upper limit was set to 10 7 response molecules per cell) and relative channel capacity corresponding to different VG-SPN configurations (see also fig. S12). The VG-SPN architecture provides ligand-specific tunability of response magnitude and robustness to noise. ( D ) Simulations (red bars) of the R max expected for IL-1 and TNF responses based on experimental values (tan bars) for CI max and Gain. Experimental R max is measured as the AUC of NEMO spot intensity in cells exposed to 1000 ng/ml of cytokine as indicated. ( E ) Simulations (red bars) of the channel capacity expected for IL-1 and TNF responses based on experimental values for CI max and Gain. Experimental channel capacities (tan bars) calculated downstream of NEMO via nuclear RelA dose-response dataset (see also fig. S13). All error bars represent SEM.

Identification of the VG-SPN architecture enabled us to conceptually characterize signal transduction properties that could not be measured from single-cell IKK data due to limited experimental throughput. We first considered the impact of the number of CI switches and their associated gains on noise propagation in the resulting VG-SPN. A mathematically controlled comparison ( 28 , 52 ) between different network configurations was enabled by assuming each configuration is capable of producing the same maximal steady-state response [maximum response (R max ); Fig. 7A ]. Simulations for different mathematically controlled configurations of the VG-SPN demonstrated that shot noise associated with signal detection and noise associated with the signal gain both fall off rapidly with increasing numbers of CI switches ( Fig. 7B ). Noise from these sources can be mitigated almost completely when cells are capable of forming ~100s to 1000s of CI per cell (CI max per cell). Here, noise mitigation benefits can be attributed to signal parallelization through quantized nodes. For example, gain noise at each CI can be positive or negative, and because each CI is independent, the distribution for noise across all CI in a cell is centered near zero. Summation of signaling through R max in a VG-SPN effectively averages the contribution of noise across all CI per cell, which is ideally zero assuming sufficient parallelization and absence of other biases. See also Materials and Methods for similar description of shot noise mitigation.

Together, our data are consistent with the architecture presented in Fig. 6A . Here, the number of CI max is receptor specific and varies between single cells, and the cell’s response to a subsaturating “S” is fractional to its CI max . Together with same-cell NEMO and RelA data, our results support a generalization of cytokine–IKK–NF-κB signaling that is evocative of an SPN, a model sensory system with noise-mitigating and information-compressing properties ( 50 , 51 ). An important difference, however, is that while binary detectors in an SPN transmit on-off measurements about an information source, each NEMO-recruiting complex performs amplification with a gain determined by the cytokine-receptor-complex identity. Henceforth, we refer to the network architecture as a VG-SPN.

The average number of NEMO spots within cells increases with cytokine concentration, yet we see significant differences between single-cell responses in all conditions ( Fig. 2F ). To determine whether some cells are predisposed to stronger responses or whether noise at the level of CI formation is predominant contributors to heterogeneity, we exposed cells to two sequential pulses of cytokines. Two-pulse experiments that compare the same cell in different conditions can be used to interpret the impact of stochastic noise on cell-to-cell heterogeneity ( 3 , 46 – 49 ). For these experiments, cells were first stimulated with a short low-concentration reference pulse of cytokine, followed by a nearly saturating concentration of cytokine to estimate the maximum number of CI complexes (CI max ) that each cell can produce in these conditions ( Fig. 6 ). In cells exposed to the same cytokine for both stimuli, the rank ordering of single-cell responses to the reference stimulus is strongly correlated with CI max for the same cell ( Fig. 6C ; Spearman ρ > 0.70, P < 10 −10 ). We then asked whether the correlation was due to shared mechanisms downstream of TNF and IL-1 receptors in cross-pulse experiments where cells are exposed sequentially to both cytokines. Same-cell correlations in cross-pulse experiments were greatly reduced with marginal to nonsignificant P values (fig. S11). These data suggest that cells are predisposed to different responses based on cell-to-cell variability in receptor-specific components of CI complexes.

( A ) The VG-SPN consists of four phases: detection, switching, amplification, and signal pooling. For detection, each CI-like complex (CI) is modeled as one of “j” independent switches with an activation probability proportionate to the signal strength. During switching and amplification, each activated switch amplifies the signal with gain “G L ,” related to the number of NEMO molecules recruited by each CI-like complex. Subsequently, NEMO from each CI-like complex is added to the cytoplasmic “pool” of activated signaling molecules. ( B ) Schematic for a two-pulse experiment. Cells are first exposed to a short- and low-concentration reference pulse of cytokine. Subsequently, the same cells are exposed to a second pulse of cytokine at saturating concentrations intended to reveal the maximum number of CI complexes that can be generated by each cell (CI max ). ( C ) Single cells are rank ordered by their response to reference cytokine stimulation and correlated with their CI max values. Responses to TNF (top) and IL-1 (bottom) are shown. Inset numbers indicate Spearman’s rank correlation (ρ) and associated P values (see also fig. S11).

When TNF and IL-1 responses are considered separately, trajectories of EGFP-NEMO spots are highly similar regardless of when they form during a cytokine response. This observation argues that each complex behaves as an independent switch, that when activated recruits a quantized amount of NEMO over its life span. To understand information transmission properties of the IKK–NF-κB signaling axis, we abstracted the system into four phases ( Fig. 6A ). Detection (i), where parallel and independent CI-like switches provide redundant measurements of the same extracellular signal (S). Switching (ii) and amplification (iii), when a CI-like complex (CI) activates in response to S, it amplifies with a ligand-specific and quantized gain (G L ) of NEMO activity. The total cellular response (R) in terms of NEMO activity and subsequent NF-κB translocation is the summation (iv) of all amplifier gains in a cell, also referred to as “signal pooling.” Although our live-cell data support independent binary switches, quantized amplification, and signal pooling of IKK into NF-κB ( Figs. 3 to 5 ), conditions for detection (i) and switching (ii) remained ambiguous because of cell-to-cell heterogeneity.

On the basis of simulations, if transcription is the predominant source of negative feedback, then spots that form later after cytokine stimulation are predicted to have lower maximum intensity, shorter track length, and smaller AUC ( Fig. 5B and fig. S8). To test the prediction, we performed a reverse time-course experiment where imaging started after a delay relative to the time of cytokine stimulation (0, 5, 10, and 15 min; Fig. 5C ). Only new spots that formed within the first 2 min were tracked for each condition, thereby enabling direct comparison of early- versus late-forming spots and minimizing effects of photobleaching. Biological replicates revealed that early- and late-forming single-spot trajectories share highly similar dynamics ( Fig. 5, B and C , and fig. S9), suggesting that transcriptional feedback is dispensable in regulation of dynamics for CI-like complexes. We therefore repeated reverse time-course experiments in the presence of cycloheximide (CHX) to prevent translation of NF-κB–regulated genes, effectively breaking the transcriptional negative feedback loop. To verify that CHX inhibited protein translation, time courses of nuclear NF-κB were also measured in the same cells during costimulation with IL-1 and CHX. NF-κB showed persistent nuclear localization in these cells, indicating a disruption of the negative feedback loop mediated by IκBα protein expression (fig. S10). Loss of transcriptional feedback did not increase features of single-spot trajectories as predicted by the transcriptional feedback model ( Fig. 5D and fig. S8). These results demonstrate that transcription is not the predominant mechanism of negative feedback on NEMO recruitment at CI-like complexes in the time scale of the primary cytokine response.

Using simulations to vary the strength for each source of negative feedback, their respective impacts on single-spot dynamics was apparent ( Fig. 5B and fig. S8). By increasing the strength of basal feedback, the decay phase for single-spot trajectories became steeper, and each spot displayed a sharp peak of intensity that was similar between spots regardless of when they form. By contrast, even though transcriptional feedback also increased steepness of the decay phase, the peak intensity of spots that form earlier was substantially greater than that of spots that formed later. For both sources, increasing strength of negative feedback reduces the overall peak height for all simulated spots (fig. S8).

( A ) Schematic of the HyDeS model of basal (red arrows) and transcription-dependent (blue arrows) feedback regulation of NEMO complexes. Each growing CI-like complex will recede with rates that depend on recruitment of DUB enzymes. ( B ) Simulations of individual EGFP-complex trajectories using the model in Fig. 5A , considering no feedback (top left), transcription-dependent feedback only (top right), basal feedback only (bottom left), or combined transcriptional and basal feedback (bottom right). ( C ) Schematic of reverse time-course experiments. Only new EGFP-NEMO complexes that form within the first 2 min are tracked in movies from each high-frequency imaging experiment. Cells were stimulated with 100 ng/ml concentrations of either IL-1 (top) or TNF (bottom). ( D ) Boxplots of maximum intensity of individual EGFP-NEMO complex trajectories after stimulation with IL-1 (top) or TNF (bottom) in cells pretreated for 20 min with media plus dimethyl sulfoxide (left) or CHX (right). Costimulation with CHX does not significantly increase values for single-spot descriptors as predicted by the transcriptional feedback model (see also fig. S5). Median and interquartile ranges are indicated, and biological replicates are shown side by side with increased transparency for each experiment.

To understand how different mechanisms of negative feedback affect trajectories of single EGFP-NEMO spots, we developed a model using ordinary differential equations. Here, cytokine stimulation induces formation of single spots at different times, and each spot becomes larger and brighter with Michaelis-Menten kinetics to approximate ubiquitin polymer growth and EGFP-NEMO recruitment. We reasoned that the basal expression of A20 and other DUBs might affect the half-life and activity of NEMO puncta. The model, therefore, considers two sources of negative feedback that act on NEMO-recruiting complexes by enhancing their disassembly rates. The first source aggregates the sum of NEMO-recruiting complexes to drive expression for an A20-like negative feedback mediator (“transcriptional feedback”; Fig. 5A ). The second source considers the impact from basal expression of A20 and other DUBs in resting cells before stimulation (“basal feedback”; Fig. 5A ). For both sources, the strength of negative feedback on each NEMO-recruiting spot increases with size to mimic DUB recruitment to ubiquitin polymers in the complex.

NF-κB–mediated expression of A20 and subsequent deubiquitinating (DUB) activity against NEMO-recruiting complexes constitute an essential negative feedback motif in inflammatory signaling ( 22 ). Within tens of minutes, it is feasible that nascent A20 contributes to the decay phase of EGFP-NEMO spot numbers observed in whole-cell measurements ( Fig. 2E ), thereby reducing IKK activation as observed in cell population assays ( 26 , 44 ). It was, therefore, unexpected that the decay phase for single EGFP-NEMO tracks is visible within several minutes of stimulation, which is fast for a feedback mechanism governed by transcription and translation ( Fig. 4B and fig. S7).

Descriptors for single-spot trajectories showed low variability, in particular for the TNF response, when compared within the same cell or when the average of single-spot descriptors was compared between cells (fig. S7, C and D). Low interspot variability combined with saturating spot numbers ( Fig. 2F ) and expectations from surface receptor expression ( Fig. 1 ) together suggest that most spots are single receptor–ligand complexes. Analysis of the “AUC spot intensity” (AUC i ) descriptor revealed a quadratic relationship between “mean AUC i ” and “AUC i variance” for spots when measured in the same cell (fig. S7E). Increased variance between large NEMO-recruiting complexes may be due to steric properties of supramolecular assemblies, for example, where growth for certain types of ubiquitin polymers is spatially limited, or intact portions of large ubiquitin chains are clipped off en bloc or through endocleavage ( 13 ), leading to greater interspot variability. In nearly all cases, the coefficient of variation values for distributions of single-spot descriptors indicate that noise is lower than Poisson processes (fig. S7). Together, our results suggest that dynamics of NEMO-recruiting protein complexes are strictly regulated for each cytokine response.

Tracking experiments revealed differences in fluorescence intensity time courses, where IL-1–induced EGFP-NEMO puncta were consistently brighter and longer lived than TNF-induced puncta ( Fig. 4B ). Both cytokines induce spots that peak within 2 to 3 min of detection, followed by a decay phase where spots decline in intensity ( Fig. 4B and fig. S7B). When stimulated with a cytokine pulse, most spots are detected only after the cytokine is removed ( Fig. 4C ), demonstrating that the formation of NEMO-recruiting complexes is variable and takes place up to 30 min following a stimulus. Step, pulse, and ramp stimulation further showed that the number of spots and the timing of spot formation are both modulated by the dynamics of cytokine presentation.

( A ) Maximum intensity projections from high-frequency 3D time-lapse imaging experiment of cells exposed to IL-1 (100 ng/ml). Colored lines in overlay indicate tracks for individual EGFP-NEMO complexes over time. On average, 60 to 75% of detected complexes are associated with a track of significant length (fig. S4A). Scale bar, 20 μm; see also movie S2. ( B ) Time courses of fluorescence intensity for single-tracked EGFP-NEMO complexes in cells stimulated with indicated conditions in a microfluidic cell culture system ( 45 ). Two representative single-complex trajectories are indicated in each condition (dark orange and dark blue lines). ( C ) Bar graphs for the average number and the time of formation for tracked EGFP-NEMO complexes. On average, 10 cells were analyzed in each condition. Error bars represent the SEM.

Asynchronous properties of EGFP-NEMO spots, such as when a spot forms or rates for NEMO recruitment and dissolution, will contribute to variability when spot intensities are measured from a snapshot image at a single time point. To understand the extent of interspot and between-cell variability, we used high-frequency imaging to enable tracking of each single spot over time ( Fig. 4 , fig. S7A, and movie S3). Cells were stimulated with a step, pulse, or ramp in cytokine concentration in a microfluidic cell culture system ( 45 ) to observe how dynamic environments modulate properties of NEMO-recruiting complexes.

Overall, same-cell measurements of EGFP-NEMO and mCh-RelA reveal that the aggregate NF-κB response in a cell is well determined by the sum of NEMO recruitment to CI-like signaling complexes near the plasma membrane. Because enrichment of NEMO at ubiquitin-rich structures is an induced-proximity mechanism for kinase activation, it is reasonable that the AUC of EGFP-NEMO intensity in puncta is a strong proxy for downstream signaling. However, our characterization for the EGFP-NEMO fluorescence intensity at cytokine-induced spots showed wide-based distributions that could indicate that the amount of NEMO recruited at each spot varies significantly ( Fig. 2C ). Nevertheless, the number of EGFP-NEMO puncta (MAX NEMO ) is almost interchangeable with AUC NEMO as determinants of same-cell responses.

Previously, we showed that the “area under the curve” (AUC) and “maximum” (MAX) are scalar descriptors of a nuclear RelA fold change time course that encode the most information about cytokine dose ( 3 ). Notably, among all descriptors for nuclear RelA fold change dynamics, AUC RelA and MAX RelA also had the strongest correlation with same-cell descriptors of NEMO-recruiting complexes ( Fig. 3D ; AUC NEMO and MAX NEMO , respectively). Both NEMO descriptors showed similar coefficients of determination for RelA descriptors, whether measured in terms of numbers for EGFP-NEMO spots or aggregate intensity of EGFP-NEMO within complexes (fig. S5).

( A ) Maximum intensity projections from 3D time-lapse images of EGFP-NEMO (top) and mCherry-RelA (bottom) expressed endogenously in the same U2OS cells. Cells were stimulated with IL-1 (10 ng/ml). ( B ) Diagrams indicating time-course descriptors for EGFP-NEMO complexes (top) and nuclear fold change of mCherry-RelA (bottom) in single cells. RelA features are similar to those used previously ( 3 , 24 ). ( C ) Example of a strong same-cell correlation between single-cell descriptors of EGFP-NEMO complexes and nuclear RelA fold change. IL-1 and TNF responses across all doses overlap and form a continuum that relates the paired descriptors. ( D ) Heatmap of Spearman’s rank correlation coefficients (ρ) for all pairs of same-cell descriptors for IL-1 (left) and TNF responses (right). See also fig. S5 and table S2 for summary of cell numbers per condition. FWHM, full weight at half maximum.

Next, we investigated the relationship between cytokine-induced dynamics of NEMO and RelA to ask whether NF-κB responses are well determined by properties of fluorescent NEMO puncta when measured in the same cell. Time courses for NEMO complexes and nuclear RelA localization were measured in response to a broad range of cytokine concentrations, and quantitative descriptors that summarize dynamic properties of EGFP-NEMO and mCh-RelA were extracted for each single cell [ Fig. 3, A and B ; see also ( 3 , 24 )]. Scatterplots of descriptors showed that cytokines and concentrations together form a continuum with a strong monotonic relationship between descriptors of NEMO and RelA ( Fig. 3C and fig. S5). We also compared coefficients for determination between NEMO and RelA descriptors, which showed stronger correlations when data are log transformed (fig. S5). Increased R 2 values in log space is because the NEMO-RelA relationship is likely governed by a power law, which is expected for a signal amplification mechanism with values that scale across orders of magnitude. Correlations between same-cell descriptors of NEMO and RelA revealed a strong effect size indicating that cell-to-cell variability in NF-κB response can be reasonably determined by descriptors of EGFP-NEMO ( Fig. 3, C and D ). Multiple linear regression for combinations of NEMO descriptors only marginally increased correlations (fig. S6).

Time courses for NEMO complexes in single cells showed a peak in spot numbers between 10 and 20 min and a rapid falloff thereafter ( Fig. 2E ), consistent with previous results from cytokine-induced IKK kinase assays ( 26 , 44 ). Numbers of NEMO spots per cell increased with cytokine concentration and showed a tendency of higher numbers in response to IL-1 at comparable molarities ( Fig. 2, E and F . Together, these data indicate that the size and intensity of NEMO complexes depend on the type of engaged receptor, whereas the number of complexes is modulated by cytokine dose.

( A ) Maximum intensity projections from 3D time-lapse images of endogenous EGFP-NEMO show rapid recruitment to CI-like complexes in cells exposed to IL-1 or TNF. See also movies S1 and S2. ( B ) Details of fluorescent complexes [orange and blue boxes in (A)] show differences between responses to IL-1 and TNF. Scale bar, 20 μm for all. ( C ) Histograms summarizing the size (left) or intensity (right) of EGFP-NEMO complexes across concentrations of IL-1 (orange) or TNF (blue). Distributions represent single-cell data in aggregate from three to five experiments. Vertical bar indicates the median of each population. Comparison between all conditions, IL-1–induced complexes are larger and brighter than TNF-induced complexes (P « 10 −20 ; t test). Numbers of cells analyzed and associated spot numbers are provided in table S1. a.u., arbitrary units. ( D ) Boxplot for estimates of the total number of EGFP-NEMO molecules in each CI-like complex. Median and interquartile ranges are indicated. See also Materials and Methods and fig. S2D. ( E ) Single-cell time courses for the number of EGFP-NEMO complexes in cells exposed to indicated concentrations of IL-1 (orange) and TNF (blue). ( F ) Dose response of maximum number of EGFP-NEMO complexes. Dark orange and blue lines represent the median, and vertical bars represent the interquartile range for cells stimulated with IL-1 or TNF, respectively.

In response to TNF or IL-1, enhanced green FP (EGFP)–NEMO transiently localizes to punctate structures near the plasma membrane ( 19 , 20 , 33 ) that are distinct from endosomal structures ( Fig. 2 and fig. S3, see also movie S1). To further characterize the role of cytokine identity and dose on NEMO recruitment at CI-like puncta, we compared descriptive features such as their size and intensity across different cytokines and concentrations. Although properties of IL-1– and TNF-induced puncta did not show a clear trend across doses, IL-1–induced spots were significantly larger and brighter than their TNF-induced CI counterparts ( Fig. 2C and movie S2; for all comparisons, P « 10 −10 , Student’s t test). To estimate the expected number of NEMO molecules in each complex, we evaluated intensity values for each fluorescent spot in terms of a reference live-cell reporter that recruits a known number of EGFP molecules into a diffraction-limited space ( 43 ). By comparing cells in identical imaging conditions, our analysis suggests that each of the larger IL-1–induced spots recruit approximately 300 NEMO molecules, whereas TNF-induced spots recruit around 80 ( Fig. 2D and fig. S4).

We set out to investigate how cytokine receptors engage NEMO as an integration hub to regulate NF-κB signaling. To counteract effects of NEMO overexpression, which can strongly inhibit NF-κB activation [fig. S2 and ( 42 )], we used CRISPR-Cas9 for targeted insertion of coding sequences for FPs into the U2OS cell line. The resulting cells coexpress N-terminal fusions of EGFP-NEMO and mCh-RelA from their endogenous loci, which can be used to monitor dynamic signaling events by live-cell imaging ( 33 ).

Although activated TNFR1 and TNFR2 both form TNF-induced homotrimeric complexes, the TNFR2 subtype binds with lower affinity to soluble TNF and shows enhanced activation by membrane-bound TNF ( 31 , 41 ). In contrast, ligand-activated IL-1 receptor (IL-1R1) forms a heterodimer with the IL-1R3 accessory protein, and dimerization can be inhibited through competitive sequestration by the IL-1R2 decoy ( Fig. 1A ) ( 32 ). Because surface expression of TNFR2 and IL-1R2 are comparably low in U2OS, receptor composition of TNF- and IL-1–induced oligomers will consist predominantly of TNFR1 trimers and IL-1R1–1R3, respectively. Together with known receptor-ligand stoichiometry ( Fig. 1A ), our results predict that single cells can simultaneously form a maximum of hundreds of IKK-recruiting complexes for saturating cytokine concentrations (approximately 400 and 700 for TNF and IL-1, respectively). Surface receptor expression is substantially lower than numbers for downstream signaling molecules such as NEMO that are expressed in orders of a million per cell ( 42 ).

( A ) Schematic of cognate receptors for TNF and IL-1 cytokines. Monomeric receptors and receptors that engage with decoy receptors (IL-1R2) are inactive and do not transmit signals into the cytoplasm (left). Activated receptor complexes (right), consisting of a TNFR homotrimer bound to TNF or an IL-1R1–IL-1R3 heterodimer bound to IL-1, are capable of seeding CI-like complexes in the cytoplasm. ( B ) Quantification of surface receptor expression on single U2OS cells. The average of three to seven biological replicates is shown for each condition. Error bars represent SEM.

IKK activity is a convergence point for proinflammatory signals that regulate NF-κB downstream of many cytokine receptors ( 26 , 30 ). Ligands that bind to multiple receptors with differing kinetics ( 31 ) and decoy receptors that sequester or antagonize signaling complexes ( 32 ) layer additional complexity to signal initiation at the plasma membrane. To establish expectations for numbers and types of IKK-activating complexes, we measured surface receptor expression in U2OS cells that were previously shown to form dynamic IKK puncta in response to TNF and IL-1 ( 19 , 33 ). Using flow cytometry with reference beads for absolute quantification, we estimated the number of surface receptors per cell for TNFR1, TNFR2, IL-1R1, IL-1R2, and IL-1R3 ( Fig. 1 and fig. S1). On average, each U2OS cell presented approximately 1300 TNFR1, 700 IL-1R1, and an abundance of IL-1R3 surface receptors. Only a small number of TNFR2 and IL-1R2 were detected on the cell surface. For reference, we measured surface receptors on HeLa and KYM1 cells (fig. S1) and found results consistent with previous reports for TNFRs ( 34 – 36 ) and agreement with surface receptor expression in other cell lines ( 37 – 40 ).

In summary, our live-cell experiments revealed CI-like complexes are independent and switch like, where each complex recruits a quantized amount of IKK over its life span. These were unexpected results that led to conceptual simplification and identification of the underlying signaling architecture. However, independence between detector nodes is not a required characteristic of a VG-SPN, and most receptor signaling systems can therefore be viewed through a similar lens. We expect the VG-SPN is common to receptor signaling systems, each with distinct pooling functions, feedbacks, and feedforwards that will reveal their individual trade-offs and information transmission benefits to the cell.

The classic model of chemoreception ( 53 ) established limitations for the signal-to-noise relationships in biological sensing systems that use time-averaged receptor occupancy to measure chemical concentration. By contrast, multivalent cytokines exhibit markedly increased apparent affinity over bimolecular interactions and, therefore, limit the likelihood of ligand dissociation from mature signaling complexes. Although these signaling systems are distinct, they are both reducible to the same motif with alternative definitions for amplification gain, and possibly with different pooling functions. A remaining question is why biological sensory systems would converge on SPN-like architectures. Although several conceptual benefits and limitations of VG-SPNs are already discussed here, other general properties of SPNs may also contribute to the answer. Stochastic resonance is an example of a counterintuitive property of SPNs that uses stochastic fluctuations to enhance information transmission ( 50 , 54 ). Although the number and timing of CI-like switches is stochastic and subject to noise during cytokine responses, they can theoretically enhance information transmission over an array of noiseless binary switches. An array of perfect switches that simultaneously activate at the same signal threshold is a binary system, either all “on” or all “off’,” whereas stochasticity in the fraction and timing of CI formation increases the number of system states through probabilistic interactions with the milieu.

Abstraction of our experimental observations revealed a VG-SPN signaling architecture ( Figs. 5 and 6 ). The resulting model enabled us to investigate in silico the emergent properties of IKK–NF-κB at the level of cytokine detection and signal amplification. With detection and parallel signal amplification at independent CI nodes, our model revealed that noise is effectively mitigated with 100s to 1000s of signaling complexes, and greater numbers have diminishing returns per complex for information transmission. Further benefits to signal transmission through a CI amplifier allow cells to fine-tune numbers of activated cytoplasmic signaling molecules, which are substantially more abundant than cytokine receptors at the cell surface. For NF-κB signaling, high-gain amplification enables a large repertoire of receptors to engage the same cytoplasmic pool of IKK with limited occupancy of space at the cell surface. Cells can therefore favor parsimony in receptor numbers or increase receptor numbers with reduced amplification gain to interface with the same size pool of signaling molecules while preserving information transmission. These trade-offs provide orthogonal controls to robustly fine-tune or diversify response sensitivity to stimuli, as shown here for different cytokines.

Induction of A20 transcription is considered a defining negative feedback in the NF-κB response ( 10 , 23 , 26 ). However, tracking experiments revealed that trajectories of early- and late-forming EGFP-NEMO spots are similar and insensitive to CHX ( Fig. 5 ). Although this does not preclude noncatalytic roles for A20 in regulation of IKK ( 29 ), it demonstrates that each CI-like complex is independent and not influenced by transcriptional feedbacks from complexes that form earlier during a response in the same cell. Our data therefore support another proposed role, where transcriptional feedback via A20 is not primarily directed at the initial immune response but instead establishes a new baseline for tolerance to subsequent stimuli ( 44 ). Together, signal amplification and negative feedback at each CI-like complex are determined predominantly by the resting cell state.

Notably, imaging requirements for EGFP-NEMO limit experimental throughput, in particular for spot tracking experiments that require both high magnification and high frequency time lapse. Although dynamical properties of NF-κB signaling are important mediators of information transmission, analyses that calculate information metrics require a large number of single-cell data points ( 1 – 3 ). Consequently, our analysis of EGFP-NEMO required simplifications through coarse-grained scalar descriptors that summarize dynamic properties of NEMO-recruiting complexes and nuclear RelA localization in the same single cells. Our analysis revealed two descriptors of EGFP-NEMO that are strong determinants for descriptors of RelA that were previously shown to carry the most information about cytokine concentrations in the milieu ( 3 ). As technologies emerge that enable data collection for calculation of information metrics between both reporters in the same cell, determination is likely to improve. However, it is also rational that the aggregate sum of NEMO in CI-like complexes during a primary cytokine response is a strong same-cell determinant of the accumulated NF-κB response. This deterministic relationship may therefore remain among the strongest in coming years, providing a readout for signaling flux and disease-associated perturbations at the level of CI.

Endogenously expressed EGFP-NEMO is a multifaceted reporter that reveals several aspects of signal transmission. At the level of detection, there is agreement between numbers of EGFP-NEMO puncta induced by saturating cytokine concentrations and average surface receptor numbers in U2OS cells ( Figs. 1 and 2 ). Similarly, the timing of formation of EGFP-NEMO puncta establishes when ligand-receptor-adapter assemblies become capable of signaling in the cytoplasm. Continuous stimulation experiments show that most EGFP-NEMO spots form within 5 to 10 min, indicating that the cytoplasmic components of CI-like complexes are not limiting. However, in cells exposed to a short pulse, new spots form up to 30 min following cytokine removal, demonstrating variability in timing for receptors to assemble into a signaling-competent stoichiometry. In the cytoplasm, spot-intensity time courses inform about biochemical interactions and feedbacks linked to signal amplification. Here, families of different ubiquitin ligases, kinases, and DUBs engaged at CI-like structures establish rates of EGFP-NEMO recruitment and dissolution. Distinct properties between puncta initiated by different receptor superfamilies ( Figs. 2 to 5 ) suggest that various adapters and ubiquitin requirements associated with different types of CI-like complexes will determine cytokine-specific signal amplification ( 10 – 13 , 19 ). We therefore expect that endogenous EGFP-NEMO will be valuable in unraveling these upstream mechanisms of inflammatory signaling, similar to reporters for NF-κB that have contributed richly for well over a decade.

MATERIALS AND METHODS

Cell culture Parental KYM1 (female), HeLa (female), and U2OS (female) cell lines were obtained from the American Type Culture Collection, and the HeLa cell line stably expressing scFv-GFP and tdPCP-tdTomato was a gift from X. Zhuang from Harvard University (43). All cell lines were cultured in RPMI (KYM1), Dulbecco’s modified Eagle’s medium (DMEM; HeLa), or McCoy’s 5A (U2OS) media at 37°C and 5% CO 2 . All media were supplemented with 10% Corning regular fetal bovine serum (FBS), penicillin (100 U/ml), streptomycin (100 U/ml), and 0.2 mM l-glutamine (Invitrogen). Cells were periodically monitored for mycoplasma contamination.

Quantitative flow cytometry Fluorescence-activated cell sorting (FACS) analysis of samples and phycoerythrin (PE) beads were performed in a BD LSRFortessa machine (University of Pittsburgh–Department of Immunology Flow cytometry facility). We used PE-conjugated beads (BD Biosciences, catalog no. 340495) and PE-conjugated antibodies against the following human proteins that were obtained from R&D Systems: TNFR1 (FAB225P), TNFR2 (FAB226P), IL-1R1 (FAB269P), IL-1R2 (FAB663P), IL-1R3 (FAB676P), polyclonal goat immunoglobulin G (IgG) PE-conjugated (IC108P), and mouse IgG1 PE-conjugated (IC002P). Cells were maintained and used between 5 to 15 passages. For deattaching cells from plates, we use 2 mM EDTA in cold phosphate-buffered saline (PBS). Forty-eight hours before staining, 1 × 105 to 2 × 105 cells were plated in six-well plates. On the day of the staining, cells were deattached from the plate, transferred to polypropylene microtiter tubes (2681377, Fisherbrand), and maintained on ice during the entire process. Fc receptors were blocked using human BD Fc block antibody (564219) for 10 to 15 min in the dark. Next, cells were washed with 2% FBS in PBS (FACS buffer) and incubated with primary antibodies for 1 hour at 4°C in the dark. Next, samples were washed twice and resuspended in DAPI (4′,6-diamidino-2-phenylindole) (1 ug/ml; catalog no. D1306, Invitrogen) in FACS buffer between 0.5 and 1.5 hours before data acquisition. Data analysis was performed using FlowJo (BD, version 10.6.1_CL). Samples stained with PE isotype controls were used to set the background fluorescence by subtracting their mean fluorescence intensity (MFI) to the MFI of samples stained with specific antibodies. Next, we used the MFI from PE-conjugated beads to estimate the number of PE molecules on the surface of each cell. All these antibodies have been previously used to estimate surface receptors, with accepted 1:1 ratios of PE to antibody molecules (55–57).

Fixed-cell immunofluorescence U2OS wild-type and U2OS cells stably overexpressing NEMO cells were seeded into plastic flat-bottom 96-well imaging plates 24 hours before cytokine treatment at a density of 7000 cells per well. On the day of the experiment, wells receiving IL-1 or TNF (10 or 100 ng/ml, respectively) were stimulated 45 min before fixation. Prewarmed 15× cytokine mixture was spiked into wells and mixed. After treatment and before fixation, the cells remained in environmentally controlled conditions (37°C and 5% CO 2 ). Cells were fixed with 4% paraformaldehyde and permeabilized using 100% methanol for 10 min each with one PBS wash in between and three PBS-T (PBS 0.1% Tween 20) washes at the end. Next, cells were incubated in primary antibody solution (3% BSA in PBS-T) with 1 μg/ml of both α-RELA (sc-8008; Santa Cruz Biotechnology) and NEMO (sc-8330; Santa Cruz Biotechnology) at 4°C overnight. The following morning, cells were washed three times with PBS-T and incubated for 1 hour at room temperature with secondary antibody solution in 3% BSA in PBS-T [4 μg/ml of both goat anti-mouse IgG Alexa Fluor 647 (A21235; Thermo Fisher Scientific) and goat anti-rabbit IgG Alexa Fluor 594 (A11012; Thermo Fisher Scientific)]. After incubation, cells were washed twice in PBS-T and incubated in Hoechst in PBS-T (200 ng/ml) for 20 min. Last, wells were washed in PBS-T and left in PBS to keep cells hydrated during imaging. Cells were imaged using a DeltaVision Elite imaging system at ×20 magnification with a LUCPLFLN objective (0.45 numerical aperture; Olympus).

Fixed-cell image analysis Using CellProfiler [www.cellprofiler.org; (58)], cell nuclei were segmented using the Hoechst channel labeling DNA. Next, secondary segmentation was performed using the EGFP-NEMO channel to define the cellular boundaries. The output of cell segmentation was compiled and analyzed using custom scripts in MATLAB (Mathworks, R2019b).

Establishing EGFP-NEMO/mCherry-RELA CRISPR double knock-in cells Single knock-in U2OS cell lines expressing EGFP-RELA/NEMO were generated previously using CRISPR-Cas9 technology as described in reference (33). To generate double knock-in cells, we assembled the RelA repair template consisting of DNA sequences for a left homology arm [−544 base pairs (bp), chromosome 11_65663376–chromosome 11_65662383], followed by an mCherry protein coding sequence with a start codon but no stop codon and a sequence encoding 3× GGSG linker in-frame with the right homology arm (+557 bp, chromosome 11_65662829–chromosome 11_65662276) from plasmids synthesized by GeneArt. Synonymous mutations were introduced to prevent interaction of the repair template and Cas9. Next, single EGFP-NEMO knock-in U2OS cells were seeded in six-well plates (2 × 105 cells) and transfected next day with pSpCas9n-(BB)-2A-Puro-RelA_gRNAs and Bgl2-linearized mCherry-RelA repair template donor plasmid. A ratio of 3.5:1 FuGENE HD (Promega) to total DNA was used for transfection of 4 μg of DNA. Plasmid generation and CRISPR modifications were all based on Ran et al. (59).

Live-cell imaging Live cells were imaged in an environmentally controlled chamber (37°C, 5% CO 2 ) on a DeltaVision Elite microscope equipped with a pco.edge sCMOS camera and an Insight solid-state illumination module (GE Healthcare). For detection of NEMO spots, U2OS cells expressing FP fusions of RELA and NEMO were seeded at a density of 15,000 cells per well 24 hours before live-cell imaging experiments on no. 1.5 glass-bottom 96-well imaging plates (Matriplate). Medium was changed to phenol red–free FluoBrite DMEM (Gibco, A18967-01) between 30 min and 2 hours before imaging. For detection of NEMO spots, live cells were stimulated with the indicated concentrations of TNF or IL-1. After cells adapted in the environmentally controlled chamber, we collected eight z-stack images of 0.5-μm separation in the Alexa Fluor 488 channel with an exposure of 0.04 s and a transmission of 32%. Using these stacks of images, we reconstructed NEMO spots in three dimensions (3D) using ImageJ (for display) and dNEMO (more details in the “Punctate structures detection and quantification” section). To test our different hypothesis, we performed two types of time-lapse imaging experiments: long-term and short-term high-frequency imaging. Images were collected over at least three fields per condition with a temporal resolution of 2 min per frame for long-term and 10 s per frame for short-term high-frequency imaging experiments. Wide-field epifluorescence and differential interference contrast images were collected using a 60× LUCPLFLN objective. For all treatments, cytokine mixtures were prepared and prewarmed so that addition of 120 μl added to wells results in the indicated final concentration. For two-pulse experiments, cells were exposed to a 1-min dilute “reference pulse” of either TNF or IL-1, selected from calibration experiments to produce a subtle response in most cells (2 and 5 ng/ml for TNF and IL-1, respectively), followed by a “saturating pulse” of 500 ng/ml for same-cell comparisons. Cells were washed once with prewarmed fresh media during media swaps. Because EGFP-NEMO foci form and disperse rapidly after treatment, a 30-min recovery period was selected between “reference” and “saturating” pulses. The short recovery period also minimized the effects of photobleaching that were apparent within 60 to 90 min of imaging. For small-molecule-inhibitor experiments to study the nature of the NEMO spots, we pretreated cells with 178 μM CHX, 50 μM Monodansylcadaverine (MDC), or 20 μM Dynasore for 20 min before adding cytokine on top of the medium with inhibitor. For dose-response experiments to calculate channel capacities using live-cell imaging data, U2OS cells were prestained with Hoechst 33342 (300 ng/ml) for 1 hour. Following the incubation with Hoechst, all growth medium contained trace amounts (60 ng/ml) of Hoechst to maintain the nuclear stain and assist with segmentation. Cells were then treated with six different doses of TNF or IL-1 (1000, 100, 10, 1, 0.1, and 0.01 ng/ml) plus a control of 0 ng/ml. After treatment, cells were imaged at the temporal resolution of 5 min per frame with a 20× LUCPLFLN objective.

Live-cell imaging in microfluidic devices Microfluidic device fabrication and operation were done as described previously (45). Briefly, dynamic stimulation devices of two inlets and one outlet were made with polydimethylsiloxane (PDMS; Sylgard), autoclaved, washed with ethanol, and incubated with a solution of 0.002 (v/v) fibronectin in PBS for 24 hours at 37°C. After incubation, excess of fibronectin was flushed with tissue culture media multiple times. Between 3 × 106 and 8 × 106 cells/ml were seeded by inserting 200-μl pipette tip in the outlet port. The outlet port was plugged with PDMS plugs after reaching a cell density for approximately 60% confluence. The device was incubated for at least 24 hours, and PDMS plugs at the outlet were then replaced with pipette tips filled with medium. On the day of the experiment, Tygon tubes (Fisher Scientific, 1471139) were attached to the device with the other end connected to basins with media or cytokine treatment. The basins were then placed on corresponding platforms on the gravity pump, and the attached Tygon tubes were clamped. The device was then fitted to a custom adapter and placed under the microscope for imaging. The clamps on the tubes were removed at the beginning of the experiment. Cells were stimulated with the desired pattern of cytokine stimulation (pulse, ramp, or continuous) by changing the heights of basins by controlling corresponding stepper motors through custom codes uploaded on the Arduino Mega 2560 Microcontroller. The cytokine treatment was prepared with Alexa Fluor 647–conjugated BSA (0.0025, v/v; Invitrogen) and imaged with CY5 filter to confirm the corresponding stimulus pattern. Images were collected under the same conditions as previously described for live-cell imaging.

Punctate structures detection and quantification EGFP-NEMO and SunTag polysome spots were detected and quantified using our application dNEMO (60). Briefly, dNEMO is a computational tool optimized for measurement of fluorescent puncta in fixed-cell and live-cell time-lapse images. The user-defined threshold for spot detection in dNEMO was set between 1.5 and 2.0 for all images. Reported pixel values for puncta were individually background corrected by averaging pixels from an annular ring surrounding each spot. The width and offset for the annular ring were both set to 1 pixel. We furthermore stipulated that puncta must appear in at least two contiguous slices of the 3D images (of eight slices) to be considered valid. The same user parameters were applied to all images, and single cells were manually segmented using dNEMO’s keyframing function. For each single cell, spot features were measured for each new spot that formed following stimulation, yielding sets of single-cell spot features over time. NEMO puncta flat-field and background-corrected images were prepared for display using ImageJ.

Tracking of individual NEMO spots The location and intensity data of EGFP-NEMO spots obtained with our dNEMO software (60) were used for single-spot tracking using the uTrack package in MATLAB (61). The hyperparameters in the uTrack package were adjusted to enhance tracking performance for TNF- and IL-1–induced spots. Specifically, the time gap window was lowered to three frames, and the gap penalty was lowered to 1 from their default values of 5 and 1.5, respectively. The minimum length of track segments used for gap closing was increased to three frames from the default of 1 frame. The merging and splitting events were considered while tracking. Properties of individual spots (intensity and size) were then associated with tracking data to generate single-spot trajectories for each spot property. Trajectories lasting for less than 3 min were excluded from further analyses. With these settings, approximately 60 and 75% of spots, respectively, for TNF and IL-1 responses were tracked and included in subsequent analysis (fig. S7A). Three main features were obtained from each single-spot trajectory: 1) Maximum intensity or Max i (peak intensity of a single-spot trajectory) 2) Integrated intensity or AUC i (sum of spot intensities at all time points in a trajectory) 3) Track length (trajectory length, i.e., time for which the spot is tracked). Features of single-spot trajectories were then compared for dynamic stimuli experiments involving different cytokines across different cells (fig. S7) and for spots formed at different time windows in reverse time-course experiments (Fig. 5D and fig. S9).

Estimating number of NEMO molecules per NEMO spot via calibration with SunTag-labeled polysomes To quantify the number of NEMO molecules within each NEMO complex (spot), the CRISPR labeling of NEMO with EGFP allows us to infer the number of NEMO molecules by counting GFP molecules and converting with 1:1 ratio. In counting GFP molecules in NEMO puncta, we used the live cell translation reporter developed by Wang et al., (43) to calibrate the relation between GFP counts and measured GFP intensity and imaged it in HeLa cells with the same imaging condition used for NEMO in U2OS EGFP-NEMO cells. During translation, one mRNA binds to multiple ribosomes simultaneously and forms a large polysome complex. We assume that the positioning of each ribosome on mRNA independently satisfies a uniform distribution. The signal intensity of each fluorescence foci in the cytoplasm, representing the translating polysomes, varied depending on the total number of ribosomes and the location of each ribosome on mRNA. When ribosomes reach the region after the coding sequence for the SunTag peptide, the intensity would be the maximum for fully assembled complex; when ribosomes are halfway toward the end, the signal would be approximately half maximal. To account for this variability, we used Monte Carlo simulations to randomly sample the positions of ribosomes and generate the distribution of numbers of SunTag peptides being produced when there are n ribosomes present on single mRNA, denoted by p 2 (m,n), where m is the number of SunTag. n, as the number of total ribosomes on each mRNA, satisfies Poisson distribution p 1 (n) = e^(−λ) λ^n/n!. The average number of ribosomes for each translation foci is measured to be 12 experimentally (43); therefore, we set λ ≈ 12. The possibility of m SunTag peptides being translated on each translation foci are calculated as follows P ( m ) = ∑ n = 1 + ∞ p 2 ( m , n ) p 1 ( n ) When n is large, p 1 (n) rapidly converges to zero. Therefore, we set n = 30 as the cutoff of the summation. We quantified the fluorescence intensity of individual translating polysome in HeLa cells using dNEMO (60) with the same threshold used for NEMO quantification in U2OS cells. As the GFP intensity is proportional to the number of GFP molecules, we determined the scaling factor between the measured fluorescence intensity and theoretical distribution of GFP molecule numbers. Here, we used the Nelder-Mead algorithm to minimize the square of difference between the “measured” GFP distribution and theoretical GFP distribution to obtain the scaling factor. To account for effects of photobleaching, we calculated the scaling factors independently for each frame of the time-lapse image. Using the scaling factor obtained from the translation reporter, the intensity of NEMO spots was converted to numbers of GFP molecules per spot.

Extracting descriptors from NF-κB dynamics Descriptors of the fold change of the FP-RelA mean intensity trajectories were extracted using custom MATLAB scripts and ImageJ movie explorer. The fold change transform is carried out by dividing the FP-RelA trajectories by the initial nuclear fluorescence at the zero time point. Those descriptors include the following: 1) Area under the fold change curve (AUC Fold ): A summary statistic that approximates the cell’s response over time, previously found to carry the most information about a cell’s response to cytokine stimulation (3). 2) Max fold change (F max /F initial ): Maximum value of the fold change trajectory, previously found to determine the transcriptional response of cells stimulated with TNF (24). 3) Maximal rate of nuclear entry (Rate in ): Slope of line fitted to three data points between time zero and time of maximum fold change. Point along fold change trajectory at which slope was determined was point where slope had maximal absolute value. 4) Maximal rate of nuclear exit (Rate out ): Slope of line fitted to three data points between time of maximum fold change and the end of the fold change trajectory. Point along fold change trajectory at which slope was determined was point where slope had maximal absolute value. 5) Time of max [t max (fold) ]: Time at which fold change is maximum. 6) Time of half up [t 50up (fold) ]: Time point at which fold change trajectory rises to half the max fold change. 7) Time of half down [t 50down (fold) ]: Time point at which fold change trajectory falls to half the max fold change.

Extracting descriptors from EGFP-NEMO dynamics Descriptors of trajectories for the NEMO spot number and NEMO spot intensity were extracted from single cell and single complex data using custom MATLAB scripts. The NEMO spot number trajectory is the number of NEMO spots detected per cell over time, while the NEMO spot intensity trajectory is the integrated intensity of the detected spots per cell over time. Those descriptors (per trajectory) include the following: 1) AUC: Area under the trajectory integrated. For a given cell’s NEMO spot number and spot intensity trajectories, this represents either the total number of spots detected or the total intensity of the spots detected, respectively. 2) Maximum value (Max): Maximum number of spots/spot intensity for a given cell trajectory. 3) Rate of entry (Rate in ): Slope of line fitted to three data points between time zero and time of max. Point along trajectory at which slope was determined was point where slope had maximal absolute value. 4) Rate of exit (Rate out ): Slope of line fitted to three data points between time of max and the end of the trajectory. Point along trajectory at which slope was determined was point where slope had maximal absolute value. 5) Time of max (t max ): Time at which the NEMO spot trajectory is maximal. 6) Time of half up (t 50up ): Time at which the NEMO spot trajectory rises to half the maximum value. 7) Time of half down (t 50down ): Time at which the NEMO spot trajectory falls to half the maximum value. 8) Full width at half max: Width of the peak of the NEMO spot trajectory between the time of half up and the time of half down. 9) Maximum over area of the cell (Max/Area): Max value of the NEMO spot trajectory divided by the area of the cell. The cell’s area is given by the polygon output from the dNEMO results. 10) Fold change peak difference [t max (fold) − t max ]: difference between the time of maximum fold change and the time of max of the NEMO spot trajectory.

Correlation and multiple linear regression of NEMO and NF-κB descriptors Descriptors collected from both the FP-RelA and EGFP-NEMO trajectories were correlated using Spearman’s rank and coefficient of determination in linear and log-log scales (see fig. S5). When correlation between a pair of descriptors could not be computed (e.g., an inability to compute some trajectory’s Rate in or t max in cells that did not respond significantly to stimulation), the single cell was omitted from the affected correlation dataset but not from the overall set of cell trajectories being analyzed. Multiple linear regression was also applied to the linear and log-transformed descriptors of the cells responding to TNF, IL-1, or the full combined experimental dataset (see fig. S6). Linear correlation and multiple linear correlation of the descriptor sets were carried out using custom MATLAB scripts.

Computational model (HyDeS) to examine relative effects of basal and transcriptional feedback on NEMO spot dynamics We modeled the intensity of individual NEMO spots with a hybrid deterministic-stochastic (HyDeS) model. The purpose of this model is to examine the relative effects of basal and transcriptional feedback via DUB molecules on NEMO spot dynamics, and not necessarily to capture precisely the absolute intensity dynamics for EGFP-NEMO at each single spot. Therefore, to directly model the relative effects of basal and transcriptional feedback, we used three main variables. First, we defined the intensity of a spot (I) to represent a continuous approximation for ubiquitin chain size and NEMO activity resulting from molecules recruited at that spot. Second, we defined a variable X_basal for each individual spot as a proxy for feedback due to activity from basal DUBs expressed in resting cells. We modeled the X_basal variable such that the basal feedback increases as the corresponding NEMO spot grows through proportionate recruitment of DUB molecules on the same complex. Third, we defined a variable X_trnsxl to capture the transcriptional feedback due to the aggregate activity of all NEMO spots. Here, transcription-induced DUBs act in the same way as basal DUB molecules. The model corresponding to N number of NEMO spots consists of 2N + 1 ODEs (Ordinary Differential Equations): N for the individual spot intensities (I i ), N for the basal feedback variables corresponding to each spot (X_basal i ), and one for the transcriptional feedback variable (X_trnsxl). For each simulation, we modeled 200 single spots that form at regular intervals to approximate the response of a single cell to cytokine stimulation d [ I i ] dt = P form * P bound K growth K limit + [ I i ] − Kd basal * X _ basal i * [ I i ] − Kd trnsxl * X _ trnsxl * [ I i ] d [ X _ basal i ] dt = Kx basal * [ I i ] − Kxd basal * [ X _ basal i ] d [ X _ trnsxl ] dt = Kx trnsxl * ∑ i [ I i ] − Kxd trnsxl * [ X _ trnsxl ] The variables and parameter values used in the model are as described in tables S3 and S4. We examined the relative effects of the two kinds of feedbacks by varying the parameters Kx basal and Kx trnsxl , which respectively control the amount of local and global feedback, by two orders of magnitudes (fig. S8).

VG-SPN model SPNs are a model sensory system in which noisy detectors are used to make independent and compressed measurements of a signal (50). These measurements are then pooled to average out the uncorrelated noise and reconstruct the original signal. Here, we characterize an SPN system with the decoration of variable gain (VG-SPN). In the VG-SPN system, each detector amplifies its binary measurement with a “Gain” before pooling. Thus, the VG-SPN can have different configurations depending on the number of detectors and the Gain per detector. In the context of cytokine signaling mediated by CI complexes, each complex acts as a binary detector of the extracellular cytokine presence, and the number of NEMO molecules recruited by a complex is considered its corresponding Gain. Therefore, we define two fundamental variables in our model of VG-SPN system. First, CI max is the number of detectors in a VG-SPN configuration that is given by the maximum number of CI complexes that can form in a cell at saturation. Second, Gain is the number of NEMO molecules recruited at a given CI complex over the time course of its activity. Depending on CI max and Gain configurations, VG-SPN systems will have different R max in terms of total number of NEMO molecules recruited. To systematically examine the noise mitigation properties of different VG-SPN configurations, we used a mathematically controlled approach by assuming each configuration is capable of producing the same steady-state response (Fig. 7A). We therefore fixed the value of R max and defined the Gain per complex to be inversely proportional to the maximum number of complexes as Gain = R max /CI max . We then varied CI max across orders of magnitude and calculated the noise associated with each configuration (Fig. 7B). The detector shot noise is the fraction of NEMO molecules recruited when a CI complex forms erroneously, i.e., in absence of the signal. We assigned the probability of a complex to erroneously form as the “shot noise level” (SNL). Therefore, the noise propagated when the SNL is less than the uniform probability (1/CI max ) will correspond to one complex forming erroneously, which, in turn, will be the Gain associated with that configuration. For example, erroneous activation of a one-receptor system via shot noise will activate the system to its fullest capabilities. Because of the inverse correlation of CI max with Gain, the shot noise will decrease with CI max (Fig. 7B, top). However, after reaching a certain value of CI max , the SNL will become higher than the uniform probability, allowing more than one complex to erroneously form at the same time. Configurations beyond such value of CI max will hit the noise floor, and the shot noise will not decrease further Fractional shot noise = { Gain R max = 1 CI max , when SNL < 1 CI max SNL * CI max * Gain R max = SNL , when SNL > 1 CI max } We observed considerable variability in intensities of NEMO complexes formed in response to a particular cytokine (Figs. 2D and 5, C and D). This variability suggests that although the mean value of Gain per complex is cytokine specific, the gain associated with a single complex can vary among different complexes within a cell, and these differences will introduce “Gain noise” in the pooled response. Specifically, Gain noise captures the noise in the pooled response corresponding to same number of active detector complexes due the variability in Gain per complex (i.e., the number of NEMO molecules recruited to a single complex over its life span). We examined the effect of different VG-SPN configurations on the Gain noise using the mathematically controlled approached as described earlier. We modeled the Gain per complex with Gaussian distribution and using the coefficient of variation given by the “Gain noise level” (GNL) as σ G = GNL*Gain. Therefore, for a given VG-SPN configuration, each cell will have a response corresponding to the sum of CI max random numbers obtained from the distribution N(Gain, σ G ). We simulated 1000 cells for every configuration with different GNL values and calculated the Gain noise (coefficient of variation) that contributes to the resulting responses (Fig. 7B, bottom). Next, we developed a computational VG-SPN model to generate simulated data for channel capacity calculations using the framework developed by (2) and codes from our previous study (3). In this model, the signal S in the range (0,1] represents the extracellular cytokine dose where S = 1 corresponds to saturating dose. The dose range is equally divided in discrete levels controlled by the hyperparameter “numDoses” S i ∈ ( 0 , 1 ] , i ∈ [ 1 , numDoses ] For each input signal S i , we simulate responses from N cells, where N is another hyperparameter. The response of the jth cell to ith level of input signal (R ij ) is a scalar quantity representing the total number of recruited NEMO molecules. The model thus produces simulated data as response vectors of length N for each of the numDoses number of input conditions. The data workflow is as follows S = [ S 1 S 2 ⋮ S i ⋮ S numDoses ] → R = [ R 1 R 2 ⋮ R i ⋮ R numDoses ] where S i → R i = [ R i 1 R i 2 … R ij … R iN ] As described earlier, CI max is the maximum number of NEMO complexes a cell can form at saturating conditions, i.e., when S i = 1. CI max will be related to the number of surface receptors (SR) per cell through the cytokine-dependent valency as CI max = SR/valency. Therefore, to model cell-to-cell variability in receptor expression (Fig. 1B), we assumed that CI max follows a Gaussian distribution over the N simulated cells with coefficient of variation given by the receptor noise level (RNL) CI max j ∼ N ( CI max , RNL * CI max ) Here, the mean CI max is set by the VG-SPN configuration, and RNL is another hyperparameter controlling the coefficient of variation in receptor expression. For instance, for N = 1000, CI max = 100, and RNL = 0.3, we sample 1000 normally distributed integer random numbers from N(100,30) to simulate responses for 1000 cells, where each cell can maximally form around 100 complexes with coefficient of variation 0.3. We used binomial distribution to model switch-like properties of complex formation in the detection and switching steps (Fig. 6A), wherein the probability of complex formation (P form ) increases with the input signal S i . We used modified Hill function to relate S i with P form (Fig. 2F) CI max j * ∼ Binomial ( CI max j , P form ) , where P form = ( 1 + 0.5 ) * S i 3 0.5 + S i 3 Here CI max j * is the number of complexes formed in response to the dose level S i in a cell having capacity to form CI max j complexes at saturation. Therefore, for a given jth cell, we generate CI max j random numbers between [0,1]. We then check how many of those fall below the P form value corresponding to the given signal S i . This number, CI max j * , will be used as number of active complexes in the jth cell. Thus, the output of detection step corresponding to signal S i is the following vector S i → [ CI max 1 * CI max 2 * … CI max j * … CI max N # * ] Note that N# ≤ N because for lower doses there might be some cells that do not cross the P form threshold, i.e., they do not form any spots. In the detection step, we simulated the number of complexes formed by a given cell ( CI max j * ) in response to a given dose (S i ). Next, in the amplification step (Fig. 6A), we model the Gain (G k ) in terms of number of NEMO molecules recruited per complex. The response of a cell, i.e., total number of recruited NEMO molecules, is then the sum of the Gains associated with all the complexes formed in that cell (signal pooling step, Fig. 6A) Response R ( S i , CI max j ) = ∑ k = 1 CI max j * G k We assume a linear relation between Gain per complex (in terms of NEMO molecules) and the experimentally observed integrated intensity (AUC i ; fig. S5) G k ( NEMO Molecules ) = AUCi _ to _ NEMO _ factor * G k ( Intensity ) Here AUCi_to_NEMO_factor is another hyperparameter. As previously stated for calculation of Gain noise, we use Gaussian distribution to model the variability in G k within a cell G k ∼ N ( Gain , GNL * Gain ) Here, GNL is the Gain noise level, i.e., the coefficient of variation of Gains corresponding to different complexes within a cell. We use the quadratic relation between intracellular mean and variance of AUC i (fig. S7E) to get GNL for different mean values of Gain corresponding to different VG-SPN configurations. To summarize, response of the jth cell with a VG-SPN configuration given by (CI max , Gain) to the input signal level S i is calculated as follows Response R ij ( S i , CI max , Gain ) = ∑ k = 1 CI max j * AUCi _ to _ NEMO _ factor * N ( Gain , GNL * Gain ) k where CI max j * = Count ( CI max j < P form ( S i ) ) and CI max j ∼ N ( CI max , RNL * CI max ) This process is repeated for N number of cells and numDoses number of input signal levels CI max j ∈ [ CI max 1 CI max 2 … CI max N ] S i ∈ [ S 1 S 2 … S numDoses ] There are two independent variables in this model, CI max and Gain, which completely define the VG-SPN configuration. In other words, the model allows independent “tuning” of maximum number of NEMO complexes (i.e., number of surface receptors) and the Gain (i.e., number of NEMO recruited by each complex) to achieve a certain channel capacity. The model has four hyperparameters: numDoses (controlling number of discrete dose levels), N (number of cells to simulate per dose level), RNL (coefficient of variation in CI max over cells), and AUCi_to_NEMO_factor (relating integrated intensity to number of NEMO molecules per complex). In addition to these four, the channel capacity calculation framework developed by (2) requires another hyperparameter “K” that controls the K-nn classification to estimate probability densities of response given input dose level. To interface the VG-SPN model with our previous code for channel capacity calculations (3), we adjusted these five hyperparameters to avoid technical biases that can artificially limit calculated channel capacity values for a VG-SPN configuration. After running preliminary computational experiments, we arrived at a set of values for the five hyperparameters. Then, for each hyperparameter, we varied its value while keeping other values in the set constant to confirm that the chosen value does not limit the maximum channel capacity value of a VG-SPN configuration. By iterating through this process several times and repeating it for multiple values of R max , we converged on the following settings (fig. S12A). First, the number of simulated input conditions can artificially limit the channel capacity calculated for a signaling system but becomes increasingly costly to compute when input conditions become excessively high. We simulated our model with 10, 30, 60, and 100 input conditions (numDoses parameter in Materials and Methods) and found that the maximum channel capacity does not increase beyond 60. Next, we observed that the maximum channel capacity decreases as the cell-to-cell variability in number of CI per cell (RNL parameter) increases. We therefore chose the intermediate value of stdR = 0.3 that is consistent with our FACS receptor analysis (Fig. 1 and fig. S1). Next, the value of K used for K-nn classification introduces minimal bias beyond a value of 5 and contributes to computational cost at higher values. We therefore chose K = 5, consistent with previous calculations (2). We then found that simulating more than 1000 cells per input condition does not increase the channel capacity but greatly affects simulation time, so we chose N = 1000. Last, to relate amplifier gain (G L ) to the quadratic noise relationship observed in data (fig. S7E), we estimated the factor to convert integrated spot intensity to number of NEMO molecules at a CI-like complex (“AUCi_to_NEMO_factor”). Values for the conversion factor did not have significant effects on max channel capacity beyond a value of 200. Once we fixed the hyperparameters of the model, we removed the constraints of the mathematically controlled framework and allowed the two model variables, CI max and Gain, to vary independently. We then calculated the channel capacities corresponding to VG-SPN configurations given by different (CI max , Gain) pairs using previous codes (3) that are based on the framework developed by (2). We also calculated the R max corresponding to each (CI max , G) configuration, which is simply given by CI max *Gain. Effectively, the value of R max relates to the number of NEMO or other response molecules in the cytoplasm that can be activated by a VG-SPN configuration. The variables and hyperparameter values used in the model are as described in tables S5 and S6.

Statistical analyses Distributions of physical properties of NEMO puncta were compared to see whether they were statistically different from one another (Fig. 2C). The distributions were compared with a Student’s t test (using native MATLAB functions), yielding P values for the null hypothesis that the two distributions were from the same (normal) distribution, with identical means and SDs. Each null hypothesis was rejected at the 5% significance level with P values «10−20 for every TNF and IL-1 pair with the same dosage condition. Spearman’s rank correlation and linear regression were performed on NEMO predictors for the RelA responses (and the log-transformed predictors/responses) for cells responding to TNF and IL-1. Native MATLAB functions were used to calculate all coefficients. Multiple linear regression was also performed on NEMO predictors for RelA AUC fold response in cells responding to TNF and IL-1 to see whether combinations of descriptors improve R2 values. R2 values shown in fig. S6 are the best R2 determined for the indicated number of predictors when considering all possible predictor combinations. Combination of predictors provided only marginal improvements, with maximum R2 values approaching 0.7 and 0.8 for the linear and log transformed (respectively).

History of Flight: Breakthroughs, Disasters and More

img ]

From hot-air balloons floating over Paris to a dirigible crashing over New Jersey, here are some of the biggest moments of aviation history.

For thousands of years, humans have dreamed of taking to the skies. The quest has led from kite flying in ancient China to hydrogen-powered hot-air balloons in 18th-century France to contemporary aircraft so sophisticated they can’t be detected by radar or the human eye.

Below is a timeline of humans’ obsession with flight, from da Vinci to drones. Fasten your seatbelt and prepare for liftoff.

1505-06: Da Vinci dreams of flight, publishes his findings

Self-portrait by Leonardo da Vinci DEA / A. DAGLI ORTI/Getty Images

Few figures in history had more detailed ideas, theories and imaginings on aviation as the Italian artist and inventor Leonardo da Vinci. His book Codex on the Flight of Birds contained thousands of notes and hundreds of sketches on the nature of flight and aerodynamic principles that would lay much of the early groundwork for—and greatly influence—the development of aviation and manmade aircraft.

November 21, 1783: First manned hot-air balloon flight

Two months after French brothers Joseph-Michel and Jacques-Étienne Montgolfier engineered a successful test flight with a duck, a sheep and a rooster as passengers, two humans ascended in a Montgolfier-designed balloon over Paris. Powered by a hand-fed fire, the paper-and-silk aircraft rose 500 vertical feet and traveled some 5.5 miles over about half an hour. But in an 18th-century version of the space race, rival balloon engineers Jacques Alexander Charles and Nicholas Louis Robert upped the ante just 10 days later. Their balloon, powered by hydrogen gas, traveled 25 miles and stayed aloft more than two hours.

READ MORE: 6 Little-Known Pioneers of Aviation

1809-1810: Sir George Cayley introduces aerodynamics

At the dawn of the 19th century, English philosopher George Cayley published “On Aerial Navigation,” a radical series of papers credited with introducing the world to the study of aerodynamics. By that time, the man who came to be known as “the father of aviation” had already been the first to identify the four forces of flight (weight, lift, drag, thrust), developed the first concept of a fixed-wing flying machine and designed the first glider reported to have carried a human aloft.

WATCH: ‘The Machines That Built America’ premieres Sunday, July 18 at 9/8c. Watch a preview now.

September 24, 1852: Giffard’s dirigible proves powered air travel is possible

Half a century before the Wright brothers took to the skies, French engineer Henri Giffard manned the first-ever powered and controllable airborne flight. Giffard, who invented the steam injector, traveled almost 17 miles from Paris to Élancourt in his “Giffard Dirigible,” a 143-foot-long, cigar-shaped airship loosely steered by a three-bladed propeller that was powered by a 250-pound, 3-horsepower engine, itself lit by a 100-pound boiler. The flight proved that a steam-powered airship could be steered and controlled.

1876: The internal combustion engine changes everything

Building on advances by French engineers, German engineer Nikolaus Otto devised a lighter, more efficient, gas-powered combustion engine, providing an alternative to the previously universal steam-powered engine. In addition to revolutionizing automobile travel, the innovation ushered in a new era of longer, more controlled aviation.

December 17, 1903: The Wright brothers become airborne—briefly

Orville and Wilbur Wright on a test flight

Flying from Kitty Hawk, North Carolina, brothers Orville and Wilbur Wright made the first controlled, sustained flight of a heavier-than-air aircraft. Each brother flew their wooden, gasoline-powered propeller biplane, the “Wright Flyer,” twice (four flights total), with the shortest lasting 12 seconds and the longest sustaining flight for about 59 seconds. Considered a historic event today, the feat was largely ignored by newspapers of the time, who believed the flights were too short to be important.

READ MORE: 10 Things You May Not Know About the Wright Brothers

1907: The first helicopter lifts off

A vintage French postcard featuring the helicopter of Paul Cornu of Lisieux, France, who piloted the first manned flight of a rotary wing aircraft in November 1907. Paul Popper/Popperfoto via Getty Images

French engineer and bicycle maker Paul Cornu became the first man to ride a rotary-wing, vertical-lift aircraft, a precursor to today’s helicopter, when he was lifted about 1.5 meters off the ground for 20 seconds near Lisieux, France. Versions of the helicopter had been toyed with in the past—Italian engineer Enrico Forlanini debuted the first rotorcraft three decades prior in 1877. And it would be improved upon in the future, with American designer Igor Sikorsky introducing a more standardized version in Stratford, Connecticut in 1939. But it was Cornu’s short flight that would land him in the history books as the definitive first.

1911-12: Harriet Quimby achieves two firsts for women pilots

Journalist Harriet Quimby became the first American woman ever awarded a pilot’s license in 1911, after just four months of flight lessons. Capitalizing on her charisma and showmanship (she became as famous for her violet satin flying suit as for her attention to safety checks), Quimby achieved another first the following year when she became the first woman to fly solo across the English channel. The feat was overshadowed, however, by the sinking of the Titanic two days earlier.

October 1911: The aircraft becomes militarized

Italy became the first country to significantly incorporate aircraft into military operations when, during the Turkish-Italian war, it employed both monoplanes and airships for bombing, reconnaissance and transportation. Within a few years, aircraft would play a decisive role in the World War I.

READ MORE: The High-Flying Fraudster of Early Aviation

January 1, 1914: First commercial passenger flight

On New Year’s Day, pilot Tony Jannus transported a single passenger, Mayor Abe Pheil of St. Petersburg, Florida across Tampa Bay via his flying airboat, the “St. Petersburg-Tampa Airboat Line.” The 23-mile flight (mostly along the Tampa Bay shore) cost $5.00 and would lay the foundation for the commercial airline industry.

READ MORE: How America’s Aviation Industry Got Its Start Transporting Mail

1914-1918: World war accelerates the militarization of aircraft

World War I became the first major conflict to use aircraft on a large-scale, expanding their use in active combat. Nations appointed high-ranking generals to oversee air strategy, and a new breed of war hero emerged: the fighter pilot or “flying ace.”

According to The Illustrated Encyclopedia of Military Aircraft, France was the war’s leading aircraft manufacturer, producing nearly 68,000 planes between 1914 and 1918. Of those, nearly 53,000 were shot down, crashed or damaged.

READ MORE: 6 Famous World War I Flying Aces

June 1919: First nonstop transatlantic flight

The first nonstop transatlantic flight ended with a nosedive into a bog in western Ireland. The pilots walked away unscathed. Topical Press Agency/Getty Images

Flying a modified ‘Vickers Vimy’ bomber from the Great War, British aviators and war veterans John Alcock and Arthur Brown made the first-ever nonstop transatlantic flight. Their perilous 16-hour journey, undertaken eight years before Charles Lindbergh crossed the Atlantic alone, started in St. John’s, Newfoundland, where they barely cleared the trees at the end of the runway. After a calamity-filled flight, they crash-landed in a peat bog in County Galway, Ireland; remarkably, neither man was injured.

READ MORE: The First Nonstop Flight Across the Atlantic Lasted 16 Harrowing Hours

1921: Bessie Coleman becomes the first Black woman to earn a pilot’s license

Pioneering aviatrix Bessie Coleman Michael Ochs Archives/Getty Images

The fact that Jim Crow-era U.S. flight schools wouldn’t accept a Black woman didn’t stop Bessie Coleman. Instead, the Texas-born sharecropper’s daughter, one of 13 siblings, learned French so she could apply to the Caudron Brothers’ School of Aviation in Le Crotoy, France. There, in 1921, she became the first African American woman to earn a pilot’s license. After performing the first public flight by a Black woman in 1922—including her soon-to-be trademark loop-the-loop and figure-8 aerial maneuvers—she became renowned for her thrilling daredevil air shows and for using her growing fame to encourage Black Americans to pursue flying. Coleman died tragically in 1926, as a passenger in a routine test flight. Thousands reportedly attended her funeral in Chicago.

READ MORE: 7 Death-Defying Historic American Daredevils

1927: Lucky Lindy makes first solo transatlantic flight

Charles Lindbergh with his historic plane ‘Spirit of St. Louis’

Nearly a decade after Alcock and Brown made their transatlantic flight together, 25-year-old Charles Lindbergh of Detroit was thrust into worldwide fame when he completed the first solo crossing, just a few days after a pair of celebrated French aviators perished in their own attempt. Flying the “Spirit of St. Louis” aircraft from New York to Paris, “Lucky Lindy” made the first transatlantic voyage between two major hubs—and the longest transatlantic flight by more than 2,000 miles. The feat instantly made Lindbergh one of the great folk heroes of his time, earned him the Medal of Honor and helped usher in a new era of interest in the possibilities of aviation.

READ MORE: 10 Fascinating Facts About Charles Lindbergh

1932: Amelia Earhart repeats Lindbergh’s feat

Aviatrix Amelia Earhart SSPL/Getty Images

Five years after Lindbergh completed his flight, “Lady Lindy” Amelia Earhart became the first woman to fly solo across the Atlantic Ocean, setting off from Harbour Grace, Newfoundland on May 20, 1932 and landing some 14 hours later in Culmore, Northern Ireland. In her career as an aviator, Earhart would become a worldwide celebrity, setting several women’s speed, domestic distance and transcontinental aviation records. Her most memorable feat, however, would prove to be her last. In 1937, while attempting to circumnavigate the globe, Earhart disappeared over the central Pacific ocean and was never seen or heard from again.

READ MORE: Tantalizing Theories About the Earhart Disappearance

1937: The Hindenburg crashes…along with the ‘Age of Airships’

The Hindenburg bursts into blames above Lakehurst, New Jersey, on May 6, 1937. (Popperfoto/Getty Images)

Between WWI and WWII, aviation pioneers and major aircraft companies like Germany’s Luftshiffbau Zeppelin tried hard to popularize bulbous, lighter-than-air airships—essentially giant flying gas bags—as a mode of commercial transportation. The promise of the steam-powered, hydrogen-filled airships quickly evaporated, however, after the infamous 1937 Hindenburg disaster. That’s when the gas inside the Zeppelin company’s flagship Hindenburg vessel exploded during a landing attempt, killing 35 passengers and crew members and badly burning the majority of the 62 remaining survivors.

READ MORE: The Hindenburg Disaster: 9 Surprising Facts

October 14, 1947: Chuck Yeager breaks the sound barrier

An ace combat fighter during WWII, Chuck Yeager earned the title “Fastest Man Alive” when he hit 700 m.p.h. while testing the experimental X-1 supersonic rocket jet for the military over the Mojave Desert in 1947. Being the first person to travel faster than the speed of sound has been hailed as one of the most epic feats in the history of aviation—not bad for someone who got sick to his stomach after his first-ever flight.

READ MORE: 6 Renowned Tuskegee Airmen

1949: The world’s first commercial jetliner takes off

Early passenger air travel was noisy, cold, uncomfortable and bumpy, as planes flew at low altitudes that brought them through, not above, the weather. But when the British-manufactured de Havilland Comet took its first flight in 1949—boasting four turbine engines, a pressurized cabin, large windows and a relatively comfortable seating area—it marked a pivotal step in modern commercial air travel. An early, flawed design however, caused the de Havilland to be grounded after a series of mid-flight disasters—but not before giving the world a glimpse of what was possible.

1954-1957: Boeing glamorizes flying

With the debut of the sleek 707 aircraft, touted for its comfort, speed and safety, Seattle-based Boeing ushered in the age of modern American jet travel. Pan American Airways became the first commercial carrier to take delivery of the elongated, swept-wing planes, launching daily flights from New York to Paris. The 707 quickly became a symbol of postwar modernity—a time when air travel would become commonplace, people dressed up to fly and flight attendants reflected the epitome of chic. The plane even inspired Frank Sinatra’s hit song “Come Fly With Me.”

READ MORE: The Chinese-Born Engineer Who Helped Launch US Commercial Aviation

March 27, 1977: Disaster at Tenerife

In the greatest aviation disaster in history, 583 people were killed and dozens more injured when two Boeing 747 jets—Pan Am 1736 and KLM 4805—collided on the Los Rodeos Airport runway in Spain’s Canary Islands. The collision occurred when the KLM jet, trying to navigate a runway shrouded in fog, initiated its takeoff run while the Pan Am jetliner was still on the runway. All aboard the KLM flight and most on the Pan Am flight were killed. Tragically, neither plane was scheduled to fly from that airport on that day, but a small bomb set off at a nearby airport caused them both to be diverted to Los Rodeos.

1978: Flight goes electronic

The U.S. Air Force developed and debuted the first fly-by-wire operating system for its F-16 Fighting Falcon fighter plane. The system, which replaced the aircraft’s manual flight control system with an electronic one, ushered in aviation’s “Information Age,” one in which navigation, communications and hundreds of other operating systems are automated with computers. This advance has led to developments like unmanned aerial vehicles and drones, more nimble missiles and the proliferation of stealth aircraft.

READ MORE: Automation of Planes Began 9 Years After the Wright Brothers Took Flight, But It Still Leads to Baffling Disasters

1986: Around the world, without landing

American pilots Dick Rutan and Jeana Yeager (no relation to Chuck) completed the first around-the-world flight without refueling or landing. Their “Rutan Model 76 Voyager,” a single-wing, twin-engine craft designed by Rutan’s brother, was built with 17 fuel tanks to accommodate long-distance flight.

Biden insists US ‘following the science’ on masks

img ]

Axios

CDC officials are concerned about a strain of the Candida auris fungus that’s resistant to all drugs and appears to have spread in small clusters in health care settings, rather than in individuals who had taken antifungals.Why it matters: “The concern is that it could spread to any of the patients who are at high risk, not just the ones who’ve been treated before — and that the population who could acquire these potentially untreatable infections could be much larger,” Meghan Lyman, medical off