Embryos assist morphogenesis of others through calcium and ATP signaling mechanisms in collective teratogen resistance

Based on the experimental data presented above, our hypothesis is that inter-embryonic signaling is a key mechanism in the collective robustness of morphogenesis. We thus simulated the development of embryos as singletons or in groups, similar to how embryos were grown in small or large groups. The results (Fig. 4) indicate that inter-embryo communication is essential for normal development in the presence of noise and that large groups are much more effective at achieving normal development than smaller groups, matching the embryo data in Fig. 1c. Overall, our agent-based simulation demonstrates that in morphogenetic systems in which emergent outcomes arise from local rules, collective dynamics can help stabilize against perturbation and noise (uncertainty) through local signaling mechanisms between agents.

Morphogenesis is the result of numerous biochemical, bioelectrical, and biomechanical signals, all of which feed into transcriptional changes that canalize physiological events into cell differentiation and other long-acting morphogenetic modules. We next asked to what extent the CEMA effect involves changes at the transcriptional level: does gene expression reflect the differential responses to challenges that we observe in small vs. large groups? We performed RNA sequencing (RNA-seq) on Xenopus samples derived from large (300-embryo) and small (100-embryo) groups of conspecifics that were either untreated or exposed to thioridazine. A randomly selected group of 15 embryos were collected from each group (N = 3 samples/group) either immediately after thioridazine treatment (developmental stage 25) or after brief recovery in MMR media (developmental stage 35). At stage 25, we observed weak separation across groups when batch-corrected data was assessed by principal component analysis, with particular overlap across the untreated group with a large number of conspecifics and both thioridazine groups (Supplementary Fig. 6). Comparison of samples between the large vs. small conspecific groups by differential expression analysis found that no genes passed FDR-correction for either the control naïve or thioridazine-treated embryos. We conclude that by stage 25, the mechanisms implementing CEMA have not had a significant impact on gene expression profiles.

However, at stage 35, transcriptional separation was observed between differently sized and treated groups (Supplementary Fig. 6). Differential expression analysis revealed 32 significantly changed genes (FDR < 0.05) between embryos housed in untreated control groups of 300 vs. 100 in control media, with 30 genes increased and 2 genes decreased in the larger (300-embryo) group compared to the smaller group (Fig. 5a and Supplementary Table 1). This reveals that even in the absence of an external morphogenetic stressor, gene expression is sensitive to the size of the developmental cohorts.Fig. 5Increasing cohort size induced changes to the transcription of different small sets of genes in control and thioridazine-treated embryos.Volcano plots of differential expression data for large (300 embryos) compared to small (100 embryos) groups of a control and b thioridazine-treated embryos. Genes undergoing significant changes in expression (FDR < 0.05) are highlighted in red (increased expression) and blue (decreased expression). Pathway analysis for comparison of large vs. small numbers of conspecifics in the c control and d thioridazine treatment conditions. Pathway analysis was performed by over-representation testing of the differentially expressed genes (FDR < 0.1) in the gene ontology (GO) terms/pathways and over-represented pathways were identified at a significance threshold of FDR < 0.1. The size of the dots reflects the gene ratios (number of significant genes associated with the GO term/total number of significant genes associated with any GO term), and the adjusted p-value (FDR) reflects the significance. Source data are provided as a Source Data file.

