Temporal signals drive the emergence of multicellular information networks

Edited by Robert Austin, Princeton University, Princeton, NJ; received February 7, 2022; accepted July 18, 2022
September 6, 2022
119 (37) e2202204119

Significance

Understanding how a group of cells cooperatively processes an environmental signal is an essential step to decoding the organizing principles of multicellular organisms. Here we demonstrate that neuronal cells form information-bearing causal networks through gap junction–mediated communication. Our experimental and theoretical results uncover the mechanism by which excitable cells self-organize in response to external stimuli and the rich collective dynamics enabled by nonsynaptic intercellular interactions. Decoding such systems will lead to a better understanding of a diverse range of physiological processes, offering insights into disease mechanisms and treatment.

Abstract

Coordinated responses to environmental stimuli are critical for multicellular organisms. To overcome the obstacles of cell-to-cell heterogeneity and noisy signaling dynamics within individual cells, cells must effectively exchange information with peers. However, the dynamics and mechanisms of collective information transfer driven by external signals are poorly understood. Here we investigate the calcium dynamics of neuronal cells that form confluent monolayers and respond to cyclic ATP stimuli in microfluidic devices. Using Granger inference to reconstruct the underlying causal relations between the cells, we find that the cells self-organize into spatially decentralized and temporally stationary networks to support information transfer via gap junction channels. The connectivity of the causal networks depends on the temporal profile of the external stimuli, where short periods, or long periods with small duty fractions, lead to reduced connectivity and fractured network topology. We build a theoretical model based on communicating excitable units that reproduces our observations. The model further predicts that connectivity of the causal network is maximal at an optimal communication strength, which is confirmed by the experiments. Together, our results show that information transfer between neuronal cells is externally regulated by the temporal profile of the stimuli and internally regulated by cell–cell communication.
Sensing and responding to chemical signals is of fundamental importance to living systems. For single cells, chemosensing is achieved by specialized receptors, which recognize molecules (ligands) in the microenvironment of cells (1, 2). Such interactions trigger a cascade of intracellular events, which regulate the functional responses of cells, such as motility (3), differentiation (4), and gene expression (4, 5). However, many chemosensing architectures determine the ligand concentration from time-integrated information, such as receptor–ligand binding and dissociation times (6, 7). As such, dynamic external stimuli can present a challenge to cell sensing. For instance, oscillatory stimuli may be misinterpreted by cells, as external and internal time scales interfere in the signaling dynamics (8).
In multicellular organisms, chemosensing is rarely accomplished by isolated single cells. Instead, collective chemosensing by communicating cells leads to rich dynamics that may be necessary to encode complex information (9). In collective chemosensing, environmental signals can induce specific single-cell dynamics that are regulated by cell–cell communication (10). For instance, we and other groups have shown that when chemosensing pathways support bifurcating signaling dynamics, cell–cell communication can shift the bifurcation boundary, so that the resulting cell response reflects both the external signal as well as the degree of communication (11, 12).
Collective chemosensing can be manifested as orchestrated multicellular dynamics, such as intercellular synchronization. Synchronized cellular dynamics have been observed in cardiac tissues (13), endothelium (14), and in the hypothalamic suprachiasmatic nucleus (15). Synchronization often requires strong external stimuli and efficient cell–cell communication to offset the intrinsic and extrinsic noise in the dynamics of individual cells (16).
Alternatively, collective chemosensing may induce a group of communicating cells to self-organize into networks that support asymmetric interactions and directed information flow. In particular, environmental stimuli facilitate the emergence of leader, follower, and pacemaker cells such as in beating cardiac tissues (17), in social amoebae that form fruiting bodies (18), and in the neuronal regulation of circadian rhythms (19). However, the hierarchical organization is often obscured by fluctuations of single-cell dynamics and requires sophisticated data analysis to reconstruct the underlying network. Information-theoretic metrics, such as Granger inference (20) and its nonparametric form of transfer entropy (21), have been instrumental in elucidating the intercellular wiring hidden from direct observations (22). Despite its biological significance, the underlying mechanisms, upstream control, and downstream function of collective chemosensing are still far from fully understood.
In this study we combine quantitative experiments and computational modeling of excitable cells to investigate the emergence of information-bearing networks when monolayers of neuronal cells sense extracellular ATP (adenosine triphosphate). We examine the calcium dynamics of KTaR cells, a neuronal cell line we derived from KNDy (Kisspeptin, neurokinin B, and dynorphin) neurons within the arcuate nucleus of an adult female mouse (23). We show that under periodic stimuli a group of interacting cells forms a directed causal network which maintains dynamic equilibrium over consecutive cycles of stimuli. The network characteristics depend not only on the level of communication between cells but also on the temporal profile of the external driving. Together, we demonstrate that temporal signals from the environment instruct the self-organization and communication dynamics of a multicellular system.

Results

Periodic Stimuli Drive the Formation of Multicellular Information Network in Dynamic Equilibrium.

