Swarming assay

All experiments were performed with a naturally competent derivative of the undomesticated B. subtilis NCIB 3610 strain, carrying the comIQ12L allele50. Bacterial cultures for swarm inoculation were grown in Luria-Broth (LB) medium (Carl Roth, X968.3) for 16 h at 37 °C with shaking at 220 r.p.m. Swarming plates containing 9 ml of 0.5% Bacto agar (BD, 214010) in LB in 9 cm Petri dishes were dried for 10 min with an open lid at room temperature and then kept upside down at 37 °C for an additional 10–30 min before inoculation. Cells were transferred onto the swarming plate by dipping a sampling tip (see description of robotic sampling tips below) into the overnight culture and then lightly touching the surface of the swarming plate with this tip. This procedure resulted in a culture spot of diameter ~1 mm containing 3,000–5,000 cells on the agar surface.

After inoculation, the swarming plate was moved into a humidity chamber whose bottom was a microscope stage insert. This humidity chamber contained water-filled basins to ensure high humidity. The humidity chamber and large parts of the inverted microscope (Nikon Ti-E) were enclosed by a microscope incubator (Okolab). During the swarming assay, the temperature inside the Okolab enclosure was kept constant at 30 °C. The humidity chamber containing the agar plate and the stage insert was mounted on a motorized xy stage (Applied Scientific Instrumentation). Brightfield images were acquired using a ×16 water-immersion objective with numerical aperture 0.8 and an Andor Zyla 4.2PLUS sCMOS camera at 100 fps.

Robotic sampling of cells for transcriptome measurements and imaging of cells during swarm development

We performed automated sampling of cells from the developing swarm using a custom-built robotic setup, which consisted of several motorized linear or rotational stages that collectively formed a robotic sampling arm, as schematically illustrated in Extended Data Fig. 1. All stages within the robotic sampling arm were controlled by the same Matlab algorithm that also controlled the motorized microscope xy stage, microscope camera and microscope focus to ensure that the robotic sampling arm is synchronized with the microscope movement and imaging.

One robotic sampling cycle consists of the following processes: Stages 1, 2 and 3.4 (stage numbers are defined in Extended Data Fig. 1) were translated such that a sampling tip was picked up from the tip box onto the tip holder attached to stage 3.4. Using rotational stage 3.3 and stage 1, the tip was then transported to the xy stage of the microscope where the swarming plate is located. The sampling tip was moved into the humidity-controlled incubation enclosure on the microscope xy stage (within which the swarming plate is located) through a hole in the cover of the enclosure, which can be opened and closed by stage 4. The sampling tip was then brought into the correct xy position for picking up cells using stages 1, 3.1 and 3.4, lowered onto the sample and brought in contact with the swarm surface by movement of stage 3.2, using live image analysis of the microscope camera’s field of view. After establishing contact with the swarm surface, the tip was kept at the surface for ~20 s before reversing the movement to retrieve the sampling tip. The duration for which the tip was in contact with the swarm surface was adapted according to the sampling position, with less contact time for dense regions. Control experiments using counting of colony forming units have shown that the sampling tip after this maneuver contains 10³–105 cells. Finally, by moving stage 3.5, the tip was moved to the sample box and ejected into a 1.5 ml Eppendorf tube containing 50 µl of lysis buffer (40 U µl−1 Ready-lyse lysozyme (Lucigen, R1804M), 0.04 U µl−1 SUPERase in RNase Inhibitor (Thermo Fisher, AM2694), 10 mM Tris, adjusted to pH 8.0 with HCl (Thermo Fisher, 15568025), 1 mM EDTA (Thermo Fisher, AM9261)). The 50 µl of lysis buffer covered the cell-containing part of the tip. The Eppendorf tube was immediately manually collected, closed and snap-frozen in liquid nitrogen, followed by storage at −80 °C until RNA isolation. The Eppendorf tube containing the lysis buffer was kept on ice until a few minutes before being inserted into the sample ejection box on the robot. To maintain the low temperature of the lysis buffer during the time that it spends in the warm environment surrounding the microscope, Eppendorf tubes were kept in an ice-cold aluminum block, which was surrounded by isolation material. The aluminum block contains slots for up to six Eppendorf tubes and the aluminum block was exchanged every 5–8 min with an identical copy previously resting on ice. These steps to ensure a constant low temperature of the lysis buffer were necessary to retain the enzyme activity and minimize transcriptome changes and RNA degradation.

