Systematic microscopical analysis reveals obligate synergy between extracellular matrix components during Bacillus subtilis colony biofilm development

Single-species bacterial colony biofilms often present recurring morphologies that are thought to be of benefit to the population of cells within and are known to be dependent on the self-produced extracellular matrix. However, much remains unknown in terms of the developmental process at the single cell level. Here, we design and implement systematic time-lapse imaging and quantitative analyses of the growth of Bacillus subtilis colony biofilms. We follow the development from the initial deposition of founding cells through to the formation of large-scale complex structures. Using the model biofilm strain NCIB 3610, we examine the movement dynamics of the growing biomass and compare them with those displayed by a suite of otherwise isogenic matrix-mutant strains. Correspondingly, we assess the impact of an incomplete matrix on biofilm morphologies and sessile growth rate. Our results indicate that radial expansion of colony biofilms results from the division of bacteria at the biofilm periphery rather than being driven by swelling due to fluid intake. Moreover, we show that lack of exopolysaccharide production has a negative impact on cell division rate, and the extracellular matrix components act synergistically to give the biomass the structural strength to produce aerial protrusions and agar substrate-deforming ability.


Introduction
Bacteria growing at an interface between two phases often form biofilms. In the laboratory, biofilms can present in many different forms including floating on the surface of a stationary liquid (pellicle biofilm), within a liquid but adherent to an object (submerged biofilm), or on a solid substrate open to the air (colony biofilm) [1]. Common across each form of biofilm is a self-produced extracellular matrix comprising polymeric proteins and saccharides, extracellular DNA, surfactants, and other substances [2]. The biofilm matrix provides cohesion within the biomass, protects cells from environmental and chemical assault, and aids in hydration and nutrient acquisition [3]. The matrix also confers "emergent" properties to the collective [4].
A key feature of laboratory-grown colony biofilms is the complex architecture that develops. Patterns of ridges and wrinkles are often seen, where the ordered and reproducible nature of the resulting architecture can be indicative of the species under investigation [5]. These complex patterns in the growing colony biofilm biomass are a response to the physical environment (for example, stiffness of the substrate), the availability of both nutrients and oxygen, and are dependent on the composition of the biofilm matrix [6].
Microscopical analyses of localised areas of colony biofilms have been performed to observe the physical movement of Gram-negative bacteria forming colony biofilms. Pseudomonas aeruginosa, in the context of a colony biofilm, divide into discrete cells and use the Type IV pilus to move across the surface and expand the biomass at the very edge of the community [7]. At higher cell densities further into the biomass, varying velocities of this pilus-driven movement generate degrees of orderly packing of cells, which has an impact on the overall biofilm morphology [7]. Additionally, a model has been proposed to simulate the growth and morphology-generating forces of Vibrio cholerae biofilms, accounting for the friction of cells on the substrate, the presence of extracellular matrix, and the diffusion of nutrients [6].
Bacillus subtilis is a Gram-positive soil-dwelling bacterium that is frequently used to study biofilms [8]. While recent systematic analyses of pellicle biofilm formation have been conducted and revealed complex dynamics in cell behaviour [9,10], the same has not been done for biofilms growing on a solid substrate. In contrast to P. aeruginosa or V. cholerae, B. subtilis grows on a semi-solid substrate as long chains of connected cells [11][12][13]. Work to elucidate the forces at play regarding movement of these chains of cells has been undertaken. As the cells grow on an agar surface and form microcolonies, deformations in the orderly chains occur at higher cell densities [14], with contributing structural effects from various matrix molecules [15][16][17]. It has been proposed, for example, that the hygroscopic nature of exopolysaccharides in the matrix results in osmotic spreading that is responsible for the lateral expansion of the colony biofilm [15]. It is well documented that colony biofilms formed by bacterial strains lacking genes that lead to matrix production display morphologies of reduced complexity [18]. Multiple end-point techniques have been employed to probe the importance of the matrix to biofilm architecture and development such as macro-scale microscopy [18]; destructive techniques including thin-section preparation for light and electron microscopy [19,20]; or biofilm disruption to single cells for flow cytometry [21]. While these approaches are informative for questions of protein localisation within the biomass or a snapshot of morphology in the presence or absence of specific matrix components, we lack a motion-based non-destructive analysis of the cell growth and movement responsible for the emergence of structural complexity in the mature community.
Here we present a systematic time-lapse microscopy-based examination of colony biofilm growth by B. subtilis NCIB 3610. We cover defined periods of growth including the point of initial deposition of the founding cells to the agar surface, through to the maturation of the colony biofilm architecture. We attain an unprecedented understanding of the different stages of colony biofilm formation and ascribe quantitative phenotypic differences to genetic mutations altering the constituents of the matrix. The matrix genes included in this study are the three genes of the tapA operon, namely tapA, encoding the secreted TapA protein which may anchor TasA to the cell exterior and can template TasA fibre formation [22]; sipW, encoding the signal peptidase SipW which processes TapA and TasA proteins for secretion [23]; and tasA, encoding the TasA polymeric fibre protein [24]. The genes also include epsH, part of the epsA-O operon whose gene products produce the matrix (WT) and derived strains carrying specific deletion mutations for matrix component genes (tapA NRS3936, sipW NRS5488, tasA NRS5267, tapA tasA NRS5748, epsH NRS5906, bslA NRS2097). The same biofilm is presented at 24 h and 48 h for each strain as indicated, incubated at 30 • C; (B) Footprints of colony biofilms formed by the indicated strains at 24 h and 48 h after cell deposition. Reflected light macroscopic images were segmented in Matlab and areas calculated from the pixel size of the images. n = 18. *p = 0.0106, **p = 0.0033, ****p < 0.0001; (C) Indentation patterns caused by colony biofilms grown at 30 • C for 72 h revealed by gentle removal of the biomass using a plastic inoculation loop. Scale bars represent 5 mm. exopolysaccharide [18]; and bslA, encoding for the colony-coating hydrophobin-like protein BslA [19,25]. For a detailed review of all these matrix components see Ref. [8]. We reveal the contribution of the polymers of the matrix to the morphological development of biofilms by analysing the behaviour of the tissue and linking the matrix components with the overall growth mechanics of the B. subtilis colony biofilm.