In order to understand the collective dynamics of communicating cells under periodic stimuli, we employ a microfluidic device as shown in Fig. 1A. A computer-interfaced flow switch alternates growth medium and ATP solution into the cell culture chamber, where a confluent monolayer of KTaR cells sense the ATP stimuli (see SI Appendix, section S1a–c for more details of device and cell characterization). To detect cellular response, we preload the cells with a calcium indicator (Calbryte, AAT Bioquest) and record the fluorescent calcium images at single-cell resolution at 1 Hz for over 15 min.
Fig. 1.
Experimental setup to uncover the underlying self-organization of KTaR cell monolayers. (A) A schematic showing the microfluidics device to deliver alternate growth medium (GM) and ATP solution (ATP) to a confluent monolayer of KTaR cells. (Inset) A fluorescent Ca2+ image of a monolayer of KTaR cells. (B) Typical calcium responses R(t) of KTaR cells to cyclic ATP stimuli at a period of 200 s. Heterogeneity among the cells leads to fluctuations in the magnitude and temporal delay of the calcium dynamics. (C and D) An example of reconstructed multicellular network via Granger inference. Direction of arrows point from a causal cell to its affected cells. Each node corresponds to the location of a cell in the field of view. The node are compositely colored by leader (red channel) and follower (blue) scores in C and authority and hub scores in D. See also SI Appendix, section S2e.
KTaR cells recognize extracellular ATP with purigenic receptors, which trigger IP3-mediated release of Ca2+ from endoplasmic reticulum stores into the cytoplasm, as well as calcium influx from extracellular space (24) (SI Appendix, section S1c and d). While overall, the relative change of intensity [Ri(t), where i is the cell index] follows the temporal profile of ATP stimuli, individual cells show variable phase delays to the global driving signal (Fig. 1B; see also SI Appendix, section S2).
To quantify if the asynchronous responses of individual cells encode information transfer, we employ Granger inference (25, 26) to construct a directed graph that represents the causal influence between cells. Qualitatively, Granger inference designates a time series C as causing a second time series E if the combined history of both {C, E} is significantly more predictive of time series E than E’s own history alone. Because the rapid flow effectively washes away secreted factors (12, 27) and because KTaR cells do not grow extended axons in our culture condition, we focus on nearest-neighbor cells where gap junctional communication is dominant (28) (SI Appendix, section S1c).
As we have shown previously, the time-derivative of fluorescent calcium intensity R˙i(t) has the benefit of being independent of the basal intensity while still measuring communication effects (11, 12). We have further confirmed that {R˙i(t)} are stationary time series (SI Appendix, section S2a) and therefore suitable for the application of Granger inference (26).
For each nearest-neighbor pair, we calculate the statistical significance of Granger difference (25) using the time series from a particular cycle. If higher than a threshold (95% confidence), an edge from the causal cell to the affected cell is drawn (see SI Appendix, section S2b for more details). Fig. 1 C and D show an example of a reconstructed causal network, where each node represents the location of a cell in the field of view, and the arrows show direction of causality.
After reconstructing the directed graph, we have calculated the leader scores (number of outgoing edges) and follower scores (number of incoming edges) for each cell. The leader/follower scores distribute randomly in space (Fig. 1C), indicating the absence of centralized organization. Indeed, we find the nodes generally have very low authority and hub scores as measured by Kleinberg centrality (29) (Fig. 1D), and the networks come with small Estrada heterogeneity index (30, 31) (0.1; SI Appendix, section S2c). Also, there are few closed loops in the networks (4% of edges form loops; SI Appendix, section S2d), indicating asuppression of feedthrough information relays, presumably by cell–cell communication noises (32). These observations suggest that cyclic external stimuli trigger information transfer between communicating KTaR cells. Although heterogeneity among the cells prevents fully synchronized responses, the cells are able to self-organize into a decentralized causal network.
Having established methods to reconstruct the underlying networks of cells performing collective chemosensing, we first examine the evolution of the network structure over consecutive cycles of ATP stimuli (see also SI Appendix, section S3a). To this end, we compute Padd, the rate (probability per cycle) of adding a new edge; Pdel, the rate of deleting an existing edge; and Pflp, the rate of flipping the direction of an existing edge (Fig. 2A). We find ~60% of edges are deleted from one stimulus cycle to the next, while a new edge would emerge from ~30% of the unconnected neighbor cell pairs (Fig. 2B). Among existing edges, less than 10% of them will flip direction in the next cycle, indicating a memory effect that stabilizes the causal relation between cell pairs.
Fig. 2.
Dynamic evolution of the multicellular networks driven by cyclic ATP stimulation. (A) An example multicellular network rewires between two consecutive cycles. The cells are exposed to 50 µM ATP at a period of 120 s. (Insets) Three types of rewiring events governed by their respective rates (probabilities per cycle): removing an edge (Pdel), adding a new edge (Padd), and flipping the direction of an existing edge (Pflp). See also SI Appendix, section S3f. (B) The rates for removing, adding, and flipping an edge at various driving periods T. Cyan indicates T= 20 s, magenta indicates T= 120 s, and yellow indicates T= 200 s. For a given period, the rates do not depend on ATP concentration (SI Appendix, section S3b). (C) Scatterplot showing the numbers of added (Eadd) and removed (Edel) edges between consecutive cycles normalized by the total number of nearest neighbors (Etot). Here colors represent the driving period as in B. Different symbols represent ATP concentration. Circle indicates 10 µM, triangle indicates 50 µM, and upside-down triangle indicates 100 µM.
Although the values of {Padd,Pdel,Pflp} do not depend on the concentrations of ATP (SI Appendix, section S3b), nor the local connectivity (SI Appendix, section S3c), we find that at a period of 20 s all three rates are smaller compared with larger periods (Fig. 2B). In all conditions, the network remains approximately stationary as the number of new edges matches the number of removed edges over consecutive cycles (Fig. 2C; see also SI Appendix, section S3d). Taken together, these observations show that cyclic ATP stimuli drive monolayers of KTaR cells into networks that maintain their dynamic equilibrium.