For embryos exposed to thioridazine, 19 unique genes were differentially expressed between embryos housed in groups of 300 vs. 100, with 3 genes increased and 16 genes decreased in the larger (300-embryo) group (Fig. 5b), revealing that the presence of large numbers of conspecifics modifies the transcriptional response to teratogen exposure. It must be emphasized that our analysis is not simply picking up conventional transcriptional responses to a teratogen: the data in Fig. 6 reflect differential responses to the same insult in large vs. small cohorts. Thus, distinct transcriptional signatures exist for both conditions: simply developing in a large group (with no exogenous stressors) induced changes in gene expression, as did specifically being in a large group responding to a teratogen.Fig. 6CEMA requires diffusion but not physical contact.a Images of embryos undergoing thioridazine treatment in a standard physical isolation device (left) and an image of a diffusion-enhanced physical isolation device (right). b Groups of n = 300 embryos were treated with thioridazine and were isolated using windowed wells, solid wells, or unseparated. ****p < 0.0001 and ***p = 0.0010. c Groups of n = 300 embryos were either put into clear or opaque separation devices or left unseparated and treated with thioridazine. **p = 0.0021 between -separated,+vision and +separated,+vision and **p = 0.0012 between -separated,+vision and +separated,-vision. b, c were both analyzed using one-way ANOVA and Tukey tests. Values are means ± SD. All treatments are normalized to the mean of non-treated controls’ spontaneous defects. Source data are provided as a Source Data file.

We computationally identified biological structures and pathways potentially impacted by these transcriptional changes through over-representation testing amongst highly expressed genes (FDR < 0.1). In the control group, pathways related to nucleosome activity, signaling receptor binding, translation, and ribosomal pathways were over-expressed (FDR < 0.1) in the 300-embryo group compared to the 100-embryo group, while those related to extracellular matrix, glycolytic, and mitochondrial processes decreased (Fig. 5c). In thioridazine-exposed embryos, the expression of genes associated with signaling pathways, including caveolae and nuclear steroid receptors, was increased in embryos housed in groups of 300 vs. 100 (Fig. 5d). In contrast, genes involved in multiple membrane transport pathways were decreased in the larger group, including proton transmembrane transport. Together, these findings suggest that embryonic transcriptional programs are sensitive to the number of individuals in a cohort. In addition, the regulation of a small, diverse set of genes involved in multiple modes of signaling, between stages 25 and 35, is a signature of group resistance to teratogens.

We next performed several functional blocking experiments to identify the modality responsible for the CEMA effect. First, to evaluate whether physical contact between embryos is required, 3D-printed separation devices were made. A clear, transparent plastic was used to create a platform with individual wells, measuring 1.80 mm tall and 1.80 mm wide, in which embryos could be placed (Fig. 6a). The entire separating device was placed inside of a dish thus allowing embryos to be physically separated, but still share media. A large (n = 300) group of embryos was separated into individual solid, clear plastic wells while treated in thioridazine. While the wells were separate, the medium could diffuse freely across the top of the wells. Separating the embryos dropped the survival rate to 1.3% compared to 96% in groups where embryos were not separated from each other (Fig. 6b) (N = 3, p < 0.001 using Welch’s t-test). To increase diffusion in the separation devices, new constructs were designed so that the walls of each well had windows that would not hinder diffusion. Animals treated with thioridazine in these new wells had an increased survival rate (63%) compared to the solid well-separated animals (1.3%) but were still lower than that for unseparated animals (96%) (Fig. 7b) (N = 3, p = 0.0005 using ANOVA). While testing for physical touch, optical cues were also tested through the use of the same separation device, but made from an opaque black plastic that prevented the sharing of optical cues. The results of blocking vision were not significantly different from blocking physical contact alone. Embryos left unseparated had an average survival of 93.7% while those treated in clear plastic had an average survival rate of 61.7% and 58.0% from those that were visually isolated (Fig. 6c).