Manufacturing of the sampling tips for the robotic arm: We designed custom sampling tips to avoid scratches in the agar surface during the sampling procedure, which would alter the swarm development. These sampling tips consisted of a standard Eppendorf 10 µl pipette tip dipped into polydimethylsiloxane (PDMS, SYLGARD 184, Dow Corning) to produce a drop at the tip edge. To ensure reproducible PDMS drop sizes on the tip, custom holders for dipping 48 pipette tips simultaneously into a levelled bath of liquid PDMS were designed. The PDMS on the tips was then cured at 50 °C for 24 h. The cured sampling tips were then sterilized by dipping the tips into 70% ethanol, followed by drying at room temperature in the closed tip rack for 2 h. The sterilized sampling tips were used within 24 h.

Robotic sampling procedure in space and time: After the inoculation of the swarming plate as described in the section ‘Swarming assay’, growth of the cells was monitored by brightfield imaging until the end of the lag phase was reached. Once the swarm started expanding, robotic sampling was initiated. For the whole duration of the swarm expansion, until a final swarm diameter of 6 cm was reached, sets of spatially separated samples were acquired in 15–20-min intervals. For each set of spatially separated samples, a line was drawn from the centre of the swarm to the outer edge of the swarm, and samples were collected from up to 9 regularly spaced positions along this line (the exact number of sampled positions depends on the diameter of the swarm at any given time, and the spacing between sampling positions was at least 1 mm). The time for each sample to be acquired by the robot was ~30 s. To avoid sampling the same spot within the swarm too many times and therefore causing disruption to this area of the swarm, the innermost sample was chosen with a distance of 1.5 mm to the point of inoculation, and the angle between the x axis of the microscope stage and the line at which samples were taken was adapted by 30° after each run such that a sampling line with the same angle was only used for every 8th set of spatially separated samples. Using this procedure, samples were acquired from only half of the swarm area, whereas the other half of the swarm was untouched by the sampling tip. Immediately before and after sampling, an overview image of the sampling area was taken to document the impact of the sampling process on the swarm (Extended Data Fig. 4). These images were also used to estimate the number of cells extracted.

In addition to acquiring brightfield microscopy videos from the sampling location before sampling, we also acquired videos at spatial locations between the sampling locations. Furthermore, as control experiments, microscopy videos were acquired at corresponding spatial locations on the half of the swarm area that was untouched by the robotic sampling tip to determine whether the half of the swarm from which samples were acquired for transcriptome measurements developed differently from the half of the swarm that was untouched by the sampling tip.

Sampling of cells during the lag phase for transcriptome measurements before swarm expansion

As the number of bacterial cells present on the agar surface during the lag phase is low and cell density plays an important role for the differentiation of cells during this period (for example, surfactin secretion, which is necessary to exit the lag phase, is controlled by quorum sensing), it was not possible to continuously acquire samples from the same swarm during the lag phase. Instead, several swarm plates were inoculated and incubated at the same time, and any given plate was only used for sampling once. To achieve a high reproducibility of lag phase development across the different plates, swarms were inoculated on single-well plates using the Echo 525 Acoustic liquid handler (Beckman Coulter, 001-10080) using an inoculation volume of 50 nl.

To collect cells from these plates, a small volume of lysis buffer previously kept on ice (3–4 µl) was pipetted onto the swarm surface and immediately retrieved by pipetting back up. This process of up-and-down pipetting was repeated 3 times to increase the number of collected cells. To increase the number of collected cells further, the process was repeated for several swarm colonies, and samples from the same timepoint were pooled in a final volume of 50 µl of lysis buffer, which was then immediately snap-frozen in liquid nitrogen, followed by storage at −80 °C until RNA isolation.

RNA isolation and sequencing