Temporal Profiles of External Stimuli Controls Multicellular Network Connectivity.

After showing the multicellular network to be stationary, we investigate whether the degree of network connectivity depends on the spatial relations between cells. We calculate the edge probability Pedge, which is defined as the number of edges divided by the number of nearest neighbors. In particular, we compare the edge probability of the original networks (directly obtained from experiments) and ones obtained by randomizing the original networks. The randomization is done by shuffling the time series of 10% of the cells with another 10% cells in the same experiment: RiRj, where ij are randomly chosen pairs. The randomized data encode identical driving signal to the original data. We find that for all experiments, even a 10% partial randomization significantly dilutes the edges (Fig. 3A), and Pedge can be reduced by as much as 20% (Fig. 3 A, Inset). The result highlights the locality of cell–cell interaction, which is consistent with gap junction–mediated information flow between nearby cells.
Fig. 3.
Multicellular network connectivity is regulated by the period of ATP stimulation. (A) The edge probability of original experimental data and the change of edge probability ΔPedge after randomization. In a particular experiment the randomization is done by switching the calcium dynamics of 10% of randomly selected cell pairs in the same field of view. (Inset) Histogram of the relative change of edge probability after randomization. (B) The dependence of edge probability with respect to driving period. (Inset) The dependence of edge probability with respect to ATP concentration. (C) The dependence of percolation degree k2/k with respect to driving period. Here k represents the degree of the (undirected) network. (Top Inset) The dependence of percolation strength with respect to ATP concentration. (Right Insets) Maps of clusters in two typical experiments (T= 20 s and T= 120 s). Nodes are colored by the size of clusters they belong to. Here the normalized cluster size is defined as the ratio between the number of cells in the cluster to the total number of cells in the field of view. See also SI Appendix, section S3f. In B and C, colors of symbols represent the driving period. Cyan indicates T= 20 s, magenta indicates T= 120 s, and yellow indicates T= 200 s. The types of symbols represent ATP concentration: circle indicates 10 µM, triangle indicates 50 µM, and upside-down triangle indicates 100 µM. Statistical comparisons are done with ANOVA. **P < 0.01, ***P < 0.001; n.s., not significant.
Having demonstrated that the short-range intercellular communication is manifested by the edge probability, we now examine what aspects of the external signal control the network connectivity. To this end we systematically vary the ATP concentration [ATP] and period of stimuli T. We find Pedge increases dramatically when the driving period increases from 20 to 120 s and plateaus without further changing when the driving period is further increased to 200 s (Fig. 3B). At long driving periods, the edge probability falls in the range of 0.3 to 0.5, corresponding to two to three edges per cell (on average, each cell has six nearest neighbors). On the other hand, the network characteristics do not depend on the ATP concentration over the physiological range of 10 to 100 µM (33) (Fig. 3 B, Inset).
To further compare the self-organized structure of multicellular networks at varying driving signals, we compute the percolation degree N=k2/k. N is the average degree (i.e., total number of edges to a node) of a node if the node is linked to another node. When N > 2, the network is above the percolation threshold and features large connected components (clusters). When N < 2, the network is below the percolation threshold and consists of many small clusters (34). We find when the driving period equals 20 s, the percolation degree is less than 2 (Fig. 3C). A typical network in this case (Fig. 3 C, Top Right Inset) indeed shows fractured topology where none of the clusters contain more than 20% of the cells in the field of view. In contrast, at larger driving period the percolation degree is greater than 2 (Fig. 3C) such that a typical network is dominated by a single large cluster (Fig. 3 C, Bottom Right Inset). These results show that there exists a time scale dependence of the external driving signal that dictates the underlying information flow of collective sensing. At a small driving period, KTaR cells form a loosely connected, fractured network. Conversely at large driving periods, highly connected and percolating networks emerge thanks to the elevated information flow between cells. These two distinct types of multicellular organization are induced by the temporal profiles of the driving signal, rather than the concentration of the stimuli.
To understand the mechanisms by which temporal signals drive the emergence of multicellular networks, we developed a mathematical model of communicating excitable cells. Each cell is modeled using a reduced form of the Hodgkin–Huxley model (35), which is widely accepted to replicate neuronal dynamics. Specifically, because the most probable experimental nearest-neighbor number is six, we model cells on a six-neighbor triangular lattice as in Fig. 4A (adding neighbor number variability to the same degree as that observed in experiments does not significantly affect our results; SI Appendix, section S4). Within each cell, two chemical species interact (Fig. 4A), which is the minimum needed for excitable dynamics (36): X, which represents calcium abundance, and Y, which represents a slower recovery variable. The following minimal reactions are chosen to produce excitations: X activates both itself (37, 38) (Fig. 4A, first two reactions) and Y (third reaction), while Y represses X (fourth reaction). Both X and Y degrade spontaneously (first and fifth reactions), and X is exchanged between neighboring cells to model the gap junction communication (sixth reaction). Transforming the rate equations of this model into a standard form (SI Appendix, section S4) makes clear that the dynamics are specified by 1) a characteristic molecule number xc, 2) a time scale separation ϵ between X and Y, and 3) an external field h that tunes the system among four regimes: stable dynamics at low molecule number, excitable dynamics, oscillatory dynamics, and stable dynamics at high molecule number (Fig. 4A). The standard form is akin to the FitzHugh–Nagumo model (39, 40), which is a reduced representation of more complex excitation models such as the Tang–Othmer (41) and Hodgkin–Huxley (35) models. Here we focus on the transition between the stable low and excitable regimes at h = hc and the effects of communication on this transition.
Fig. 4.
Mathematical model of collective excitable dynamics. (A) Each cell in a triangular lattice contains a calcium variable subject to positive and negative feedback reactions (via a slower recovery variable) to enable excitations and exchange reactions to enable nearest-neighbor communication. (B) An excitable cell (cell 1, h>hc) can induce a nonexcitable cell (cell 2, h>hc) to excite with a delay via cell–cell communication. (C) Fraction of nearest neighbors with causal edges (edge probability) increases and saturates with stimulus period in model. Here the presence of an edge is determined by the delay between neighbors’ excitations (SI Appendix, section S4). (D) Edge probability decreases for small or large duty fraction (fraction of period for which stimulus is on) in (Left) model and (Right) experiments. Experiments are done with 50 µM ATP at a period of 120 s. In D and E, Left, error bars are SD over 100 simulations of eight cycles.
Modeling ATP as setting the value of the field h, we find that communication between an excitable cell (h>hc) and a nonexcitable cell (h<hc) can induce an excitation in the nonexcitable cell (Fig. 4B), albeit with a delay. Indeed, applying Granger inference to stochastic simulations (42) of the model, we find that the excitable cell “Granger-causes” the nonexcitable cell in this case. For simplicity, in larger networks we take the peak delay between neighboring cells as a proxy for the Granger metric when determining edges, as we find that the two are correlated in the simulations (SI Appendix, section S4). We then investigate networks of similar size to the experimental viewing window (we find that network size does not significantly affect our results; SI Appendix, section S4c). Field strengths h are drawn from a normal distribution centered just above hc, such that slightly more than half of cells are excitable from the stimulus alone.
In these model networks we find that the edge probability increases and then saturates with the driving period (Fig. 4C). This finding is consistent with the experiments (Fig. 3B), which validates the model. The intuitive reason is that when the driving period is shorter than the excitation time scale, cells are still in the recovery phase and cannot respond, which reduces causal information. This intuition continues to hold if the on- and off-times are not constant but instead drawn from an exponential distribution (SI Appendix, section S4d). This intuition also holds when fixing the period but varying the duty fraction: if either the on- or off-portion of the cycle is too brief, the edge probability is reduced (Fig. 4 D, Left). Varying the duty fraction in the experiments, we see that this prediction is upheld (Fig. 4 D, Right), in further support of the model.