Gross colony biofilm morphology varies with matrix composition
Through macro-scale imaging of colony biofilms, the inability of B. subtilis to produce the full array of biofilm matrix components has been linked with i) gross morphological differences in colony biofilm architecture [18], ii) the footprint occupied by the biomass [26], and iii) the community's ability to deform the substratum on which it is grown [27]. Here, reflected light images of colony biofilms formed by NCIB 3610, and a suite of otherwise isogenic mutant strains carrying deletions in the tapA, sipW, tasA, tapA tasA, epsH, and bslA coding regions, demonstrate these impacts (Fig. 1A) (see Table 1). Qualitatively, as expected, each of the deletion strains showed a gross reduction in structural complexity as compared to the parental strain colony biofilm (e.g. Refs. [18,24]), although it should be noted that the exact colony biofilm morphology varies with each genotype. The mutant strains also showed differences in the area occupied by their colony biofilms (Fig. 1B). Specifically, at 24 h of growth, the footprint (the surface area occupied by the colony biomass) of the epsH and tapA tasA mutant colonies were smaller than those formed by NCIB 3610 (p < 0.0001), and by 48 h all the mutant strains occupied significantly smaller footprints (Fig. 1B).
Colony biofilms formed by NCIB 3610 force indentations to develop in the agar substratum, thought to be caused by osmotic pressure exerted by the presence of EPS in the biofilm [27]. As expected, multiple indentations (impressions made in the agar by the colony biofilm) in the agar under NCIB 3610 colony biofilms could be seen below the initial deposit of founding cells (within the 'coffee ring' of high-density cells produced by inoculum evaporation), with further indentations in the agar extending in a radial pattern to the colony biofilm periphery (Fig. 1C). Moreover, as previously reported only low levels of substratum deformation could be induced by the epsH strain ( Fig. 1C) [27]. Our new analysis of the other matrix deletion strains revealed that the only matrix mutant to show a peripheral pattern of substratum indentation was the bslA strain, although visually the indentations formed appeared to be of lower frequency and less pronounced than those formed by NCIB 3610. The remaining matrix mutant strains (sipW, tapA, tasA, tapA tasA) deformed the substratum to varying degrees in the agar directly under the initial inoculum, but deformations were mainly absent under the colony biofilm periphery. We therefore conclude that the absence of any single one of the matrix components alters how the biofilm interacts with the substratum.

Growth rate in the initial colonisation phase is not linked with the mature biofilm footprint
We proposed that novel insight into the microscale processes that influence the macroscale presentation of colony biofilms may be gained by in-depth systematic imaging that quantifies the dynamics of cell movement during the formation process. Therefore, a series of confocal imaging experiments were designed and conducted (Fig. S1, Fig. S2B). The initial experiments involved collecting data from the earliest point of founder cell deposition onto the agar substrate (0-10 h) (Fig. S1A). The starting planktonic NCIB 3610 cultures were normalised to an optical density of 1 at OD 600 and a 1 μl inoculum was spotted on the agar surface. The founding cells contained a ratio of 20% GFP-positive and 80% GFP-negative cells to aid visualisation, by adding contrast to the biomass that quickly becomes optically dense. The inoculum was dried onto the surface of the agar forming a patch of founding cells in a characteristic 'coffee ring' pattern [28] (Movie 1). Over time, the bacteria increased in number and started to grow as chains ( Fig. 2A), forming coherent biomass where the cells filled the agar surface within the coffee-ring (Movie 1).

Movie 1.
We systematically assessed the behaviour of the suite of biofilm matrix mutants over the same timeframe (see Table 1) (Movies 2-7). As with NCIB 3610, each of the matrix mutants contained populations of cells that were physically constrained within a coffee-ring formed upon deposition of the inoculum (Fig. 1A). To determine if differences in mature colony biofilm morphology and footprint were due to changes in growth rate during colonisation of the agar surface, the normalised mean intensity of GFP signal was measured every 10 min over 10 h, with the GFP intensity value serving as a proxy for cell density measurements (see materials and methods) ( Fig. 2A) (Fig. 2B) (Fig. 2C). There was a difference in the mean doubling times for each of the matrix mutant strains compared to NCIB 3610, with the tapA strain being slightly slower and bslA strain slightly faster (p = 0.018 and p = 0.002 respectively, Fig. 2C). Although statistically different, most of the calculated doubling times for these two strains were within the standard deviation of NCIB 3610 (Fig. 2C). Therefore, it is unlikely that differences in the ability of the biofilm matrix mutant strains to initially colonise the agar surface of the inoculation zone would lead to consistently smaller colony biofilm footprint or architecture.