Our experiments yielded two types of sample, which were snap-frozen in liquid nitrogen and then stored at −80 °C until RNA isolation. Sample type A: cells from the spatiotemporal sampling experiments during the swarm expansion phase, for which the cells are located on a sampling tip submerged in 50 µl of lysis buffer within an Eppendorf tube. Sample type B: cells from the temporal sampling experiments during the lag phase, for which the cells were directly collected in 50 µl of lysis buffer within an Eppendorf tube. In both cases, samples were thawed at room temperature, followed by incubation at room temperature for 5 min, interrupted every 1 min by vortexing to ensure complete lysis of cells. Then, the tip was removed and total RNA was extracted from this lysate using the hot SDS/hot phenol method51 with some modifications as follows. To the lysate, we added 6 µl of 1 M sodium acetate (pH 5.5, Sigma, S7899) and 62.5 µl of Roti-Aqua-Phenol (Carl Roth, A980) and incubated the mixture at 65 °C for 8 min. The whole mixture was transferred to a phase lock gel tube (VWR, 733-2478), followed by the addition of 62.5 µl chloroform (Sigma, C2432). The mixture was centrifuged at 21,130 × g for 15 min at 12 °C. The aqueous phase (~65–70 µl) was transferred to a new 0.2 ml PCR tube. RNA in the solution was purified by adding 120 µl of Agencourt RNAClean XP kit (Beckman Coulter, A63987). Samples were then treated with TURBO DNase (Thermo Fisher, AM2238) and the total RNA quality was analysed with a TapeStation 4150 (Agilent, G2992AA). According to the TapeStation results, several RNA samples with relatively high concentrations were diluted to ensure that a consistent amount of RNA was used for further processing to minimize sample-to-sample bias. For ribosomal RNA (rRNA) depletion in the total RNA, we used the ‘do-it-yourself’ method52 with reduced reaction volume, followed by purification of rRNA-depleted RNA using Agencourt RNAClean XP kit (Beckman Coulter, A63987). Sequencing library preparation was performed using NEBNext Ultra II Directional RNA Library Prep with sample purification beads (NEB, E7765S). Sequencing for swarming cells was carried out at the Max Planck Genome Centre (Cologne, Germany) using an Illumina HiSeq 3000 with 150 bp single reads (aiming for >5 M reads per library). Sequencing for lag-phase cells before swarm expansion was carried out at the Basel Genomics Facility (Basel, Switzerland) using an Illumina NovaSeq 6000 with 101-bp single reads (aiming for >5 M reads per library).

Extracellular metabolite measurements during swarm development

To measure extracellular metabolites in the agar underneath the cells, swarms were allowed to develop to the desired diameter. Then, cells were removed from the agar by gently scratching the surface with a razor blade. At specific radii, a biopsy of the agar was acquired (the full depth of the agar) to acquire ~20 mg of agar, which was placed in an Eppendorf tube. Samples were weighed using a fine scale, and the appropriate amount of metabolite extraction solution (10 µl mg−1) containing 50% TE buffer (10 mM Tris, adjusted to pH 7.0 with HCl (Thermo Fisher, AM9850G), 1 mM EDTA (Thermo Fisher, AM9261)) and 50% methanol was added. Eppendorf tubes were then kept shaking at 4 °C for 2 h and afterwards centrifuged at 4 °C for 10 min at 21,130 × g. The liquid phase was then filtered using a 0.22 µm filter (regenerated cellulose, 4 mm diameter, Phenomenex, AF0-3203-52) and stored at −20 °C until further processing by mass spectrometry. In addition to the dilution of the samples by 1:10 resulting from adding the extraction solution, the extracts were diluted further by 1:20 for amino acid measurements to optimize the detection by the instruments. For organic acids, the extracts were used without further dilution.

Amino acids

Quantitative determination of amino acids was performed using HRES LC–MS. The chromatographic separation was performed on an Agilent Infinity II 1260 HPLC system using a ZicHILIC SeQuant column (150 × 2.1 mm, 5 μm particle size, 100 Å pore size) connected to a ZicHILIC guard column (20 × 2.1 mm, 5 μm particle size) (Merck KgAA). We used a constant flow rate of 0.3 ml min−1, with mobile phase A being 0.1% formic acid in a 99:1 mixture of water:acetonitrile (Honeywell) and phase B being 0.1% formic acid in a 99:1 mixture of acetonitrile:water (Honeywell) at 25 °C. The injection volume was 1 µl. The profile of the mobile phase consisted of the following steps and linear gradients: 0–8 min from 80% to 60% B; 8–10 min from 60% to 10% B; 10–12 min constant at 10% B; 12–12.1 min from 10% to 80% B; 12.1–14 min constant at 80% B. An Agilent 6470A mass spectrometer was used in positive mode with an electrospray ionization source and the following conditions: ESI spray voltage 4,500 V, nozzle voltage 1,500 V, sheath gas 400 °C at 12 l min−1, nebulizer pressure 30 psig and drying gas 250 °C at 11 l min−1. Compounds were identified on the basis of their mass transition and retention time compared to standards. Chromatograms were integrated using MassHunter software (Agilent). Absolute concentrations were calculated on the basis of an external calibration curve prepared in the sample matrix. To mimic the sample matrix, an aliquot of freshly prepared agar was treated similar to the extraction performed on samples for exometabolome determination. Mass transitions, collision energies, cell accelerator voltages and dwell times were optimized using chemically pure standards. Parameter settings of all targets are given in Supplementary Table 1.