External and Internal Time Scales Jointly Regulate the Structure of Multicellular Information Networks.

Having investigated the relationship between the external driving and intracellular excitation time scales, we now use the model to investigate the effects of changing the strength of intercellular communication. Upon varying the cell–cell coupling constant g over 4 decades, we find that stronger communication leads to a higher fraction of cells exhibiting an excitable response, x>xc (Fig. 5A). Evidently, for sufficiently strong communication, all nonexcitable cells can be induced to excitation by communication alone. This result is consistent with our previous report in fibroblast monolayers that communication can augment the effects of external stimuli by modulating the bifurcation threshold of excitable cell dynamics (12).
Fig. 5.
Collective sensory responses at varying levels of communication. (A) Model prediction of the fraction of firing cells as the coupling constant g changes over 4 decades. (B) The fraction of KTaR cells with NWS greater than 0.5 for three cases: cells are treated by 10 µM palmitoleic acid to inhibit gap junction and exposed to 50µM ATP at a period of 120 s (CX), untreated cells are exposed to 50µM ATP at a period of 120 s (control), and cells are exposed to 30 mM KCl at a period of 120 s (KCl). Error bars represent SD from 500 bootstrap samples. (Top) The NWS values shown as heat maps. Horizontal axis of the heat maps represents time, where the black dashed lines indicate t=2T= 240 s. NWS values of individual cells are stacked vertically. See SI Appendix, section S3 for more details. (C) Model prediction of the edge probability as the coupling constant g changes over 4 decades. (D) Experimentally measured edge probability for the same conditions as in B. See also SI Appendix, section S3f. Statistical comparisons are done with ANOVA. **P < 0.01, ***P < 0.001; n.s., not significant. N = 3. Error bars are SD.
To test the model prediction, we devised methods to either reduce or enhance intercellular communication experimentally (SI Appendix, section S1). To reduce communication, we treat the KTaR cells with 10 µM palmitoleic acid, a broad spectrum gap junction inhibitor (43). To enhance communication, we replace ATP with KCl as the external stimulus. KCl depolarizes the KTaR cell membrane, triggering an action potential as well as intracellular calcium responses (44). As a result, the cells couple electrically through gap junctions, which is faster compared with the diffusion-limited molecular exchange.
Unlike in the model where the criterion for excitation is self-evident (x>xc), in the experiments the criterion must be defined from the response itself. To this end, we define a normalized wavelet score (NWS) to quantify the cellular dynamics in the frequency domain using a time-resolved wavelet transformation (SI Appendix, section S3). If the calcium dynamics of a cell perfectly follows the driving frequency, its NWS equals 1 at all times (except for boundary effects that affect the beginning and end of the time series). Otherwise, the NWS will fluctuate between 0 and 1 when irregular response occurs.
We find that the NWS of cells treated with palmitoleic acid (abbreviated as CX-) show significantly stronger fluctuations compared with untreated cells, whereas the NWS of cells stimulated with KCl quickly reach and stay close to 1 (Fig. 5 B, Top). To compare with the model prediction, we calculate the fraction of cells with an NWS greater than 0.5 at time point t=2T to avoid boundary effects, where T is the driving period set to be 120 s. Consistent with the model, we find the fraction of cells with NWS greater than 0.5 is highest for KCl excited cells and lowest for gap junction–inhibited cells (Fig. 5 B, Bottom).
Interestingly, our model also predicts that the network connectivity reaches a maximum at an optimal coupling constant g and decreases in either direction from the optimal value (Fig. 5C). The intuitive reason is that with weak communication, only the inherently excitable cells are responding, such that causal edges do not form with nonexcitable cells. Intermediate communication induces nonexcitable cells to excite with a delay, introducing new edges. Strong communication synchronizes cells, reducing causality and removing edges.
Experiments confirm that the edge probabilities for gap junction–inhibited monolayers and for KCl excited monolayers are both lower than the untreated KTaR cells exposed to ATP stimuli (Fig. 5D). Consistently, typical networks of untreated cells show characteristics of percolation, while networks under the other two conditions are evidently fractured (Fig. 5 D, Inset). These observations suggest that under the control condition the KTaR monolayers are posed close to the optimal coupling strength for causal information flow. Inhibiting gap junction curtails cell–cell communication, reducing network connectivity, whereas accelerating cell–cell communication leads to rapid synchronization between neighboring cells, also reducing information flow. Together, our results demonstrate that the self-organization of multicellular networks is modulated by the level of cell–cell communication.