Biofilm edge expansion is linked with biofilm matrix molecules
Concurrent with the growth of the cell population that was constrained within the central coffee-ring zone, chains of NCIB 3610 cells emanated from the exterior boundary of the coffee-ring and expanded radially across the agar surface increasing the footprint (Movie 1) [11]. To record the expansion of the colony biofilm from this point forward, a separate imaging experiment following the growth of the monolayer cell chains was initiated. At 10.5 h after the cells were initially deposited, imaging data were collected for the next 12-15 h (Fig. S1B) (Movie 8, Fig. 3A). By comparing images collected over time we observed that the NCIB 3610 chains travelled >700 μm radially across the agar before leaving the field of imaging. The chains typically traversed the field of view at a rate of 4.25 μm/min and with complex motion of the constituent loops and chains in the expanding biomass (Movie 8). We also examined the biofilm matrix mutant strains during this phase of colony biofilm formation (Movies 9-14). All the strains grew chains of cells emanating from the coffee ring biomass and we calculated the rate of expansion of the biomass across the agar surface for each genotype (see materials and methods). All bar one of the matrix mutant strains displayed a significantly slower mean expansion rate than NCIB 3610, ranging from 0.74 μm/min for the epsH strain (p < 0.0001) to 2.54 μm/min for the bslA mutant strain (p = 0.001). The sipW strain mean expansion rate was not significantly different from NCIB 3610 at 3.72 μm/min, but the data showed two populationssome at the upper bounds of NCIB 3610 rates and the others as slow as the tapA and bslA strains (Fig. 3C). Overall, these data link specific biofilm matrix molecules to the expansion of the biomass from the coffee-ring and are consistent with observed differences in the final colony biofilm footprint

The absence of the biofilm exopolysaccharide increases cell cycle times
Slower expansion across the agar surface by chains of cells produced by the biofilm matrix mutants could be the consequence of slower cell division rates at the single cell level culminating in reduced biomass production. Therefore, to calculate the time taken for a single cell to divide during biofilm expansion, strains were constructed in NCIB 3610 and matrix-gene deletion backgrounds that harboured a gene encoding an FtsZ-GFP fusion protein inserted at the native ftsZ locus. The strains also constitutively expressed mKate2 (Table 1). FtsZ is one of the first proteins to localise at the mid-point of dividing cells, forming a scaffold ('z-ring') for the machinery that will carry out cytokinesis (Movie 15) [29]. After cytokinesis, the z-ring disassembles and re-forms at the mid-point of the daughter cells. The interval between the formation of z-rings of parent and daughter cells can therefore be monitored to calculate the generation time (Fig. 4B). The presence of the FtsZ-GFP fusion protein did not impact growth, as no phenotypic difference in the colony biofilm morphology compared with the parental strain lacking the reporter fusion was detected (Fig. S2A).

Movie 15.
To measure the generation time at the single cell level in conditions where the biofilm matrix is produced, the suite of strains was grown under biofilm-forming conditions for 12 h before high-resolution imaging for multiple generations (Fig. 4A, Fig. S2B). These analyses provided data to calculate the mean generation time for each strain growing on a solid surface (Fig. 4B). NCIB 3610 showed a distribution of doubling times, with a mean of 54 min (±11 min SD) (Fig. 4B). Analysis of the biofilm matrix mutants yielded comparable values ( Fig. 4B) with two exceptions: the sipW strain (p = 0.04), with a slight increase in mean division time of 4 min compared to NCIB 3610 but with similar standard deviation; and more significantly the epsH strain, which took 79 min (±22 min SD) to divide (Fig. 4B). These findings reveal that the absence of exopolysaccharides during sessile growth impacts cell division, and this is likely to contribute to the reduction of the footprint of the epsH mutant compared to NCIB 3610.

A reduction in collective biomass movement is associated with an incomplete matrix
Alongside differences in the expansion rate attained by the biofilm matrix mutants, variations in the form of motion displayed by the cell chains expanding across the agar surface were observed. In NCIB 3610, the density of the biomass to the inside of the expanding edge of the colony biofilm increased to the point where the lateral movement became restricted. At this point, the chains of cells in the biomass (due to the definition afforded by the GFP-positive cells contrasting with the GFP-negative) became bundled and appeared twisted, ultimately leading to slower-moving biomass above the plane of the agar (Movie 16, arrows) [14]. The formation of vertical structures in the peripheral zone of the biomass is consistent with the behaviour of cells within the coffee-ring zone which also appear to rise because of physical constraints (Movie 1). Likewise, at this point, as the cells pile up vertically, they continue to increase the overall biomass. In contrast to NCIB 3610, the expanding chains of the matrix mutant strains, except for the bslA mutant, formed fewer 3-dimensional structures and maintained a monolayer appearance in a larger zone of the colony biofilm periphery (Fig. S3 WT vs sipW). These observations are consistent with the static single-time point perspective of the biomass obtained by confocal imaging of later-stage colony biofilms at both the periphery and interior (Movie 9-14, Fig. 7).

