Controlling periodic longrange signalling to drive a morphogenetic transition
Abstract
Cells use signal relay to transmit information across tissue scales. However, the production of information carried by signal relay remains poorly characterised. To determine how the coding features of signal relay are generated, we used the classic system for longrange signalling: the periodic cAMP waves that drive Dictyostelium collective migration. Combining imaging and optogenetic perturbation of cell signalling states, we find that migration is triggered by an increase in wave frequency generated at the signalling centre. Wave frequency is regulated by cAMP wave circulation, which organises the longrange signal. To determine the mechanisms modulating wave circulation, we combined mathematical modelling, the general theory of excitable media, and mechanical perturbations to test competing models. Models in which cell density and spatial patterning modulate the wave frequency cannot explain the temporal evolution of signalling waves. Instead, our evidence leads to a model where wave circulation increases the ability for cells to relay the signal, causing further increase in the circulation rate. This positive feedback between cell state and signalling pattern regulates the longrange signal coding that drives morphogenesis.
Editor's evaluation
This fundamental work substantially advances our understanding of how multicellular structures transmit information over long ranges. Compelling approaches combining experiments and theory unravel the mechanism by which amoeba form migrating cellular waves by chemotaxis. The work will be of broad interest to cell and developmental biologists.
https://doi.org/10.7554/eLife.83796.sa0Introduction
Tissues require longrange signalling strategies to coordinate cell behaviour. As signal diffusion is only efficient over a few cell diameters, cells relay signals to rapidly disseminate information over tissue scales (Gelens et al., 2014). Signals are relayed by the activation of cells in response to active neighbours, which in turn activate their own neighbours (Dieterle et al., 2020). In this way, longrange signalling through tissues is propagated as waves of cell activation that travel outwardly from signalling centres (Deneke and Di Talia, 2018). Prominent examples of biological activation waves include ERK waves in the developing skin, breast, and colon and during wound closure (Aoki et al., 2017; Ender et al., 2022; Pond et al., 2022), periodic activation waves in the brain (Huang et al., 2010), chemoattractant waves in neutrophil swarms (Afonso et al., 2012), morphogen waves during gastrulation (Liu et al., 2021), waves of cell stress in bacterial biofilms (Chou et al., 2022), and cAMP waves in Dictyostelium aggregation (Tomchik and Devreotes, 1981). Waves can also be pathological, such as in cardiac fibrillation, where waves drive muscle contraction independently of the pacemaker (Pandit and Jalife, 2013). In these examples, the behaviour of the cells is organised by a selfsustaining signalling pattern; however, it is unclear how cell populations regulate these patterns.
One of the most intensively studied paradigms for signal relay is the Dictyostelium transition from a singlecell state to a multicellular state (Weijer, 2004). Upon starvation, single cells begin to communicate via cAMP (Tomchik and Devreotes, 1981; Gregor et al., 2010). Upon receipt of cAMP, cells secrete their own cAMP, which then activates the next set of cells, hence propagating the signal. As cAMP is a chemoattractant, cells move towards signalling centres, resulting in aggregation into the multicellular phase of their lifecycle. At the population level, the cAMP signal propagates in waves with circular and spiral patterns characteristic of other relay systems (Durston, 1973; Lee et al., 1996).
Despite decades of research into cAMP signalling in Dictyostelium, the mechanism of information production and control remains poorly characterised. Several lines of evidence suggest information is carried by cAMP wave gradient, frequency, and speed and decoded by cells in terms of changes in motility (Siegert and Weijer, 1989). One view is that cells interpret an increase in the frequency of the cAMP waves as a trigger to begin collective chemotaxis (Skoge et al., 2014). Alternative views are that other changing features of the waves, such as their speed and gradient, convey information instructing chemotaxis (Fisher et al., 1989; Wessels et al., 1992; Song et al., 2006; Nakajima et al., 2014). The current paradigm explaining the changes in wave behaviour over developmental time is that the excitability (the ability to relay the signal) of cell groups increases (individual cell excitability remains constant) at higher cell densities, which then manifests in the changing wave dynamics coupled to the onset of aggregation (Höfer et al., 1995; van Oss et al., 1996; Dallon and Othmer, 1997). One caveat with this view is that it requires some degree of aggregation to initiate the signal that drives the aggregation process. This paradigm was informed by extensive theory and studies on excitable chemical systems. In these systems, changes in wave dynamics emerge from rotating spiral waves whose rotation rate and curvature specifies the frequency and speed of activation waves. The curvature and rotation rate are specified by the circulation of the signal at the signalling centre (Tyson and Keener, 1988; Winfree, 2001). However, the major limitation in uniting this theory with experiment to decipher wave generation and interpretation is that studies have looked either at macroscopic (population) or microscopic (cell) levels but not both at the same time. As a result, it is not clear how cell signalling, movement, and spatial patterning at the signalling centre modulate the periodic cAMP waves that coordinate this morphogenetic transition.
To determine how populationlevel signalling dynamics are modulated by the signalling centre, we have leveraged a sensitive cAMP reporter to characterise collective cell signalling and motion during Dictyostelium aggregation, with singlecell resolution, over developmental time and length scales. We find that the onset of cell aggregation coincides with an increase in wave frequency and decrease in wave speed. Optogenetic activation of cAMP signalling with defined frequencies reveals how an increase in the frequency of wave production from the signalling centre is required for efficient aggregation. These changing wave dynamics result from the increasing circulation rate of cAMP waves at the signalling centre. We show that the circulation of cAMP waves spatially organises the signal as a spiral pattern, and the rotation rate and curvature of the spiral sets wave frequency and speed. We evaluate contrasting established and new hypotheses by which the circulation dynamics are regulated. We challenge the current paradigm, that increasing cell density drives an increase in collective cell excitability, and show that a model in which wave dynamics result from cell rearrangement at the signalling centre also cannot explain the observed signalling patterns. Instead, we show that the ability for individual cells to relay signalling waves increases as the wave frequency increases. This creates positive feedback that gradually increases the rate of wave circulation and the ability for cells to relay the signal. Together these observations imply the dynamic behaviour of signalling driving morphogenesis results from a mutual reinforcement between cell state and signal pattern.
Results
The onset of collective cell signalling
To track the signalling state of individual cells during morphogenesis, we imaged cells expressing the intensiometric cAMP signalling reporter Flamindo2 (Odaka et al., 2014; Hashimura et al., 2019) at an endogenous genetic locus. For simultaneous analysis of cell motion, the cells also expressed a nuclear marker (H2BmCherry) for cell tracking (Corrigan and Chubb, 2014). The fluorescence of Flamindo2 has an inverse relationship to the intracellular [cAMP]. This allows visualisation of the onset and progression of cAMP signalling wave patterns over cell populations during early development (Figure 1A and Figure 1—video 1). To monitor the motion and signalling of individual cells in their population context, we tracked the central position of the nuclei and measured the internal [cAMP] from the mean Flamindo2 intensity of the cell body (Figure 1B, Figure 1—figure supplement 1A and Figure 1—video 2). Analysing all cell tracks simultaneously provided the means to study the onset and adaptation of signalling behaviour and collective motion of single cells in their population context (Figure 1C–E).
Early during development, cells predominantly fired asynchronously. The proportion of spontaneously activating cells declined as the signalling behaviour became more correlated, with periodic waves of signalling including an increasing proportion of the cells (Figure 1D). Prior to the onset of cAMP waves (around 4 hr after starvation), cells produced cAMP pulses spontaneously at a mean rate of around 1 pulse per 12 min, with a minimum of around 1 pulse per 4 min (Figure 1E). During this time, the population distribution of interactivation intervals follows a gamma distribution (Figure 1E and Figure 1—figure supplement 1B). We interpret this fit as there being a minimum waiting time between pulses (a refractory phase), with the exponential tail of the distribution suggesting the probability of firing is independent of the previous firing event.
The initial cAMP waves (around 4–5 hr) are broad bands of broken waves amalgamating from numerous circular patterns that expand from spontaneously activating cells (Figure 1A). Activation events prior to synchronisation are spontaneous, spatially variable, and mostly do not initiate a cAMP wave. Together with the independence of successive firing events, these observations are inconsistent with the model that individual ‘pacemaker’ cells act as the organising source of global cAMP waves. However, it may be possible that cell clusters (arising for example from random motion) may collectively act as a pacemaker (Mirollo and Strogatz, 1990; Gregor et al., 2010). For example, if the cell cluster is sufficiently dense such that at least one cell in the cluster is likely to activate just after the refractory time. Nevertheless, any given cell will experience circular waves from several different sources with a variable periodicity (5–9 min).
Circular waves produce cAMP signals at irregular intervals which are typically greater than the interpules interval (Figure 1D). Consequently, not all cells will relay a circular wave because cells that fire synchronously with one wave can potentially reactivate and so be refractory for the next wave. This is highlighted by a multipeaked population distribution for interactivation intervals where the smallest peak is associated with spontaneous cAMP pulses, while the other peaks are associated with synchronised firing during a wave (Figure 1E and Figure 1—figure supplement 1B). Rotating spiral cAMP waves form from broken waves around 5–6 hr following starvation (Figure 1A). To interrogate how waves break and give rise to spiral waves, we engineered cells to express bluephotoactivable adenylyl cyclase (bPAC; Stierl et al., 2011) to control the time and place of cAMP waves by exposing a subset of cells to blue or green light (Figure 1—figure supplement 2 and Figure 1—video 3). In the absence of light pulses, we observed a steady train of periodic waves (Figure 1—figure supplement 2B). Photoactivation of a small patch of cells generated a circular wave which fused with oncoming waves (Figure 1—figure supplement 2A first and second row and Figure 1—figure supplement 2C). Wave breaks (Figure 1—figure supplement 2A forth – sixth row) occurred only when the oncoming wave collided with an asymmetric wave generated by the activated patch (Figure 1—figure supplement 2A third row). We interpret the asymmetry arising when cells activate at the border of two regions; one where cells are in a refractory state (preventing wave propagation) and the other in which cells are recovered (allowing wave propagation; Figure 1—figure supplement 2D). This experiment suggests that spiral waves might naturally arise from spontaneous cell activation behind a circular wave, where neighbouring cells are transitioning from a refractory to rest state.
Rotating spiral patterns generate waves that are more frequent and regular than circular patterns (Figure 1 56 hr onwards) and hence become the dominant pattern over time, overriding the less frequent circular waves, as reported in previous studies (Lee et al., 1996). The mean wave frequency steadily increased towards the maximum firing frequency of cells, around once every 3–4 min (Figure 1E). Consequently, the proportion of spontaneously activating cells tends to zero as there is no time to activate between waves (Figure 1D), and the population distribution of interpulse intervals tends to a normal distribution (Figure 1E). To understand the basis of the heterogeneity in interpulse intervals, we quantified the singlecell differences in activation times relative to the wave centre (Figure 1—figure supplement 3A, B). Aided by Markov modelling, we observed a shortterm bias in firing tendencies: for example, a cell that fires ahead of its neighbours has a 60% chance of repeating this early firing, 35% chance of switching to late firing, and a 5% chance of skipping the subsequent wave (Figure 1—figure supplement 3C, D). Whether a cell fires before or after the wave does not affect its probability of skipping a wave. Overall, these data suggest no intrinsic difference between cells beyond a slight history dependence.
The onset of collective migration
To relate how information encoded by spiral waves coordinates collective motion, we compared changes in cAMP wave dynamics to cell chemotaxis (Figure 2 and Figure 1—video 2). Chemotaxis is directed by the external cAMP concentration field, which is determined by the release of internal cAMP (reported by Flamindo2) from individual cells. The spatial organisation and motility of cells changed rapidly, occurring around 7 hr following starvation, coincident with the dominance of spiral waves (Figure 2A and Figure 2—figure supplement 1). During this process, cAMP waves steadily became more frequent, slower and thinner, and cells became faster, more directional, and persistent in their movement (Figure 2B and C). The profile of the mean cAMP pulse in single cells remained constant with each wave. Since the duration of a cAMP pulse ($\tau $) remains constant, and the wave speed ($v$) decreases, then the wave width ($w=v\times \tau $) also decreases. A thinner wave generates a steeper external cAMP gradient, which may enhance chemotaxis (Song et al., 2006).
To quantify collective cell motion, we used the following two measures: (i) the mean chemotaxis index, the proportion of cell movement in the direction of the oncoming wave (Skoge et al., 2014; Figure 2B and D), and (ii) the local correlation of cell velocity (De Palo et al., 2017; Figure 2E). Cells consistently moved towards the wave source at wave fronts but remained spatially uniform before 6 hr of starvation due to random cell motion at the wave rear (Figure 2B bottom left panel). This is highlighted by a peak in the local correlation of cell motion (Figure 2E) and the mean chemotaxis index (Figure 2D bottom panel) at the wave front, both of which then decay towards zero (indicating random cell motion). For infrequent, fast and wide waves (i.e. circular patterns), the chemotaxis index and local correlation of cell motion reached their baseline values at the following wave front, indicating cells move randomly. This random cell motion disrupts the order in collective cell motion initially induced by the front of the wave. For frequent, slow and thin waves (i.e. rotating spiral waves), the chemotaxis index and velocity correlation significantly change as the cAMP wave frequency increases between 6 and 8 hr following starvation. From 6.5 hr onwards, when the wave period is less than 4.5 min, the measures of collective cell motion do not relax back to their previous values by the next wave and hence gradually increase over time to constant, elevated values across the entire wave (Figure 2D and E). In parallel, the cell speed at the front of the wave increases incrementally with each wave until it peaks at around 7 µm min^{−1} (Figure 2D second panel). These data are consistent with models where decreasing the wave period causes cells to gradually become more persistent towards the direction of the wave source and hence gain consistent collective motion (Skoge et al., 2014).
To test the relationship between wave frequency and collective migration in the full population context, we photoactivated bPAC in a small patch of cells at different intervals and monitored the resultant activation waves (Figure 3 and Figure 3—video 1). Consistent with our observations shown in Figure 2, cells collectively migrated towards the signalling centre only when the induced wave period was less than 4.5 min (Figure 3A). However, photoactivation of cells with a period less than 3.5 min did not result in aggregation. We infer that this is because cells cannot relay a steady periodic signal close to their maximum firing rate initially, resulting in a missed wave every second photoactivation event (Figure 3C), leading to the signalling centre becoming ‘taken over’ by another set of period waves. Motivated by the observed change in wave frequency during aggregation, we sought to test if it were possible to induce large cell aggregates by photoactivating cells at a faster rate than what we observe experimentally by reducing the wave period gradually from 5 min to 3.5 min from 5 to 6 hr following starvation, rather than 6–8 hr (Figure 3C second panel). This was indeed achievable, resulting in the accumulation of a significantly larger number of cells around the signalling centre and conglomerates of several aggregates (Figure 3B). These artificial signalling centres induce large aggregates by entraining the territories of neighbouring signalling centres. Our data support the notion that the increase in wave frequency is a crucial determinant of the transition to collective motion. We now turn to the question of how cells generate this increasing wave frequency.
Spiral wave progression
To characterise the regulation of spiral wave dynamics, we imaged and analysed the motion and signalling of individual cells around the spiral core following the onset of spiral formation. The spiral core is the region about which the spiral tip (the end of the activation wave) moves. Based on the general theory of excitable media, the spiral tip will sustain itself and the entire spiral wave by virtue of circular motion (Tyson and Keener, 1988; Winfree, 2001). In other words, the circulation of the spiral tip is the organising source of the spiral wave. To better understand the dynamics of spiral tip circulation, and how it relates to the changes in wave dynamics, we tracked the position of the spiral tip and the shape of the spiral wave from spiral formation to aggregation (Figure 4, Figure 4—figure supplement 1, Figure 4—videos 1–3). To track tip motion, we mapped the [cAMP] of individual cells to a phase variable measuring the sequence of changes in the internal [cAMP] during a cAMP pulse (Figure 4B, Figure 4—figure supplement 1B, C). During a cAMP pulse, this variable increases from 0 (baseline [cAMP]) to ±π (maximum [cAMP]) and then back to 0. This approach reveals the spiral singularity — or topological defect — which is the point in space around which cells are distributed across the full spectrum of phase states. At each time point, we estimated the spiral tip location as the mean position of maximally activated cells closest to the singularity (Figure 4B). Additionally, we fitted a spiral curve to the spatial distribution of maximally activated cells (Figure 4A and C first panel).
Spiral waves significantly change in rotation rate and curvature during aggregation, which we refer to as spiral wave progression (Figure 4C and Figure 4—videos 1–3). These changes in spiral wave structure and dynamics are coupled to changes in spiral tip circulation (Figure 4B–D). Each rotation of the spiral tip generates a wave: the rate at which the activation waves circulate around the signalling centre sets the frequency of global activation waves disseminated from the signalling centre. The spiral tip rotational period decreased from 5.5 min to 3.5 min. These changes in the rotational period of the spiral tip translate to the changes in the rotation rate of the spiral wave: the decrease in the spiral tip rotation period (Figure 4D first panel) is consistent with the decrease in the wave period away from the spiral core that coincides with the onset of collective migration (Figures 2C and 3). The acceleration of spiral tip circulation is associated with spiral core contraction; the radius of circulation steadily reduced from around 100 µm to 20 µm (Figure 4D second panel). Meanwhile, during this contraction, the spiral tip speed decreased from 100 to 120 to 40–50 µm min^{−1} (Figure 4D third panel).
How do cells at the signalling centre drive spiral wave progression? To measure the fundamental signalling parameters and how they change over time, we utilised the general theory of excitable media (Tyson and Keener, 1988). Dictyostelium cAMP waves can be shown to satisfy this theory (Tyson et al., 1989; Foerster et al., 1990). Our highresolution datasets now make it possible to apply this theory to understand how spiral waves change during aggregation, independently of the specific equations that model cAMP release, which are still disputed (Sgro et al., 2015; Kamino et al., 2017). According to this theory, the shape of the spiral wave must satisfy a curve that is approximately Archimedean (Materials and methods, Equation 2). The spiral curve depends on a single parameter $k={N}_{0}/\epsilon D$ (referred to here as the curve parameter) which relates the signalling parameters ${N}_{0}$ (the planar wave speed), $D$ (the signal diffusion constant), and $\epsilon $ (the ratio of the rates of cell inactivation and activation, which we call the cell activation parameter; Materials and methods, Equation 3). Crucially, the inverse of $\epsilon $ is the measure of cell excitability and specifies the rate at which cells activate. By fitting the spiral curve to the spatial distribution of maximally activated cells at each time point (Figure 4A; Figure 5A, Figure 4—videos 1–3), throughout spiral wave progression, we obtain a measurement of how the cell excitability and planar wave speed change over time (Figure 5B second column). Specifically, both measures of the cell signalling phenotype can be expressed in terms of spiral wave properties that we can measure (Materials and methods, Equation 4), such as the normal wave speed and wave curvature (Figure 5B first column).
This analysis reveals that the activation parameter and planar wave speed both decrease during spiral wave progression (Figure 5B). Firstly, we find that the planar wave speed halves in value as the wave period decreases towards its minimum (Figure 5D). The theory of excitable media provides the following interpretation: as the wave period decreases, cells have less time to reset between waves, causing slower cell activation and wave speed (the dispersion relation), which gives rise to the deceleration of wave speed during core contraction. This deceleration quantitatively matches the changing wave dynamics away from the spiral core (Figures 2C and 5B), meaning the circulation of the spiral tip determines the wave speed away from the core, in addition to the frequency. Secondly, the theory interprets the increase in spiral wave rotation rate and curvature to be a consequence of cells becoming more excitable over time (Figure 5C). In other words, the spiral waves progress because cells can activate faster. Indeed, although there was temporal variation between experimental replicates, the activation parameter collapses onto a highly consistent curve when plotted against the wave period (Figure 5C), showing how the wave dynamics are strongly coupled to the excitability of the cells. A simple quadratic provides a good fit regardless of derivation used for quantifying the activation parameter (Materials and methods, Equations 4 and 5). In summary, this analysis implies the change in cell excitability is a key element driving spiral wave progression. We now evaluate different hypotheses for mechanisms driving spiral wave progression.
Hypotheses
Hypothesis 1: increased collective excitability by cell reorganisation with constant individual excitability
The existing paradigm, based on mathematical modelling, is that spiral wave progression is driven by cell aggregation (Höfer et al., 1995; van Oss et al., 1996). In this view, an increase in cell density elevates cAMP production and degradation rates per unit volume, increasing the collective excitability around the aggregate. In this way, the spiral wave rotation rate and curvature increase as cells accumulate at the spiral core, even though the intrinsic excitability of individual cells remains constant. Consistent with this paradigm, we measured a substantial decrease in the cell activation parameter (inverse of cell excitability) and increase in cell density during spiral wave progression (Figure 6). The relationship between density and the activation parameter is incredibly sensitive, meaning slight fluctuations in density (such as occurring through random motility) would correspond to dramatic changes in the activation parameter (Figure 6C first panel). However, unlike the strong relationship between wave period and cell excitability (Figures 5C and 6C third panel), a causative relationship between density and excitability is not supported by reliable fits (Figure 6C second panel). In addition, our data show that spiral wave progression occurs prior to clear changes in cell density (Figure 6A). Specifically, spiral core contraction and the increase in rotation rate are immediate and steady for around 1.5 hr following spiral wave formation (Figure 6B second panel), while the cell density sharply increases only after around 1 hr (Figure 6B first panel) once the wave period is below 4.5 min. Indeed, the rotational period is close to the minimum value before substantial changes in density occur (Figure 6B third panel). Altogether, these observations imply that an increase in cell density is not sufficient to account for spiral wave progression, at least at the onset of collective migration.
A more recent model for cAMP release (Sgro et al., 2015), which is based on the FitzHughNagumo model rather than the MartielGoldbeter model (Martiel and Goldbeter, 1987), shows improved quantitative behaviour of cAMP release from single cells (Huyan et al., 2021). It is conceivable that incorporating this more recent model into the earlier models of spiral waves (Höfer et al., 1995; van Oss et al., 1996; Dallon and Othmer, 1997) could provide spiral wave progression consistent with our data  an increase in the spiral wave rotation rate prior to cell aggregation. To test this, we developed a hybrid partial differential equation (PDE) agentbased model (ABM; see Materials and methods and Supplementary Information) in the same vein as van Oss et al., 1996, but replacing the MartielGoldbeter model with the modified FitzHughNagumo model. To explore how cellular reorganisation influences the signalling pattern, we fixed all signalling parameters such that the intrinsic cell excitability does not change. This model reproduced a rotating spiral wave and experimentally observed spatial patterning of cells (Figure 7, Figure 7—videos 1 and 2). In the model, the rotating spiral wave creates a limit cycle in cell motion (see the circular streamline in Figure 6A and Figure 6—figure supplement 1). Cells collect along the limit cycle via chemotaxis from both inside and outside the path of tip circulation, forming a ring pattern. The emergence of this limit cycle is matched in our data (Figures 6A and 7, Figure 6—videos 1 and 2). In data and model, both the radius of spiral tip circulation and the ring pattern of cells contract over time. In our model, this contraction occurred in the absence of cellcell physical interactions (i.e. no crowding effects) but not in the absence of cell motion (Figure 7D and E first panel). These tests indicate that spiral core contraction can be attributed to the spatial reorganisation of cells via chemotaxis towards the circulating wave, but not necessarily due to crowding effects from cells funnelling towards the spiral core. This effect is also observed in data (Figure 6A): the limit cycle is ‘blunted’ at the location of the spiral tip, suggesting that cell chemotaxis is deflected around the external [cAMP] field of the spiral tip, inwards from the edge of the cell ring. Put simply, the tip is deflecting cells inwards.
Since our data point to the contraction of the radius of spiral tip circulation as the driver of spiral wave progression (Figure 4B–D), our model suggests the following scenario for the progression of the spiral wave: as the external [cAMP] field of the spiral tip rotates about the spiral core, it draws cells from the cellular ring inwards, causing a decrease in the circulation radius hence an increase in the rotation rate. In this view, the spiral progression associated with the onset of collective chemotaxis is driven by spiral core contraction. This contrasts the original formulation where cell reorganisation increases the cell density, which increases the collective excitability.
This contraction scenario sufficiently explains spiral core contraction but not spiral wave progression. Specifically, the spiral wave does not increase in rotation rate during spiral core contraction (Figure 7E second and third panel). The speed of spiral tip motion decreases too much to offset the contraction in the radius of circulation (Figure 7E fourth panel). This is a flaw in the model since the increase in wave frequency is the key element of change coinciding with the onset of collective cell migration. By adding an additional component to the model, which allows the excitability of individual cells to increase over time, this flaw can be bypassed, with the model outputting a strong increase in wave frequency (Figure 7—figure supplement 1 and Figure 7—video 3).
Hypothesis 2: increased individual cell excitability driven by wave dynamics
An increase in cell density is not sufficient to explain spiral wave progression. This suggests that spiral wave progression does not arise because of the increase in collective excitability resulting from aggregation, but cells become more excitable individually. How would this increase in individual cell excitability occur? Based on the relationship in Figure 5C, we propose that individual cell excitability (i.e. the rate of cell activation) increases as the wave period decreases. That is, cell excitability is coupled to the wave dynamics. This hypothesis states that spiral wave progression is due to positive feedback between cell excitability and the activation wave period: in the context of the spiral wave, which generates periodic waves, the excitability of cells increases the rotation rate of the spiral wave, thereby decreasing the wave period, which raises cell excitability and so forth.
Due to the conservation of the rotational period at the signalling centre and wave period away from the signalling centre, one consequence of Hypothesis 2 is that cell state should be similar across the field. To test this property, we used the quadratic relationship between the cell activation parameter and wave period (Figure 5C) to rederive the relationship between the planar wave speed and wave period (the dispersion relation) in a different context – away from the spiral core. To measure the wave dynamics and structure away from the spiral core (Figure 6B), we exploited that the spiral curve approaches a circle distant from the core (Figure 6A) and fitted an equation for a cone to the time evolution of the spatial distribution of maximally activated cells (Materials and methods). Using this approach, we find a consistent dispersion relation both at and away from the spiral core (Figure 5D). In other words, the changes in cell state due to wave dynamics are the same across the population. This match was critically dependent on the relationship between cell excitability and wave period that is stated in Figure 5C.
To test whether the coupling of excitability to the wave dynamics sufficiently describes spiral wave progression, we modelled cAMP spiral waves at constant cell density, with the intrinsic excitability of cells changing in response to external cAMP (Figure 8 and Figure 8—video 1). The earlier models analysed by Höfer et al., 1995 and van Oss et al., 1996 were based on the spatial continuum extension of the MartielGoldbeter model analysed by Tyson et al., 1989. The key extension to the Tyson model by Hofer/van Oss is that the cell density can change. By using the Tyson model, we return to the original formulation, but instead of allowing density to change, we couple cell excitability to the cAMP wave dynamics. We set the cell activation variable to be a simple timedependent function where the variable tends to a baseline value and decreases in response to external [cAMP] (Materials and methods). Figure 8A shows the time evolution of the simulated spiral wave. This solution reproduces features of spiral wave progression such as the rate of spiral core contraction (Figure 8B), increase in spiral wave curvature, increase in the spiral rotation rate, and the coupling between the cell activation parameter and the wave period (Figure 8C). The similarity between our experimental data and this model (even with parameter values rounded to the nearest order of magnitude) suggests that the positive feedback between cell excitability and the wave period is a sufficient mechanism driving spiral wave progression (Figure 8—figure supplement 1).
Lastly, the general theory allows us to predict the radius of spiral tip circulation from the measured cell excitability (the curvature relation, Materials and methods Equation 5). Is the change in cell excitability driving the change in circulation radius (Hypothesis 2) or is the radius driving the change in excitability (Hypothesis 2)? If Hypothesis 1 was true, then a change in the measured radius of circulation should precede the change in excitability and hence the predicted radius. In Hypothesis 2, the circulation radius adapts to changes in cell excitability, in line with previous experiments using caffeine, which inhibits cAMP relay (i.e. reduces excitability), caused a larger radius of circulation (Siegert and Weijer, 1989). Further supporting Hypothesis 2, our data show the change in predicted radius precedes the change in the measured radius, with convergence at the end of spiral wave progression (Figure 5—figure supplement 1). That is, the radius of spiral tip circulation is a consequence, not driver, of spiral wave progression.
Exploring the properties of the radius contraction in the adapted Tyson model shows that the inclusion of a hole in the cell density field does not hinder spiral wave progression. The spiral wave becomes pinned to the hole, which sets the minimum radius of spiral tip circulation (Figure 8B and Figure 8—video 2). The wave period reduces (Figure 8C third panel) because the speed of the spiral tip speed remains elevated (Figure 8C second panel) by virtue of an increase in excitability (Figure 8C fourth and fifth panel), compensating for a large radius of circulation (Figure 8C first panel). This simulated result contrasts Hypothesis 1, suggesting that preventing core contraction does not impede spiral wave progression.
Hypothesis testing
We have investigated two hypotheses for spiral wave progression. The prevailing paradigm (Hypothesis 1) is that excitability increases due to an increasing cell density; however, this model is not supported by our data comparing density with wave period and excitability (Figure 6). Our reformulation of this hypothesis predicted that cell rearrangements drive the increase in the spiral wave rotation rate via contraction of spiral tip circulation (Figure 7). In contrast, Hypothesis 2 proposes that the wave dynamics change cell excitability (Figure 8). The key difference between these hypotheses is that in Hypothesis 1, spiral tip contraction is the driver of spiral progression (Figure 7F), whereas in Hypothesis 2, contraction is a consequence (Figure 8—figure supplement 1). Exploring the properties of the models underlying these hypotheses provides a clear experimental test: to constrain the radius of spiral tip circulation.
To experimentally test the competing hypotheses, we imaged Dictyostelium aggregation in the presence of physical blocks (low melt agarose particles with radii ranging between 25 and 150 µm; Figure 9 and Figure 9—video 1). Spiral waves can become pinned to the blocks (Figure 9A), with their minimum radius of circulation set by the block radius (Figure 9B). Spiral waves that are pinned to a block cannot fully contract (Figure 9C first panel) but nevertheless undergo a steady increase in rotation rate and hence decrease in wave rotational period (Figure 9C third panel). The circulation rate increases because the spiral tip speed does not decrease to the same extent as the unconstrained spirals (Figure 9C second panel). This observation is in agreement with the model simulations shown in Figure 8, revealing the same effects as in the data – wave speed increased to compensate for larger radius of circulation, generating a reduced rotational period. This observation also reflects our simulations predicting that spiral waves can progress in the presence of physical blocks due to a coupling between the cell excitability and wave dynamics. In addition, this coupling between wave frequency and excitability is also apparent in our photoactivation experiments, in which a gradually increasing pulse frequency enhanced aggregate formation (Figure 3). Together, these observations are consistent with Hypothesis 2, but not Hypothesis 1: the contraction of the spiral tip circulation is a consequence, not a driver, of spiral wave progression.
Overall, our approach provides a clear interpretation of our data that spiral wave progression results from positive feedback between cell excitability and wave frequency, arising only in the context of the spiral wave pattern (Figure 8—figure supplement 1). That is, from the initial symmetry breaking event (Figure 1—figure supplement 2), the formation of a rotating spiral wave produces stable periodic activation waves that increase cell excitability, then increases the spiral rotation rate and hence wave frequency, which then feedback onto further increases in cell excitability.
Destabilisation of the spiral wave
To conclude, we now describe the cAMP wave dynamics following aggregation. As development proceeds, more and more cells move into the aggregate. Previous studies show the activation waves are organised as a multiarmed spiral in aggregates (Siegert and Weijer, 1995). How does this continued cell accumulation affect the cAMP dynamics at the signalling core and, consequently, the encoding of information by cAMP waves?
To characterise the single to multiarmed spiral transition, we monitored the structure and dynamics of the activation waves following ring collapse (Figure 10, Figure 10—figure supplement 1 and Figure 10—video 1). Singlearmed spirals destabilise only once the cellular ring fully contracts. Contraction results in a phase singularity (convergence of all signalling phases) that is fixed at the aggregate centre (Figure 10—figure supplement 1A). The phase singularity and spiral wave are unstable and disintegrate over time. This is likely because the cells are not infinitely small and can only be in one signalling phase at a time. In other words, the cells at the centre cannot activate at a rate necessary to support the full spectrum of signalling phases. The destabilisation process begins by cells signalling outofphase with the spiral wave, resulting in the emergence of transient patches of cell activation that gradually increase in number, size, and lifetime (Figure 10, around 2 hr following spiral formation). These activation patches rotate about the aggregate while fusing with, and breaking apart from, other activation patches and the original spiral wave, which progressively disintegrates (Figures 10A and 2–3 hr following spiral formation). This process concludes with the appearance of three stable, slowmoving and uncurved activation waves that rotate about the same fixed point at the aggregate centre (Figure 10A, from 3 hr following spiral formation). At this point, the signalling regime can be classed as a multiarmed spiral. Although the ring has collapsed, the centre of the aggregate contains relatively few cells (Figure 10—figure supplement 2). At this point, there is no phase singularity: the few cells that occupy this central region appear to have signalling phases that are uncorrelated with each other and the surrounding multiarmed spiral wave (Figure 10—figure supplement 1 – last panel). This central region of uncorrelated signalling phases is characteristic of theoretical chimera spiral waves (Totz et al., 2018).
To analyse the transition from the single to multiarmed spiral, we generated a set of kymographs of cell activation states across the aggregate: one that is crosssectional and the others that are circular (Figure 10B). The crosssectional kymograph shows a disordered signalling regime between 2 and 3 hr following spiral formation, from around 45 min after cellular ring closure. In the circular kymographs, both the single and multiarmed spiral waves appear as diagonal lines across space and time (<2 and >3 hr), meaning a constant rotation rate. These data show that the singlearmed spiral (<2 hr) has a rotational period of 3.4±0.2 min (the minimum cell firing period, see Figures 1, 3 and 4), which is three times faster than the threearmed spiral whose rotational period is 10.1±0.7 min and wave period is 3.3±0.5 min. This means that the overall frequency with which a cell experiences a wave is the same because the threearmed spiral provides three waves in its longer rotation period. As such, the frequency of cAMP waves disseminated from the signalling centre, along the streams, remains the same for both signalling patterns. That is, both patterns encode the same temporal information.
Slow and fastwaves coexist at the outside of the aggregate midway through the transition to the multiarm spiral. This is visualised as a mixture of lines of two different gradients in the kymograph (Figure 10B, fourth panel, 2–3 hr). This observation suggests that both signalling regimes coexist by the mixing of cells that relay either the single or multiarmed spiral wave. The multiarmed spiral wave eventually becomes dominant as the proportion of cells that relay the singlearmed spiral wave declines. Altogether, our data suggest that the signalling pattern adapts to the change in aggregate geometry while maintaining an active phase wave that generates waves at the maximum frequency.
Discussion
To determine the properties of signal coding during signal relay, we have carried out a characterisation of collective signalling and cell motion during Dictyostelium morphogenesis. Our analysis shows how cell aggregation is coordinated by the dynamics of circulating waves of cell activation at the signalling centre. The circulation rate modulates the speed, frequency, and width of longrange periodic activation waves, spatially organised as a spiral. Leading up to aggregation, the circulation rate of the spiral tip gradually increases, resulting in an increase in the global wave frequency and a decrease in wave speed by virtue of the spiral wave rotating faster and becoming more curved. Using the general theory of excitable media, mathematical modelling, and physical perturbations of population structure, our data imply the increase in the circulation rate arises due to positive feedback between the cell signalling pattern and cell signalling state. Specifically, the wave circulation rate is determined by the ability of cells to relay the signal (cell excitability), whilst the ability to relay increases as the circulation rate increases. As such, the information encoding a morphogenetic transition unfolds from positive feedback between the signalling pattern and the cell signalling properties at the signalling centre. Altogether, our results describe how the collective behaviour of 10^{5}–10^{6} cells across several millimetres is organised by circulation of a selfsustaining signal over 10^{2}–10^{3} cells.
This study distinguishes biological activation waves from those produced by idealised mathematical models and chemical reactions (Winfree, 1972). The difference is that biological activation waves can encode information that can be modulated and translated to changes in the spatial organisation and gene expression of individuals beyond the spatial scale of the individual, which can feedback to changes in activation wave dynamics and structure. For Dictyostelium aggregation, information is generated locally by the rate of circulation of the spiral tip and decoded globally by cell chemotaxis. Here, by information, we mean the properties of the periodic waves, such as frequency and speed, which are translated by cells into a chemotactic and/or gene regulatory response. Information decoding via cell chemotaxis has been studied at both the single cell and population levels. There are several contrasting views regarding the wave properties which cells interpret for chemotaxis. One view is that the ability for cells to move collectively increases with an increase in wave frequency. Skoge et al., 2014 showed, by controlling the frequency of cAMP waves, that the chemotactic index increased as the frequency became maximal (Skoge et al., 2014). This observation is consistent with both our correlative and optogenetic data showing onset of collective migration as the wave frequency approaches the maximum cell activation rate. An alternative view is that cells ‘count’ waves: artificial rotating waves revealed that 4–5 waves at high frequency can induce cell polarity migration (Nakajima et al., 2016). While we cannot rule out a counting mechanism, based on our data here and early data on population measurements of aggregation, in which cells would be counting upwards of 30 waves before making the decision to collective migrate, any counting mechanism would be contextdependent, for example, using counts above a certain frequency of pulses. An alternative view is based on the wave gradient, where cells move maximally towards steep stationary (Fisher et al., 1989; Song et al., 2006) and/or increasing (Varnum et al., 1985; Wessels et al., 1992; Nakajima et al., 2014) gradients over time. Again, this view is consistent with our data, where both the wave speed and width decrease as the wave frequency increases, meaning that cells are exposed to a steeper gradient for a longer time. Overall, these different mechanisms are looking at different properties of the wave. What our study shows is that all the wave properties: speed, frequency, and width are specified by the circulation of the spiral tip. Specifically, once the frequency is set by the tip, all other features are set by default. Altogether, our study implies that the circulation rate of the signal at the signalling centre is the mechanism of information production, regardless of the specific informationdecoding mechanism.
The mechanism by which cell state changes in response to the changing wave dynamics is explained, at least in part, because information contained within periodic waves is also decoded in terms of gene expression. In Dictyostelium, and many other systems, gene expression is responsive to the frequency of signalling. Expression of pdsA, which encodes a cAMP phosphodiaterase, and csaA, which encodes a cellcell adhesion protein, is sensitive to the changing frequency of cAMP pulses (Masaki et al., 2013; Corrigan and Chubb, 2014). The csaA gene is induced during infrequent cAMP signalling, then inactivated as the pulse frequency increases as cells aggregate. The transcription factor GataC displays rapid shuttling between the nucleus and cytoplasm in response to external cAMP oscillations (Cai et al., 2014). This process is lowpass filtered, with nucleocytoplasmic shuttling occurring only when external cAMP oscillations reach a threshold frequency. Here, we describe a mechanism controlling a gradual increase in the global wave frequency. In the light of these findings, it is conceivable that Dictyostelium use high and/or lowpass filtered transcription factors to exploit the gradual increasing wave frequency to control the timing of gene activation. In this way, the sequential timing of activation of several transcription factors across the cell population can be coordinated by signalling centres that disseminate periodic waves that steadily change frequency over time. A change in cell excitability could also arise independent of changes in gene expression. For example, we can imagine a scenario in which the degradation of cAMP cannot keep up with the increased wave frequency. In this situation, the baseline cAMP increases over time.
In Dictyostelium, the temporal coding of cAMP waves allows cells to sense the extent of cell aggregation through the frequency of activation waves produced by the aggregate, which increases as aggregation progresses. In other words, the temporal code allows for longrange sensing of aggregate progression, which can be interpreted as a nonlocal quorum sensing mechanism. A potential advantage of this adaptive information coding system is that by moving only towards highly frequent and consistent cAMP waves, effort is then concentrated on migration directed only towards an established signalling centre. If cells were to move maximally towards each initial wave, they would meander and unnecessarily expend energy while cAMP waves settle to a stable pattern, in addition to missing out on the opportunity to chemotaxis towards sources of nutrition.
Our data reveals the transition from the single to multiarmed spiral wave during aggregate progression (Agladze and Krinsky, 1982). In the established aggregate, the singlearmed spiral is associated with an unstable signalling phase singularity at the aggregate centre. The multiarmed spiral lacks a clear singularity and instead is associated with relatively small numbers of desynchronised cells at the spiral core. This is the hallmark of the chimera spiral state, a phenomenon yet to be observed in a biological system (Totz et al., 2018). Overall, we view the transition to the multiarmed spiral as an adaption to the change in aggregate geometry to maintain an active phase wave at the signalling centre that generates activation waves at the maximum frequency.
In their physiological context, most cell populations, including Dictyostelium, are structured in three spatial dimensions. It is not trivial to extrapolate our findings (in the standard 2D plane context) to understand threedimensional activation waves, owing to the limitless number of selfsustaining wave patterns—the scroll waves—afforded by the extra spatial dimension (Siegert and Weijer, 1992; Winfree, 2001; Durston, 2013). To better understand how biological signals are assembled, sustained, and controlled in a more general class of tissues, it would be important to extend the experimental, computational, and mathematical approaches outlined in this study to characterise the structure, dynamics, and biological effects of cAMP wave patterns that coordinate Dictyostelium morphogenesis in 3D.
Materials and methods
Generation of cell lines
Request a detailed protocolDictyostelium AX3 cells were engineered to express Flamindo2 and H2BmCherry from safehaven genomic locations for stable and uniform expression. Flamindo2 is a cAMP biosensor (Odaka et al., 2014; Hashimura et al., 2019) with a fluorescence intensity that is inversely proportional to internal cAMP levels. A codonoptimised sequence encoding the Flamindo2 reporter was knockedin to the act5 locus for stable and uniform expression (Paschke et al., 2018; Tunnacliffe et al., 2018). H2BmCherry, knocked into the rps30 locus, labels the nucleus for cell tracking at highcell densities (Corrigan and Chubb, 2014). For optogenetic activation of cAMP signalling, we used the lightactivatable adenylyl cyclase, bPAC (Stierl et al., 2011). A codonoptimised bPAC cDNA was cloned into the extrachromosomal expression vector, pDM1203, then transfected into the act5Flamindo2 cells.
Cell handling
Request a detailed protocolCells were grown at 22°C in 12 ml of Formedium HL5 medium supplemented with penicillin and streptomycin in 10 cm tissue culture dishes. Cultures were diluted to 25% of their density each day to stop cell entry into stationary phase. Cells were not used beyond 9 days of culture. For development, logphase cells at 80% confluency were washed and resuspended in KK2 (20 mM KPO_{4} pH 6.0) and then plated on 1.5% agar at between 4×10³ and 5×10³ cells mm^{−}². After 1 hr in a humidified chamber, a 1.2 cm × 1.2 cm slab of agar was cut and gently placed cellfacedown onto a glass imaging dish (Bioptechs deltaT) and then immediately submerged in silicon oil to prevent desiccation. To incorporate physical blocks into the cell aggregation field, we sprinkled lowmelting point agarose powder onto the KK2 agar bed before the agar had set.
Cell imaging
Request a detailed protocolCells were imaged at 22°C on a custombuilt inverted wide field microscope (Cairn Research) with 470 nm and 572 nm LED light sources (Miermont et al., 2019). This was equipped with a Prime 95B CMOS camera (Photometrics) and 10× UplanFL N objective (Olympus). A 1.2 mm × 1.2 mm area was imaged at resolution 1 µm² per pixel every 4 s for 8 hr from 2 hr after starvation. Frame stitching was used to expand the effective fieldofview size. To image cell motion and signalling about the spiral core, including the experiments with the physical block, we imaged a 1.2 mm × 1.2 mm area for 2 hr following from 5 to 6 hr after starvation. For this, we monitored cAMP wave patterns in real time and identified waves that had started to rotate about an axis.
Photoactivation
Request a detailed protocolFor light activation of cAMP signalling, we used bPACexpressing act5Flamindo2 cells and a 3i spinningdisk confocal microscope with a 10× objective and Prime 95B CMOS camera (Photometrics). To initiate periodic cAMP waves, cells contained within a circular area with a ≈100 µm diameter were pulsed with 488 nm light for 1 s, at intervals to mimic physiological cAMP pulsing. Flamindo2 images (excitation 515 nm) were collected throughout the pulse sequence, every 5 s. To initiate a rotating spiral wave, we repeated the photoactivation protocol above, except that cells were photoactivated at operatorcontrolled points in time and space. We monitored waves in real time and identified periodic waves away from the signalling centre (around one per 5 min, 5 hr following starvation). To induce circular wave patterns, we photoactivated cells around 3 min after the previous wave. To induce asymmetric wave patterns (resulting in the broken and hence spiral wave), we photoactivated around 2 min after the previous wave.
Cell tracking and measuring cell signalling states
Request a detailed protocolThe following image analysis pipeline was implemented in Matlab 2020b. To track the position of individual cells from image sequences, the IDL particle tracking algorithm (Crocker and Grier, 1996) was applied to the cell nuclei channel (H2BmCherry) of the image time series. The velocity of each cell at each time point was estimated by the first derivative, second order, central finite difference of cell position. To measure the mean fluorescence intensity of each tracked cell, watershed segmentation (the Fernand Meyer algorithm) seeded from the tracked cell positions was applied to the cell body channel (Flamindo2). To infer internal cAMP levels (nondimensional), the tracked mean fluorescent intensity of each cell was normalised and smoothed using Matlab’s SavitzkyGolay finite impulse response smoothing filter. Two measures were used to identify cAMP pulses: (i) the similarity between the normalised intensity and a square pulse of length 2 min and (ii) the signalling phase $\phi $. Variable $\phi $ was estimated by stratifying the time series of the internal cAMP level $c$ by its time derivate $\dot{c}$ (estimated by the first derivative, second order, central finite difference). The value of $\phi $ is set as the angle around the point $c={c}^{*}$ and $\dot{c}=0$:
In addition to the time series of the inferred cAMP levels, we recorded the time points of the maximal cAMP level during each pulse, and the time points of cell activation ($\dot{c}>0$) and deactivation ($\dot{c}<0$), measured by the time point where the cAMP levels are halfmaximal.
To quantify the internal [cAMP] signalling and phase field following aggregation, we tracked the fluorescence intensity of each individual pixel as a function of time, over the entire 2D image sequence. The fluorescence intensity time series was oscillatory and typically increased over time due to cell crowding. We used Matlab’s SavitzkyGolay finite impulse response smoothing filter to separate the oscillatory and background component of the signal. The background component was used to estimate the cell density of the aggregate. The oscillatory component was used to estimate the spatial distribution of signalling states and phases across the aggregate, using the same technique described above.
Activation wave curve fitting
Request a detailed protocolTo quantify the evolving dynamics of the spiral wave from cell track data, we first estimated the location of the spiral phase singularity at each time point as the central point of a circular region of radius 100 µm with the largest variance in internal [cAMP] across each cell. The location of the spiral tip was estimated by the mean position of five maximally activated cells closest to the location of the phase singularity. We estimate the axis of rotation at each time point by fitting a circle to the track of the spiral tip, half a revolution backwards and forwards in time. Using the Matlab function ‘fit’, we fit the spiral curve given by Equation 2 (see below) to the spatial distribution of maximally activated cells, starting from the axis of spiral tip circulation (Tyson and Keener, 1988). To achieve this, we determined the magnitude ($r$ in Equation 2) and orientation ($\theta $ in Equation 2) of the vector extending from the axis of spiral rotation to the location of each maximally active cell. This fit provides the spiral orientation (${\theta}_{0}$ in Equation 2) and curve parameter ($k$ in Equation 2).
Away from the spiral core, the time and place of activation for each cell were clustered in three dimensions (two spatial and one temporal) using the DBSCAN (densitybased clustering of applications with noise) algorithm as a means to identify cells with correlated signalling states (e.g. during the propagation of cAMP waves) and also, by exclusion, the cells that activate independently. This approach reveals the clusters of cell activation associated with each cAMP wave. For each wave, we estimated the dynamic curve of the wave by fitting an equation for a cone ($t{t}_{0}=Nr,r=\sqrt{{\left(x{x}_{0}\right)}^{2}+{\left(y{y}_{0}\right)}^{2}}\mathrm{}$, where ${t}_{0}$ and ${(x}_{0},{y}_{0})$ is spatiotemporal origin of the wave) to the time evolution of the spatial distribution of maximally activated cells. The slope of each of the cones fitted to the spatial distribution of maximally activated cells at each time point was used to deduce the speed of each wave.
Quantifying cell motion
Request a detailed protocolTo estimate the mean cell velocity field and density around the spiral wave, we measured the mean velocity and total number of cells at regular lattice areas (step size 10 µm) in the reference frame that rotates with the fitted spiral parameter ${\theta}_{0}$ (Figure 6—video 2). The streamlines across this field were analysed using the Matlab function ‘streamlines’.
To analyse the behaviour of single cells during stable periodic waves, the fitted cone equations were used to estimate the mean internal [cAMP] and velocity of each cell as a function of distance from the wave centre. This allows an estimation of measures such as the angle difference between the direction of wave travel and cell motion, for each cell:
where $\widehat{{\mathit{v}}_{i}}$ is the normal velocity of cell $i=1,\dots ,N$, and $\widehat{\mathit{w}}\left({\mathit{X}}_{i}\right)$ is the normal velocity of the cAMP wave at the location of cell $i$, ${\mathit{X}}_{i}$ . The following two formulae were used to measure collective cell motion: (i) the mean chemotaxis index (the proportion of cell movement in the direction of the oncoming wave) as a function of distance from the wave:
where ${\Omega}_{j}$ represents the subset of cells (subset size ${N}_{{\Omega}_{j}}$) located at a distance $j\times 10$ µm ($j=\dots \mathrm{1,0},1,\dots $) from the wave centre, and (ii) the local correlation of cell velocity:
where ${\Phi}_{i}$ denotes the subset cells (subset size ${N}_{{\Phi}_{i}}$) within a 50 µm radius of cell $i$ (De Palo et al., 2017). Local cell density was measured by the mean density of cells within a 60 µm radius of each cell.
General theory of excitable media
Request a detailed protocolThis section outlines the elements of the theory we have used, and we used these to interpret spiral wave progression (Tyson and Keener, 1988; Winfree, 2001).
The eikonal relation dictates how the normal speed along the wave front $N$ changes due to wave curvature $K$:
where ${N}_{0}$ is the planar (no curvature) wave speed, $D$ is the signal diffusion coefficient (which we assume remains constant), and $\epsilon $ is the ratio of the rates of cell inactivation and activation (which we call the cell activation parameter; Keener, 1986). For a wave extending to the pivotal axis, the solution to Equation 1 is accurately approximated by the following curve (in polar coordinates, $x=r\mathrm{cos}\theta (r,t),y=r\mathrm{sin}\theta (r,t)$):
where $m$ is a constant. Parameter $k$ (the curve parameter) defines the curvature of the wave, and ${\theta}_{0}$ is the spiral angle of rotation, with rotation rate, ${\omega =\dot{\theta}}_{0}$ (Tyson and Keener, 1988). We fit Equation 2 to the spatial distribution of maximally activated cells at each time point, throughout the entire progression of the spiral wave. The fit is accurate for all time points, which validates the use of the theory of excitable media even though cells are discretely distributed across space. This fit provides a measurement of the time evolution of the curve parameter $k$, rotation rate $\omega $, and hence the wave period $T=2\pi /\omega $ (Figure 4D). The value of $k$ relates the planar wave speed to the cell excitability:
Substituting Equation 3 into Equation 1 for ${N}_{0}$ allows the cell activation parameter $\epsilon $ and planar wave speed ${N}_{0}$ to be expressed in terms of variables that can be directly measured from data:
Consistent with the traditional paradigm of spiral wave progression (Höfer et al., 1995; van Oss et al., 1996), we find that the cell activation parameter, which is the inverse of cell excitability, decreases during spiral wave progression (Figure 5). This result illustrates how we can identify key signalling properties of the cells from their population signalling patterns, at each time point during aggregation. Plotting ${N}_{0}$ from Equation 4 against the wave period $T$ provides a prediction of the dispersion relation, which quantifies how a decrease in wave period reduces the wave speed because cells do not fully reset between waves, causing them to activate slower. The curvature relation describes the relationship between radius and rate of circulation of the spiral tip:
where ${r}_{0}$ is the radius of spiral tip circulation, which can measure from our data (Tyson and Keener, 1988). In support of our approach, the curvature and eikonal relations produce similar values of the cell activation parameter at all time points. We use the curvature relation to predict the circulation radius of the spiral tip based on our spiral measurements.
Mathematical model (ABM)
Request a detailed protocolA hybrid PDE ABM was developed to study the interplay between cell motion, spiral wave dynamics, and spatial patterning of cells around the spiral core.
External cAMP is secreted by active cells. We define active cells as having a positive activator value ${A}_{i}>0$ (Sgro et al., 2015). Once secreted, we assume external cAMP diffuses and decays. The 2D external [cAMP] field $c$ was modelled by the following reactiondiffusion equation:
where $\beta $ and $\sigma $ are the rates of external cAMP decay and release from individual cells, $D$ is the external cAMP diffusion coefficient, $N$ is the number of cells, $H$ is the Heaviside function, and $\delta $ is the Dirac delta function, modelling the point source release of cAMP from cell $i$ ($1\le i\le N$), with position ${\mathit{X}}_{i}$ and activation state ${A}_{i}$. The spatial domain is $\mathrm{\Omega}={\left[\frac{L}{2},\frac{L}{2}\right]}^{2}\subset {\mathbb{R}}^{2}$ where $L>0$ is the domain size. We used same parameter values of VidalHenriquez and Gholami (VidalHenriquez and Gholami, 2019). The following modified FitzHughNagumo model was used to describe the excitable release of cAMP from each cell:
where ${I}_{i}$ is the inhibitor level of cell $i$, and the value of each parameter was the same as Sgro et al., 2015, with exception of the time scale $\tau $ which was adjusted to match experimental data.
We tested several mathematical models of cell chemotaxis against the data for the velocity of individual cells as a function of the external cAMP field. The external cAMP field was inferred in our experimental data by substituting the measured position and discretised activation state of each cell in the reactiondiffusion equation above, which provides a prediction for the external [cAMP] at the location of each cell ${c}_{i}={c}_{i}\left({\mathit{X}}_{i}\right(t),t)$ and the spatial gradient $\left{c}_{i}\right$ (Figure 7—videos 1 and 2). The model that best fit the experimental data was a simple system of linear ordinary differential equations that encapsulate chemotactic memory:
where $0\le {\theta}_{i}<2\pi$ and ${v}_{i}\ge 0$ are the direction and speed of cell movement, and ${s}_{i}\ge 0$ is the sensation variable to external cAMP as used in Höfer et al., 1995. This variable causes cell chemotaxis to be sensitive to the external cAMP gradient at the front of the wave, but not the rear. Variable $\left\nabla {c}_{i}\right$ and $0\le {\theta}_{i}^{c}<2\pi$ are, respectively, the magnitude and direction of the external [cAMP] spatial gradient evaluated at ${\mathit{X}}_{i}$. Model parameters ${\alpha}_{m}$, ${\beta}_{m}$, ${\gamma}_{m}$, ${\sigma}_{m}$, and ${\lambda}_{m}$ were determined during the model identification using the Matlab function ‘fminsearch’ to obtain the best fit to experimental data. The first equation models the turning of the cell towards the local gradient of external cAMP. The turning speed depends on a cell’s sensitivity and the size of the cAMP gradient. The second equation describes adaptation of the cell speed. It drives cell speed to be proportional to ${s}_{i}\left\nabla {c}_{i}\right$ , i.e., higher sensitivity and larger local gradients will lead to faster speeds. Finally, the third equation models sensitivity: If a cell experiences a high concentration of external cAMP ${c}_{i}$, it will become desensitised, and it takes some time to become sensitive again. The variable ${s}_{i}$ provides ‘inertia’ or ‘memory’ that makes cell chemotaxis to be dependent on the history of external [cAMP] dynamics. When comparing the model predictions to the data, this effect was crucial to ensure that cells move towards the wave source at the wave front while also neglecting the wave rear. Specifically, between waves, the value of ${s}_{i}$ recovers to a state where the cells are sensitive to the external cAMP gradient at the front of the wave but subsequently ${s}_{i}$ reduces to a state where the cells are insensitive to the external cAMP gradient at the rear of the wave. For the final cell movement model, we also consider size exclusion effects between cells by including a cellcell repulsion term. This means ${\theta}_{i}$ and ${v}_{i}$ are understood to mean the orientation and speed a cell would have if not in contact with other cells. We define as ${d}_{ij}={\mathit{X}}_{i}{\mathit{X}}_{j}$ the distance between the cell centres of the $i$th and $j$th cell. A cell’s movement now follows:
i.e., cells move with speed ${v}_{i}$ in direction ${\theta}_{i}$ in the absence of cellcell repulsion. The second term describes a cellcell repulsion with maximal repulsion magnitude ${\rho}_{R}$ where only cells within a distance of ${d}_{R}$ of each other exert a pushing force. See Table 1 for parameter values.
We use noflux boundary conditions for the external cAMP concentration, i.e., cAMP cannot leave the domain. Similarly, cells cannot leave the domain. Since the dynamics cause a concentration of cells towards the spiral centre, cells are being depleted from the very edge of the simulation domain, and hence, we only considered the middle part of the simulation domain to represent the biological situation. For all numerical experiments, we initialise cell positions ${\mathit{X}}_{i}$ and orientations ${\theta}_{i}$ using a uniform random distribution on $\mathrm{\Omega}$ and $[0,\mathrm{}2\mathrm{\pi})$, respectively. Cell sensitivities ${s}_{i}$ were initially set to 1 and cell speeds ${v}_{i}$ to 0. The initial external cAMP concentration, activator, and inhibitor values were initialised as shown in Figure 7—video 3. The initial state was created by disabling cell movement and: (1) letting a vertical cAMP wave moving from left to right across a periodic domain until it is equilibrated, (2) setting the external cAMP in the upperhalf of the domain to zero and setting the inhibitor value in the upper half on the domain to a high value (we used 2.5), and (3) setting the boundary conditions to noflux and restarting the simulation.
Model equations were solved numerically using Matlab. The external [cAMP] field was solved across a 50 grid points per mm with the finite difference, implicit Euler method in time and central difference in space. Cell velocities and activation states were solved using the Matlab function ‘ode45’.
Mathematical model (PDE)
Request a detailed protocolTo test whether the coupling of cell excitability to the external [cAMP] dynamics sufficiently describes spiral wave progression, we extended a mathematical model of cAMP spiral waves to simulate the time evolution of the external [cAMP] field over time. We adopt the model by Tyson and Murray, 1989, which is of the general form of continuum excitable media with the dynamics of cAMP release governed by the MartielGoldbeter model:
where $u$ and $v$ (nondimensional) are the spatial distribution of external [cAMP] (which can diffuse) and an internal inactivator (which cannot diffuse). We adopted the following nondimensional parameter values from Tyson and Murray, 1989, most of which were rounded to the nearest order of magnitude: ${\alpha}_{u}={10}^{3},{\beta}_{u}={10}^{2},\kappa ={10}^{1},{L}_{1}={10}^{1},{L}_{2}={10}^{1},{\lambda}_{1}={10}^{4},{\lambda}_{2}=0.5,c={10}^{1}$. Equation 6 is dimensionalised by the time scale $\tau =0.075$ min and diffusion coefficient $D=3\times {10}^{4}$ µm^{2} min^{−1}. We note that solutions to this equation rapidly tend to an equilibrium state. Motivated by our finding that the cell excitability is coupled to the cAMP wave dynamics, we extended Equation 6 by letting the cell activation variable, which we previously called the cell activation parameter (this inverse being a measure of cell excitability), to be a simple timedependent function of the external [cAMP]:
Where ${\alpha}_{\epsilon}=\frac{1}{50},{\beta}_{\epsilon}=\frac{35}{50},{\epsilon}_{0}=10$. Equation 7 states that the cell activation variable tends to a baseline value ${\epsilon}_{0}$ and decreases in response to external cAMP. The simulations were performed on a 4.5 mm × 4.5 mm grid with 0.025 mm spacing and 0.05 min time step. The initial condition was the state where $=0,v=0$ , $\epsilon =\frac{{\epsilon}_{0}}{2}$, except for a rectangular region representing a broken wave, with $u=10$ and $v=0$ and to one side $u=0$ and $v=1$. To numerically solve Equations 6; 7, we used the implicit finite difference method for the diffusion term and the RungeKutta (Matlab ‘ode45’ function) method for the reaction terms. The Neumman boundary condition (no flux) was used.
Data availability
The high resolution movies of cAMP signalling and optogenetic treatments can be accessed at UCL's institutional repository using the DOI: 10.5522/04/21360975. All image analysis and mathematical modelling methodology described in full in methods section. Matlab code for the agentbased model is provided.
References