Discussion

A group of interacting cells encodes environmental information in different forms than single cells do. Revealing the underlying principles of collective chemosensing is an essential step to understanding the rules of life. Here we study the external ATP-triggered calcium dynamics of neuronal cell monolayers. We employ microfluidics to deliver alternate ATP solution and pure growth medium to KTaR-1 cells, a neuronal cell line we derived from KNDy neurons within the arcuate nucleus of an adult female mouse. KTaR cells express connexin proteins in vitro which constitute gap junction channels between adjacent cells. Using Granger inference, we show that during each ATP–growth medium cycle, there is asymmetric information flow between adjacent cells manifested as causal relations between their intracellular calcium dynamics. As a result, the external stimuli drive the neuronal cell monolayers to establish directed networks. These networks display hierarchical structure where leader and follower cells distribute spatially without any apparent centralized organization (Fig. 1).
The information networks are highly dynamic from one cycle of stimuli to the next, while the overall connectivity remains stationary. For all conditions tested, most structural fluctuations of the networks manifest as adding or removing edges, whereas less than 10% of the edges flip directions over consecutive cycles (Fig. 2B). This suggests that the network reconfiguration is due to stochastic disappearance and reappearance of deterministic causal relationships that presumably arise from cell-to-cell heterogeneity (45). The time evolution of the networks shows characteristics of detailed balance. For instance, the number of edges remains approximately constant (Fig. 2C). The probability flux in the cellular state space defined by the leader/follower scores also vanishes (SI Appendix, section S2). This is in contrast to the nonequilibrium stationary states observed in other living systems, especially in the macroscopic brain dynamics (46). It is conceivable that higher-order organization in the brain leads to the emergence of entropy production that is absent at the scale of locally communicating neuronal cells.
Many neuronal systems demonstrate characteristics of learning and reinforcement (47). In contrast, we find that under repeated stimuli, gap junctions mediate a Markovian evolution of KTaR networks that keeps the system stationary. This is expected as gap junctions alone have rapid turnover time (48). It will be interesting for future studies to elucidate the mechanisms by which neuronal cells stabilize their information exchange dynamics.
We find that the edge probability of the multicellular network primarily depends on the time scale and is impervious to the magnitude of external stimuli (Fig. 3). Interestingly, both the experiments and the theoretical model show that the effective time scale of the external signal is determined by the lesser of the on- and off-duty cycles (Fig. 4). It makes sense that short on-times may be insufficient to trigger excitations (or sustain neighbor-induced excitations), but this result implies that short off-times are also insufficient to do so. This is likely a result of the need for a postexcitation recovery time, which is a generic property of excitable systems. Indeed, the minimal nature of the model suggests that our findings on network responses to temporal signals may be generalized to other multicellular excitable systems.
Our finding that an intermediate communication strength maximizes causal connectivity (Fig. 5) has implications for information propagation in multicellular systems. A metaanalysis which compares the coupling constant in the model and effective diffusivity between cells further validates this idea (SI Appendix, section S4e). In systems unlike ours, where a stimulus is localized or the medium itself is spatially directed, one expects that causal information should increase indefinitely with the communication strength between units. However, in systems like ours, where neither the stimulus nor the medium break symmetry, our results highlight an interesting regime where intermediate communication amplifies heterogeneity to create random but reproducible pathways of information flow. Such a regime may be important in systems where these initially spontaneous pathways are reinforced and built upon to break symmetry permanently, facilitating the formation of differentiated structures (49).
In vivo, KNDy neurons are critical for pubertal progression and sex steroid feedback involved in the neuroendocrine regulation of reproduction, expressing and secreting Kisspeptin, a peptide stimulatory to gonadotropin-releasing hormone (GnRH) secretion (50). It is also interesting for future investigations to determine the physiological effects of temporal sensitivity for KNDy neurons in vivo, which receive episodic afferent signals including glutamate, GABA, and neurokinin B in addition to ATP.
Broadly speaking, collective chemosensing in a multicellular system may behave as isolated noninteracting units, may form hierarchical information flow networks, or may achieve synchronization. Here we show that the transition between these scenarios is controlled by the interplay of two time scales: one that is set internally by the communication channels between the cells and one that is set externally by the driving signal. Our results highlight the rich dynamics exhibited by spatially coupled excitable units. Decoding such systems will lead to a better understanding of a diverse range of physiological processes from development to neuronal dynamics to tissue organization (51), offering insights into disease mechanisms and treatment (52, 53).