Movie 16.
To compare the dynamics of motion across the array of strains, we marked and tracked the 2D movement of 36 visibly distinct features in the biomasses (six features were followed per image, six images were in each dataset, and the features selected were distributed across the initial biomass in the field of view imaged) (see Materials and Methods) (Movie 17). The field of view imaged is positioned at the edge of each biofilm such that movement in a straight line along the x-axis from left to right is equivalent to, in the context of a polar coordinate system of a whole biofilm, increasing in radius without deviating in angle theta from the centre (radial growth). We analysed the path, rate of motion, and change in distance from the biofilm edge for each feature for the maximum duration of the fastest-moving features, to the point at which the edge of the biofilm left the field of view (which was 1 h and 40 min). Plotting of feature trajectories from NCIB 3610 biofilms revealed the direction of movement is not strictly aligned with the x-axis, with features deviating to different degrees along the y-axis (Fig. S4). To compare the types of tracks of each genotype, the starting point for each trajectory was normalised to the origin and overlaid in separate plots (Fig. 5A). Comparing all the feature trajectories of the matrix mutant strains to those of NCIB 3610 revealed differences in the motion profiles -for example, the path taken by features within the sipW strain colony biofilm tended to extend outwards along the x-axis to a similar degree to, but less prone to motion along the y-axis than, those of NCIB 3610 (Fig. 5A sipW). Isolating the xaxis component of these trajectories revealed the sipW strain was the only mutant strain that did not produce feature tracks that were significantly shorter than NCIB 3610 in this direction (Fig. 5B, p = 0.1617 compared to p < 0.0001 for all other mutant strains). Isolating the y-axis component of the trajectories shows that all the mutant strains had a significantly reduced range of motion away from an x-axis trajectory compared to NCIB 3610 (Fig. 5C, p = 0.0008 for tapA mutant strain, all other mutant strains p < 0.0001).

Movie 17.
Analysis of the rate of movement for the features tracked revealed further insight into the collective biomass motion (Fig. 5D, see annotation marked 's' on Movie 16). Calculating all the distances moved between each time-point, for NCIB 3610 three distinct populations of stepsizes are evident as shoulders in the density histogram (Fig. 5D, WT, black brackets), implying that there are different categories of motion. There was a similar pattern of shoulders in the distribution of tapA mutant strain step sizes as NCIB 3610, albeit with an increase in density of smaller step sizes (Fig. 5D tapA). The sipW mutant had a broad distribution of step-size with some features moving as far as NCIB 3610 with step-sizes up to almost 60 μm per 10 min interval but lost the three categories of step-sizes evident from NCIB 3610 and the tapA mutant strain (Fig. 5D sipW). The remainder of the matrix mutant strains displayed a unique profile of step-size distribution, generally showing much smaller step sizes than NCIB 3610 or the tapA and sipW mutant strains (Fig. 5D).
Finally, we tracked each landmark relative to the outer leading edge of its colony biofilm over time to determine if step-size correlated with the expansion rate of the colony biofilms (Fig. S5, see annotation marked 'd' on Movie 16). Across each of the mutant strains, individual features either maintained pace with the expansion rate (approx. horizontal lines), travelled more slowly, increasing the distance from the expanding edge (positive gradient lines), or had an increased pace, catching up with the edge (negative gradient lines). Re-plotting these data to represent the change in distance to the edge over time, rather than the absolute distance, highlights each strain's propensity to the movement categories described here (Fig. 6A). The predominant characteristic of NCIB 3610 feature trajectories is to fall behind the leading edge, as evidenced by the very few negative gradients when plotted (Fig. 6A). The trends of the sipW and bslA strains could be described similarly, but those of tapA, tasA and epsH strains include a higher proportion of negative gradients (Fig. 6A). The variability in the trajectory travel profiles and landmark motion across the space implies that there are multiple forces acting on the biofilm biomass and that during growth these forces push or pull the biofilm biomass in different directions in a matrix-dependent manner, revealing the impact of the matrix components on the material properties of the colony biofilm as a collective biomass.

Colony biofilm expansion is a result of peripheral mono-layer growth
To address whether the collective movement of biofilm mass, from the interior of the community to the exterior, contributes to the expansive capability of B. subtilis colony biofilms or establish if expansion is driven by the outer edge of the community, we performed confocal imaging at regions immediately adjacent and interior to the areas of mono-layer chain expansion (Figs. S1C and D). We hypothesised that if movement of biofilm from the interior of the community contributed to expansion of the biofilm footprint, motion of the material moving towards the exterior would be at a rate that exceeds the rate of expansion directly measured for the biofilm periphery (Fig. 3C). Image collection started after 12 h of biofilm growth and was sustained for 15 h. Data were acquired in two separate experimental set-ups capturing different perspectives of the colony biofilm growth: the biofilm-agar interface (Movie 18 (WT); Movie 19-24 (matrix mutants)) and biofilmair interface (Movie 25 (WT); Movie 26-31 (matrix mutants)).
Supplementary Imaging of the biofilm-agar interface from below showed that NCIB 3610 and the bslA strain formed channels that hosted clusters of motile planktonic cells close to the expanding edge of the biofilm (Movie 18, WT and bslA; Fig. S1D). Both NCIB 3610 and the bslA mutant also had sectors of material in the growing colony biofilm material that showed concerted motion, first moving slowly in the outward direction before reversing towards the interior. The reverse motion is evident in biomass on the interior side of the field of view first, before appearing to stretch and drag the rest of the material towards the interior (Movie 18, WT and bslA). The NCIB 3610 colony biofilm material can be seen to press into the agar substratum, while the bslA biomass folds in on itself. This concerted motion within the biomass was not evident in the material formed by any of the other matrix deletion strains (Movie 18, compare NCIB 3610 and bslA to tasA and see Movies 19-23 tapA, sipW, tasA, tapAtasA, epsH). Given that the indentations under biofilm peripheries were only observed in the cases of NCIB 3610 and the bslA mutant strain (Fig. 1C), it seems likely that, as movement at the agar interface was only observed in only these two strains, it is this movement that leads to formation of the indentations. Overall, no gross outward expansion of the biomass was apparent for any of the strains imaged at the agar interface (Movies 18-24).
Imaging of the biofilm-air interface (Fig. S1C) of NCIB 3610 showed transition of the peripheral mono-layer chains into a thicker tissue before the interior biomass entered and traversed the imaging field, moving towards the periphery of the colony biofilm at a rate of 0.87 μm/ min (Fig. 6B), (left to right, Movie 25, WT). In contrast, although each of the matrix deletion strains showed a transition of the peripheral monolayer chains into a thicker tissue, the material displayed no coordinated translational motion for the 15-h duration of imaging (Movies 26-31).
When the mean outward movement of NCIB 3610 biomass at 0.87 μm/ min distant from the edge is compared to the mean peripheral expansion rate of 4.25 μm/min (Fig. 3C) and is combined with the lack of any outward movement of the same biomass in the matrix deletion strains, we conclude that expansion of B. subtilis colony biofilms is primarily the consequence of growth at the biomass periphery.