Based on these data, we concluded that the most likely medium of inter-embryo influence is a short-range chemical signal.Fig. 7ATP/P2X receptors may mediate CEMA.a Suramin eliminates the protective effect of large group size against thioridazine-induced craniofacial defects and death. Tadpoles in groups of 300 were treated with 100 µM suramin, 90 µM thioridazine, 90 µM thioridazine and 100 µM suramin together, and compared to 90 µM thioridazine exposure as singletons. Blue arrows showing normal eyes, pigmentation, and head shape. Purple arrows show defects in those regions. b Percent defect of each group (n = 300 for each group) was then quantified. ****p < 0.0001 and ***p = 0.0002. c Thioridazine concentration was constant at 90 µM across all treatments and calcium or P2 receptor blockers were added. Calcium was depleted using 5 µM BAPTA. P2 receptors were inhibited with 100 µM PPADS. Each treatment group had a total of n = 300 embryos. ****p < 0.0001. d Quantification of extracellular ATP concentration. *p = 0.0443. Each dot represents the average of a group of three replicates. All animal images are dorsal views, and are oriented with the anterior end at the top of the image. For all data: values are plotted as mean ± SD, one-way ANOVA, and Tukey tests were conducted, and ***p ≤ 0.001, **p ≤ 0.01, *p ≤ 0.05. All treatments are normalized to the mean of non-treated controls’ spontaneous defects. Scale bars represent 2 mm. Source data are provided as a Source Data file.

To help identify the molecular nature of the diffusible signal, we focused on ATP, which is a known signaling molecule implicated in morphogenetic coordination in a range of model systems86–88, as well as mediating long-range contraction responses after injury51. It is also rapidly diffusible, which makes it a plausible target given all of the above data. Studies have shown calcium and ATP acting in a range of different communication pathways, including injury response and intercellular communication89–94. Calcium waves in particular have been shown to act as communicators for damage and regeneration93,95–97, suggesting these as potential candidates for inter-embryo signaling mechanisms for CEMA. We made use of a loss-of-function reagent known to block ATP signaling via the P2X/P2Y receptors: the purinergic receptor type 2 antagonist suramin98–101 to see whether it would prevent the CEMA effect and render large groups as susceptible to defects as small ones93,94,102.

Consistent with our results above, animals treated with thioridazine as singletons had a 96.2% incidence of defects; treating embryos in groups of 300 reduced this to 38.8%. Including suramin with thioridazine in a large group (300) treatments resulted in an incidence of 86.4%, not significantly different from the rate for embryos treated as singletons, suggesting suramin blocked the protective effect of the large group. Suramin by itself induced some defects in embryos treated in groups of 300, but only in an average of 7% of embryos (Fig. 7a, b) (N = 3). In an effort to identify a more precise molecular target, we focused on suramin’s two major targets: calcium and P2 receptors. Previous studies have shown that PPADS blocks a number of different metabotropic and ionotropic receptors, including the P2 receptors103,104, and BAPTA was used as a calcium chelator to investigate calcium’s role. When thioridazine was used in conjunction with BAPTA, the average survival rate of a group of 100 embryos was 80.3% which was not significantly different from thioridazine alone, 86.7%. Using PPADS with thioridazine also caused a drop of average survival to 20.7% (Fig. 7c). These data indicate that calcium stores and P2 receptors are involved in the CEMA effect. To further investigate whether ATP is involved, we measured extracellular ATP levels with the hypothesis that larger groups expel more ATP into the media that signal to other embryos. Media from treated and control animals were harvested, and an ATP determination kit was used to quantify concentrations in the two conditions. Ratios of ATP concentration of the large group over the small group revealed that the average extracellular ATP concentration in treated groups was 1.073 μM while control groups had an average of 0.989 μM (Fig. 7d). This increase further supports a potential ATP-based mechanism for CEMA.