Materials and Methods

See SI Appendix for details of cell culture, microscopy, and image analysis. The statistical analysis and computer simulations are performed with MATLAB (MathWorks). See SI Appendix for details of the theory and the simulations.

Data, Materials, and Software Availability

Experimental recordings of calcium dynamics are available at https://figshare.com/s/a681463e1e69134f7320 (54). Simulation codes are available at https://github.com/rwl23/mutlicellular_information_networks (55).

Acknowledgments

G.L. and B.S. are supported by NSF grant PHY-1844627. R.L., A.M., G.L., A.S., and P.C. are supported by National Institute of General Medical Studies grant R01GM140466. B.S. is supported by National Institute of General Medical Sciences grant R35GM138179.

Supporting Information

Appendix 01 (PDF)

References

1
P. G. Data, H. Acker, S. Lahiri, Neurobiology and Cell Physiology of Chemoreception (Springer Nature, 1993).
2
B. Alberts et al., Molecular Biology of the Cell (Garland Science, ed. 4, 2002), chap. 15.
3
J. d’Alessandro et al., Cell migration guided by long-lived spatial memory. Nat. Commun. 12, 4118 (2021).
4
M. A. Basson, Signaling in cell differentiation and morphogenesis. Cold Spring Harb. Perspect. Biol. 4, a008151 (2012).
5
W. Tang et al., BSKs mediate signal transduction from the receptor kinase BRI1 in Arabidopsis. Science 321, 557–560 (2008).
6
W. Bialek, S. Setayeshgar, Physical limits to biochemical signaling. Proc. Natl. Acad. Sci. U.S.A. 102, 10040–10045 (2005).
7
R. G. Endres, N. S. Wingreen, Maximum likelihood and the single receptor. Phys. Rev. Lett. 103, 158101 (2009).
8
A. Mitchell, P. Wei, W. A. Lim, Oscillatory stress stimulation uncovers an Achilles’ heel of the yeast MAPK signaling network. Science 350, 1379–1383 (2015).
9
D. Ellison et al., Cell-cell communication enhances the capacity of cell ensembles to sense shallow gradients during morphogenesis. Proc. Natl. Acad. Sci. U.S.A. 113, E679–E688 (2016).
10
Y. Fujii, S. Maekawa, M. Morita, Astrocyte calcium waves propagate proximally by gap junction and distally by extracellular diffusion of ATP released from volume-regulated anion channels. Sci. Rep. 7, 13115 (2017).
11
B. Sun, G. Duclos, H. A. Stone, Network characteristics of collective chemosensing. Phys. Rev. Lett. 110, 158103 (2013).
12
G. D. Potter, T. A. Byrd, A. Mugler, B. Sun, Communication shapes sensory response in multicellular networks. Proc. Natl. Acad. Sci. U.S.A. 113, 10334–10339 (2016).
13
N. N. Agladze et al., Synchronization of excitable cardiac cultures of different origin. Biomater. Sci. 5, 1777–1785 (2017).
14
B. Ubezio et al., Synchronization of endothelial Dll4-Notch dynamics switch blood vessels from branching to expansion. eLife 5, e12167 (2016).
15
S. Bernard, D. Gonze, B. Cajavec, H. Herzel, A. Kramer, Synchronization-induced rhythmicity of circadian oscillators in the suprachiasmatic nucleus. PLOS Comput. Biol. 3, e68 (2007).
16
D. Wilson, S. Faramarzi, J. Moehlis, M. R. Tinsley, K. Showalter, Synchronization of heterogeneous oscillator populations in response to weak and strong coupling. Chaos 28, 123114 (2018).
17
V. M. Christoffels, G. J. Smits, A. Kispert, A. F. M. Moorman, Development of the pacemaker tissues of the heart. Circ. Res. 106, 240–254 (2010).
18
T. Gregor, K. Fujimoto, N. Masaki, S. Sawai, The onset of collective behavior in social amoebae. Science 328, 1021–1025 (2010).
19
E. D. Herzog, Neurons and networks in daily rhythms. Nat. Rev. Neurosci. 8, 790–802 (2007).
20
C. W. J. Granger, Investigating causal relations by econometric models and cross-spectral methods. Econometrica 37, 424–438 (1969).
21
L. Barnett, A. B. Barrett, A. K. Seth, Granger causality and transfer entropy are equivalent for Gaussian variables. Phys. Rev. Lett. 103, 238701 (2009).
22
C. Zou, C. Ladroue, S. Guo, J. Feng, Identifying interactions in the time and frequency domains in local and global networks—A Granger causality approach. BMC Bioinformatics 11, 337 (2010).
23
D. C. Jacobs, R. E. Veitch, P. E. Chappell, Evaluation of immortalized AVPV- and arcuate-specific neuronal Kisspeptin cell lines to elucidate potential mechanisms of estrogen responsiveness and temporal gene expression in females. Endocrinology 157, 3410–3419 (2016).
24
A. W. Lohman, B. E. Isakson, Differentiating connexin hemichannels and pannexin channels in cellular ATP release. FEBS Lett. 588, 1379–1388 (2014).
25
F. H. Lin et al., Increasing fMRI sampling rate improves Granger causality estimates. PLoS One 9, e100319 (2014).
26
A. K. Seth, A. B. Barrett, L. Barnett, Granger causality analysis in neuroscience and neuroimaging. J. Neurosci. 35, 3293–3297 (2015).
27
B. Sun, J. Lembong, V. Normand, M. Rogers, H. A. Stone, Spatial-temporal dynamics of collective chemosensing. Proc. Natl. Acad. Sci. U.S.A. 109, 7753–7758 (2012).
28
R. A. Rossello, D. H. Kohn, Gap junction intercellular communication: A review of a potential platform to modulate craniofacial tissue engineering. J. Biomed. Mater. Res. B Appl. Biomater. 88, 509–518 (2009).
29
J. M. Kleinberg, Authoritative sources in a hyperlinked environment. J. Assoc. Comput. Mach. 46, 604–632 (1999).
30
E. Estrada, Quantifying network heterogeneity. Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 82, 066102 (2010).
31
E. Estrada, Degree heterogeneity of graphs and networks. II. Comparison with other indices. J. Interdiscip. Math. 22, 711 (2019).
32
A. Erez, T. A. Byrd, M. Vennettilli, A. Mugler, Cell-to-cell information at a feedback-induced bifurcation point. Phys. Rev. Lett. 125, 048103 (2020).
33
M. W. Gorman, E. O. Feigl, C. W. Buffington, Human plasma ATP concentration. Clin. Chem. 53, 318–325 (2007).
34
R. Cohen, S. Havlin, Complex Networks: Structure, Robustness and Function (Cambridge University Press, 2010).
35
A. L. Hodgkin, A. F. Huxley, A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 117, 500–544 (1952).
36
S. H. Strogatz, Nonlinear Dynamics and Chaos (CRC Press, 2018).
37
F. Schlögl, Chemical reaction models for non-equilibrium phase transitions. Z. Phys. 253, 147–161 (1972).
38
A. Erez, T. A. Byrd, R. M. Vogel, G. Altan-Bonnet, A. Mugler, Universality of biochemical feedback and its application to immune cells. Phys. Rev. E 99, 022422 (2019).
39
R. Fitzhugh, Impulses and physiological states in theoretical models of nerve membrane. Biophys. J. 1, 445–466 (1961).
40
J. Nagumo, S. Arimoto, S. Yoshizawa, An active pulse transmission line simulating nerve axon. Proc. IRE 50, 2061–2070 (1962).
41
Y. Tang, H. G. Othmer, Frequency encoding in excitable systems with applications to calcium oscillations. Proc. Natl. Acad. Sci. U.S.A. 92, 7869–7873 (1995).
42
D. T. Gillespie, Exact stochastic simulation of coupled chemical reactions. J. Phys. Chem. 81, 2340–2361 (1977).
43
A. Salameh, S. Dhein, Pharmacology of gap junctions. New pharmacological targets for treatment of arrhythmia, seizure and cancer? Biochim. Biophys. Acta 1719, 36–58 (2005).
44
R. A. de Melo Reis, H. R. Freitas, F. G. de Mello, Cell calcium imaging as a reliable method to study neuron–glial circuits. Front. Neurosci. 14, 569361 (2020).
45
R. Wollman, Robustness, accuracy, and cell state heterogeneity in biological systems. Curr. Opin. Syst. Biol. 8, 46–50 (2018).
46
C. W. Lynn, E. J. Cornblath, L. Papadopoulos, M. A. Bertolero, D. S. Bassett, Broken detailed balance and entropy production in the human brain. Proc. Natl. Acad. Sci. U.S.A. 118, e2109889118 (2021).
47
Y. Niv, Reinforcement learning in the brain. J. Math. Psychol. 53, 139–154 (2008).
48
C. E. Flores et al., Trafficking of gap junction channels at a vertebrate electrical synapse in vivo. Proc. Natl. Acad. Sci. U.S.A. 109, E573–E582 (2012).
49
A. B. Goryachev, Symmetry breaking as an interdisciplinary concept unifying cell and developmental biology. Cells 10, 86 (2021).
50
A. M. Moore, L. M. Coolen, D. T. Porter, R. L. Goodman, M. N. Lehman, Kndy cells revisited. Endocrinology 159, 3219–3234 (2018).
51
J. I. Alsous, J. Rozman, R. A. Marmion, A. Košmrlj, S. Y. Shvartsman, Clonal dominance in excitable cell networks. Nat. Phys. 17, 1391–1395 (2021).
52
M. F. Arnsdorf, Cardiac excitability, the electrophysiologic matrix and electrically induced ventricular arrhythmias: Order and reproducibility in seeming electrophysiologic chaos. J. Am. Coll. Cardiol. 17, 139–142 (1991).
53
S. F. Kazim et al., Neuronal network excitability in Alzheimer’s disease: The puzzle of similar versus divergent roles of amyloid β and tau. eNeuro 8, ENEURO.0418-20.2020 (2021).
54
G. Li, A. Starman, P. Chappell, B. Sun, video for KTaR-1 cells calcium signaling respond. Figshare. https://figshare.com/s/a681463e1e69134f7320. Deposited 18 January 2022.
55
R. LeFebre, A. Mugler, Mutlicellular-Information-Networks. GitHub. https://github.com/rwl23/mutlicellular_information_networks. Deposited 2 February 2022.