LTB4 is a signalrelay molecule during neutrophil chemotaxisDevelopmental Cell 22:1079–1091.https://doi.org/10.1016/j.devcel.2012.02.003

Propagating wave of ERK activation orients collective cell migrationDevelopmental Cell 43:305–317.https://doi.org/10.1016/j.devcel.2017.10.016

Methods of digital video microscopy for colloidal studiesJ Colloid Interface Sci 179:298–310.https://doi.org/10.1006/jcis.1996.0217

A discrete cell model with adaptive signalling for aggregation of Dictyostelium discoideumPhilosophical Transactions of the Royal Society of London. Series B, Biological Sciences 352:391–417.https://doi.org/10.1098/rstb.1997.0029

Chemical waves in cell and developmental biologyThe Journal of Cell Biology 217:1193–1204.https://doi.org/10.1083/jcb.201701158

Dictyostelium discoideum aggregation fields as excitable mediaJournal of Theoretical Biology 42:483–504.https://doi.org/10.1016/00225193(73)902427

Spatial trigger waves: Positive feedback gets you a long wayMolecular Biology of the Cell 25:3486–3493.https://doi.org/10.1091/mbc.E14081306

Dictyostelium discoideum: Cellular selforganization in an excitable biological mediumProceedings. Biological Sciences 259:249–257.https://doi.org/10.1098/rspb.1995.0037