We next sought to directly visualize the physiological dynamics that could underlie coordination between embryos using an acute mechanical injury assay. We reasoned that we may be able to observe such dynamics in the embryonic collective, and made use of a genetically encoded fluorescent calcium reporter, GCAMP6S. To directly visualize if Xenopus are sensitive to the damage of conspecifics, we first measured calcium signaling dynamics in GCAMP6S-expressing embryos exposed to mechanical injury as well as calcium activity in distant, uninjured embryos that share the same pool of 0.1x MMR media. We chose mechanical injury because it enabled a precise, known time point relative to which we could examine embryonic neighbors’ activity. Experiments were conducted in custom channels to ensure repeatable media volume and embryo spacing paradigms (Fig. 8a). Mechanical damage is well-known to induce calcium waves within a single organism that aid in cell-cell coordination and permit rapid wound closure in Xenopus embryos in response to superficial injury49,50.Fig. 8Mechanical damage in Xenopus induces calcium waves in distant, uninjured conspecifics.a Dimensions for holder and embryo arrangement for injury experiments. b Fluorescent images of GCAMP6S-expressing embryos at 0-, 10-, and 20-min post-injury for untreated stage 10–12 embryos. The injured embryo is on the left and the neighboring uninjured embryo is on the right. c The speeds of GCAMP6S signal propagation between and within embryos. Speeds were compared using the Kruskal–Wallis test, assessing differences in speed between embryos (n = 7) and within injured (n = 8) and receiver (n = 9) embryos; *p < 0.05. d Kymographs of receiver embryo calcium dynamics before injury (spontaneous activity) and post-injury under control, suramin, and ATP conditions.

Maximum GCAMP6S signal for e injured and f receiver embryos from the control (n = 10) and suramin-treated (n = 5) groups. Intensities are normalized to the mean of the pre-injury spontaneous activity. Groups were compared separately for injured and receiver embryos using the two-tailed Welch’s t-test; **p = 0.009. g Maximum GCAMP6S signal for embryos exposed to a bolus of MMR media (control; n = 3) or ATP (n = 4). Groups were compared using the two-tailed Welch’s t-test; *p = 0.01. All error bars represent the standard deviation. Source data are provided as a Source Data file.

Prior to injury, a low level of spontaneous flashes of calcium activity was detected. These spontaneous activity levels were quantified over 10 min of recording and used as a normalization factor for injury experiment quantification. Mechanical injury by a pulled glass needle induced a robust cell-to-cell calcium wave within the injured embryo (Supplementary Movie 1) and injury dynamics were recorded for 20 min post-injury. We observed that the injury wave was transferred to the second embryo sharing the same channel (Fig. 8b and Supplementary Movie 2), which was not itself manipulated in any way. Furthermore, we observed multiple injury wave transfers when an individual was injured within a group of 10 embryos (Supplementary Movie 3).

We measured cell-to-cell injury wave speeds of 6.7±3.48 μm/s in the injured embryo (N = 7), 5.28±1.89 μm/s between embryos (N = 8), and 2.36±1.66 μm/s in the uninjured embryo (N = 9). The wave speed in the uninjured embryo was slower than both the inter-embryo (p = 0.023) and injured embryo (p = 0.045) speeds (Fig. 8c). To more easily visualize calcium dynamics over time, the microscopy videos of the uninjured embryo were transformed into kymographs depicting the pre-injury spontaneous activity and post-injury activity (Fig. 8d) and intensity vs. time curves (Supplementary Fig. 7). For untreated receiver embryos, post-injury calcium activity was detected within 5–10 min following mechanical damage in the adjacent injured embryo. Calcium dynamics in the receiver embryo included a combination of sporadic bursting activity as well as lower, background levels of consistent calcium activity. Across replicates (N = 10), the median peak of calcium activity in the injured embryo was 1.77-fold higher than spontaneous activity levels (p < 0.0001) and 1.08-fold change over spontaneous activity in the uninjured embryo (p = 0.002). Furthermore, the injury wave crossed distances of >1 mm in the medium, suggesting that diffusible molecules may allow for long-range signaling and detection of tissue damage between conspecifics.