Information & Authors

Information

Published in

The cover image for PNAS Vol.119; No.37
Proceedings of the National Academy of Sciences
Vol. 119 | No. 37
September 13, 2022
PubMed: 36067282

Classifications

Data, Materials, and Software Availability

Experimental recordings of calcium dynamics are available at https://figshare.com/s/a681463e1e69134f7320 (54). Simulation codes are available at https://github.com/rwl23/mutlicellular_information_networks (55).

Submission history

Received: February 7, 2022
Accepted: July 18, 2022
Published online: September 6, 2022
Published in issue: September 13, 2022

Keywords

  1. calcium dynamics
  2. neuron dynamics
  3. intercellular communication
  4. temporal stimuli

Acknowledgments

G.L. and B.S. are supported by NSF grant PHY-1844627. R.L., A.M., G.L., A.S., and P.C. are supported by National Institute of General Medical Studies grant R01GM140466. B.S. is supported by National Institute of General Medical Sciences grant R35GM138179.

Notes

This article is a PNAS Direct Submission.

Authors

Affiliations

Department of Physics, Oregon State University, Corvallis, OR 97331
Department of Physics and Astronomy, University of Pittsburgh, Pittsburgh, PA 15260
Alia Starman
Department of Biomedical Sciences, Carlson College of Veterinary Medicine, Oregon State University, Corvallis, OR 97331
Department of Biomedical Sciences, Carlson College of Veterinary Medicine, Oregon State University, Corvallis, OR 97331
Andrew Mugler2 [email protected]
Department of Physics and Astronomy, University of Pittsburgh, Pittsburgh, PA 15260
Department of Physics, Oregon State University, Corvallis, OR 97331