A geometrical theory for spiral waves in excitable mediaSIAM Journal on Applied Mathematics 46:1039–1056.https://doi.org/10.1137/0146062

Competing patterns of signaling activity in Dictyostelium discoideumPhysical Review Letters 76:1174–1177.https://doi.org/10.1103/PhysRevLett.76.1174

Synchronization of pulsecoupled biological oscillatorsSIAM Journal on Applied Mathematics 50:1645–1662.https://doi.org/10.1137/0150098

Rectified directional sensing in longrange cell migrationNature Communications 5:5367.https://doi.org/10.1038/ncomms6367

The microfluidic lighthouse: An omnidirectional gradient generatorLab on a Chip 16:4382–4394.https://doi.org/10.1039/c6lc00898d

Rotors and the dynamics of cardiac fibrillationCirculation Research 112:849–862.https://doi.org/10.1161/CIRCRESAHA.111.300158

Cellular memory in eukaryotic chemotaxisPNAS 111:14448–14453.https://doi.org/10.1073/pnas.1412197111

Dictyostelium discoideum chemotaxis: Threshold for directed motionEuropean Journal of Cell Biology 85:981–989.https://doi.org/10.1016/j.ejcb.2006.01.012

Light modulation of cellular camp by a small bacterial photoactivated adenylyl cyclase, bpac, of the soil bacterium beggiatoaThe Journal of Biological Chemistry 286:1181–1188.https://doi.org/10.1074/jbc.M110.185496