We then tested if inter-embryo injury wave propagation uses signaling mechanisms necessary for CEMA by testing the P2-purinoceptor antagonist suramin, which interferes with ATP signaling. Embryos treated with suramin (N = 5) were tested in the injury wave assay and compared to a compiled set of controls (N = 10). No treatment suppressed GCAMP6S intensity in injured embryos below that of controls (Fig. 8e). Calcium waves traveled from the epicenter of injury to the edge of the embryo as observed in controls. However, injury wave transfer to adjacent embryos was attenuated compared to untreated embryos in the case of suramin treatment (p = 0.009; Fig. 8e, f), suggesting that suramin suppresses only the propagation of the wave from the injured embryos to adjacent embryos. The mean peak of calcium activity was 1.08-fold change over spontaneous activity levels in the uninjured embryos in control media and 1.02-fold change over spontaneous activity in uninjured embryos treated with suramin (Fig. 8f). As suramin is defined as a P2 receptor antagonist, this finding suggests that injury signaling between conspecifics could involve ATP signaling. Indeed, when we apply a bolus of ATP, GCAMP6S is increased (p = 0.01) in uninjured embryos (2.11-fold change from spontaneous) over embryos exposed to a bolus of control media (1.05-fold change from spontaneous) (Fig. 8d, g and Supplementary Movie 4). However, the representation of the GCAMP6S signal in a kymograph (Fig. 8d) shows that the calcium dynamics are substantially different between the mechanical injury scenario and the application of a bolus of ATP. This difference in dynamics may occur because the mechanically injured embryo continuously releases diffusible molecules, including ATP, into the shared media pool, and the receiver embryo detects these molecules, resulting in prolonged calcium bursting behavior.

When a bolus of ATP is applied in the media, the receiver embryo detects all the molecules in bulk, resulting in a strong, but short-lived, signal. Together, these differences suggest that simply applying ATP to media may not mimic the more complex signaling dynamics provided by nearby conspecifics that could provide an additional layer of information apart from the molecule itself. Additionally, because suramin has been shown to exhibit multiple off-target effects beyond its defined P2-purinoceptor antagonism105, it is possible that additional molecules beyond ATP contribute to inter-embryo injury wave propagation and ongoing work in our lab aims to fully characterize the diffusible molecules involved in the injury wave effect as well as CEMA.

It is clear that although the hardware components of embryogenesis are encoded by the genome, phenotypic plasticity integrates inputs from the environment11,21,22,106. Here, we addressed the hypothesis that embryonic outcomes result not only from the activities of a single embryo’s developmental mechanisms but also from a kind of collective computation82,107–112, with stability properties arising from inter-embryo influence. Our results identify a source of lateral influence on development and morphology which resembles molecular- and cell-level collective repair14,113.

Increasing cohort size (while scaling the amount of media and drug to be constant per embryo) was sufficient to mitigate the effects of different teratogens, including thioridazine, nicotine, and forskolin (Figs. 1 and 2). As cohort size increased, survival increased and the incidence of defects decreased (Fig. 1c, e), an effect we term Cross-Embryo Morphogenetic Assistance (CEMA). The observation of this protective effect of cohort size against defects induced by thioridazine, forskolin, and nicotine57,114 (Fig. 2b, d), indicates that the effect is not specific to disruptors of a specific pathway. To further understand the boundaries of this assistive effect, we also looked at the effect of cohort size on defects known to be induced by overexpression of a dominant mutant of Kir6.1, a transmembrane potassium channel critical for bioelectric signaling during morphogenesis56. The protective effect of cohort size held for this disruptor as well (Fig. 2f). Finally, we also show that genetic diversity among cohorts does not affect the magnitude of CEMA and that CEMA is evident in perturbed rearing conditions, where larger treated cohorts have reduced baseline incidence of defects (Fig. 3). We conclude that individual Xenopus laevis embryos communicate with each other via instructive signals that protect against teratogenic insults.

These results, particularly the broad effect of CEMA across diverse disruptors and in control conditions, support the idea that robustness in morphogenesis is partly the result of lateral interactions between conspecifics. This has precedence in several related biological systems. Perhaps the earliest (pre-metazoan) involves quorum sensing and coordination in bacteria115,116. Others include clutch hatching synchronization in response to perceived predation in birds31,117, body color and brain changes in locusts depending on presence of others118,119, and growth pattern of plant roots in order to avoid other plants or forage120–122. The boundary between the body and the outside world shifts on both ontogenic (developmental) and phylogenic (evolutionary) timescales16. Thus, it is perhaps not surprising that a kind of multi-scale dynamic is in play in morphogenetic homeostasis: interactions between cells, known to be crucial for normal embryogenesis, have a parallel in the interaction between embryos in a cohort14.