Thickening of mature biofilms requires the full suite of matrix components
Our final imaging dataset provided a more detailed perspective of the gross colony biofilm architecture by the application of confocal imaging to both the outer edge and the interior of NCIB 3610 colony biofilms. Both the fine structures seen in the periphery of the colony biofilm and the cylindrical structures of the colony biofilm interior at 48 h of growth appeared to become thicker and more pronounced when equivalent regions were imaged at 72 h of growth (Fig. 7A). To quantify this visual change, the images from the periphery were analysed for the total volume of material in the visible biomass, using the position of the agar as a reference (Fig. 7B). The analysis showed continued development of complex structures, alongside an ~1.76-fold thickening of the biomass between 48 and 72 h (Fig. 7C). In sharp contrast, similar imaging of the biofilm matrix mutants revealed few complex structures and no increase in biomass volume over a comparable time frame, mostly reporting a volume ratio very close to 1 (Fig. 7C, dashed line). The bslA mutant strain showed a slight reduction in volume, with a mean ratio of 0.87. Therefore, only wild-type NCIB 3610 colony biofilms have the propensity to increase in thickness (Fig. 7B, C), revealing the dependency on the full suite of biofilm matrix molecules for vertical expansion of the colony biofilm.

Discussion
By developing and implementing time-lapse imaging and quantitation methods we observe and describe the dynamics of Bacillus subtilis colony biofilm formation, revealing the contribution of the extracellular matrix molecules. Through our analyses, we provide new insight to the microscale processes underpinning the formation of the resulting macroscale community. Understanding the specific function and impact of each of the polymeric matrix molecules secreted by cells in a colony biofilm is a challenging goal. By deleting the genes needed to produce each matrix component in turn (tapA, sipW, tasA, epsH, bslA), and in combination (tapA tasA), it was hoped that these functions would present themselves in the physical dynamics of the biomass. Our data show that BslA is not required to produce internal architecture leading to channel formation, but this feature is lost in the absence of epsH and any of the tapA operon gene products. Moreover, all the matrix component deletion strains showed a loss of coherence in the movement that accompanies growth.

Matrix influence on cell division and physiology
Acknowledging that cells undergoing sessile growth experience physiological conditions that are different from those in planktonic culture [30], we developed a new method to calculate the cell division time of bacteria growing in the context of a biofilm with single-cell resolution. Each of the matrix-deficient strains, except the epsH mutant, divide at comparable rates to NCIB 3610 during expansion of cell chains at the periphery of the colony biofilm. That the absence of exopolysaccharide results in increased cell-cycle times helps to explain the eps mutant colony biofilms' diminutive form. The increase in division time shown by the eps mutant cells was surprising given the metabolic burden of synthesising the exopolysaccharide [31]. It could be that the eps mutation has a global impact on metabolism under biofilm growth conditions, but why cell division is slowed is currently unexplained and warrants further investigation.
There is a disparity in observed doubling times at different stages and physical locations in biofilm growth. During initial colonisation of the agar, in the period after deposition of the cells, NCIB 3610 had a mean doubling time of around 4.3 h (Fig. 2C), while after 12 h of growth the chains of cells at the expanding periphery double every 54 min, on average (Fig. 4B). We speculate that the physiology of the cells is very different in these two cases: during colonisation the cells are switching from planktonic to sessile growth mode and are ultimately constrained within the area of inoculation; at the periphery at later times the cells have adapted to sessile growth and have immediate access to free space in which to grow. Given the range of morphological structures in NCIB 3610 colony biofilms, from dense regions of cells, aerial projections, and fluid-filled cavities, we surmise that there are many physiological pressures present that will lead to regulatory adaptations within a local area of cells.
Altered gene regulation in micro-environments of the biofilm may also account for the bimodality of expansion rates measured for the sipW strain, where two distinct clusters of speeds were detected in the fields of view that were selected for imaging (Fig. 3C). SipW is a signal peptidase responsible for processing and releasing nascent TapA and TasA from the cell [23], but also has a matrix gene regulatory function in the context of biofilms [32]. It is conceivable that disturbing the balance of matrix gene regulatory signals in situations of heterologous physiological pressures could result in inconsistent responses of cells in a local area. However, on a whole-biofilm scale, sipW gene deletion leads to reduced morphological complexity and smaller biomass footprints ( Fig. 1A and B) compared to NCIB 3610.