Spatial pattern formation during aggregation of the slime mould Dictyostelium discoideumJournal of Theoretical Biology 181:203–213.https://doi.org/10.1006/jtbi.1996.0126

Spontaneous center formation in Dictyostelium discoideumScientific Reports 9:3935.https://doi.org/10.1038/s41598019403734

Dictyostelium morphogenesisCurrent Opinion in Genetics & Development 14:392–398.https://doi.org/10.1016/j.gde.2004.06.006

Behavior of Dictyostelium amoebae is regulated primarily by the temporal dynamic of the natural camp waveCell Motility and the Cytoskeleton 23:145–156.https://doi.org/10.1002/cm.970230207
Article and author information
Author details
Funding
Wellcome Trust (202867/Z/16/Z)
 Jonathan R Chubb
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication. For the purpose of Open Access, the authors have applied a CC BY public copyright license to any Author Accepted Manuscript version arising from this submission.
Acknowledgements
This study was supported by a Wellcome Trust Senior Fellowship (202867/Z/16/Z) to JRC. The funder had no role in study design, data collection and interpretation, or the decision to submit the work for publication. We thank Elizabeth Westbrook, Ricardo Barrientos, Courtney Lancaster, Kees Weijer, and Allyson Sgro for comments on the manuscript.
Copyright
© 2023, Ford et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics

 4,435
 views

 345
 downloads

 6
 citations