We also gained insight into the mechanism of inter-embryo communication. Evidence of inter-individual communication has been noted in a variety of cases, including quorum sensing123, viral lysis-lysogeny124, and regeneration125. In an effort to understand the modality of communication in CEMA, special culture devices were used to physically separate embryos. Results indicated that CEMA required adjacency but not physical contact (Fig. 6), identifying diffusible chemical(s) as a likely inter-embryo signaling mechanism underlying CEMA. When embryos were physically distanced and prevented from touching, the average survival rate decreased significantly from non-separated, densely situated individuals. While separating embryos influenced survival, visual cues did not have any significant impact on the outcome.

Previous studies have found that following damage, ATP is released in high concentrations into extracellular space and can induce calcium waves126–129. In our study, we have shown that increasing the number of embryos in a group has a protective effect when faced with a stressor (chemical or physical) and that either depleting calcium or blocking P2 receptors can eliminate this shielding effect (Fig. 7).

The involvement of P2 receptors in craniofacial development is not surprising as purinergic receptors are found in nearly every tissue in adult animals130. In zebrafish, P2X3.1, a paralog of P2X3 mammalian receptors, is expressed in the embryonic head, specifically neural and ectodermal cells. Knockdown of P2X3.1 was shown to result in malformed lower jaws, malformed branchial arch regions, and smaller heads due to defective pharyngeal skeletal development131. It was also noted that the loss of signaling also disrupts normal neural crest behavior within the branchial arches. In relation to Xenopus laevis, P2 receptors have been detected in the central and peripheral nervous system132. The authors note that while each subunit has a distinct expression profile, there are several that are expressed in the head region. At stages 27 and 33/34, the p2rx1.L subunit is localized to a specific region of the hindbrain as well as in the head region, in the branchial arches, and in the cement gland. Expression of the p2rx1.S starts at the cleavage stage and then localizes to the eye field during neurulation till stage 41. The expression of p2rx2.L is mostly found in developing mesoderm derivatives and the nervous system. Expression can be seen in the head region, the brain, and sensory organs like the eyes at stage 41. Finally, a third member of the p2x receptor family is weakly found at the somites (stage 33/34), p2rx6.L. At stage 41, this receptor expression is strong in the head and dorsal regions. Temporally, the start of expression may be unique, but the level of expression for all of them increases during development and reaches a peak at stage 45.

When embryos are physically injured, they elicit a calcium wave that starts at the site of injury, propagates across the embryo, and spreads to a neighboring unharmed embryo (Fig. 8). It should be noted that, as always, the identification of the molecular components underlying the effect is just part of a full explanation. A simple concentration of ATP or any other molecule does not have the bandwidth to specify the information needed for complex morphogenesis as observed here. While it’s formally possible that CEMA involves a counter-teratogen effect rather than a pro-morphogenetic effect (e.g., a signal to degrade teratogenic drugs or mRNA, rather than a signal to improve morphogenesis), we think it unlikely based on the specificity of the effect (failure to cross-protect) and because evolutionarily, it would be much easier to develop mechanisms that support one particular target morphology vs. ones that try to anticipate the huge variety of possible teratogens. Future work will address the dynamic encoding of patterning cues via these molecular implementations, which could include spatial or temporal patterns (e.g., pulse-coding as has been observed with calcium133–136).

