Predicting global dynamics of spatial microbial communities from local interaction rules

Interactions between cells drive biological processes across all of life, from microbes in the environment to cells in multicellular organisms. Interactions often arise in spatially structured settings, where cells mostly interact with their neighbors. A central question is how the properties of biological systems emerge from local interactions. This question is very relevant in the context of microbial communities, such as biofilms, where cells live close by in space and are connected via a dense network of biochemical interactions. To understand and control the functioning of these communities, it is essential to uncover how community-level properties, such as the community composition, spatial arrangement, and growth rate, arise from these interactions. Here, we develop a mathematical framework that can predict community-level properties from the molecular mechanisms underlying the cell-cell interactions for systems consisting of two cell types. Our predictions match quantitative measurements from an experimental cross-feeding community. For these cross-feeding communities, the community growth rate is reduced when cells interact only with few neighbors; as a result, some communities can co-exist in a well-mixed system, but not in a spatial one. In general, our framework shows that key molecular parameters underlying the cell-cell interactions (e.g. the uptake and leakage rates of molecules) determine community level properties. Our framework can be extended to a variety of systems of two interacting cell types, within and beyond the microbial world, and contributes to our understanding of how biological functions arise from interactions between single cells.


Introduction
Biological interactions pervade all of life. Interactions at lower levels of organization can confer new functionality at higher levels. For example, interactions between different cell types determine the functioning of organs and tissues in multicellular organisms and interactions between different species determine the processes an ecosystem performs [1,2]. In natural systems, interactions often arise in spatially structured settings, where individual entities interact mostly with others that are close by in space [3]. When interactions are local, the spatial organization of the different entities defines their network of interaction. A central question is how the properties of biological systems emerge from this network of interactions. This has primarily been studied