Views, downloads and citations are aggregated across all versions of this paper published by eLife.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading

 Chromosomes and Gene Expression
 Developmental Biology
Most teleost fishes exhibit a biphasic life history with a larval oceanic phase that is transformed into morphologically and physiologically different demersal, benthic, or pelagic juveniles. This process of transformation is characterized by a myriad of hormoneinduced changes, during the often abrupt transition between larval and juvenile phases called metamorphosis. Thyroid hormones (TH) are known to be instrumental in triggering and coordinating this transformation but other hormonal systems such as corticoids, might be also involved as it is the case in amphibians. In order to investigate the potential involvement of these two hormonal pathways in marine fish postembryonic development, we used the Malabar grouper (Epinephelus malabaricus) as a model system. We assembled a chromosomescale genome sequence and conducted a transcriptomic analysis of nine larval developmental stages. We studied the expression patterns of genes involved in TH and corticoid pathways, as well as four biological processes known to be regulated by TH in other teleost species: ossification, pigmentation, visual perception, and metabolism. Surprisingly, we observed an activation of many of the same pathways involved in metamorphosis also at an early stage of the larval development, suggesting an additional implication of these pathways in the formation of early larval features. Overall, our data brings new evidence to the controversial interplay between corticoids and thyroid hormones during metamorphosis as well as, surprisingly, during the early larval development. Further experiments will be needed to investigate the precise role of both pathways during these two distinct periods and whether an early activation of both corticoid and TH pathways occurs in other teleost species.

 Developmental Biology
 Neuroscience