Consistent with previous studies, we found that an intra-embryo calcium wave is induced when mechanical damage occurs50,51. Furthermore, we observed for the first time that the calcium wave not only travels across the injured embryo but also into and across adjacent, uninjured embryos (Fig. 8 and Supplementary Movies 1–4). This occurs at a velocity in the same order as that observed for cell-to-cell signaling in developing Xenopus laevis embryos (5 μm/s) and injured Xenopus neuroepithelium (9 μm/s)93,137, but much slower than that observed for cardiac electrical conduction (2–5 m/s), neuronal action potentials (0.5–100 m/s), or cranial blood flow (24–100 cm/s)138,139. It is also slower than the long-range signaling observed between injured and uninjured limbs within 5 s of amputation in a model of hindlimb regeneration in Xenopus125. Further, targeting ATP signaling with suramin, an antagonist of purinergic receptors P2X and P2Y, abrogated both the inter-embryo calcium wave and the CEMA effect against morphogenesis disruptors (Figs. 7 and 8), suggesting a common ATP signaling mechanism. While these results do not directly implicate inter-embryo calcium waves in CEMA, they suggest that embryo damage can cause diffusible signaling molecules to be detected by nearby conspecifics. While a mechanical injury and chemical treatment both induce stress or injury onto the embryo, the methods and scale of doing so are quite different. The mechanical injury serves as an example of a short, discrete stressor, compared to a longer-term, continuous, teratogen-based stressor. The mechanical damage observed in the injury wave assay may induce different signaling responses than that induced by chemical injuries, like thioridazine and nicotine. By probing both modalities of trauma, we were able to examine whether there is a shared mechanism despite their differences.

Additional modes of injury should be investigated in the injury wave assay, including possibly chemical and laser-induced damage. The above-implicated mechanisms are not meant to be exclusive, and it’s entirely possible that other modalities are also involved.

Communication between cells and embryos that affects morphogenetic outcomes is likely to involve downstream changes in the expression of developmental genes. In an effort to identify unique transcriptome features between large and small conspecific groups, an RNAseq analysis (using rRNA depletion, to capture more unique features) was done on samples from the thioridazine experiment. Our analyses showed that changes are first detectable between stages 25 and 35, and involve a small number of genes that specifically respond to cohort size and not simply to the action of a teratogen. There were 16 genes that were down-regulated and 3 that were up-regulated (Supplementary Table 1). Many of these are currently of unknown function, but the list includes the V-ATPase, which is well-known to drive morphogenetic events140–146, and interesting new components linked to DNA damage and endopeptidase activity. Gene set enrichment analysis showed that genes associated with signaling pathways, including caveolae and nuclear steroid receptors, were increased in large groups, but genes involved in multiple membrane transport pathways were decreased in the larger group (Fig. 5). These changes also distinguish the context of a cohort responding to a challenge vs. simply to the size of a developmental clutch (Fig. 5).

Taken together, the data suggest the following signaling elements (Fig. 9). We propose that in normal conditions, the extracellular concentration of ATP is low. Injury induces an internal calcium wave and increased ATP release from the injured embryo. The fact that suramin or PPADS both block the calcium wave in neighboring intact embryos and, blocking P2 receptors prevents a calcium response, suggests that ATP released from injured embryos binds to P2 receptors on the surface of neighboring embryos, triggering the secondary calcium response in uninjured embryos (Fig. 9).Fig. 9Mechanistic model of CEMA effects.a Naive condition in which embryos are undisturbed. Damaged conditions in which embryos are subjected to injury and ATP expulsion increases. Damaged embryo + PPADS/Suramin showing the environment in which P-type receptors are blocked and inhibiting the binding of ATP. Damaged + BAPTA shows that internal calcium stores are depleted and there is a reduction of calcium response. b Enlarged view of internal events during injury. Top: The embryo on the left is injured and elicits a calcium response and increased ATP release. ATP binds to P2X receptors on neighboring, uninjured embryos and elicits its own calcium response. Middle: In the presence of PPADS/Suramin, the injured embryo has a calcium and ATP response but P2X receptor blockade prevents ATP binding and subsequent calcium wave in the uninjured embryo. Bottom: In the presence of BAPTA, the injured embryo no longer has a calcium response, but is still able to increase the release of ATP. ATP can still bind to P2X receptors in the neighboring embryo but there is no calcium wave, indicating that the calcium responses in both the injured and uninjured neighbor are dependent on calcium. Created with BioRender.com.