Organic acids

Quantitative determination of organic acids was performed using LC–MS/MS. The chromatographic separation was performed on an Agilent Infinity II 1290 HPLC system using a Kinetex EVO C18 column (150 × 2.1 mm, 3 μm particle size, 100 Å pore size, Phenomenex) connected to a guard column of similar specificity (20 × 2.1 mm, 3 μm particle size, Phenomoenex). We used a constant flow rate of 0.2 ml min−1, with mobile phase A being 0.1% formic acid in water and phase B being 0.1% formic acid in methanol (Honeywell) at 25 °C. The injection volume was 0.5 µl. The profile of the mobile phase consisted of the following steps and linear gradients: 0–2.5 min constant at 0% B; 2.5–6 min from 0% to 100% B; 6–8 min constant at 100% B; 8–8.1 min from 100% to 0% B; 8.1–12 min constant at 0% B. An Agilent 6495 ion funnel mass spectrometer was used in negative mode with an electrospray ionization source and the following conditions: ESI spray voltage 2,000 V, nozzle voltage 500 V, sheath gas 260 °C at 10 l min−1, nebulizer pressure 35 psig and drying gas 100 °C at 13 l min−1. Compounds were identified on the basis of their mass transition and retention time compared to standards. Chromatograms were integrated using MassHunter software (Agilent). Absolute concentrations were calculated on the basis of an external calibration curve prepared in the sample matrix. To mimic the sample matrix, an aliquot of freshly prepared agar was treated similar to the extraction performed for exometabolome determination. Mass transitions, collision energies, cell accelerator voltages and dwell times were optimized using chemically pure standards. Parameter settings of all targets are given in Supplementary Table 1.

Transcriptome data analysis

For each of the FASTQ files from the 284 samples of the swarm expansion phase (3 replicates: replicate 1 contained 96 samples, replicate 2 contained 92 samples, replicate 3 contained 96 samples), we performed read trimming using Trimmomatic (v.0.39)53 and mapped the trimmed reads to the B. subtilis NCIB 3610 reference genome and plasmid (NCBI accession number: NZ_CP020102 and NZ_CP020103) using HISAT2 (v.2.2.1)54 in single-end mode. Reads mapped to all functional genes were counted using featureCounts (v.2.0.1)55 with fractional counting for multimapping and multi-overlapping reads (−OM−−fraction) in a strand-specific manner (−s 2). Any short transcripts such as small RNA or transfer RNA should not be detected because the above-described RNA purification steps using RNAClean XP kit can only retain transcripts bigger than 200 nt. Reads that were mapped to protein-coding sequences, non-coding RNAs, transfer-messenger RNA, signal recognition particle RNA and ribonuclease P RNA (in total 4,342 genes) were used for downstream analyses.

Raw read counts of the expansion phase samples were loaded into R and filtered by keeping all genes for which there were at least 2 samples with a read count of at least 10, leaving 3,932 genes out of 4,342 genes. After further removing samples with total read counts of <106 mapped to those genes, there were 278 samples left: 95 samples in replicate 1, 92 samples in replicate 2, and 91 samples in replicate 3, which were analysed further. We then generated a DGEList object containing all 278 transcriptomes and applied the TMM normalization method implemented in edgeR (v.3.26.8)56,57 for samples pooled from all three experiments, using the calcNormFactors function to enable sample-to-sample comparison of the data. All subsequent analysis was then performed with normalized log2 values.

For Supplementary Fig. 21, we performed the same processing for the FASTQ files of the 33 lag phase samples and then performed the normalization using both the 278 expansion phase and 33 lag phase samples, following the same criteria. This resulted in a final dataset comprising 4,083 genes across 311 samples. The normalization in Supplementary Fig. 21 therefore differs from that used for all other figures.

Spectral representation of spatiotemporal gene expression