Results
We present a mathematical framework that predicts community-level properties from local interaction rules for a variety of systems consisting of two interacting cell types. First we apply it to microbial cross-feeding communities and compare its predictions to experimental observations. Subsequently, we derive general predictions for other systems.
To address our first goal, we build on previous work were we measured the local interaction rules for a cross-feeding community [11]. This community consists of two bacterial cell types that exchange essential amino-acids: one cell type (∆P ) is unable to produce proline, while the other (∆T ) is unable to produce tryptophan. Amino acids are passively released into the extracellular environment and cells in this community can thus grow by exchanging proline and tryptophan [17,[26][27][28][29]. We found that cells receive amino acids from partner cells residing within a finite interaction range; the growth rate of a cell increases linearly with the frequency of the partner type within this interaction range [11]. The local interaction rules can be fully described by the interaction range and the maximum growth rate that a cells achieves when all its neighbors are of the partner type (Fig. 1B).  Fig. 1: A mathematical framework to scale up from molecular-level properties, to individual-level properties, to community-level properties. We previously measured the local-interaction rules for a crossfeeding community (B) and showed that these can be derived (arrow 1) from the molecular mechanisms of the interaction (A). Here we developed a mathematical model that can predict community-level properties (C) either from measured local rules (arrow 2) or directly from the underlying molecular mechanisms (arrow 3). (A) The community consists of two types of Escherichia coli: ∆P is unable to produce the amino acid proline and ∆T is unable to produce the amino acid tryptophan. Cells exchange amino acids with the environment through active uptake (with rate r u ) and passive leakage (with rate r l ). Amino acids are exchanged between cells through diffusion in the environment (with rate D). All rates differ between the two amino acids. (B) Local interaction rules can be fully described by two fundamental quantities: the size of the interaction neighborhood (r, left); and the growth function of a cell (charecterized byμ, right). Each dot corresponds to the measured growth rate of a single cell, n = 2610 for ∆P and n = 2162 for ∆T, the line shows the result of a linear regression, data reproduced from [11]. (C) Our framework makes predictions for the community-level properties, such as the equilibrium frequency of the two types, their spatial arrangement, and growth rate.
We also showed that the local interaction rules can be derived from key biophysical parameters of the underlying amino acid exchange (arrow 1 in Fig. 1). The interaction range primarily depends on the cell density and on the ratio between the uptake and diffusion rate of the molecules (Methods, Eq. 4). The maximum growth rate primarily depends on the leakage rate of the molecules (Eq. 6). Uptake, leakage, and diffusion rates are molecule specific; as a result the interaction range and maximum growth rate will generally differ between cell types, as we also observed for our specific community (Fig. 1B). Our goal here is to use these molecular parameters, which determine the local interaction rules, to predict community-level properties (arrow 3 in Fig. 1).

A mathematical framework for cross-feeding communities
Our framework builds on previous work using evolutionary graph theory [30][31][32]. We consider a spatial system of two cell types, A and B, that exchange essential cellular building blocks. We assume that the total number of cells is constant and that each cell interacts with a fixed number of neighbors. In contrast to previous work, we allow the two cell types to interact with a different number of neighbors, as suggested by our experimental observations [11]. We assume that type A interacts with r A neighbors and that type B interacts with r B neighbors. We call this quantity the neighborhood size and set r B ≥ r A (Fig. 2B). To allow for different neighborhood sizes, we need to use directed graphs, instead of the undirected graphs used in previous work. As a result, several simplifications no longer hold and we discuss in the SI the required extensions to previous work.
Each cell (of type A or B) interacts with all (r A or r B ) cells in its interaction neighborhood and these interactions determine its growth rate. We assume that the growth rate of a cell increases linearly with the frequency of the partner type within its interaction neighborhood Pair approximation allows for a quantitative description of spatial systems. (A) Pair approximation assumes that a spatial system (left side) can be fully described by tracking the number of all pairwise links (right side). NX←Y indicates that focal cell X interacts with its neighbor Y. (B-D) A system is fully described by the local rules, as specified by: (B) the interaction neighborhoods of both types, characterized by the neighborhood sizes rA and rB; (C) the replication neighborhood, assumed to be identical to the smallest interaction neighborhood; and (D) the growth functions of both types, characterized by the maximum growth ratesμA andμB. (E) Pair approximation has three dynamical variables that describe the global composition, P (A), and local composition P (B|A, rA) and P (A|B, rB) of the system.
( Fig. 2D), as observed in the experimental communities [11]. Each type achieves its maximum growth rate,μ, when it is completely surrounded by the partner type, and its value is different for the two types (Fig. 2D). When a cell divides, it replaces a random neighbor within the replication neighborhood. We assume that the replication neighborhood is identical to the smaller interaction neighborhood (r R = r A , Fig. 2C). This assumption has its limitations: in reality neighboring cells are not replaced, but rather pushed aside. Moreover, the replication neighborhood could be different from the interaction neighborhood. However, this assumption does capture one of the most essential features of real biological systems: that cells place their offspring close by in space and thus become surrounded by their own type (see SI for a more detailed discussion). We implemented our model in two complementary ways: we used a cellular automaton to simulate the dynamics numerically and we used pair approximation to derive analytical predictions. Pair approximation assumes that the dynamics of a spatial system can be described by only specifying how often the different cell types are found in each others neighborhoods ( Fig.  2A) [30,33,34]. In our case, this means that the outcome of an interaction between a focal individual and one of its neighbors depends only on the identity of these two cells and not on the wider context in which these two cells are found. This approximation allows us to describe the dynamics of a spatial community (Fig 2A left) by only tracking how often cells of the different types are found within each others interaction neighborhood, i.e. by tracking the number of A←A, A←B, B←A, and B←B links (Fig 2A right). In our convention, X←Y means the focal cell X interacts with Y, i.e. that the focal cell X receives building blocks from Y. These links are directional: because the two types have different interaction ranges, cell X might receive building block from cell Y, even when the inverse is not true. This is the primary difference with previous work, where interactions were assumed to be symmetric (and links undirected) [30].
With these assumptions, we obtain the dynamical equations for the number of pairwise links in time (see SI). Because the number of cells is constant, we only have three independent variables (e.g. the number of A←A, A←B, B←A links) and we can express the fourth one (the number of B←B links) as a function of the other three. Instead of tracking the number of links directly, it is convenient to change variables and track the following dynamical variables: P (A), the global frequency of A, P (B|A, r A ), the local frequency of type B within the interaction neighborhood of type A, and P (A|B, r B ), the local frequency of type A within the interaction neighborhood of type B (Fig. 2E). The first quantity describes the global composition of the community, while the last two quantities describe the average local composition of a cells' neighborhood. The average growth rate of a cell in a spatial system can be directly calculated when we know the average local composition, i.e. the average frequency of the partner type within its interaction neighborhood.

Predicted community properties at steady state
Pair approximation allows us to predict global properties of the community from the local interaction rules between cells. We input the neighborhood size of each type (r A and r B ) and the growth functions characterized by the maximum growth rate obtained when a cell is completely surrounded by the other type (μ A andμ B , Fig. 2 B and D). As output we obtain a system of dynamical equations describing the global, P (A), and local, P (B|A, r A ) and P (A|B, r B ), composition of the system (Fig. 2E). By solving the dynamical equation for steady state, we can find analytical expressions for the global and local compositions at equilibrium (see SI).
Pair approximation predicts that the global composition of the community reaches a steady state in which the frequency of type A is given by: When the neighborhood size is large (r A , r B ≫ 1), this equation simplifies to the predicted frequency in a well-mixed system (see SI): The same global community composition can correspond to many different spatial arrangements of the two cell types. For example, the two cell types could be highly mixed in space or completely segregated. The spatial arrangement matters because cells interact only locally: the growth of a cell depends on the local frequency of the partner type, P (B|A, r A ) and P (A|B, r B ), which can be different from the global frequency, P (B) and P (A), when the cell types are clustered in space. A relevant quantity to describe a spatial system is thus the ratio between the The difference between the model prediction and data is less than 13%. (F) The average growth rate of the community is reduced due to cell clustering. An in-silico (see Methods) shows that the growth rate in clustered communities, with experimentally observed spatial arrangements, is reduced by a factor of 0.87 (CI: 0.85-0.90) compared to randomized communities, where cell clusters were disrupted (Data reproduced from [11]). Pair approximation predicts a decrease by a factor of 0.94 (SI Eq. [22][23], the difference between the model prediction and data is less than 6%.
local and global frequency of the partner type, and pair approximation predicts that at steady state: where P (B) = 1 − P (A). This equation shows that the frequency of the partner type in the interaction neighborhood is lower than one would expect from the global composition of the community. This happens because cells place offspring close to themselves when they divide. When the neighborhood size is small (e.g., r = 3) this reduction in the local frequency is as much as 50%, but when it is large (r ≫ 1), the local frequency approaches the global frequency.

Predicting global properties of an experimental microbial community
We validate our mathematical framework by comparing its predictions to experimental observations of our cross-feeding bacterial community. We parameterize our model directly from the biophysical parameters describing the amino acid exchange (arrow 3 in Fig. 1, Table S2). The model can be fully parameterized by specifying the size of the interaction neighborhood (r), which mostly depends on the uptake and diffusion rates (Eq. 4), and the maximum growth rate of a cell (μ), which mostly depends on the leakage rate (Eq. 6, see Methods). Alternatively, we can parameterize the model directly from the measured interaction range and measured maximum growth rate (arrow 2 in Fig. 1), the predictions obtained from both parameterization match closely (Table S3). Our model recapitulates the experimentally observed community-level dynamics. In our experiments, the frequency of the ∆T type converges over time to a stable frequency across 22 communities (Fig. 3A). Our model predicts similar dynamics, as we demonstrate by numerically solving the dynamical equations given by pair approximation and by simulating the system using a cellular automaton. Both pair approximation (Fig. 3B) and simulations (Fig. 3C) show that the community reaches a stable equilibrium. From Eq. 1 we find a predicted equilibrium frequency of ∆T of 0.20, which matches well with both the experimentally observed value of 0.19 and cellular automaton simulations value of 0.22 (Fig. 3D).
Cells in microbial communities tend to be surrounded by cells of their own type because when a cell divides the two new cells remain close in space. Therefore, the local frequency of the partner type found in the interaction neighborhood of a focal cell is typically lower than the global frequency in the whole community. Our model predicts (Eq. 3) that the local frequency of the partner type is reduced by a factor of 0.99 for ∆P and 0.89 for ∆T , compared to the global frequency. We measured this reduction in the frequency for the experimental communities and found values of 0.87 for ∆P and 0.85 for ∆T . Qualitatively the experimental observations thus match our model predictions: cells tend to be surrounded by their own type, thus reducing the local partner frequency.
For cross-feeding communities, cells in a spatially structured community grow slower than cells in an equivalent well-mixed community, because cells become surrounded by their own type. Our model predicts that spatial structure decreases the average growth rate in our experimental community by a factor of 0.92 (SI Eq. 22-23). Experimentally we cannot directly compare the growth in spatially structured communities with those in well-mixed communities. However, we previously used an individual-based biophysical model to perform an equivalent in-silico experiment: we compared the predicted growth rate for clustered communities, in which the spatial arrangement of cells was based on experimental measurements, with randomized communities, in which the spatial arrangement of cells was randomized [11]. We found that the growth rate in the clustered communities was reduced by a factor of 0.87 compared to the randomized communities, a value that is comparable with the value of 0.94 predicted by pair approximation.

Global predictions from local interactions
We now investigate what general predictions we can make for cross-feeding communities if we know the local interaction rules. Our model predicts that the equilibrium composition of the community is mostly set by the maximum growth rates of the two types. In general, the type with highest maximum growth rate constitutes the majority (Fig. 4A). When the neighborhood size is large (r ≫ 1), the equilibrium frequency approaches that of a well-mixed system, while if the neighborhood size is small (r = 3), the equilibrium frequency is shifted more to the type that grows faster (Fig. 4A).
When both cell types have the same maximum growth rate, the ratio of neighborhood sizes ( r B r A ) does not affect the equilibrium state of the community (Fig. 4B, black line). When the the two cell types have different growth rates, increasing the neighborhood size of even a single type moves the equilibrium frequency closer to the expected frequency in a well-mixed system (Fig.  4B, colored lines). , asymmetric communities, where the cell types have different maximum growth rates, grow poorly. When the asymmetry is too large communities cannot grow in a spatially structured environment even though they could grow in a well-mixed one. When the neighborhood size is large (r = 30), the growth rate of the spatial community is close to that of the well-mixed one even when there is an asymmetry. (B) Increasing the neighborhood size of even one of the types can increase the growth rate of spatial communities and prevent collapse. rA is held constant at 3, while rB is increased.
The neighborhood size strongly alters the local composition of a cells' neighborhood; as the neighborhood size becomes smaller, the partner fraction in each cells' interaction neighborhood decreases compared to the global fraction (Fig. 4C). As a result the average growth rate of the community decreases relative to an equivalent well-mixed community (Fig. 4D).

Small interaction neighborhoods can lead to community collapse
We have shown that having a small interaction neighborhood reduces the growth rate of mutualistic communities. Here we show that this growth reduction is more pronounced when the two cell types have different maximum growth rates. We call such communities asymmetric. Our model predicts that the more asymmetric a mutualistic community is, the more their growth is hindered by small neighborhood sizes, to the point that the community collapses when neighborhood sizes are very small (Fig. 5A). Increasing the neighborhood size of either one (Fig. 5B), or both (Fig. 5A) of the types can prevent the collapse of the community. This finding shows that there is a limit to the stability of mutualistic interactions in spatially structured communities: mutualistic types with different maximum growth rates that can stably coexist in well-mixed environments might not be able to survive in spatially structured environments. For cross-feeding communities, spatial structure can thus have a detrimental effect.

Spatial effects for other types of communities
So far, we have focussed on cross-feeding systems were the growth rate increases with the frequency of the other type, however for other systems we expect different growth rate functions to be of relevance. As long as the growth rate depends linearly on the neighborhood composition, we can find an analytical solution for the steady state community properties (for non-linear functions the model can be solved numerically, see SI). The global composition of the community at steady state (P (A)) depends strongly on the chosen growth functions (SI, Eq. 28). However, the relation between the local and global composition of the community is always the same: i.e. for all growth functions we find that equation 3 holds (see SI).
Here we illustrate results for two different types of system. First, we consider a system where a type's growth rate increases with the number, rather than the frequency, of partner type within the interaction neighborhood (Fig. 6). In this case, the type with the highest product ofμ · r reaches the higher global frequency (see SI). Therefore, the global composition of the community depends strongly both on the maximum growth rate and on the neighborhood size (Fig. 6).
Second, we consider systems where two cell types inhibit each other's growth. In such systems, a cell's growth rate decreases (linearly) with the frequency of the other type within the interaction neighborhood. We find that cells that inhibit each other, grow slower when they are well-mixed compared to when they interact within a small neighborhood (see SI). This happens because cells typically are surrounded by their own type, and thus have reduced interactions with the growth inhibiting cells.

Discussion
We developed a mathematical framework that can make quantitative predictions for the dynamics of spatially structured communities consisting of two interacting cell types (Fig. 1). We primarily focussed on microbial cross-feeding systems and used pair approximation to derive analytical predictions for the community level properties from knowledge of the local interaction rules (Fig. 2). These rules are unequivocally defined by two fundamental quantities: the size of the interaction neighborhood and the maximum growth rate that cells achieve when they are completely surrounded by the partner type; both quantities typically differ for the two cell types. We furthermore showed that these local interaction rules can be directly derived from key biophysical parameters of the underlying molecular mechanisms (Fig. 1). We validated this framework using an experimental community of two cross-feeding bacteria and showed that our model reproduces the experimentally observed community level properties (Fig 3). We thus elucidated how the local interaction rules arise from the molecular exchanges, and how the local interaction rules scale up to determine community properties. Generally, we are thus able to scale up from molecular-level properties, to individual-level properties, to community level properties.
For spatial cross-feeding communities we found that the local and global properties can change independently. The global composition of the community (i.e. the equilibrium frequency of the two types) is set by the ratio of the maximum growth rates, which mostly depends on the leakage rates of the exchanged metabolites. In contrast, the composition of the local neighborhood, the one that matters most for the cells, is set by the neighborhood size, which mostly depends on the ratio of the uptake and diffusion rates. Different biophysical parameters thus control the global and local properties of the community.
Small neighborhood sizes reduce the frequency of the partner type around each cell. This happens because cells place their offspring close-by in space. In our model we assumed that cells cannot actively move. However, cells could overcome the negative effect of having a small neighborhood size if they could actively move to locations with a higher frequency of the partner type. Without such active movement, cross-feeding cells in a spatial system grow slower than in equivalent well-mixed system (Fig. 4). This reduction in growth rate can even lead to the collapse of the community: cells in asymmetric communities (i.e. communities where the two type have different maximum growth rates) can stably coexist in well-mixed system, yet they might not coexist in a spatial system, if they interact at small ranges (Fig. 5).
Previous work suggested that short-range interactions provide an advantage to cooperative interactions, because they separate cooperative types from non-cooperators [24,[35][36][37]. In general, when it is good to be surrounded by your own type, we expect short-range interactions to be beneficial, as for example when two cell types inhibit each other. When it is good to be surrounded by the other type, we expect short-range interactions to be detrimental, as for example when two cell types exchange beneficial resources. When a community consists of more than two cell types the situation can be more complex. Previous studies have shown that short-range interaction can be beneficial for cross-feeding communities that contain non-producing cells: the short-range interactions harm the producing cells, but they harm the non-producers even more, and thus prevent them from taking over [20,38,39].
Our mathematical framework can model a variety of systems of two interacting cell types. Here we considered systems where a type's growth rate increases with the number, or the frequency, of the partner type within the interaction neighborhood, and systems where the two cell types inhibit each other's growth rate. Our framework can model systems beyond the bacterial world, as long as they are composed of two interacting types. For example, it could be adapted to model homeostasis in tissues composed of two different cell types that control each other's growth via growth factors [40,41].
Many biological systems are spatially structured and multi-scale: they consist of individual entities (e.g. cells or species) that interact with each other in space. The sum of these interactions determines the global properties of the system (e.g. multicellular organism or microbial community). To understand such biological systems, it is essential to scale between levels of organization. Our work provides a contribution to this effort by creating a mathematical framework that can scale from molecular mechanisms, to local interaction rules, to global system properties.

Model assumptions
• Two cell types, A and B, fully occupy a regular directed graph. Type A interacts with r A neighbors, type B interact with r B neighbors, where r B ≥ r A .
• A cell's growth rate increases linearly with the frequency of the other type within its interaction neighborhood.
• Individuals reproduce with a probability proportional to their growth rate, and their offspring replaces a random neighbor within the replication neighborhood.
• The replication neighborhood is identical to the smallest interaction neighborhood r R = r A .
We used pair approximation to derive analytical predictions and a cellular automaton to run simulations of the model.

Pair approximation
We track the number of all pairwise links N X←Y , where X, Y ∈ A, B. The total number of links changes through time, because the neighborhood sizes are different for the two cell types (see SI). This is in contrast to previous work were the two neighborhood sizes are the same and the total number of links remains constant. Two events can change the number of links: a type A cell reproduces and replaces a B neighbor, with rate T + , or a type B cell reproduces and replaces an A neighbor, with rate T − . During a T + event the number of A cells thus increases by one, and during a T − event it decreases by one. The rate T + is given by: where the first factor gives the probability of choosing a type A focal cell, the second the probability that this A cells has a type B neighbor, and the third the average growth rate of a type A cell that has at least 1 type B neighbor, relative to the average population growth rate ⟨µ⟩, which is given by: Similarly, for T − we have: The extra factor at the end is the probability that a type A neighbor that is in the (larger) interaction neighborhood is also part of the smaller replication neighborhood. We can express all probabilities as function of the number of links: where N is the total number of cells in the systems. When a T + or T − event happens, the number of X ← Y link changes by ∆ + XY and ∆ − XY , respectively. These changes can be fully expressed as function of N X←Y , as we show in the SI. We can then write the dynamical equation for N X←Y as: We can solve these equation numerically to obtain the temporal dynamics, and we can solve them analytically for the steady state, see SI for full derivation.

Model parametrization
To parametrize the model, we need to specify the size of the interaction neighborhood r and the maximum growth rate of a cellμ. We previously showed that that interaction-range (R) of a cell (the distance over which amino-acids can be exchanged) is given by [11]: where β and γ are constants, ρ is the volume fraction occupied by cells, and D, r u , and r l are the rates of diffusion, uptake, and leakage of the exchanged metabolite (see SI for details and parameter values). The number of neighbors can be directly calculated from the interaction range: where l and w are the average length and width of cells, and ρ 2D is the number of cells per area. We previously showed that that maximum growth rate of a cellμ is given by [11]:

Cellular automaton
In the cellular automaton cells are placed on a square grid of size 100 x 100 with periodic boundary conditions. Cells interact within an extended Moore neighborhood with range d X , and thus have a neighborhood size of (2d X + 1) 2 − 1. For each cell we calculate the growth rate (using linear growth function) and we randomly pick a cell to reproduce with a probability that is proportional to its growth rate. We then randomly replace one of its neighbors. The grid was initialized with a random arrangement with a given frequency of the two types. We parameterized the experimental community by choosing d ∆T = 1 (8 neighbors) and d ∆P = 5 (120 neighbors). Simulations are implemented in C++.

Experimental communities
The experimental methods are described in detail in reference [11], here we summarize the most relevant details. The community consists of two strains of Escherichia coli: ∆T: MG1655 trpC::frt, PR-sfGFP and ∆P: MG1655 proC::frt, PR-mCherry. The deletion of trpC and proC prevents the production of tryptophan and proline, respectively. Cells are labeled with constitutively expressed fluorescent proteins. Cells were grown in microfluidic chambers of 60x60x0.76µm, the small height forces cells to grow in a monolayer. The chambers open on one side into a feeding channel of 22µm high and 100µm wide. The microfluidic devices were fabricated from Polydimethylsiloxane (PDMS) using SU8 photoresist molds. Overnight cultures of the two strains (started from single colonies) were mixed in a 1:1 volume ratio and loaded into the microfluidic devices by pipette. Cells were grown on M9 medium (47.8mM Na 2 HPO 4 , 22.0 mM KH 2 PO 4 , 8.6 mM NaCl and 18.7 mM NH 4 Cl) supplemented with 1 mM MgSO 4 , 0.1 mM CaCl 2 , 0.2% glucose, and 0.1% Tween-20. For the first 10h the medium was supplemented with 434 mM of L-proline and 98 mM of L-tryptophan to allow cells to grow independently and fill the chambers. Subsequently, cells were grown without externally supplied amino acids.
The growth of the communities was followed using time-lapse microscopy. Phase contrast and fluorescent images were taken every 10min using fully automated Olympus IX81 inverted microscopes, equipped with a 100x NA1.3 oil objective, a 1.6x auxiliary magnification, and a Hamamatsu ORCA flash 4.0 v2 sCMOS camera. The sample was maintained at 37 • C with a microscope incubator.

Image analysis
Microscope images were analyzed using Vanellus (version v1.0 [42]). Images were registered, cropped to the area of the growth chambers, and deconvoluted. Images were segmented on the fluorescent channels using custom build Matlab routines [11] or on a combination of the phase contrast and fluorescent channels using Ilastik (version v1.3.3 [43]). Cell tracking was done using custom build Matlab routines [11] followed by manual correction.
Cell length and width were measured by fitting an ellipse to the cell shape and taking the major and minor axis length, respectively. Cell growth rates were estimated by fitting a linear regression to the log transformed cell lengths over a 40min time window [11].

Data analysis
For each chamber we estimated P (∆T ) as the number of pixels occupied by ∆T cells divided by the total number of pixels occupied by cells. We estimated P (∆T |∆P, r ∆P ) by first calculating the local ∆T frequency for each ∆P cell, and averaging this over all ∆P cells within the chamber.
The local frequency was calculated as the number of pixels occupied by ∆T cells divided by the total number of pixels occupied by cells, considering only pixels within 12.1µm (the interaction range of ∆P ) of the cell perimeter. P (∆P |∆T, r ∆T ) was estimated in similar way, using the interaction range of 3.2µm for ∆T.
The difference between the average growth rates of clustered and randomized communities was estimated using an in-silico experiment. For the clustered communities experimentally observed spatial arrangements of 22 communities were converted to a 40x40 square grid. Each grid site was assigned the cell type that occupied the majority of the corresponding pixels in the real image. Grid sites that remained empty after this procedure were randomly assigned one of the types. This grid was used as input for a previously described individual based model that can predict single cell growth rates [11]. The resulting growth rates were averaged over all cells to obtain the average growth rate for the clustered communities. Subsequently the 40x40 grids were randomized, keeping the frequency of the two cell types constant but permuting their locations. We calculated the average community growth rate for 20 independent randomizations, and we assigned the average over all these permutations as the growth rate of the randomized communities.