Matrix influence on tissue motion and structure
Described as an active matter system, chains of B. subtilis cells on an agar substrate have been observed to behave nematically and produce topological defects to alleviate stresses incurred by growth, creating 3dimensional structures [14]. When considering why mutant strains deficient in matrix components form biofilms with smaller footprints and display less complex morphology, we find that it is not the initial growth rate of the cells, as they transition into sessile growth from planktonic culture, that is the root cause of macroscale impact. Rather, we observe 2D chains growing into 3D structures a short distance from the edge of expanding biofilms, and we propose that the dynamism of movements of these ordered structures away from the purely lateral direction is dependent on strength provided by the combination of matrix molecules present. This tensile stress is transmitted over a larger length-scale when the full suite of matrix molecules is produced, confirmed by the deviations of feature motions from the radial direction of biofilm growth, and the formation of structures that extend both aerially from the biomass and embed into the agar substratum. This process produces and maintains internal cavities hosting fluid and pockets of swimming planktonic cells that contribute to heterogeneity in the isogenic community [33]. The matrix-mutant biofilms show a lack of 3-dimensional organised structures. If this is due to a loss of scaffolding in the spacing between chains of cells, it follows that there might be tighter cell packing, that the overall surface area is reduced and total cell density greater in these cases. Tighter cell packing in this way has negative implications with regards to nutrient and oxygen permeation [34], and would go some way to explaining why matrix mutant biofilms with cell division rates similar to NCIB 3610 at the periphery produce consistently smaller biofilms.

Mode of biofilm growth
Osmotic spreading, or the swelling of the B. subtilis colony biofilm by uptake of water from the substratum, has been proposed [15]. However, by directly observing the motion of biomass at various locations throughout biofilm expansion we see that instead of swelling from within, NCIB 3610 colony biofilms have the highest rate of lateral growth at the periphery, and this growth results from cell division in mono-layer chains rather than biomass accumulating and pushing from the interior. Moreover, for NCIB 3610 the slow biomass movement inside the biofilm towards the exterior was only observed at the biofilm-air interface and not at the biofilm-agar interface. This model of colony biofilm expansion is consistent with waves of biofilm matrix gene expression at the periphery of the biofilm in B. subtilis [35] and zones of active growth in similar colony biofilm communities formed by Escherichia coli [36]. Furthermore, motion of material at the underside (agar proximal) of the colony biofilm, that was initially at a short distance from the edge, indicates that as growth continues this foundational biomass not only lacks expansive sliding, but also becomes embedded into the agar substratum, forming peripheral indentations in the substrate [27].

Overarching conclusion
In toto, this study has provided unprecedented views of a complex macroscale microbial community forming from an initial dispersed population of cells. The synergy of the matrix molecules in providing the tensile strength needed to form the complex architecture becomes apparent, and collective motion of the community reveals how the forces generated by cell growth are dissipated in the material. Whilst the specific functions of each matrix component remain elusive, in total, the matrix components function with dependence on each other to confer structure and strength, as indicated by the mutant strains being unable to maintain even small-scale morphological features or to increase in thickness as the biofilm matures. Combining imaging analysis of this nature with microscale selection of distinct micro-pockets of cells would allow the molecular processes underpinning the transition from the motile state to the embedded social community to be followed. In the future it will be intriguing to apply these combined techniques to identify the benefits of matrix molecules on a broad range of bacterial species in their varied natural environments.

Construction of strains
Strains generated for this study that carried genes encoding GFPMut2 or mKate2 under constitutive promoters were produced by SSP1 bacteriophage transduction of genetic material from previously published strains (see Table 1). Briefly, bacteriophages were generated from the donor strain carrying the desired genetic element (sacA:Pspachy-gfp-kan or amyE:Pspachy-mKate2-spc) using previously published methods and introduced into the recipient strain [37].

Construction of FtsZ-GFP fusion protein strains
A marker-less ftsZ-gfp fusion gene was created at the native ftsZ locus by homologous recombination using vector pMAD2 [38] containing the gene for GFPMut2 flanked by homology arms (HA) of ftsZ at the 5 ′ end and bpr at the 3 ′ end. The two HAs were created by PCR of NCIB 3610 genomic DNA using primers NSW3027 and NSW3033, and NSW3037 and NSW3032 respectively (see Table 2). A gfpmut2 fragment was generated by PCR from genomic DNA of strain NRS2388 (see Table 1) using primers NSW3034 and NSW3035. These PCR fragments contained flanking homology for assembly using a Gibson Assembly kit (New England Biolabs) to yield the final vector pMAD2-5 ′ HA-gfpmut2-3 ′ HA (pNW2410) which was verified by sequencing using primers NSW3023 and NSW3024. This vector was introduced into B. subtilis strain 168 cells, then into final recipient strains (NCIB 3610, NRS3936, NRS5488, NRS5267, NRS5748, NRS5906 and NRS2097) (see Table 1) via phage transduction [37]. Integration of ftsZ-gfpmut2 into the genomes of these strains was achieved by a method published previously [38] and positive colonies were screened for GFP fluorescence on a stereomicroscope and PCR analysis. The resultant strains had the addition of constitutively expressing mKate2 at the amyE locus via phage transduction [37] from strain NRS5841, selecting on spectinomycin (100 μg/ml) plates to yield final strains NRS6781, NRS6793, NRS6794, NRS6795, NRS7525, NRS7526 and NRS7527.