Much as in traditional developmental studies of single embryo morphogenesis, it is not enough to identify the molecular components necessary for CEMA. Previous work indicates that the robustness of normal craniofacial development is due to the computational capabilities of cell collectives, which allow them to reach appropriate target anatomy despite perturbations in the starting state or environment. To begin to identify system-level collective dynamics sufficient to produce group robustness across embryos, computational modeling can help interrogate the ways in which instructive lateral signals might protect embryos in large groups against disruptive signals. In addition, it can motivate intervention strategies. There are many examples of stability in groups, for example, coupled oscillators—from Huygens’ clocks to networks of networks147–150 or of robot swarms151; the “wisdom of crowds”152 is now a well-known phenomenon.

We modeled this effect using agent-based, cellular automata, with an emphasis on anatomical homeostasis—the ability of morphogenetic systems to resist noise. In this model, embryos that are affected proportionally seek necessary morphogenetic signaling information from others in the cohort. Numerous examples exist of embryos sensing and progressively correcting for powerful disruptors to early developmental processes, such as those that affect left-right asymmetry153, but to our knowledge this has not before been studied at the group level. This aspect of the model is a formalization of the idea that stress propagation is an important aspect of collective problem-solving because it leads to multiple subunits taking action to address the same unmet need. It is easy to imagine how this could have evolved: autocrine signaling can readily become paracrine signaling, at whatever scale of organization. Selection could progressively reward the leakiness of signals that represent the current error state as broadcasts that enable integrated activity to reduce the error.

Our model reinforced the conclusion that inter-embryonic signaling is a key factor in CEMA. A key component is that we modeled both embryonic development and CEMA using local interactions, not requiring a central controller mechanism. These types of local interactions are a hallmark of complex collectives that are found across biology and are hypothesized to support a number of functions154,155. In our simulations, local interactions were sufficient to solve the majority problem, overcoming noise and supporting ‘normal development’, with the assumption that the local interactions are further supported by the local cell’s neighbors (neighbor’s neighbors). We see that when these interactions are taken away, or reduced to a significant degree, the population becomes much less robust to perturbation (Supplementary Fig. 5). However, one open question in this study, as well as in the field of complex systems, is how many agents (i.e., embryos) do you need to become a ‘collective’ that functions differently than smaller groups? Here, the experimental and simulation data support a number around 300 for 100% survival in the presence of a teratogen but we do not have a clear picture as to why this specific number and if this number is unique to this species or can be generalized to other systems. In future studies, this “behavior vs. cohort size” phenomenon could be investigated if the relationship between them is similar. For example, in both the experiment and model we see a sharp, seemingly non-linear relationship between cohort size and survival.

While our results do reflect what is seen in Fig. 1c, they do not match perfectly. In our simulation groups of fewer than 25 conspecifics demonstrate small robustness against disruption by noise (<10% survival), whereas in the experimental data in Fig. 1c, groups of 5 and 10 embryos show no significant resistance to the effects of teratogen. Our model could be extended in several ways in future work. First, while there is embryonic death during the simulation, there is no ‘death-signaling’, wherein a dying embryo may emit a destructive signal that affects other embryos and could greatly affect development. Second, all ECAs are given a random initial configuration of 0’s and 1’s, with a majority (60%) set to 1. However, it may be that initial embryonic stages are not so random but have some amount of reproducible structure. Further experimental and computational studies are needed to address this issue. Third, while our inter-embryonic signaling is indeed local, it does not necessarily account for inter-embryonic traveling waves. If signals can propagate across a larger space over time in a relay fashion, this may change how the collective develops. We are currently developing a more in-depth model to address these issues in future work.