Notes

2
To whom correspondence may be addressed. Email: [email protected], [email protected], or [email protected].
Author contributions: A.M. and B.S. initiated the project; G.L., A.S., and P.C. conducted the experiments; R.L. and A.M. developed the theoretical model; G.L., R.L., A.S., P.C., A.M., and B.S. analyzed data; and G.L., R.L., A.S., P.C., A.M., and B.S. wrote the paper.
1
G.L. and R.L. contributed equally to this work.

Competing Interests

The authors declare no competing interest.

Metrics & Citations

Metrics

Note: The article usage is presented with a three- to four-day delay and will update daily once available. Due to ths delay, usage data will not appear immediately following publication. Citation information is sourced from Crossref Cited-by service.


Citation statements




Altmetrics

Citations

Export the article citation data by selecting a format from the list below and clicking Export.

Cited by

    Loading...

    View Options

    View options

    PDF format

    Download this article as a PDF file

    DOWNLOAD PDF

    Login options

    Check if you have access through your login credentials or your institution to get full access on this article.

    Personal login Institutional Login

    Recommend to a librarian

    Recommend PNAS to a Librarian

    Purchase options

    Purchase this article to access the full text.

    Single Article Purchase

    Temporal signals drive the emergence of multicellular information networks
    Proceedings of the National Academy of Sciences
    • Vol. 119
    • No. 37

    Figures

    Tables

    Media

    Share

    Share

    Share article link

    Share on social media