The transcriptome measurements and the microscopy-based measurement of phenotypic properties were sampled at a set of radial space–time points \(\}}_}\left(_,_\right)}_^\), where \(}}_\) is the time and \(_\) is the radial position from the centre of the swarm at which the sample was acquired. For each gene, we have a sample vector \(}}_\) with length \(L\), where the \(l\) th entry is the gene expression at the point \(}}_}}\). Similarly, for each phenotypic property, we have a sample vector \(}}_\) with length \(L\).

To form a spectral representation across the three replicates, we fit a common domain to the three replicates. Experimentally, the radial position of the boundary \(_\) of the swarm at each time \(_\) was determined automatically by detecting the presence and location of bacteria in the microscopy field of view and moving the microscope stage until the field of view was split between colonized agar containing bacteria and uncolonized agar in approximately equal proportions. We simultaneously fitted a boundary of the form \(b\left(t\right)=_\exp \left(t/\tau \right)\) to all replicates by first minimizing the loss function,

$$}\left(},_^,_^,_^\right)=\mathop\limits_^\mathop\limits_^^}_^-_^\exp \left(_^/}\right)\right)}^$$

(1)

where superscript (n) denotes the index of the three different replicates. This exponential fit approximates the experimental data very well (Extended Data Fig. 9a).

To obtain non-dimensional data, we rescaled data as follows. Let \(r\) be the index corresponding to the largest \(_^\); we defined the time shift \(_^=_^\) and scaled the initial value \(_=_^\exp \left(_^/\tau \right)\). The other time shifts are given by

$$_^=}\left[\log \left(_^\right)-\log \left(_^\right)\right]+_^.$$

(2)

We then non-dimensionalized the data (\(\widetilde\) variables) using \(}_=\left(}}_-_^\right)/}\) and \(}_=_/}}_\). The domain boundary is then given by \(0\le \widetilde\le T\) and \(0\le \widetilde\le \exp \left(\widetilde\right)\) where \(T\) is maximum non-dimensional time present in all three replicates. Data points that lay outside the domain after rescaling were not used in the spectral representation. The non-dimensionalized domain and how the individual sampling data points are distributed within this domain are shown in Extended Data Fig. 9b.

We built a domain-specific orthogonal polynomial basis \(}}_\left(\widetilde,\widetilde\right)}_^\) by applying Gram–Schmidt orthogonalization58 to the monomial set,

$$\,\widetilde,}^,\widetilde\widetilde,}^,\ldots \},$$

(3)

under the inner product

$$\left\langle f,g\right\rangle =_^}\text\widetilde^\right)}}_\text\widetilde^}.$$

(4)

The space–time dependence of each gene was compressed by expanding each sample vector,

$$}}_=\mathop\limits_^}}_}}_}}$$

(5)

where \(}}_}}\) is the length \(}\) vector formed by evaluating \(}}_}}\left(\widetilde,\widetilde\right)\) at each non-dimensional space–time point \((}_,}_)\). The coefficients were fitted using least squares on the matrix equation,

$$}}_=\left[\begin}}_ & }}_ & \cdots & }}_\end\right]\left[\begin}}_}}\\ }}_}}\\ \vdots \\ }}_},}}\end\right]$$

(6)

for each gene and property.

To combine information from all replicates, we fitted a single coefficient vector for the replicates by stacking the least square problems on top of each other to form a single linear regression problem,

$$\left[\begin}}_^}}}}\\ }}_^}}}}\\ }}_^}}}}\end\right]=\left[\begin\begin}}_^ & }}_^ & \cdots & }}_^\end\\ \begin}}_^ & }}_^ & \cdots & }}_^\end\\ \begin}}_^ & }}_^ & \cdots & }}_^\end\end\right]\left[\begin}}}_}}\\ }}}_}}\\ \vdots \\ }}}_},}}\end\right],$$

(7)

from which a smooth average spatiotemporal gene expression was formed,

$$}}}_(\widetilde,\widetilde)=\mathop\limits_^}}}_}}_}}(\widetilde,\widetilde).$$

(8)

The same procedure was also used to produce average properties \(}_\).

Spatiotemporal pattern identification

The coefficients \(_\) for the basis functions \(}_}}\) encode information about the spatiotemporal expression pattern of the genes and phenotypic properties, which allowed us to cluster genes on the basis of their spatiotemporal expression pattern using the spectral coefficients \(_\) directly. With the clustering analysis, we intended to identify the underlying patterns independent of global shifts and scaling. Therefore, we defined coefficients with the mean subtracted and scaled by the standard deviation,