Culture of cells and biofilm growth
Colony biofilms were founded from a pre-culture grown in 3 ml liquid Lysogeny Broth (LB, per litre: 10 g yeast extract, 5 g NaCl, 10 g tryptone) from a single colony streaked from frozen strain stocks on an LB agar (1.5% w/v SelectAgar, Invitrogen) plate. The liquid cultures were incubated for three to 4 h at 37 • C shaking at 200 rpm and adjusted to the same optical density (OD 600 ) as the lowest culture in a set, in the range of OD 600 of 0.9-1. Strains expressing fluorescent molecules were then mixed with their parental strains (see Table 1) after normalisation of cell density (80 μl non-fluorescent strain with 20 μl fluorescent strain). Colony biofilms were inoculated from this mixture by spotting 1 μl of the pre-culture onto MSgg (5 mM potassium phosphate and 100 mM MOPS at pH 7.0 autoclaved and cooled to 55 • C then supplemented with 2 mM MgCl 2 , 700 μM CaCl 2 , 50 μM MnCl 2 , 50 μM FeCl 3 , 1 μM NaCl 2 , 2 μM thiamine, 0.5% (v/v) glycerol, 0.5% (v/v) glutamic acid) agar (1.5% w/v). The pre-culture spot was air dried under a flame for 5 min before incubating the covered colony biofilms at 30 • C in a humidified box or transferring them immediately to a pre-warmed microscope chamber.

Colony biofilm and indentation imaging
Macro imaging of whole biofilms and areas of agar indentation were performed on a Leica stereoscope (M205FCA) using reflected light from a ring illuminator and either a 0.5 × 0.2 NA or 1×0.2 NA objective.

Confocal imaging of constitutive GFP colony biofilms
LabTek II 2-chamber microscope vessels (Thermo Scientific) were used to grow and image biofilms 30 min, 10.5 h, 48 h, and 72 h after cell deposition (as indicated), with 2 biofilms per chamber. Two vessels could be imaged concurrently, totalling 8 biofilms. For imaging the biofilm-agar interface, the chambers were filled with 750 μl MSgg agar to put the biofilm into the working distance of the objective without excessively scattering light with the agar, in all other cases 4.4 ml MSgg agar was used, leaving a small headspace for the biofilms to grow, and be imaged from the biofilm-air interface. After solidifying, the MSgg agar was air-dried in a laminar flow hood for 45 min.
Non-fluorescent and constitutive GFP-expressing equivalent strains of wild-type and matrix gene mutants (see Table 1) were grown, normalised for optical density and mixed, as detailed above. The cells were placed onto 4 ml MSgg agar (1.5% w/v) set into each chamber of LabTek II 2-chamber microscope vessels (Thermo Scientific) set into various containers, as detailed in specific methods. Two vessels were prepared per experiment, giving four chambers, and the MSgg agar was dried in a flow hood for 45 min before spotting. Two biofilms per strain analysed were spotted in a single chamber, spaced evenly apart.
A Leica SP8 upright confocal was used with a chamber pre-warmed to 30 • C. The 488 nm argon laser was set to 2% power and rapid imaging was performed with a resonant mirror at 8 kHz with bi-directional scanning and line averaging of 16x. Gain of the PMT detector was set to between 550 v and 650 v to avoid saturating pixel values. The pinhole was set to 1 AU for a wavelength of 525 nm. Nyquist acquisition was achieved using a 10x, 0.3 NA long working distance objective for zstacks with 3.87 μm steps and a scan of 1024 × 1024 pixels over an 890 μm area. Time-lapse imaging was set with an interval of 10 min. For imaging the biofilm-agar interface the vessels had their lids removed and were inverted onto the microscope stage. For imaging the biofilm- air interface the vessels were placed on the stage and had their lids replaced with a 60 mm × 24 mm #1.5 coverglass to prevent biofilm exposure to the flowing heated air of the microscope chamber. The biomass proximal to the Petri dish edge was imaged.