Adolescence is characterized by changes in rewardrelated behaviors, social behaviors, and decisionmaking. These behavioral changes are necessary for the transition into adulthood, but they also increase vulnerability to the development of a range of psychiatric disorders. Major reorganization of the dopamine system during adolescence is thought to underlie, in part, the associated behavioral changes and increased vulnerability. Here, we utilized fast scan cyclic voltammetry and microdialysis to examine differences in dopamine release as well as mechanisms that underlie differential dopamine signaling in the nucleus accumbens (NAc) core of adolescent (P2835) and adult (P7090) male rats. We show baseline differences between adult and adolescentstimulated dopamine release in male rats, as well as opposite effects of the α6 nicotinic acetylcholine receptor (nAChR) on modulating dopamine release. The α6selective blocker, αconotoxin, increased dopamine release in early adolescent rats, but decreased dopamine release in rats beginning in middle adolescence and extending through adulthood. Strikingly, blockade of GABA_{A} and GABA_{B} receptors revealed that this α6mediated increase in adolescent dopamine release requires NAc GABA signaling to occur. We confirm the role of α6 nAChRs and GABA in mediating this effect in vivo using microdialysis. Results herein suggest a multisynaptic mechanism potentially unique to the period of development that includes early adolescence, involving acetylcholine acting at α6containing nAChRs to drive inhibitory GABA tone on dopamine release.