$$_=\frac_-_\mu }_/_}_}$$

(9)

where \(_\) is the standard deviation of the expression of gene \(n\), \(_\) is the mean of the expression of gene \(n\), \(_\) is the Kronecker delta that is 0 when \(n\ne m\) and 1 when \(n=m\), and \(_\) is constant since it is a degree 0 polynomial. Note that under this rescaling, \(_\) is no longer an independent parameter since it is fully determined by the means of the higher-order polynomials and their respective coefficients. Using these rescaled coefficients \(_\), we defined a score for how strongly patterned a gene expression profile is. A spatiotemporally patterned expression profile should have two components:

The pattern should vary smoothly so that the pattern should be well approximated by the spectral representation, which means that the scaled representation error

$$}}_}}=}}_}}-_}_}-\mathop\limits_^_}}_}}\right\Vert}^$$

(10)

should be small.

The pattern should not just be constant, meaning that the higher coefficients should be important, which also means that the pattern score

$$}}_}}=\mathop\limits_^_^$$

(11)

should be large.

We therefore defined a space–time ranking,

$$}}_}}=\frac}}_}}}}}_}}}$$

(12)

which is large when a gene expression has a strong spatiotemporal pattern, and low when there is no spatiotemporal pattern. For identifying different types of spatiotemporal pattern on the basis of the spectral coefficients \(_\), we kept only those genes with a high ranking \(}}_}}\). We defined the cut-off for which genes were designated as displaying a spatiotemporal pattern by ordering the genes by their ranking (smallest to largest) and then finding the largest integer \(_\) such that,

$$\frac\nolimits_^_}}_}}}\nolimits_^}_}}}\le 0.5.$$

(13)

For the genes that displayed a spatiotemporal pattern on the basis of the above criterion, we calculated the cosine similarity between the coefficients,

$$_}=\frac\nolimits_^__}\limits_^_^}\sqrt\limits_^_^}}.$$

(14)

We used cosine distances since we needed a metric that is independent of a global scaling to the expression level. The distance matrix \(D=\left(_}\right)\) was then separated into clusters using the k-medoids algorithm59 to produce k distinct patterns. The choice of the number of clusters was based on a plot of the total cost versus number of clusters and choosing a value in the elbow of the curve (Extended Data Fig. 10). From each cluster, we chose the highest ranked gene as the most representative spatiotemporal pattern.

Image analysis

To perform single-cell segmentation on the brightfield microscopy images of the swarm, we used Stardist60. Bacterial movement fields were calculated using the Horn–Schunck optical flow method61, applied to consecutive images of a video. From the segmentation and bacterial movement fields, all phenotypic properties listed in Supplementary Table 2 were calculated. For properties calculated for each individual cell, median values were used across one image field of view, and mean values across all 48 frames in a video were used for visualization in heat maps. The following paragraphs describe the definitions of emergent property parameters.

Nematic order

The nematic order parameter \(S\left(\vartheta \right)\) for two cells with a relative angle \(\vartheta\) was quantified as

$$S\left(\vartheta \right)=1.5\cdot ^-0.5.$$

(15)

The local nematic order parameter for a specific cell in a swarm was then defined as the mean of the nematic order parameter of the cell with each of its neighbours within a centroid–centroid distance of ≤10 µm.

Non-motile clusters

Non-motile cells were identified by thresholding of the bacterial movement field with a cut-off of 8 µm s−1 and by hierarchical clustering based on centroid–centroid distance with a cut-off of 2 µm. Only clusters with 10 or more cells were considered.

Rafts

The neighbourhood of each cell was defined as all cells within a centroid–centroid distance of 30 µm to the cell under investigation. To quantify rafting behaviour, we counted the number of motile cells (speed of 10 µm s−1 or more) in this neighbourhood that share the same orientation as the cell under investigation up to a tolerance of 15°. This number was then divided by the total number of cells in the neighbourhood. This ratio, called the ‘local rafting ratio’, was calculated for each cell and used as a measure for local rafting activity.

Density fluctuations

To calculate density fluctuations, the image was split into sub-images of size 120 × 120 pixels (~48 × 48 µm2) and for each sub-image, the density was calculated as the fraction of the area covered by cells to the entire area. For the density fluctuations in space, the standard deviation between the local densities in sub-images was calculated for each timepoint. The final value was taken to be the mean across all timepoints.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.