Confocal imaging of FtsZ-GFP fusion strains
Vessels were prepared by adhering a single Gene Frame (65 μl) (Thermo Scientific) onto a glass microscope slide, leaving the plastic aperture cover in place. A stepped flat mould was created by gluing 22 mm × 22 mm glass coverslips to a second microscope slide at a distance apart that would align to the inside edges of the Gene Frame. 150 μl of MSgg agar (at 55 • C) was pipetted into the aperture of the Gene Frame, avoiding forming bubbles on the surface, and the stepped mould was immediately inverted over the aperture and held in place for 30 s to allow the agar to set, slightly raised above the height of the Gene Frame, and removed by gently sliding off. This was done to ensure that the sample could still be mounted after agar shrinking during incubation before imaging. After mixing cultures (see above), three biofilms were inoculated onto the agar, spaced equally. The slide was incubated in a humidified box at 30 • C for 12 h before mounting. To maximise oxygen availability during imaging, most of the agar in the Gene Frame was cut away with a sterile scalpel, leaving approximately one quarter of each biofilm (Fig. S2Biii) before peeling off the aperture cover from the Gene Frame and applying half a gas-permeable polymer coverslip (Ibidi) that had been cut in two equal pieces. Imaging was performed on a Leica SP8 inverted confocal microscope with a chamber pre-warmed to 30 • C. Excitation of GFP and mKate2 was achieved with 488 nm and 568 nm lasers, respectively, set to 2% and 0.1% power and read by hybrid detectors with a gain of 10 v and 2 v. The pinhole was set to 1 AU at 525 nm for optimal resolution of GFP and both channels were acquired simultaneously. Fast imaging was employed with a resonant mirror at 8 kHz and line averaging of 16x and the image field was 752 × 752 pixels. Three image fields were set around the periphery of each of the three biofilms on the slide, generating nine movies per acquisition, with a z-stack specified to capture at Nyquist sampling with the 63×1.4 NA objective. Images were collected every 10 min for 5 h.

Image handling and figure preparation
All microscopy data were imported into an OMERO server [39] for organisation, annotation, and analysis. Maximum intensity image projection and annotation with ROIs were performed with OMERO.insight. Static image figures were created using OMERO.figure while movie figures were prepared using bespoke code in Matlab (MathWorks) (See Image Analysis). Raw, projected and ROI-annotated image data are available at https://www.ebi.ac.uk/biostudies/studies/S-BIAD474 which contains a link to live data on the OMERO server.

Image analysis
Analysis of intensity images and ROI positions were performed in Matlab (R2019a) using the OMERO.matlab toolkit (v5.5.4) to read data from the OMERO server (v5.6.3). All code is available at (https://doi. org/10.5281/zenodo.6563834). A summary of analysis workflows is as follows: Initial doubling timethree identical square ROIs were placed on images in the internal region of the newly-spotted biofilm and propagated through every time-point of the images. The ROIs were used to define patches of the image to be downloaded, and the mean intensity of each patch was calculated. The means were normalised to the highest value over all time-points, producing a sigmoid curve. The highest gradient was detected and used to define the doubling time.
Expansion rate -Time-lapse images of chains of cells, oriented to move left to right across the image space, were segmented using Otsu thresholding [40] until the chains reached the far edge of the image. The right-most pixel segmented had its x-axis location recorded for each time point, converted from pixel units to μm using the pixel size metadata.
These position data were plotted, and a straight line was fitted with high confidence. Expansion rate was calculated from the gradient of the fitted line.
Tracking of feature landmarks -For each time-lapse image, six visibly recognisable features were marked with a point-ROI, which was propagated to the time-point where the chains reached the far edge of the image. Features were chosen to be spread across the height and width of the observed chains in the image at the first time-point. The locations of the points were then updated manually for each time-point. For fair comparison, the number of time-points used for downstream analysis was limited to that of the shortest tracks detected. The XY location of each point was used to generate tracks across the image space, from which the step-sizes (distance moved between time-points), distance to the edge (using the edge location from Expansion rate, above) and speed were calculated.
FtsZ cycle time -Time-lapse imaging was examined for the appearance of new Z-rings. Upon first sighting, a Z-ring was marked with a single point ROI and designated as the 'root' of a pedigree with a '0' entered as the ROI's comment. Progeny were followed over time and when new daughter Z-rings became visible they were marked with a point and given the id of their parent's ROI as a comment. The pedigrees were followed over multiple generations, using the constitutive mKate2 signal of this strain to visually delineate daughter cells. Each image's ROIs were then analysed in Matlab, sorted into pedigrees using the root ('0') and parent ROI id comments to form links between them, calculating the interval from the difference in time-point of the ROIs.
Footprint size analysis -The area occupied by whole biofilms was calculated in Matlab by segmentation using Otsu thresholding [40] of reflected-light images acquired on a stereomicroscope. Pixel units were converted to μm using pixel size metadata from the image and reported as μm 2 .
Biofilm thickness ratio -Volumetric images from the edge of biofilms were segmented using Otsu thresholding [40], and for each XY coordinate the upper-most pixel in the Z dimension containing segmented data was recorded. The total volume of biomass was then calculated using the z-section size metadata and a ratio was generated between biomass volumes at 48 h and 72 h of growth.
Internal biomass movement -Visibly recognisable structures moving across the image field were identified and marked with a point ROI that was used to manually track progress across the images for a minimum duration of 11 h. The Euclidean distance given by the change in coordinates of the point ROIs during each 18-min time-interval was then used to calculate movement rate.
Movie figure preparation -Time-lapse images were downloaded into Matlab plane by plane from OMERO to generate frames onto which ROIs could be used to highlight specific features. Scale bar lengths were chosen and drawn onto each frame using the pixel size metadata of each image. The frames were then collated and saved as.avi files with a framerate of 2 per second.

Graphing and statistical analysis
Histogram ridge plot of step-sizes was generated in R from data exported from Matlab, code for which is available at (https://doi. org/10.5281/zenodo.6563834). Other graphs were drawn in Matlab and Prism (Graphpad). Statistical analyses were performed in Prism and p-values reported are the results of 1-way ANOVA using Dunnett's multiple comparisons test with significance set to p = 0.05.

Data availability statement
Computational code is available on Github and has been archived by Zenodo -doi:10.5281/zenodo.6563834.