Copyright © 2018 Logsdon and Aldridge. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# The Empirical Fluctuation Pattern of *E. coli* Division Control

Jacopo Grilli 1 , Clotilde Cadart 2,3, Gabriele Micali 4,5, Matteo Osella6,7 and Marco Cosentino Lagomarsino8,9,10 \*

1 Santa Fe Institute, Santa Fe, NM, United States, 2 Centre National de la Recherche Scientifique, Institut Curie, PSL Research University, UMR 144, Paris, France, 3 Institut Pierre-Gilles de Gennes, PSL Research University, Paris, France, 4 Department of Environmental Microbiology, Eawag, Dübendorf, Switzerland, 5 Department of Environmental Systems Science, ETH Zurich, Zurich, Switzerland, 6 Physics Department, University of Turin, Turin, Italy, 7 Istituto Nazionale di Fisica Nucleare Sezione di Torino, Turin, Italy, 8 Sorbonne Université, Paris, France, 9 Centre National de la Recherche Scientifique, UMR 7238, Paris, France, 10 IFOM, FIRC Institute of Molecular Oncology, Milan, Italy

In physics, it is customary to represent the fluctuations of a stochastic system at steady state in terms of linear response to small random perturbations. Previous work has shown that the same framework describes effectively the trade-off between cell-to-cell variability and correction in the control of cell division of single E. coli cells. However, previous analyses were motivated by specific models and limited to a subset of the measured variables. For example, most analyses neglected the role of growth rate variability. Here, we take a comprehensive approach and consider several sets of available data from both microcolonies and microfluidic devices in different growth conditions. We evaluate all the coupling coefficients between the three main measured variables (interdivision times, cell sizes and individual-cell growth rates). The linear-response framework correctly predicts consistency relations between a priori independent experimental measurements, which confirms its validity. Additionally, the couplings between the cell-specific growth rate and the other variables are typically non zero. Finally, we use the framework to detect signatures of mechanisms in experimental data involving growth rate fluctuations, finding that (1) noise-generating coupling between size and growth rate is a consequence of inter-generation growth rate correlations and (2) the correlation patterns agree with a near-adder model where the added size has a dependence on the single-cell growth rate. Our findings define relevant constraints that any theoretical description should reproduce, and will help future studies aiming to falsify some of the competing models of the cell cycle existing today in the literature.

Keywords: linear response theory, single-cell growth and division, fluctuation patterns, control of cell division, models, theoretical, data interpretation, statistical

### 1. INTRODUCTION

Today, dynamically tracked data of many dividing cells offer the possibility to analyze with great precision and detail the decision process leading to cell division (Osella et al., 2017). While such data are starting to be abundant, measurements of many single cells yield a complex tangle of correlation patterns between growth-related variables, which is often difficult to grasp. Consequently, new theoretical and data-analysis tools need to be developed in order to extract from such data the relevant information to understand the control of cell-cycle progression, cell division and their impact on cell proliferation.

### *Edited by:*

Ariel Amir, Harvard University, United States

### *Reviewed by:*

Hanna Salman, University of Pittsburgh, United States Sattar Taheri-Araghi, California State University, Northridge, United States

*\*Correspondence:*

Marco Cosentino Lagomarsino<EMAIL_ADDRESS>
### *Specialty section:*

This article was submitted to Microbial Physiology and Metabolism, a section of the journal Frontiers in Microbiology

> *Received:* 15 March 2018 *Accepted:* 20 June 2018 *Published:* 30 July 2018

### *Citation:*

Grilli J, Cadart C, Micali G, Osella M and Cosentino Lagomarsino M (2018) The Empirical Fluctuation Pattern of E. coli Division Control. Front. Microbiol. 9:1541. doi: 10.3389/fmicb.2018.01541

Specifically, one can restrict the question to the control of cell size, where some general principles are emerging from empirical data. E. coli cells grow with a cell-specific growth rate, they are born with a cell-specific size, and they divide with a cell-specific cell-cycle time. A cell that deviates from average behavior in size or growth rate can correct by a compensatory deviation in doubling time. A "near-adder" control of cell division (Campos et al., 2014; Taheri-Araghi et al., 2015), whereby the differential size extension in a single cell cycle is independent of the initial size of the cell, is an effective principle governing cell division. Additionally, the distributions of interdivision times and cell sizes across conditions (Iyer-Biswas et al., 2014b; Kennard et al., 2016) show a clear link between average values and the variability of these variables from cell to cell, thus suggesting a "universal" shape of the size distribution. However, we still understand relatively little regarding how these simple principles emerge (Ho and Amir, 2015; Harris and Theriot, 2016; Wallden et al., 2016; Osella et al., 2017). For example, there is a debate on whether near-adder behavior could emerge from specific mechanisms or molecular circuits (Taheri-Araghi et al., 2015; Harris and Theriot, 2016), as a byproduct of chekpoint control of cell-cycle progression (Adiciptaningrum et al., 2015; Wallden et al., 2016), or as a consequence of external constraints or trade-offs (Osella et al., 2017).

Comprehensive and precise quantitative methods are necessary in both data analysis and theory to deal with the new data. Broadly, the open question is how much a mechanism can be isolated and specified with available data. The simplest theoretical framework for cell-cycle control describes the cell cycle as a discrete-time process, relating the measured variables (size, interdivision times, etc.) across generations (Amir, 2014; Campos et al., 2014; Taheri-Araghi et al., 2015; Grilli et al., 2017). We have recently shown the general existing equivalence between this simple discrete-time formalism to more detailed continuous-time models describing the division process by a division rate varying with cell state and time (Grilli et al., 2017). In the limit where deviations around the mean initial size (or interdivision time) are small (as first proposed in Amir, 2014), we have shown explicitly how the resulting "linearresponse" framework describes a wide range of division control mechanisms and characterizes with remarkable precision the available experimental data. These results make the linearresponse framework perfectly suitable to model cell-cycle control in E. coli given the available datasets.

However, studies are usually restricted to a subset of the measured variables. Most importantly, while single-cell growth rate clearly fluctuates (Wang et al., 2010), very few studies have addressed the consequences of these fluctuations on cell division control, and the studies that did started from very specific assumptions (Osella et al., 2014; Wallden et al., 2016; Logsdon et al., 2017). More specifically, previous studies adopting the linear-response framework typically neglected cell-to-cell variability in growth rates.

Here, we extend this linear-response framework to incorporate growth rate fluctuations, and we show how it can be used in general to evaluate exhaustively all the possible correlations and fluctuations in the data. We use this systematic approach to evaluate jointly all the homeostatic and noisegenerating couplings measured in different experimental studies, and to connect correlation patterns with possible mechanisms underlying cell division.

# 2. DATA SETS

We tested our theoretical considerations and their implications on the analysis of empirical data with E. coli single-cell growth division data from Kiviet et al. (2014), Taheri-Araghi et al. (2015), Kennard et al. (2016), and Wallden et al. (2016). Scripts and formatted data are available with the authors. A detailed list of the growth conditions is available in the Supplementary Information of Cadart et al. (2017).

### 3. LINEAR-RESPONSE FRAMEWORK FOR CELL GROWTH-DIVISION FLUCTUATIONS

### 3.1. General Features of Fluctuations of Cell Size and Growth Rate

This section describes the main model assumptions and definitions. We assume exponential growth (Campos et al., 2014; Iyer-Biswas et al., 2014b; Osella et al., 2014; Taheri-Araghi et al., 2015) of cell size V(t) = V(0) exp(αt), where α is the single-cell growth rate. A cell divides at size Vf = V0 exp(ατd), where V0 is the size at birth and τd is the division time. We neglect the fluctuations around binary fission or the process of filamentation and recovery (Osella et al., 2014), thus cell size at division Vf is equally partitioned between the two daughter cells in our description. Given the assumption of exponential growth, it is useful to introduce the logarithmic size q(t) = log(V(t)/V ∗ ). With this notation, the exponential growth translates into a linear relationship q(t) = q(0) + αt. In the definition of q, we introduced a size scale V ∗ , which in principle can be an arbitrary scale, to make the argument of the logarithm dimensionless. A particularly useful choice is to define V ∗ = hV0i , so that, at least in the small noise limit, hq(0)i = 0 independently of the condition (Grilli et al., 2017).

The distribution of cell size at birth gives a static picture of cellsize fluctuations. It is well accepted (Taheri-Araghi et al., 2015; Kennard et al., 2016; Grilli et al., 2017) that the distributions of size at birth (or, equivalently, at division) obtained for different strains and conditions collapse when rescaled by their average. This collapse is the consequence of the existence of a unique size-scale parameter of cell size control that varies across conditions (Kennard et al., 2016; Grilli et al., 2017). A consequence of the collapse is that the expectation of any function of the rescaled size at birth V0/hV0i is constant across conditions. **Figure 1** shows that the variance of logV0/hV0i (which is by definition equal to σ 2 q ) is constant across conditions, as predicted by the collapse of the size distribution.

Growth rates fluctuations are usually neglected. Contrarily to the well accepted collapse of the size distributions (Giometto et al., 2013; Taheri-Araghi et al., 2015; Kennard et al., 2016; Grilli et al., 2017), there is no consensus on whether the distribution of the single-cell growth rates collapse when rescaled by the average

Each symbol represents an average over all single cells available in the dataset (defined by shape and color) for a specific strain in a specific condition. (A) Shows that the standard deviation of the rescaled log-initial size log(V0/hV0i), is constant across conditions and experiments, and does not depend on the average log-initial size hlog V0i. This fact is a consequence of the collapse of initial size distribution (Taheri-Araghi et al., 2015; Kennard et al., 2016; Grilli et al., 2017), since the standard deviation log(V0/hV0i) is a function of the coefficient of variation of the initial size only (Grilli et al., 2017). (B) Shows the standard deviation of the rescaled growth rate α/hαi (which corresponds to the coefficient of variation of α). This quantity is approximately constant for some datasets (e.g., Taheri-Araghi et al., 2015, blue squares), suggesting a collapse of the growth rate distribution (Taheri-Araghi et al., 2015), while it shows a decreasing trend in others (e.g., Kennard et al., 2016, green circles).

growth rate, with some studies reporting the collapse (Taheri-Araghi et al., 2015) and others showing the opposite (Kennard et al., 2016). **Figure 1** shows that the coefficient of variation of the growth rate appears to be constant, independently of the condition and the average growth rate (consistently with the collapse hypothesis) for some datasets (Taheri-Araghi et al., 2015), while it shows a decreasing trend for others (Kiviet et al., 2014; Kennard et al., 2016; Wallden et al., 2016). A possible explanation of these differences could lay in the growth mode of the different experiments (microcolonies in agar or polyacrylamide, wide microfluidics channels, single-file microfluidics channels, etc). Theoretical efforts to provide a rationale of growth rate distributions are available (Iyer-Biswas et al., 2014b; Pugatch, 2015; De Martino et al., 2016; De Martino et al., 2017).

Independently of the controversy about the collapse of the growth rate distributions, it is clear from **Figure 1**, that the fluctuations of the growth rates are not negligible. The coefficients of variation are in the range 0.1–0.3, which is comparable with the typical fluctuations of the sizes at birth, whose CVs are of the order of 0.2 (Amir, 2014; Grilli et al., 2017).

### 3.2. Linear-Response Framework in Presence of Growth Rate Fluctuations

Since the coefficient of variation of both sizes and growth rates are around 0.2, it is reasonable to assume that the sizecontrol mechanism can be expanded around the average size and average growth rate (Amir, 2014; Grilli et al., 2017). This maps the problem into a discrete-time linear response framework. This section generalizes the usual linear response approach to include single-cell growth rates fluctuations, with the goal of disentangling how size and growth rate fluctuations are connected.

We use a discrete-time description, focusing on a single-cell at a given generation i. This cell has initial size V (i) 0 , which corresponds to a log-size at birth q (i) 0 , a growth rate α (i) and divides at a log-size q (i) f after a time τ (i) d . These variables are connected by the relation

$$q\_f^{(i)} = q\_0^{(i)} + \alpha^{(i)} \mathbf{r}\_d^{(i)}.\tag{1}$$

We also indicate the net logarithmic multiplicative growth (sometimes referred to as elongation) by

$$G^{(i)} = \alpha^{(i)} \mathfrak{r}\_d^{(i)} = \log(V\_f^{(i)} / V\_b^{(i)}) = q\_f^{(i)} - q\_0^{(i)} \dots$$

Here, α (i) is a random variable with mean hαi (which depends on the condition) and variance σ 2 α . We assume that individual cells maintain the same growth rate for a whole cell cycle (this approximation corresponds to how growth rate distributions are often evaluated in empirical data). The growth rates of subsequent generations are drawn from a common (Gaussian) distribution, but may be correlated across generations, with Pearson correlation coefficient ρ. An alternative approach to model growth fluctuations would consist in modeling the stochastic fluctuations of α (i) in continuous time (Iyer-Biswas et al., 2014a,b). The more minimalistic discrete-time approach that we employ can be seen as a coarse-graining of an appropriate continuous-time model. We also note that the discrete-time modeling approach is equivalent to a continuoustime description for single-exponential growth (Grilli et al., 2017; Ho et al., 2018), but in presence of a growth rate fluctuating in continuous time, this equivalence does not hold anymore. For the small fluctuations observed in the data, the two approaches are qualitatively equivalent and likely quantitatively not distinguishable given the available statistical power.

We introduce the variable fluctuations

$$\begin{aligned} \delta G^{(i)} &:= G^{(i)} - \langle G \rangle \\ \delta \mathfrak{r}\_d^{(i)} &:= \mathfrak{r}\_d^{(i)} - \langle \mathfrak{r} \rangle \\ \delta \alpha^{(i)} &:= \alpha^{(i)} - \langle \alpha \rangle \\ \delta q\_0^{(i)} &:= q\_0^{(i)} - \langle q\_0 \rangle \ . \end{aligned} \tag{2}$$

Under the above assumptions, the evolution equations for the fluctuations of the three main variables can be written as

$$\begin{split} \frac{\delta G^{(i)}}{\sigma\_{G}} &= -\lambda\_{Gq} \frac{\delta q\_{0}^{(i)}}{\sigma\_{q}} - \lambda\_{G\alpha} \frac{\delta \alpha^{(i)}}{\sigma\_{\alpha}} + \nu\_{G}^{(i)} \\ \frac{\delta \tau\_{d}^{(i)}}{\sigma\_{\sharp}} &= -\lambda\_{\tau q} \frac{\delta q\_{0}^{(i)}}{\sigma\_{q}} - \lambda\_{\tau \alpha} \frac{\delta \alpha^{(i)}}{\sigma\_{\alpha}} + \nu\_{\tau}^{(i)} \end{split} \tag{3}$$

where λXY are linear coupling parameters (analogous to "susceptibilities") between measured variables. Such linear couplings are central in comparing this theory to experimental data. ν (i) τ and ν (i) G are zero-mean random variables representing cell-to-cell variability ("noise"). Previous approaches neglecting cell size fluctuations assume δα = 0. In order to be solved, Equation 3 requires the specification of how the fluctuations δα(i) are related to the fluctuations of logarithmic size δq (i) 0 . We assume that generally the individual growth rate may depend on the initial size, i.e., within a linear-response framework,

$$\frac{\delta\alpha^{(i)}}{\sigma\_{\alpha}} = -\lambda\_{\alpha q}\frac{\delta q\_0^{(i)}}{\sigma\_q} + \upsilon\_{\alpha}^{(i)}.\tag{4}$$

The values of these coupling constants can be associated to specific mechanisms of size control (Amir, 2014; Grilli et al., 2017). For instance a timer corresponds to λτq = 0 and λτ α = 0, while a sizer mechanism implies λGqσG/σq = 1 and λGα = 0. A standard adder corresponds instead to λGqσG/σq = 1/2 (Amir, 2014; Grilli et al., 2017).

This general framework can be used to evaluate the direct couplings between fluctuations of birth size, growth rate and division time. In particular, all the couplings λs can be evaluated from the Pearson correlations of these variables or from the slope of conditional averages. The two methods are equivalent if the slopes are estimated using least-squares method. Indeed, assuming a relationship y = β ∗ x + noise , and performing a linear regression, the least-square estimate of beta is β = cov(x, y)/var(x). Any discrepancy of the values of the coupling parameters obtained in these two different ways should be considered as a signature of the violation of the linear-response regime (Cadart et al., 2017). For instance, by multiplying both sides of Equation (4) by δq0/σq and averaging over the fluctuations, it is easy to obtain (see Supplementary Information)

$$
\kappa\_{\alpha q} = -\lambda\_{\alpha q} \,, \tag{5}
$$

where cαq is the Pearson correlation between the growth rate α and the log-size q. Similarly, multiplying by δq/σq Equation (3) and averaging gives

$$
\omega\_{\tau q} = -\lambda\_{\tau q} - \lambda\_{\tau a} c\_{aq} = -\lambda\_{\tau q} + \lambda\_{\tau a} \lambda\_{aq} \,. \tag{6}
$$

Finally, multiplying the same equation by δα/σα and averaging over fluctuations gives

$$
\omega\_{\tau a} = -\lambda\_{\tau q} c\_{aq} - \lambda\_{\tau a} = \left. -\lambda\_{\tau a} + \lambda\_{\tau q} \lambda\_{aq} \right. \tag{7}
$$

A detailed comparison between the two ways to evaluate the couplings in empirical data is shown in the Supplementary Information of Cadart et al. (2017). The very small discrepancies found reinforce the ideas that the linear-response framework is effective.

The last two equations clearly show that the presence of a non-zero coupling between single-cell growth rate and birth size induces a non-trivial dependence of the correlations from all the coupling parameters. In other words, when the coupling between division time and cell size is estimated uniquely from cτq (Amir, 2014; Grilli et al., 2017), one is measuring a combination of two different effects: the direct coupling between division time and size at birth and an indirect coupling mediated by growth rate fluctuations.

Direct and indirect effects can be disentangled by solving Equations (6) and (7) to calculate λτq and λτ α in terms of measurable covariances,

$$
\lambda\_{\tau q} = -\frac{c\_{\tau q} - c\_{\alpha q} c\_{\tau \alpha}}{1 - c\_{\alpha q}^2} \,, \tag{8}
$$

and

$$
\lambda\_{\tau\alpha} = -\frac{c\_{\tau\alpha} - c\_{\alpha q}c\_{\tau q}}{1 - c\_{\alpha q}^2} \,. \tag{9}
$$

**Figure 2** illustrates this procedure over one example dataset.

The main assumption of the above modeling framework is that all the relations between the variables are linear. Deviations from linear response are not captured, but in the data these are fairly small (see **Figure 2**). Additionally, it is difficult to evaluate precisely the non-linear behavior from data because the statistics decreases radically for large fluctuations where such non-linear trends are stronger. Given the statistical power available from existing datasets, it is typically not possible to discriminate nonlinear terms from stochastic fluctuations (Grilli et al., 2017).

It has previously been noted that measured Pearson correlations between observed variables in such experiments are very sensitive to (experimental and biological) noise (Eun et al., 2018), so that it is not simple to reconstruct mechanisms from correlations. Our approach, which directly tackles this issue, has two main advantages. First, it subtracts the contributions of indirect correlations and only measures direct couplings. Second, the coupling constants are insensitive to noise, because they are defined as linear slopes of conditional averages, and therefore robust to experimental noise. This is supported by the analysis of Cadart et al. (2017) mentioned above, verifying the equivalence between theory coupling constants measured directly and their expressions in terms of covariances. Indeed, this equivalence not only supports the validity of the linear regime, but also reinforces the idea that scatter in the data (due e.g., to experimental noise) does not bias too much the coupling constants, when they are estimated from the covariances.

We also note that the choice of parameters in Equations (3) and (4) is not unique. There are several alternative (equivalent) choices of fluctuating variables, linked by the condition that G =

fluctuations are not independent (C), the slopes observed in (A,B) are determined by a combination of the coupling between size fluctuation and interdivision time (with strength λτq, see Equation 2) and the effect of growth rate fluctuation on interdivision time (with strength λτ α, see Equation 2). The method described in the text disentangles direct and indirect effect on the slopes of (A–C) to obtain the direct couplings between τ , α, and q, which are shown in (D–F), for different datasets (symbols and colors) and conditions. Interestingly, all the possible couplings are non-zero. Current models do not describe these nontrivial correlations, which poses a challenge to future models.

ατ : α, τ , q; α,G, q; G, τ , q. Each triplet carries a set of linearresponse coupling constants between pairs of variables that can be mapped to the ones defined above.

### 4. RESULTS

### 4.1. All the Possible Couplings Between Measured Variables Are Non-zero in Empirical Datasets

One might expect that not all the possible couplings are different from zero, supporting the simplifying assumptions made by most studies. **Figure 2** shows instead that all the independent couplings λτq, λτ α, and λαq are in fact not negligible.

The existence of nonzero couplings between division time and growth rate λτ α or, equivalently, elongation and growth rate λGα (see next section and **Figure 3**) suggests that cell size control depends in a non-trivial way on growth rate fluctuations. Moreover, the existence of a coupling between growth rates and cell size (measured by the coefficient λαq), affects the observed correlation between division time and size at birth cτq (see Equation 6), or, equivalently, the slope of the sizegrowth plot cGq (Skotheim, 2013), which is normally used to evaluate the strength of the size control. Different values of these correlations are associated with different strengths of cell size control (Skotheim, 2013; Osella et al., 2014), and, more indirectly, to different cell size control mechanisms (Amir, 2014; Grilli et al., 2017). The presence of non-zero couplings between interdivision time and growth rate and between growth rate and size effectively reduces or increases the observed strength of homeostasis, and might be a signature of the mechanisms effecting such control.

### 4.2. The Linear-Response Framework Allows to Define Consistency Relations

As discussed above, the framework provides equivalent alternative descriptions for the interdivision time τ or the net growth G as state variables. These two parameter settings are equivalent under the linear-control assumption. This section deals with the mapping between these parametrizations, and reports consistency tests defined by the expected relations using empirical data.

The consistency relations can be obtained from the relation G (i) = α (i) τ (i) , which, under the linearization assumption reduces to (see Supplementary Information)

$$
\delta G^{(i)} = \langle \mathfrak{r} \rangle \delta \alpha^{(i)} + \langle \alpha \rangle \delta \mathfrak{r}^{(i)}.\tag{10}
$$

Using this relation together with Equation (3), one obtains (see Supplementary Information)

$$
\lambda\_{Gq} = \langle \alpha \rangle \lambda\_{\tau q} \frac{\sigma\_t}{\sigma\_G} \tag{11}
$$

and

$$
\lambda\_{G\alpha} = \langle \alpha \rangle \lambda\_{\tau\alpha} \frac{\sigma\_\sharp}{\sigma\_G} + \langle \tau \rangle \frac{\sigma\_\alpha}{\sigma\_G} \,. \tag{12}
$$

An important consequence of Equation (12) is that a nonzero λGα (or λτ α) could appear even if the corresponding coupling λτ α (or λGα) was null. For instance, if the coupling λτ α is negligible, indicating that the interdivision time is independent from the single-cell growth rate, one would still observe a nonzero λGα. In fact, in that case, the division time would be independent of the growth rate and therefore, in the same amount of time, cells growing faster would grow more than slow-growing cells. Interestingly, this is not the case. Both λτ α (see **Figure 2**) and λGα

FIGURE 3 | Alternative parameterizations of cell-size control linear response are consistent across different experimental datasets and conditions (symbols and colors are the same as in Figure 1). One can parametrize the effect of size and growth rate fluctuations on interdivision time, or, alternatively, their effect on the final size (or the elongation G). (A,C) Show the effect of log-size fluctuation on elongation (with strength λGq, see Equation 3) and the effect of growth rate fluctuation on interdivision time (with strength λGα ). The linear framework predicts consistency relations between the alternative parameterizations (B,D), which are verified across datasets and conditions in the plots in the right panel.

are in fact different from zero (see **Figure 2**), suggesting that some specific mechanism must be at play in determining their value.

**Figure 3** shows that these consistency relations are verified in the data. This result shows that the generalization of the linear-response framework to include growth rate fluctuations can correctly describe the data. As previously discussed, small deviations from these relations could be due to the presence of non-linearities that cannot be captured by the theoretical framework. However, non-linearities can be appreciated only in presence of large fluctuations (Grilli et al., 2017). Also note that while these deviations could be the signature of not-yetunderstood biological mechanisms, they could also come from measurement, segmentation or tracking errors associated to specific experimental and image-analysis protocols.

### 4.3. Growth Rate Fluctuations Can Confound Standard Tests of Division Control Mechanisms, the Case of "Grower" vs. "Timer"

A standard simple analysis to understand what kind of cell division control is in place is based on size-growth plots (Skotheim, 2013) and more generally on measures of correlation between cell elongation and initial size. Since growth rate fluctuations are typically neglected, a value of λGq = 0, i.e., no correlation between elongation and size at birth (**Figure 4A**), is usually interpreted as the signature of a "timer." However, in presence of growth rate fluctuations the value of this coupling parameter alone is not enough to pinpoint a single mechanism. In fact, λGq = 0 can be the actual result of a timer, i.e., λτq = 0 which then requires that λτ α = 0 (**Figure 4B**). However, the relevant variable that is actually uncoupled to size could be the net growth (elongation) itself, rather than the doubling time and the growth rate. In this case, one has a "grower" (**Figure 4C**), for which again λGq = 0 and λGα = 0. For a grower, λτq = 0, but λτ α = −hτ iσα/(hαiστ ). Importantly, if λαq 6= 0, then the correlation between τ and q (cτq) is also different from zero in a grower. This illustrative example emphasizes the importance of keeping into account growth rate fluctuations in order to correctly interpret empirical correlations patterns. It also shows how the linear-response framework proposed here can be used to delineate the correct interpretation. A more general decomposition of contribution of growth vs. timing coupling to size control is discussed in Cadart et al. (2017).

We now proceed to describe two situations where our analysis points directly to specific mechanisms contributing to the complexity of the observed correlation patterns.

### 4.4. Mother-Daughter Growth Rate Correlations Explain the Prevalent Negative Coupling Between Growth Rates and Cell Size Observed in Data

Since the framework describes exhaustively the measured correlations in the datasets, it can be used to explore signatures of mechanisms by which division is coupled to size. This section provides a rationale for the emergence of a non-zero coupling between single-cell growth rate and size at birth λαq observed in most datasets. In particular, we show that a negative λαq may emerge from the presence of a correlation between the growth rates across generations.

**Figure 5** shows the Pearson correlation ρ between the mother and the daughter(s) growth rates in several data sets. We computed the correlation using both the daughters (when available in the dataset). Some previous studies have reported, or assumed, small mother-daughter correlations in growth rate (Wang et al., 2010; Eun et al., 2018), while others have found larger ones (Wallden et al., 2016). An important result of the present analysis is that this correlation is never negligible in the analyzed data sets and appears very conserved across conditions and experimental setup. Specifically, all the datasets considered here show a significantly positive correlation, which does not seem to depend strongly on the condition and takes values around 0.3. Variations of single-cell growth rates can be both due to an external origin (e.g., local variation of nutrients in an agar experiment) or inherent to the mechanisms of growth (e.g., from fluctuations in gene expression Iyer-Biswas et al., 2014a). In both cases, one may expect the presence of positive correlation across generations.

We can model mother-daughter correlations in growth rate as a linear constraint on the single-cell growth rates over generations

$$
\delta \alpha^{(i+1)} = \rho \delta \alpha^{(i)} + \sigma\_{\alpha} \upsilon\_{\alpha}^{(i)} \,. \tag{13}
$$

This equation substitutes Equation (4), and must be considered together with Equation 3. The question we address is whether

2 q = 0.2.

inter-generational correlations give rise to a non-zero coupling λαq between single-cell growth rate and size.

from numerical simulations of the corresponding linearized model with σα/hαi = 0.2 and σ

Since we are considering a fully linear framework, a non-zero coupling λαq can be directly mapped into a non-zero correlation between single-cell growth rate and size cαq via Equation (5). We find that a non-zero correlation cαq emerges in presence of a nonzero correlation of growth rate across generation and a non-zero coupling λGα. In the Supplementary Information we obtain the relation

$$
\lambda\_{\alpha q} = -c\_{\alpha q} = \rho \frac{\lambda\_{G\alpha}}{1 - \rho(1 - \lambda\_{Gq}\sqrt{\frac{\sigma\_G^2}{\sigma\_q^2}})} \sqrt{\frac{\sigma\_G^2}{\sigma\_q^2}} \,. \tag{14}
$$

From the above equation we observe that, for positive ρ and negative λGα (the case of empirical data) only negative couplings λαq can emerge. Hence, mother-daughter growth rate correlations provide a rationale for this widespread negative coupling.

**Figure 5** shows that the prediction of Equation (14) quantitatively reproduces the values of λαq for the datasets that display negative values of this coupling. On the other hand, some datasets have weakly positive couplings λαq, which cannot be reproduced via this mechanism. This suggests that some other unknown mechanisms might play a role in coupling growth rates and sizes.

### 4.5. An Adder Model Captures the Experimental Correlation Patterns Only If Added Size Depends on Single-Cell Growth Rate

This section explores how the observed negative values for the coupling λGα emerge from a simple extension of the adder model that explicitly includes growth rate fluctuations. According to the adder model (Amir, 2014; Taheri-Araghi et al., 2015), the size at division Vf is given by

$$V\_f = (V\_0 + \Delta)e^v,\tag{15}$$

where 1 is the added size and ν is a noise term. We also know that the average size hV0i has an exponential dependence on the average growth rate hαi, which is often referred to as Schaechter's law (Schaechter et al., 1958; Taheri-Araghi et al., 2015; Kennard et al., 2016; Si et al., 2017). Since for steadily dividing cells the average size at birth is equal to the added size, the exponential dependence of the average size on the average growth rate also implies that the average added size must have an exponential dependence on the average growth rate.

Equation (15) does not take into account fluctuations of the growth rate, and there are several possibilities to extend it. Schaechter's law can be seen as a constraint, as the added size averaged over fluctuations of the growth rate has to scale nearly exponentially with the average growth rate. Two possible extreme scenarios can be disentangled with the data. The first scenario assumes that in a given growth condition the added size does not depend on the fluctuation of the growth rate, i.e.,

$$
\Delta = \text{S}\_0 \exp(\langle \alpha \rangle T) \,. \tag{16}
$$

In such case, two cells with different individual growth rates, but growing in the same condition, will add (on average) the same size. In the opposite extreme scenario, the added size depends on the single-cell growth rate, following Schaecter's law even for small fluctuations

$$
\Delta = S\_0 \exp(\alpha T) \,. \tag{17}
$$

In this second scenario, cells in the same growth conditions growing at different growth rates will add on average different sizes. Conversely, cells growing under different conditions but having the same growth rate due to single-cell fluctuations will B).

mechanism correctly predicts negative values of λαq, but the prediction cannot capture the datasets with positive values of λαq, which are shown in (C) (not visible in

add (on average) the same size. Since fluctuations of α are small, in both cases the average added size h1i will scale exponentially with the average growth rate.

This question was approached in Kennard et al. (2016), but the linear-response theoretical framework gives us the tools to reformulate it more precisely. **Figure 6** shows how we proceed to test the two scenarios. For each dataset we estimated the parameter T (which can be related to the duration of the replication/segregation ("C+D") period Donachie, 1968; Zheng et al., 2016; Si et al., 2017) from hlog 1i averaged across cells within each condition. Note that for both Equations (16) and (17), the slope of hlog 1i vs hαi is T. For each condition, one can then compare the relation between the (logarithmic) added size and the individual growth rate. Under the first scenario, the relation should show no dependency (a flat slope). Under the second scenario, one should instead expect a linear relation with slope T (where T is the mean timing from replication initiation to division, assumed to be nearly constant across conditions).

Under the assumption that fluctuations are small, one can linearize the models and compare the coupling coefficients. In the first case, one obtains λGα = 0, while in the second case a non-trivial coupling arises with value

$$
\lambda\_{G\alpha} = \frac{\sigma\_{\alpha}}{\sigma\_{G}} - \frac{T}{2}.\tag{18}
$$

**Figure 6** compares this prediction with the empirical values of this coupling parameter measured across datasets, obtaining a remarkable agreement. It is important to stress that the parameter T, which in turn determines the value of λGα, has been measured only using the dependence of average added size on the average growth rate across different conditions, without considering single-cell data. This procedure assumes that the parameters S0 and T are constant within each dataset across different conditions (Donachie, 1968; Si et al., 2017). Variation around Schaecter's law, which could arise for many reasons (Zheng et al., 2016; Si et al., 2017), can affect the values of these parameters. Precise estimates, or measurements of T can improve the prediction of Equation (18). Moreover, the variability around

FIGURE 6 | Dependency of added size on single-cell growth rates and effect on the values of λGα . (A) Shows that the average logarithmic added size is linear when plotted against the average growth rate (for all the conditions from Wallden et al., 2016). We consider two extreme scenarios for the logaritmic added size of individual cells: it could depend only on the average growth rate (gold line) or it could linearly grow with the single-cell growth rate (purple line), with the same slope observed in (A). (B) Shows the results of the two scenarios for one condition (intermediate growth rate data set from Wallden et al., 2016). These two scenarios translate into two different predictions for the value of λGα , which are reported in panel (C) and tested in panel (D). The non-zero values of λGα are therefore a consequences of the dependence of the added size on the the single-cell growth rate.

the adder mechanism (visible in the values of λGq plotted in **Figure 3**) also affects the predictions.

### 5. DISCUSSION AND CONCLUSIONS

With a long list of recent published studies, the current literature on cell division control remains fragmented in terms of conceptual tools and conclusions on the data (Osella et al., 2017).

This work proposes an extension of the linear-response framework to include growth rate fluctuations, and at the same time provides a comparative meta-analysis of different datasets. As such, we hope it could be useful to the community. The results lead us to propose the linear-response model as a basic common data-analysis framework for the recent (and future) growth division data obtained on bacteria.

Analyzing the experimental data, we found two important results. First, the linear response framework including growth rate fluctuations provides a remarkably good description of the data. Couplings are in all cases only moderately nonlinear (Grilli et al., 2017), and all the predicted self-consistency relations are verified to a satisfactory degree. This carries the important lesson that models where nonlinear couplings are involved might be true, but they would be hard to test with the currently available data. Second, we show that the linear correlation patterns store a wealth of information that is not yet fully appreciated. In particular, single-cell growth rates seem to have a complex effect on cell cycle progression and cell-division control.

The study of the linear coupling parameters between growth rate and size reveals an unexpected variability across experiments. This variability depends on the growth conditions and on the measured mother-daughter correlations in growth rate. Thus, the effect of growth rate fluctuations on cell size varies from one experiment to the other and appears, in some cases, to contribute to size correction. In contrast, the nearadder correlation for the net growth G and the dependency of interdivision time on initial size are robustly observed in all datasets, although quantitatively, they also vary across conditions.

Since all these linear-response effects are clearly testable, we believe that future efforts should be funneled into explaining these patterns comprehensively. Each particular model for a size control mechanism carries a set of constraints for the linearcoupling parameters, which are easily testable in the data. Such testing involves all the measured variables and correlations. Hence, our framework provides a comprehensive and stringent tool to test proposed mechanisms.

The question of specifying the mechanisms for cell division and cell-cycle progression remains largely open, and therefore could benefit from such a systematic tool. One main limitation of current studies is that they discard relevant information and only consider a subset of correlation patterns. Conversely, linear-coupling parameters allow to control all the measured correlations at once. In this work, we provided a proof of principle that this plan of action can be effective by studying two simple examples: (i) the role of mother-daughter correlations in size control, and (ii) the role of growth rate fluctuation in determining the added size.

In particular, we find that growth rate fluctuations are not consistent with an adder model where the added size does not depend on the single-cell growth rate. This could be an important clue on the mechanism underlying the adder behavior and could

### REFERENCES

Adiciptaningrum, A., Osella, M., Moolman, M. C., Cosentino Lagomarsino, M., and Tans, S. J. (2015). Stochasticity and homeostasis in the e. coli replication and division cycle. Sci. Rep. 5:18261. doi: 10.1038/srep18261

help selecting one of the currently proposed models (Ho and Amir, 2015; Taheri-Araghi et al., 2015; Harris and Theriot, 2016; Wallden et al., 2016).

Finally, we address the coherence between different experimental data sets. The discrepancies we found could be due to biological factors, but also to possible sources of experimental noise and systematic errors in the data. Current experimental and data-analysis pipelines are different in both the segmentation and tracking steps. This problem is largely disregarded and in particular we know very little about possible induced biases in the measured correlations, which would lead to different conclusions on the mechanistic aspects. As end-users of the data, we propose our pipeline, i.e., comparing all measurable linear correlations as a useful testing ground for different datasets. In our analysis, we have noticed different trends depending on the experimental laboratory, and growth device, as well as outlier points where possibly the growth conditions were not steady or the analysis had problems. Interestingly, while some observations are robust, others are more erratic across datasets. Beyond this, we can just encourage the experimentalists to adopt common protocols and shared comparison pipelines in future studies. It would be sufficient to define a well defined "reference experiment" with fixed strain, growth device, and growth conditions that each data set needs to share, in order to make precise comparisons possible.

### AUTHOR CONTRIBUTIONS

MCL and JG planned the study. JG, MCL, and MO developed the modeling framework, JG, CC, and GM, analyzed data. JG and MCL outlined and wrote the paper. MO, CC, and GM contributed with scientific discussion and critical reading.

### ACKNOWLEDGMENTS

This work was supported by the International Human Frontier Science Program Organization, grant RGY0070/2014. We acknowledge useful discussions with A. Amir. MO is supported by the Departments of Excellence Grant awarded by the Italian Ministry of Education, University and Research (MIUR) (L. 232/2016). GM was supported by grant no. 31003A 169978 from the Swiss National Science Foundation to Martin Ackermann. JG was supported by an Omidyar Postdoctoral Fellowship at the Santa Fe Institute.

### SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb. 2018.01541/full#supplementary-material

growth and organization. Curr. Biol. 27, 3367–3374.e7. doi: 10.1016/j.cub.2017. 09.046

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Copyright © 2018 Grilli, Cadart, Micali, Osella and Cosentino Lagomarsino. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Phenotypic Heterogeneity in Sugar Utilization by *E. coli* Is Generated by Stochastic Dispersal of the General PTS Protein EI from Polar Clusters

Sutharsan Govindarajan† , Nitsan Albocher † , Tamar Szoke, Anat Nussbaum-Shochat and Orna Amster-Choder\*

*Department of Microbiology and Molecular Genetics, IMRIC, Faculty of Medicine, The Hebrew University of Jerusalem, Jerusalem, Israel*

### *Edited by:*

*Ariel Amir, Harvard University, United States*

### *Reviewed by:*

*Simon Ringgaard, Max Planck Society (MPG), Germany Thomas Shimizu, Fundamental Research on Matter Institute for Atomic and Molecular Physics (NWO), Netherlands*

> *\*Correspondence: Orna Amster-Choder<EMAIL_ADDRESS>† These authors have contributed equally to this work.*

### *Specialty section:*

*This article was submitted to Microbial Physiology and Metabolism, a section of the journal Frontiers in Microbiology*

> *Received: 03 August 2017 Accepted: 26 December 2017 Published: 17 January 2018*

### *Citation:*

*Govindarajan S, Albocher N, Szoke T, Nussbaum-Shochat A and Amster-Choder O (2018) Phenotypic Heterogeneity in Sugar Utilization by E. coli Is Generated by Stochastic Dispersal of the General PTS Protein EI from Polar Clusters. Front. Microbiol. 8:2695. doi: 10.3389/fmicb.2017.02695* Although the list of proteins that localize to the bacterial cell poles is constantly growing, little is known about their temporal behavior. EI, a major protein of the phosphotransferase system (PTS) that regulates sugar uptake and metabolism in bacteria, was shown to form clusters at the *Escherichia coli* cell poles. We monitored the localization of EI clusters, as well as diffuse molecules, in space and time during the lifetime of *E. coli* cells. We show that EI distribution and cluster dynamics varies among cells in a population, and that the cluster speed inversely correlates with cluster size. In growing cells, EI is not assembled into clusters in almost 40% of the cells, and the clusters in most remaining cells dynamically relocate within the pole region or between the poles. In non-growing cells, the fraction of cells that contain EI clusters is significantly higher, and dispersal of these clusters is often observed shortly after exiting quiescence. Later, during growth, EI clusters stochastically re-form by assembly of pre-existing dispersed molecules at random time points. Using a fluorescent glucose analog, we found that EI function inversely correlates with clustering and with cluster size. Thus, activity is exerted by dispersed EI molecules, whereas the polar clusters serve as a reservoir of molecules ready to act when needed. Taken together our findings highlight the spatiotemporal distribution of EI as a novel layer of regulation that contributes to the population phenotypic heterogeneity with regard to sugar metabolism, seemingly conferring a survival benefit.

Keywords: bacterial polarity, cell poles, proteins localization, dynamic localization, PTS system, general PTS proteins, phenotypic heterogeneity

# INTRODUCTION

The poles of rod-shaped bacterial cells play an important role in various molecular processes, including DNA segregation, metabolic regulation and aggregate clearance (Bowman et al., 2011; Govindarajan et al., 2012; Laloux and Jacobs-Wagner, 2014). The poles are also emerging as hubs for localization of many proteins, including sensory systems (Alley et al., 1992; Janakiraman and Goldberg, 2004; Briegel et al., 2009; Rudner and Losick, 2010; Amster-Choder, 2011). Currently, a handful of mechanisms have been suggested to underlie targeting of proteins to the poles. These include interaction with the polar anionic lipid cardiolipin, recognition of membrane domains with strong negative curvature, present in the poles and septa, and self-assembly in nucleoid-free spaces (reviewed in Govindarajan et al., 2012; Laloux and Jacobs-Wagner, 2014; Treuner-Lange and Sogaard-Andersen, 2014).

Notably, while the list of proteins that localize to the poles is constantly growing, less is known about their dynamics within or between the poles. Several polar complexes were shown to exhibit dynamic behavior, relocating from pole to pole, from pole to midcell, or from lateral sites to the poles, with their spatiotemporal dynamics often linked to developmental pathways (Jensen et al., 2002). A classic example for proteins that exhibit pole-topole dynamics is the Escherichia coli MinCDE complex, which negatively regulates FtsZ polymerization at the poles and restricts Z-ring formation to mid-cell (Lutkenhaus, 2007). In Caulobacter crescentus, the polarity-establishing proteins PopZ and TipN were shown to relocate from one pole to the other in a cell cycledependent manner, such that they mark the old and new pole, respectively (Huitema et al., 2006; Lam et al., 2006; Bowman et al., 2008; Ebersbach et al., 2008). FtsZ, which assembles into a ring in mid-cell during division, was recently shown to relocate to the cell poles of non-dividing cells and re-assemble there (Yu et al., 2017). Dynamic localization from lateral sites to the poles was observed for chemoreceptor complexes in E. coli (Thiem and Sourjik, 2008).

The bacterial phosphotransferase (PTS) system controls hierarchical uptake and utilization of preferred carbohydrates from complex environments (Deutscher et al., 2006). Additionally, the PTS controls other pathways, including carbon catabolite repression and inducer exclusion (Deutscher et al., 2006). The spatial organization of the PTS system has been studied in our lab. We found that the control center of the PTS system, which comprises the major PTS proteins EI and HPr, localizes to the cell poles of E. coli. Although EI and HPr interact with each other, polar localization of each of them occurs independently of the other (Lopian et al., 2010). Moreover, activation of the PTS system was shown to affect the localization of HPr but not of EI. Our recent study, which focused on identifying the mechanism of EI localization revealed that a geometric cue is important for EI targeting (Govindarajan et al., 2013). Thus, similar to DivIVA, a B. subtilis negative membrane curvature sensor protein, EI localizes to regions of strong negative curvature in the membrane, which are usually present in the poles and septa (Govindarajan et al., 2013). However, unlike DivIVA, which can directly sense the membrane curvature through its membrane-binding α-helical domain (Lenarcic et al., 2009; Ramamurthi and Losick, 2009; Oliva et al., 2010), the soluble EI protein was suggested to localize to these regions via other, yet unknown, proteins that sense membrane curvature.

In this study, we employed time-lapse fluorescence microscopy in live cells to explore the temporal organization of EI in growing and quiescent E. coli cells. We show that polar EI clusters are often dynamic and that their dynamic range differs among cells in the population, with their speed negatively correlating with cluster size. EI dynamics is energy-dependent, since it is negatively affected by inhibition of cell metabolism. EI cluster dynamics does not depend on the type of sugar, whether it is PTS or non-PTS. However, regardless of the sugar source, EI clustering inversely correlates with its function, that is, EI has a higher capacity to be active in cells in which it is uniformly distributed, and its higher-order assembly into clusters prevents its activity. In line with this, during transition from inactive to active state of growth, EI molecules disperse out of the cluster in a significant number of cells in a population. Intriguingly, EI cluster formation is an event that is stochastic in time, which generates phenotypic heterogeneity within a population.

# RESULTS

## EI Clusters Exhibit Several Dynamic Localization Patterns that Are Energy Dependent

The general PTS protein EI has been shown to form clusters that localize mainly to the poles or to mid-cell (Lopian et al., 2010; Govindarajan et al., 2013). However, the temporal behavior of these clusters, as well as that of EI molecules that are not associated with clusters has not been characterized. We first set out to monitor the spatiotemporal localization of EI in actively growing cells. For this purpose, we monitored exponentially growing cells, which express EI fused to mCherry, as well as ZapA—a marker for septal location (Galli and Gerdes, 2010) fused to GFP, both expressed from the native chromosomal loci under the control of their respective promoters, by timelapse microscopy. First, we verified that the activity of the chromosome-encoded EI-mCherry protein is comparable to that of the wild-type protein by comparing the growth rate of the strains expressing EI-mCherry or EI in minimal medium supplemented with PTS sugars (glucose or sorbitol) or with a non-PTS sugar (lactose). The results in Figure S1 show that the growth of the two strains on all sugar sources is alike, indicating that EI-mCherry has a similar activity to EI. Next, cells expressing EI-mCherry and ZapA-GFP were grown in minimal medium supplemented with glucose till they reached OD600 = 0.2, and images were acquired every minute for 10 min Based on the presence and the spatiotemporal localization patterns exhibited by EI, we divided the cells to four groups (see **Figures 1A,B** and Movie S1): (a) Non-dynamic (non-Dyn), that is, the EI cluster was observed in one of the poles without detectable movement (9%); (b) Dynamic within one pole (Dyn-1P), that is, the EI cluster was dynamic within the zone of one pole (43%); (c) Dynamic outside the zone of one pole (Dyn-OP), that is, the EI cluster switched from one pole to another or between the pole and mid-cell (11%); (d) EI clusters were not detected (UC), but, rather, EI molecules appeared as freely diffuse in the cytoplasm (37%) (**Figures 1A,B**). Two-dimensional trajectories of the movement tracks, presented in **Figure 1C**, show that, except for cells that exhibited a non-dynamic pattern (non-Dyn), in all cells exhibiting dynamic localization patterns (Dyn-1P and Dyn-OP), EI moved in a zigzag manner, that is, in jagged lines with sharp turns.

Heat map presentation of EI-mCherry fluorescence intensity, made for representative cells at time 0, shows an inverse correlation between the EI cluster intensity and its dynamics (**Figure 1A**, last panel). Thus, the intensity of non-Dyn clusters is

EI-mCherry in each of the cells. The contour of each cell is outlined. (B) Pie chart showing the fractions of cells that exhibit the different patterns of EI dynamics, i.e., non-dynamic (non-Dyn, green), dynamic within one pole (Dyn-1P, red), dynamic from pole to midcell or from pole to pole (Dyn-OP, orange), and undetectable cluster (UC, blue), in populations grown to early log in minimal medium supplemented with glucose as carbon source. The Standard division was between 1 and 2%. (C) Two-dimensional trajectories of EI-mCherry cluster movement overlaid on the corresponding phase contrast (gray) images. (D) Distribution diagram indicating the distribution of EI-mCherry clusters area (µm2 ) in cells (referred to as size) exhibiting the different patterns of EI dynamics (40 cells from each group). *P*-values were obtained by ordinary one-way ANOVA: \*\*\*\**p*-value < 0.0001. (E) Scatter plot of cluster speed (µm/s) vs. EI-mCherry clusters area (µm2 ), drawn for cells treated (red) or not treated (purple) with DNP. See Materials and Methods for experimental details. Spearman correlation for the untreated ρ = −0.78, *p*-value < 10−8 .

much higher compared to Dyn-1P, and the intensity of the latter is significantly higher than that of Dyn-OP. To quantitatively analyze the relationship between EI dynamics and cluster size at the population level, we plotted the area of EI-mCherry clusters (will be referred to as cluster size) for each localization pattern in the corresponding population (**Figure 1D**, n = 120). The results demonstrate that non-Dyn clusters have the largest area and are non-mobile, Dyn-1P clusters are smaller and they move within a limited zone, and Dyn-OP clusters are the smallest and the most mobile. These results indicate that the dynamics of EI clusters inversely correlates with their size, that is, the bigger the cluster, the more static it is.

Next, we quantitatively analyzed the correlation between the rate of EI clusters movement and their size. To answer this question, we performed time-lapse experiments, during which we acquired images of EI-mCherry localization every 10 s for a total of 5 min. From these images, we calculated the average speed (µm/s) of individual EI-mCherry clusters and plotted it as a function of average EI-mCherry cluster size. The results (**Figure 1E**, purple circles) show that the average speed of EI cluster inversely correlates with its size (Spearman correlation ρ = −0.78, p-value < 10−8 ), that is, clusters with smaller size moved faster than clusters with bigger size.

To figure out whether active cellular metabolism is required for EI motion, we repeated the time-lapse experiment described above with cells treated with Di-nitro phenol (DNP), an inhibitor of cell metabolism (Parry et al., 2014). The results (**Figure 1E**, red circles) show that movement of EI was restrained in DNP-treated cells, as evidenced by the dramatic reduction in cluster speed, and the lack of correlation between the speed and the cluster size.

Together, the data presented in this section suggest that the motion range and speed of EI clusters inversely correlate with their size and that this motion requires energy.

### EI Clusters Are Formed Stochastically in Time by Assembly of Pre-existing Dispersed Molecules

To shed light on the process of EI cluster formation, we asked whether EI clusters form by de novo protein synthesis or by assembly of pre-existing molecules. We assumed that if new clusters nucleate as a result of expression burst, then the amount of EI before cluster formation would be lower than its amount after cluster formation. On the other hand, if new clusters are formed by assembly of pre-existing molecules, the amount of EI before and after cluster formation are expected to be comparable. To distinguish between these possibilities, we followed, by time-lapse microscopy, the process of EI cluster formation in 78 cells, and quantified the mean intensity (MI) and the standard deviation intensity (SDI) of the EI-mCherry fluorescent signal before and after cluster formation (**Figure 2A**). Of note, whereas the MI estimates the amount of EI-mCherry in the cell, the SDI provides information on the distribution of EI-mCherry within the cell i.e., whether or not a cluster has been formed. Moreover, Figure S2 shows that the size of EImCherry clusters directly correlates with cellular EI mCherry SDI (Pearson correlation ρ = 0.83 P-value = 7.89E-34). This implies that the cluster size and the cellular distribution of the EI-mCherry signal directly correlate and, hence, we can use the SDI of the signal as an estimate for cluster size. Thus, cells with low SDI value have diffuse EI molecules, while cells with high SDI value have highly clustered EI. As shown in **Figure 2Aa**, the mean fluorescence intensity of EImCherry of individual cells, before and after cluster formation, are comparable, indicating that the new clusters are assembled from pre-existing molecules. We then plotted the SDI of EImCherry before and after cluster formation, one against the other. The results in **Figure 2Ab** show that the SDI fluorescence intensity values of EI-mCherry after cluster formation were higher than the values before cluster formation in the majority of cells, indicating that new clusters are indeed formed in these cells.

Together, the results in **Figures 2Aa,b** imply that the new clusters are formed from pre-existing molecules. To validate this observation, we tracked the distribution and the concentration of EI-mCherry molecules in 36 individual cells, in which this protein appeared diffuse at time 0, for 60 min. To determine the cell cycle status of the cells, we also measured the area of these cells. The results presented in **Figure 2Ba** show that the growth rate was similar in all the tracked cells. Moreover, the EI-mCherry MI did not change much during cell growth, i.e., was similar in all time points (**Figure 2Bb**). However, in 14 cells, the EI-mCherry SDI abruptly increased at different time points, indicating that new clusters are formed in these cells during the time course of the experiment (**Figure 2Bc**).

The results in **Figures 2A,B** suggest that the freely diffuse EI molecules assemble rapidly and haphazardly to form higher order clusters within the cytoplasm in a manner that is independent of EI-mCherry level, which does not fluctuate. Of note, the results in **Figure 2B** were not obtained with synchronized cells. To test whether there is any linkage between cluster formation and time or cell cycle, we followed newly divided cells (n = 80), in which EI-mCherry molecules appeared diffuse, and monitored the time at which the clusters formed (**Figure 2Ca**). Additionally, the length of the cell, which provides information on the cell cycle status, at the time of cluster formation was also measured (**Figure 2Cb**). As shown in **Figures 2Ca,b**, the time of cluster formation and the length of the cell during cluster formation were very different in these cells, which were synchronized according to their cell cycle, indicating that new EI clusters are formed stochastically during the cell cycle. Of note, the stochasticity is deduced from the fact that we did not find any cellular parameter that is linked to EI clustering. The possibility that EI clusters are formed in response to some unknown intracellular change exist, but this putative change, in turn, is seemingly not linked to cell growth and cell cycle and shows stochasticity in time. Taken together, our results show that E. coli populations are heterogonous with respect to whether and when EI clusters are formed during growth.

### Diffuse EI Molecules Dynamically Incorporate into Polar Clusters

After showing that diffuse EI molecules can assemble to form a new cluster, we asked whether these molecules can join a pre-formed EI cluster. To answer this question we performed Fluorescence Recovery After Photobleaching (FRAP). For this purpose, we photobleached the entire pole area, where the clusters are localized, and the recovery of fluorescence intensity at the pole was measured as a function of time and plotted relative

to the time after bleaching. For the analysis, we chose cells with big EI clusters that were non-dynamic before photobleaching. Furthermore, in order to distinguish between recovery by assembly of pre-existing diffuse molecules and recovery by assembly of newly synthesized molecules, we performed another set of FRAP with cells spotted on an agar pad that contained the translation inhibitor chloramphenicol (CM). The results in **Figure 3A** show that nearly 6% of EI-mCherry fluorescent signal in the pole area was recovered after 600 s in untreated cells (n = 14), whereas the recovery of the fluorescent signal in the CM-treated cells during this time was slightly lower i.e., nearly 4% (n = 12). That is, about one third of the recovered EI-mCherry at the pole were de novo synthesized. Of note, the time it took the polar signal to appear and the pace it developed in the untreated and CM-treated samples were very similar (**Figure 3B**). These results indicate that the majority of

FRAP microscopy images of wild-type *E. coli* cells expressing EI-mCherry in the absence of treatment (upper panel) or in the presence of chloramphenicol (lower panel). Images at selective time points are shown. The contour of the bleached area in the fluorescent images is outlined.

molecules that assemble with the clusters are pre-existing diffuse EI-mCherry molecules, which were present outside the bleached area.

### Growth Phase, but not the Carbon Source in the Growth Medium, Affects Formation of EI Clusters

To determine whether dynamic localization of EI depends on the presence of PTS sugars, we followed EI cluster dynamics in minimal medium supplemented with either glucose or sorbitol, which are PTS sugars, or with lactose, which is a non-PTS sugar, as a sole carbon source, and quantified the fraction of cells in each of the groups that exhibit the different patterns of EI dynamic localization. The results in **Figure 4A** show that EI clusters were similarly distributed between the groups when cells were grown in glucose or sorbitol containing media, and were only slightly and not significantly different in cells grown in lactose containing media. Hence, the spatiotemporal behavior of EI does not depend on the type of sugar in the growth medium and does not really change between cells grown in the presence of PTS and non-PTS sugars.

Next, we asked whether formation of EI cluster depends on the growth state of the cells. To this end, we grew EI-mCherryexpressing cells in non-PTS medium (minimal medium with lactose) or in PTS medium (minimal medium with glucose or sorbitol) and calculated the percentage of cells with EI clusters at early logarithmic (EL), mid logarithmic (ML), early stationary phase (SP), and late stationary phase (after overnight growth, ON). As shown in **Figure 4B**, in all three growth media, nearly 60% of cells in EL growth phase formed EI clusters, and this percentage increased a bit upon transition from EL to ML (65– 70%) and remained similar during early SP. Conversely, the percentage of cells with EI clusters increased significantly in ON cultures (80–95%), with cells grown in lactose medium containing the highest percentage of clusters. The difference in the percentage of cells with EI clusters in ON cells is not due to an increase in the cellular amount of EI, since Western blot analysis, using α-mCherry antibodies, show that the amount of EI-mCherry present in cells at the different growth phases and in the different growth media, is largely comparable, similar to the housekeeping control protein GroEL (**Figure 4C**). Therefore, the percentage of cells with EI clusters in the population is growth phase-dependent and is the highest in quiescent cells.

FIGURE 4 | Growth phase, but not the carbon source in the growth medium, affects formation of EI clusters. (A) Stacked column chart showing the fractions of cells that exhibit a certain pattern of EI dynamics, i.e., non-dynamic (non-Dyn, green), dynamic within one pole (Dyn-1P, red) and dynamic from pole to midcell or from pole to pole (Dyn-OP, orange), in populations grown to early logarithmic phase in minimal medium supplemented with glucose (Glu, *n* = 328), sorbitol (Sor, *n* = 311) or lactose (Lac, *n* = 286). Cells from three independent experiments were analyzed. Ordinary one-way ANOVA indicated no significant variations among the samples. (B) Column chart of the fraction of cells with detectable clusters in cultures grown on different carbon sources glucose (Glu), sorbitol (Sor) or lactose (Lac) at different growth phases: early log (EL; OD600 0.2–0.25), mid-log (ML; OD600 0.5–0.55), early stationary phase (SP; OD600 1.00–1.1) and after overnight growth (ON). More than 400 cells were analyzed from three independent experiments. *P*-values were obtained by ordinary one-way ANOVA: \*\*\**p*-value < 0.001; \*\*\*\**p*-value < 0.0001. (C) Western blot analysis showing the level of EI-mCherry and GroEL proteins, detected by α-mCherry and α-GroEL antibodies, respectively. Cells were grown on different carbon sources to different growth phases, as in (B). Equal amount of samples were collected at specific growth phases, separated on 10% SDS polyacrylamide gel, blotted onto nitrocellulose membrane and detected by the two antibodies.

# The Process of EI Cluster Dispersal Occurs after Exiting Quiescence

The finding that more ON cells contain EI clusters than growing cells (**Figure 4B**) raised the question of what happens during the transition from stationary to logarithmic phase. Do clusters disperse or a fraction of new cells are born without clusters? To address this question, we spotted diluted ON culture on an agar pad with fresh medium and observed the cells for 3 h by timelapse microscopy. At the time that the cells were spotted on the agar pad (time 0), we monitored the distribution (SDI) of the EImCherry fluorescent signal in 269 cells (**Figure 5A**, left panel). Of note, as explained above, cluster size can be inferred from the SDI (Figure S2). From this time point and on, we classified the cells to cells having or not having detectable clusters at each time point till they divided (**Figure 5A**, middle panels). Together, the results in the left and middle panels of **Figure 5A** show that the population contained the following sub-populations: (i) A relatively small sub-population with no detectable clusters at time 0 (13%), which formed clusters during cell growth in the vast majority of cases. (ii) A very big sub-population of cells with relatively small clusters at time 0 (64%) that dispersed during growth and reformed (35%, see such a cell in **Figure 5B**) or did not disperse (29%) till they divided. (iii) A medium size subpopulation of cells with big clusters at time 0 (23%) that remained throughout the experiment.

To negate the possibility of other factors, such as mechanical perturbations upon shifting cells from liquid media to solid media, contribute to cluster disassembly, we performed a control experiment in which ON cells were spotted on an agar pad made from overnight filtered medium (n = 255). In this setup, none of the EI clusters dispersed and only 4% of the cells divided during the 3 h time-lapse microscopy, suggesting that dispersal of EI clusters is specific to transition from quiescence to active growth (Figure S3).

To better understand the nature of the transition from ON to EL phase regarding EI distribution, we followed EI cluster size by tracking the SDI of 10 cells from each of the three sub-populations in **Figure 5A**. The results of this analysis are presented in Figure S4. The left panel shows the SDI at time 0, using the same scale for all sub-populations. The middle panel shows whether and when the cells had a detectable cluster at each time point till they divided. The right panel shows the SDI at each time point till division, using a different scale for each sub population. By and large, the results show a common trend for each of the first two groups. The SDI of cells in sub-population (i), which had no detectable clusters when they exited quiescence (upper panels; purple shades), remains low till cluster formation, but, upon cluster formation, the SDI doubles, suggesting that cluster formation is a rapid event, although cluster size continues to increase. In cells from sub-population (ii), which had relatively small clusters when exiting quiescence (middle panels; blue shades), there is a moderate decrease in SDI values, even in cells with clusters that remained detectable throughout the experiment, suggesting that cluster dispersal is a continuous process, rather than an abrupt event. With time, the SDI value start to increases. The SDI in cells from sub-population (iii), which had big clusters after exiting quiescence (lower panels; green shades), showed fluctuations with no consistent trend, maybe because big clusters are old and do not disperse completely, or because there is more noise in our measurements in this range. Together, these results suggest that, unlike the abrupt event of cluster formation that occurs at random time

FIGURE 5 | Assembly of EI into clusters affects EI function and dynamics (A) Diluted ON cultures were spotted on an agar pad with fresh medium and the cells were observed at time 0, and every 15 min till they divided or for 3 h, by time-lapse microscopy. Left panels: heat map of standard deviation intensity (SDI) of EI-mCherry at time 0. Middle panels: heat map of cells with (turquoise) or without (pink) detectable EI-mCherry clusters over time, from time 0 till cell division (white, time after cell division). Right panels: Heat map of dividing cells with a cluster in each sister cell (blue) or only in one of the sister cells (gray) after 3 h; white rows denote cells that did not divide within the 3 h. (B) The images show a representative cell from the ON population, which has an EI-mCherry cluster at time 0 that completely disperses during growth (see image after 45 min) and later reforms (see image after 100 min). The contour of the cell is outlined. Scale bar corresponds to 1µm. The plot shows the EI-mCherry intensity in this cell vs. the relative cell length at the different time points at which the images were taken (0, 45, and 100 min). (C) The images show a representative cells from the ON population, which express EI-mCherry, that were incubated with the glucose analog 2-NBDG. The mCherry fusion protein (red) and the 2-NBDG (green) were observed by fluorescence microscopy and are shown over the phase contrast images (gray) that were observed with phase microscopy. Scale bar 2µm. The scatter plot presents the 2-NBDG mean intensity (AU) vs. the EI-mCherry SDI (AU) (*n* = 194). Pearson correlation ρ = 0.16; *p*-value = 0.02. (D) The scatter plot shows the mean intensity (MI) of mCerulean expressed from Pnag, vs. the cellular EI-mCherry SDI (AU) at the moment of transition to NAG-continuing medium (0 min; green cycles; *n* = 222) and after an hour (60 min; purple cycles; *n* = 274). Pearson correlation for the "60 min" plot ρ = −0.16; *p*-value = 0.007. Pearson correlation for the "0 min" plot ρ = 0.06; *p*-value = 0.366.

points during the cell cycle, dispersal of EI from the cluster is a gradual process, which most cells are going through when exiting quiescence, although it might not be completed, since in some cells EI molecules start to join the cluster before its complete dispersal. Therefore, the process of cluster assembly seems stochastic in time, whereas disassembly is probably due to active sensing, together generating a heterogeneous population.

Because all cells except for two had detectable EI-mCherry clusters during their lifetime, we tested the possibility that cluster formation promote division by measuring the time of cluster formation relative to the cell cycle duration. The results in Figure S5A show no linkage between cluster formation and the cell cycle stage (n = 138). Moreover, we did not observe a significant difference in cell cycle duration between cells which were born with or without a cluster (Figure S5B). Finally, 15% of the cells which were born without cluster divide without forming a cluster. Therefore, cluster formation does not seem to promote division. The reason that 2 cells out of 269 cells recovering from ON did not form a cluster (**Figure 5A**) is probably due to long cell cycle duration. Still, the possibility that EI cluster formation is linked to cell division can be revisited if EI mutants that do not form clusters become available.

Since most cells divided within the timespan of the experiment (3 h), we also addressed the question of what is the fraction of cells that are born without clusters. Hence, we asked if both sister cells had clusters after cell division or only one. The results in **Figure 5A** show that all cells but two had clusters before division, and that 72% of the sister-cell pairs had only one cell with a cluster (**Figure 5A**, right panel, gray) and in 28% of the cases both sister cells had a cluster (**Figure 5A**, right panel, blue) at the end of the time-lapse experiment, implying that 64% of the cells after doubling of the population had clusters. Taking into consideration that at the beginning of the experiment 87% of the cells had clusters, the percentage of clusters in the population dropped after division. Summarily, two phenomena were observed when closely following cells that emerge from stationary phase: the majority of newly added cells lack detectable clusters, suggesting that all EI molecules are dispersed and active at the beginning of growth, and EI clusters in most cells form and disperse at random times during growth, yielding a phenotypically heterogeneous population.

### Clustering of EI Inversely Correlates with EI Function

The process of EI cluster dispersal and the random events of cluster formation, characterized above, are probably general rather than sugar-specific, as there were more EI clusters in ON cultures than in the EL, ML, and SP cultures, independent of the carbon source (**Figure 4B**). We hypothesized that when the carbon source in the medium is consumed and EI activity is not required, it tends to form clusters that can disperse when a fresh carbon source is supplied and EI needs to function. To explore the possible correlation between EI distribution and function, we incubated overnight cells with 2-NBDG, a fluorescent glucose analog, for 10 min. After washing the cells, images were acquired to quantify the clusters of EI-mCherry and the amount of 2- NBDG in each cell. The scatter plot of EI-mCherry SDI against 2- NBDG mean fluorescence, shown in **Figure 5C**, suggests that EI has a higher capacity to be active (higher 2-NBDG signal) in cells with lower SDI of EI-mCherry, i.e., cells in which EI-mCherry is distributed uniformly or form small clusters. These results imply that EI is less active when it is clustered in higher-order assemblies and, thus, that its function inversely correlates with clustering and with cluster size (Pearson correlation ρ = −0.16; p value = 0.02).

Since 2-NBDG is metabolized to a non-fluorescent product (Yoshioka et al., 1996), we decided to test the correlation between EI distribution and function by another approach. For this purpose, we constructed a strain that contains an EImCherry fusion, expressed from the native chromosomal locus, as well as mCerulean gene under the regulation of N-acetyl glucosamine (NAG) promoter (Pnag). To induce expression from Pnag promoter, phosphorylated NAG needs to be present inside the cell (Yamada and Saier, 1988; Westermayer et al., 2016). Since EI is the first protein in the phosphorylation cascade that enables uptake and phosphorylation of PTS sugars, including NAG (Yamada and Saier, 1988), expression of mCerulean implies that EI is active. The constructed strain, which was initially grown in minimal medium containing glycerol until EL phase, was resuspended in fresh minimal medium containing NAG. Fluorescent images were acquired in order to quantify the expression of EI-mCherry and mCerulean proteins, at time 0 and 60 min. The expression of mCerulean from Pnag (mean intensity of Pnag-mCerulean) as well as the distribution of EI-mCherry in the cell (SDI of EI-mCherry) were plotted against each other. As shown in **Figure 5D**, at time 0 min, at which point cells which were grown only in the presence of glycerol, induction of PnagmCerulean could not be observed. In contrast, expression of Pnag-mCerulean was clearly induced when cells were grown in the presence of NAG for 60 min. More importantly, the scatter plots in **Figure 5D** suggest that EI has a higher capacity to be active (higher mCerulean production) in cells with lower SDI of EI-mCherry i.e., those cells in which EI-mCherry is distributed uniformly or form small clusters. This suggests that EI is relatively less active when it clusters to form higher-order assemblies. Although the possibility of other factors, for example, changes in the uptake or degradation of the sugars, contributing to the observed heterogeneity cannot be completely ruled out, the results obtained with the two assays which tested EI activity in different manners, suggest that EI activity inversely correlates with clustering.

# DISCUSSION

Dynamic organization of proteins in higher-order assemblies, provides a versatile and attractive mechanism for metabolic regulation. In this study, by using time-lapse microscopy as a primary tool, we demonstrate for the first time that the temporal distribution of a central metabolic enzyme, the E. coli general PTS protein EI, is multimodal, and that EI clustering inversely correlates with its enzymatic activity. We also show significant

repositioning of EI molecules from being in a cluster to being in a diffuse state during transition from inactive to active state of growth, highlighting the importance of EI cluster formation for future need of sugar metabolism.

Our data show that EI clusters, observed as mostly polar in snap-shots in previous studies (Lopian et al., 2010; Govindarajan et al., 2013), exhibit rapid and continuous relocation within the pole region, between the two poles or between the pole and midcell. Thus, the EI clusters spend most of the time in the poles or pole-to-be and a small fraction of their time in relocating among or within these sites, explaining their observation as merely polar in snap-shots. Still, a significant fraction of cells (nearly 11%) contain an EI cluster that remains static at one pole. The demonstration that the mobility of newly formed EI clusters is not restricted to a single pole region, whereas that of bigger clusters is restricted only to the pole region till they become nondynamic is in line with our finding that the motion and speed of EI clusters is influenced by their size. It has recently been suggested that dynamics of macromolecular complexes within the bacterial cytoplasm inversely correlates with their size (Parry et al., 2014). Our observation of inverse correlation between the size of the EI clusters and their dynamics is in agreement with this.

We have previously shown that EI clusters are located at or near negatively curved regions, although the EI soluble protein, which cannot be anchored in the membrane, localizes to these regions via other, yet unknown, proteins that sense membrane curvature (Govindarajan et al., 2013). Of note, EI is not the first polar protein reported to form at negatively curved membrane regions and to exhibit dynamics. The membraneanchored chemoreceptors, recently shown to also recognize negative curvature (Strahl et al., 2015), exhibit dynamic cluster positions within the polar zone, which can be more subtle (Studdert and Parkinson, 2005) or less (Thiem et al., 2007). Another protein that senses membrane curvature, but shows dynamics is MreB (Izoré and van den Ent, 2017). Together, these examples highlight the fact that capture of proteins at specific sites, including the poles, is transient, with the degree of transience depending of the nature of the interactions. In light of these findings, dynamics of a soluble protein, such as EI, despite its association with a yet unknown curvature-sensing protein, is less surprising.

EI is the first protein in the PTS system cascade, which controls preferential use of sugars, i.e., preferred use of PTS sugars over non-PTS sugar. Still, EI is important also for uptake of non-PTS sugars, since phosphorylated IIAglc, which receives the phosphate from EI via HPr, is required for the activation of adenylate cyclase, whose product is required for the transcription of non-PTS catabolic operons (Deutscher et al., 2014). Indeed, we show that the nature of the sugar, whether it is PTS or non-PTS, does not affect the mobility of EI clusters. Moreover, our results suggest that EI assembly into clusters and disassembly is a mechanism to control its activity in sugar metabolism, independent of the nature of the sugar. This is deduced from the finding that the majority of cells in overnight cultures that grew in the presence of PTS or non-PTS sugars, which apparently exhausted the sugar source and are not actively growing, contained more EI clusters compared to cells in the other growth phases, suggesting that assembly of EI into clusters reflects a nonfunctional state. We provide a support for this hypothesis by showing that in a significant percentage of cells from overnight cultures (35%), the EI clusters break down and disperse in the cytoplasm shortly after inoculation into a fresh medium. However, at later time points within the cell cycle, new clusters are formed at random time points, in what seem like rapid and non-synchronized events. The negative relationship between EI clustering and sugar uptake is deduced from our demonstration that the intracellular level of a fluorescent glucose analog, as well as the expression of a PTS-dependent promoter inversely correlate with EI clustering. Having said that, disassembly and reassembly of EI clusters were not observed in all cells. It is possible that some cells have an existing pool of disperse active cytoplasmic EI, not detectable in our experiments, whereas EI clusters below our detection capability may exist in other cells, serving as "nucleation centers" for EI higher-order assemblies.

Taken together, our data suggests that EI molecules that are not needed at a given moment join the cluster, from which they can easily disperse when needed. The finding that 63% of the cells in exponential phase have clusters, in addition to diffuse molecules, imply that most cells have more EI than needed under the conditions tested. Since EI is involved in various processes beyond PTS sugar uptake, such as chemotaxis, catabolite repression and utilization of nutrients such as nitrogen and potassium (Lüttmann et al., 2015), our data suggest that the cells anticipate higher need of EI. Hence, the capability of EI clusters to assemble and disassemble contributes to dynamic regulation of the current and expected cellular metabolism.

Many cases of bacterial proteins, whose higher-order assembly hinder their active site or prevent a conformational change required for their activity, have been documented. One example is the CTP synthetase, whose large scale polymerization keeps a subpopulation of the protein in a conformationally restricted form, which can be readily activated (Ingerson-Mahar et al., 2010; Barry et al., 2014). Another example is MurG, an enzyme involved in peptidoglycan subunit synthesis that, when produced in excess, is sequestered to the poles in an inactive form, which, upon need, can be released to the cytoplasm in an active form (Michaelis and Gitai, 2010). More recently, FtsZ was also demonstrated to assemble at the cell poles of non-dividing cells, in the form of a quiescent body, and only upon sensing growth-supporting conditions, the pole-localized FtsZ structures disassemble to form active molecules that engage in cell division (Yu et al., 2017). In all of these cases, the proteins are not degraded when they are not required, but, rather, are kept as a dynamic storage depots, which can be used in the future. Of note, we have previously reported that cluster formation of B. subtilis EI occurs only upon growth in a medium supplemented with mannitol as the carbon source and not in a rich medium, such as LB, in which it appeared diffuse (Govindarajan et al., 2013). Hence, it appears that the dynamic clustering property of EI is a general phenomenon, which holds for Gram positive bacteria as well.

Populations of genetically uniform microorganisms, bacterial and yeast cells alike, were shown to exhibit phenotypic heterogeneity with regard to fitness-determining traits (e.g., Gefen and Balaban, 2009; Holland et al., 2014; Wang et al., 2014). This phenomenon, vastly studied with respect to variability among bacterial cells in their resistance to antibiotics, is considered a population-based strategy that can be beneficial to a growing population. The mechanisms that underlie adaptation to stress by heteroresistance are not fully understood at the single cell level. Population heterogeneity was also studied with respect to metabolism, which continuously adapts to unpredictable environmental changes. Upon nutrient change, including that of carbon source, a homogeneous E. coli population was shown to split into a growing and a non-growing persisters, which was suggested to stem from stochastic variation in metabolic flux (Kotte et al., 2014). We show that bacterial populations are composed of phenotypic subsets with respect to cluster mobility of a major carbon metabolic enzyme, EI, which can be explained by difference in cluster size. However, the level of heterogeneity is revealed when observing single cells in the population. Thus, despite our general observation that more cells in a quiescent population contain EI clusters than in a growing population, we show that cluster formation is a stochastic event that happen at random timings in individual cells and is not tied up to cell cycle. This stochasticity gains further importance in light of our finding that cluster formation inversely correlates with enzyme activity. Recent studies have shown that biochemical processes are inherently stochastic and cause molecule abundances to fluctuate (Elowitz et al., 2002), a phenomenon that can be used to generate distinct phenotypes (Balaban et al., 2004; Losick and Desplan, 2008). Our findings provide a new organizational explanation to enzymatic heterogeneity, which can hold also for other cellular pathways. Moreover, different reported cases of phenotypic heterogeneity might be linked in nature, as suggested by the observation that uptake of sugars and amino acids increases uptake of various antibiotics (Franklin and Godfrey, 1965; Peng et al., 2015).

# EXPERIMENTAL PROCEDURES

### Bacterial Strains and Growth Media

Unless otherwise indicated, overnight cultures were grown in LB supplemented with appropriate antibiotics at 30◦C. M9 medium, supplemented with glucose, sorbitol or lactose, were used for time-lapse microscopy experiments. For inhibition of cell metabolism, Dinitrophenol (DNP) at a concentration of 2 mM, was added to the growth medium. When appropriate, antibiotics were added at the following concentrations: kanamycin (30µg/ml) and chloramphenicol (25µg/ml) (Sigma-Aldrich). MG16558(ptsI-mCherry)strain expressing EI-mCherry from the native chromosomal locus was described previously (Lopian et al., 2010). SUT201, which expresses EI-mCherry and ZapA-GFP from the chromosome, was constructed by P1 transduction of zapA-gfp(::cat) from HC261 (Peters et al., 2011) to MG16558(ptsI-mCherry). NTS102, which expresses EI-mCherry, Psrl-mVenus and PnagmCerulean from the chromosome, was constructed by P1 transduction of MG16558(ptsI-mCherry), linked to KanR, to T1683 (Westermayer et al., 2016).

# Growth Conditions Used for Snap-Shot and Time-Lapse Imaging

AriekacellsR coverslip cell chamber (SC15012) was used for time-lapse imaging. Overnight cultures were diluted 1:100 in fresh M9 medium and cells were allowed to grow at 30◦C, unless otherwise indicated, till they reached OD600 of 0.2. Samples were then spotted on 1% agarose pads, prepared with the respective minimal medium and supplemented with the appropriate sugar which had been pre-equilibrated to 30◦C, and imaged immediately by time-lapse microscopy. For calculating the fraction of each dynamic group, time-lapse images were acquired every 3 min for a total of 1 h. For FRAP microscopy, cells were grown in M9 glucose media as described above, and samples were spotted on 1% M9 glucose agarose pads without or with chloramphenicol (25µg/ml). To measure the correlation of cluster speed vs. EI-mCherry cluster size, cultures grown in fresh M9 glucose medium until OD600 of 0.2 were spotted on 1% M9 glucose agar pad and time-lapse images were acquired every 10 s for a total of 5 min. For calculating the fraction of cells with detectable EI clusters, MG16558(ptsI-mCherry) cells were grown overnight in M9 minimal medium supplemented with casamino acids, vitamin B1 and 0.4% sugar (glucose, sorbitol or lactose). The media was supplemented with kanamycin. Overnight cultures were diluted 1:100 in the respective fresh M9 medium and cultures were grown until specific growth phase. Snap shot were taken from the cultures at OD600 0.2 to 0.25 for early log (EL), 0.5 to 0.55 for mid log (ML), 1 to 1.1 for stationary phase (SP) and after overnight growth (ON). To follow EI localization in overnight-grown cells that were inoculated into fresh medium, MG16558(ptsI-mCherry) cells, which were grown overnight in M9 minimal media supplemented with CAA, vitamin B1, and glucose, were diluted and spotted onto 1% agar pads in the same medium, which had been pre-equilibrated to 30◦C, and imaged immediately by time-lapse microscopy. To measure the correlation of 2-NBDG [(2-(N-(7-Nitrobenz-2-oxa-1,3-diazol-4-yl)Amino)-2-Deoxyglucose] uptake vs. EI-mCherry spatiotemporal organization, MG16558(ptsI-mCherry) cells were grown overnight in M9 minimal medium supplemented with casamino acids, vitamin B1 and 0.4% glucose. Cells were pelleted, washed thrice and incubated for 10 min in M9 minimal medium supplemented with casamino acids, vitamin B1 and 10µM 2-NBDG. After 10 min, cells were washed thrice and imaged in the microscope. To measure the correlation of Pnag mCerulean expression vs. EI-mCherry spatiotemporal organization, NTS102 cells were grown overnight in M63 minimal media supplemented with CCA, vitamin B1 and Glycerol (M63 glycerol) diluted 1:100 into M63 glycerol and grown till OD600 = 0.3–0.35. Cells were pelleted and resuspended in M63 supplemented with vitamin B1 and NAG. Snap shot in phase contrast, mCherry and mCerulean (CFP) channels were taken at time 0 and after 1 h.

# Fluorescence Microscopy

Fluorescence microscopy was carried out as described previously (Govindarajan et al., 2013). For snap-shot imaging, 0.5 ml cells were centrifuged, washed with 1X phosphate buffered saline (PBS) and finally resuspended in 10–100 µl of PBS. Light microscopy was performed on a Nikon Eclipse Ti-E inverted microscope equipped with Perfect Focus System (PFS), CFI PLAN Apochromat DM 100X (numeric aperture 1.45) oil Ph3 objective and ORCA Flash 4 camera (Hamamatsu photonics). Chroma filter cubes were used as follows: ET-GFP for GFP and 2-NBDG, and ET-mCherry for mCherry. Time-lapse imaging was performed with the aid of OKOLAB cage incubator. Unless otherwise indicated, cells were spotted on 1% agar pads containing M9 media supplemented with glucose, sorbitol or lactose. Agar pads were pre-equilibrated to the appropriate temperature and cells were imaged by time-lapse microscopy at the respective temperature. For each experiment, all images were collected using uniform parameters of magnification and exposure. Images were processed using NIS Elements-AR (Nikon). Two-dimensional trajectories of dynamical clusters were obtained by processing the time-lapse images using Fiji in combination with its Manual Tracking plugin (Schindelin et al., 2012). Correlations were calculated using MATLAB. Statistical tests were performed using GraphPad Prism.

For FRAP microscopy, Nikon A1R confocal microscope equipped with Apochromat 60X objective (numeric aperture 1.4) was used. Photobleaching was done over the area in which the EI cluster is localized. Recovery was measured every 10 s for a total period of 10 min. Mean fluorescence intensity was normalized to the total fluorescence intensity for each ROI after bleaching. Images were analyzed using NIS Elements AR module.

### Image Analysis

Unless otherwise indicated, image analysis was done by NIS-Elements Advanced Research (AR) version 4.4 software (Nikon).

### Analyzing of EI Cluster Size

Regions of interest (ROIs) were selected to mark the cluster region of the first image in each time-lapse series as following. First the details of the mCherry-acquired image were enhanced using the "Homogenization" tool. Next, binary layers were created using the "Define Threshold" tool. ROI were then drawn over the binary layers and copied onto the original mCherry image and the area of each ROI was exported. Finally, cells with marked ROI were manually sorted according to their EI dynamic pattern (non-Dyn, Dyn-1P or Dyn-OP). Scatter plot of the cluster area value for each pattern was obtained using GraphPad Prism. The statistical significance of the differences among cells with detectable clusters was determined by ordinary one-way ANOVA test using GraphPad Prism.

### Calculating the Fraction of Each Dynamic Group

Time-lapse of individual cells were manually classified based on their EI dynamic pattern (non-Dyn, Dyn-1P, Dyn-OP or UC). The results of this analysis are presented as a pie chart or a stacked column chart. The statistical significance of the differences among cells with detectable clusters were determined by ordinary one-way ANOVA test using GraphPad Prism.

### Measuring the Correlation between Cluster Speed and Cluster Size

EI-mCherry clusters were marked throughout each time-lapse series by binary layer by using the "Define Threshold" tool. Clusters marked by binary layer were tracked using the "Track Binaries" tool of the NIS-Elements AR "Tracking" module. The binary center was tracked from the end to the beginning without any skip. For other features, we used the default parameters of the program. Values for average cluster speed and cluster area, which were obtained from the Tracking analysis, were plotted and displayed as a scatter plot. Spearman correlation was calculated using MATLAB corr function.

### Measuring the Correlation between Cluster Area and Cell SDI

First, ROIs were selected to mark the cells outline as following. Binary layers were created using the "Define Threshold" tool over the phase contrast images. ROI were then drawn over the binary layers. Next, EI-mCherry clusters were marked by binary layer by using "Define Threshold" tool. Then the binary lyres were copied onto the original mCherry images. Finally, mCherry images data were exported. The SDI of the ROIs, which were drawn over the cell outline, were plotted and correlated with the MI of the binary layers, which were drawn over the cluster, using Microsoft Excel and MATLAB corr function.

### Measuring Cells Arrea EI-mCerry Concentration and Distribution

ROIs were selected to mark the cells outline throughout timelapse series of cells with undetectable clusters as following. ROI were selected on the first phase contrast image in each time-lapse series, using the "Define New ROI" tool of the "Tracking" module. Then the tool "Track Autodetected ROIs" of the "Tracking" module was used to outline the selected cells by appropriate ROI for each time point throughout the time-lapse series. The data were exported from the mCherry-acquired image and plotted using MATLAB custom made script.

### Calculating the Distribution of EI Clusters from Cell Cycle-Synchronized Population

SUT201, which expresses EI-mCherry and ZapA-GFP from the chromosome, was used.  From these images, newly born daughter cells, in which EI clusters were not detected, were manually monitored. The time at which a clear Z-ring was observed was noted and this was considered time 0 for the newly born daughter cells. The time-lapse series of individual cells were monitored and the time at which new clusters were formed was noted. Using the NIS Elements AR module, the cell length at the time of new cluster formation was measured.

### Calculating the Fraction of Cells with Detectable EI Clusters

For each growth condition, 400 or more cells from three independent experiments were manually classified as having or not having a detectable EI-mCherry cluster. The statistical significance of the differences among cells with detectable clusters were determined by ordinary one-way ANOVA test using GraphPad Prism.

### Following EI Localization in Overnight-Grown Cells Inoculated into Fresh

medium: ROIs were selected to mark the cells outline of the first image in each time-lapse series as following. Binary layers were created using the "Define Threshold" tool over the phase contrast acquired images. ROI were then drawn over the binary layers. From each ROI, the value of SDI from the mCherry channel, taken at time 0, for each cell were obtained and presented as a heat map. Each cell outlined by ROI was manually monitored and noted for each time point of the time-lapse series as having or not having a detectable EI-mCherry cluster until cell division was completed. At the final time point, the sister cells were classified depending on whether both had detectable clusters or not. Finally, the cells were aligned according to their EI-mCherry SDI values.

### The Relation between 2-NBDG Uptake and EI-mCherry Spatiotemporal Organization

ROIs were selected to mark the cells outline of the first image in each time-lapse series as following. Binary layers were created using the "Define Threshold" tool over the phase contrast images. ROI were then drawn over the binary layers. From each ROI, the 2-NBDG MI value and the SDI values of the mCherry signal were exported and presented as a scatter plot.

### Western Blotting

Equal amount of samples were collected, washed and their proteins were separated on 10% SDS–polyacrylamide gels. Gels

### REFERENCES

were subjected to Western blot analysis as described previously (Lopian et al., 2010). α-GroEL (Abcam) and α-mCherry (Abcam) were used for detection of GroEL and EI-mCherry proteins, respectively.

# AUTHOR CONTRIBUTIONS

SG, NA, and OA-C contributed to the design of the study and wrote the manuscript; SG, NA, and TS contributed to data acquisition; SG, NA, TS, AN-S, and OA-C contributed to data interpretation.

### ACKNOWLEDGMENTS

This research was supported by the Israel Science Foundation founded by the Israel Academy of Sciences and Humanities and the Deutsch-Israeli Project Cooperation (DIP). We thank Zakhariya Manevitch from the Core Research Facility of the Faculty of Medicine at The Hebrew University for help with FRAP microscopy. We acknowledge Thomas G. Bernhardt (Harvard Medical School) for gifts of strains. We thank members of Orna Amster-Choder's lab for fruitful discussions.

### SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb. 2017.02695/full#supplementary-material

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Copyright © 2018 Govindarajan, Albocher, Szoke, Nussbaum-Shochat and Amster-Choder. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Division-Based, Growth Rate Diversity in Bacteria

Ghislain Y. Gangwe Nana1 , Camille Ripoll 2 , Armelle Cabin-Flaman3 , David Gibouin3 , Anthony Delaune3 , Laurent Janniere4†, Gerard Grancher 5 , Gaelle Chagny 5 , Corinne Loutelier-Bourhis 6 , Esther Lentzen7 , Patrick Grysan7 , Jean-Nicolas Audinot 7 and Vic Norris 1 \*

1 Laboratory of Microbiology Signals and Microenvironment, Department of Biology, University of Rouen, Mont Saint Aignan, France, 2 Department of Biology, University of Rouen, Mont Saint Aignan, France, 3 Groupe de Physique des Matériaux, Centre National de la Recherche Scientifique, Département de Biologie, Université de Rouen Normandie, Saint-Etienne du Rouvray, France, 4 Laboratoire iSSB, Evry, France, 5 R. Salem Laboratory of Maths, UMR 6085 Centre National de la Recherche Scientifique-University of Rouen, Saint Etienne du Rouvray, France, 6 UMR Centre National de la Recherche Scientifique, 6014 COBRA, University of Rouen, Mont Saint Aignan, France, 7 Material Research & Technology Department, Luxembourg Institute of Science and Technology, Belvaux, Luxembourg

### *Edited by:*

Conrad L. Woldringh, University of Amsterdam, Netherlands

### *Reviewed by:*

Tanneke Den Blaauwen, University of Amsterdam, Netherlands Niculina Musat, Helmholtz-Zentrum für Umweltforschung (UFZ), Germany

> *\*Correspondence:* Vic Norris<EMAIL_ADDRESS>
### *†Present Address:*

Laurent Janniere, Génomique Métabolique, Génoscope, Institut François Jacob, CEA, CNRS, Univ Evry, Université Paris-Saclay, Evry, France

### *Specialty section:*

This article was submitted to Microbial Physiology and Metabolism, a section of the journal Frontiers in Microbiology

> *Received:* 27 November 2017 *Accepted:* 12 April 2018 *Published:* 11 May 2018

### *Citation:*

Gangwe Nana GY, Ripoll C, Cabin-Flaman A, Gibouin D, Delaune A, Janniere L, Grancher G, Chagny G, Loutelier-Bourhis C, Lentzen E, Grysan P, Audinot J-N and Norris V (2018) Division-Based, Growth Rate Diversity in Bacteria. Front. Microbiol. 9:849. doi: 10.3389/fmicb.2018.00849 To investigate the nature and origins of growth rate diversity in bacteria, we grew Escherichia coli and Bacillus subtilis in liquid minimal media and, after different periods of 15N-labeling, analyzed and imaged isotope distributions in individual cells with Secondary Ion Mass Spectrometry. We find a striking inter- and intra-cellular diversity, even in steady state growth. This is consistent with the strand-dependent, hyperstructure-based hypothesis that a major function of the cell cycle is to generate coherent, growth rate diversity via the semi-conservative pattern of inheritance of strands of DNA and associated macromolecular assemblies. We also propose quantitative, general, measures of growth rate diversity for studies of cell physiology that include antibiotic resistance.

Keywords: heterogeneity, asymmetry, bacteria, DNA segregation, cell cycle, secondary Ion mass spectrometry, NanoSIMS 50, isotope-labeling

# INTRODUCTION

Although phenotypic diversity is fundamental to the way populations of bacteria deal with the opportunities and risks presented by their environment including exposure to antibiotics (Smits et al., 2006; Maisonneuve and Gerdes, 2014), the extent of this diversity and the mechanisms responsible for generating it remain to be fully elucidated. Phenotypic diversity is reflected in the diversity of growth rates. The conclusion drawn from early, influential studies on nondifferentiating bacteria growing on nutrient plates (Schaechter et al., 1962) or in liquid media (Ecker and Kokaisl, 1969) was that individual cells had the same growth rate. In recent years, however, the results from microfluidics and other studies have pointed to the radically different conclusion that individual bacteria grow in apparently the same conditions with different rates (Godin et al., 2010; Campos et al., 2014; Taheri-Araghi et al., 2015; Wallden et al., 2016). That said, there is no disagreement that major changes in the growth rate at the level of the population due to growth in nutritionally different media are accompanied by major changes in the size and composition of the cells (Schaechter et al., 1958).

Growth rates affect the key steps of the cell cycle (Weart et al., 2007) and, reciprocally, the cell cycle affects growth rates (Wang et al., 2010; Bowman et al., 2013; Osella et al., 2014). If an evolutionarily useful phenotypic diversity at the level of growth rates is to be generated—for example, to allow some bacteria to escape the action of an antibiotic (Balaban et al., 2004; Kim and Wood, 2016)—each phenotype must be coherent with respect to the set of genes expressed (Norris and Amar, 2012). To achieve such coherent diversity, we have proposed that one of the parental strands of DNA can be physically associated with proteins appropriate for a survival strategy whilst the other strand can be physically associated with proteins appropriate for a growth strategy, so allowing division to generate daughters with different, coherent phenotypes (Rocha et al., 2003). To investigate growth rate diversity, we grew the model organisms Escherichia coli and Bacillus subtilis in liquid minimal media, labeled them with the rare, stable isotope, 15N, and analyzed them using the sensitive, quantitative imaging technique of Secondary Ion Mass Spectrometry (SIMS) (Musat et al., 2008; Boxer et al., 2009; Petroff et al., 2011).

### MATERIALS AND METHODS

### Cell Culture

E. coli BL21 (B F− ompT lon hsdS(rB− mB−) gal dcm λ(DE3) was grown continuously at 37◦C in M9 medium containing per liter 0.1 mmol CaCl2, 8.498 g Na2HPO4-2H2O, 3 g KH2PO4, 1 g NH4Cl, 2 mmol MgSO4, 0.5 g NaCl, 4 g D-glucose. Twentyfive milliliters of culture was shaken at 240 r.p.m. in a 250 ml Ehrlenmeyer flask in either a Buhler incubator or a New Brunswick G76 shaker. The mass doubling time was 64 min (OD600). After 15 generations at an OD600 that never exceeded 0.1, bacteria were inoculated, via preheated pipettes, at the dilution of 1:10 in new M9 medium with 15NH4Cl (98% 15N, ISOTEC, USA) as the only nitrogen source (note that an OD600 of 0.1 corresponds to early exponential growth and, in the conditions we used, plenty of nutrients are available and the cells could have continued growing exponentially for several more generations). One milliliter samples were then taken at 2, 4, 8, 16, 32, 64, and 128 min; the OD600 was <0.05 after labeling for 128 min. Growth was stopped by adding 1 mL of M9 14N medium at 0◦C. All subsequent manipulations were then performed at or below 4◦C. To exclude the artifactual incorporation of isotopes during this cold treatment, cells were grown in 14N medium and growth was stopped by adding M9 15N medium at 0◦C and cells were prepared and analyzed as described below; no incorporation of 15N was detected (Supplementary Figure 1). To exclude the presence of contaminant bacteria, a standard metabolic test, API, was performed, which confirmed that the bacteria were E. coli (not shown); moreover, the bacteria produced identical colonies on agar plates (not shown) and the mass doubling time as measured by optical density corresponded to that previously reported for this strain in the same growth conditions. Finally, if it is supposed we started with equal numbers of two different species, for example, one with a mass doubling time of 36 min and the other 72 min; after 960 min of steady-state growth—round 15 generations—the slower species would be present at a frequency of 2960/72/2960/36 so around 1/10000, effectively leaving a single species growing at the faster rate.

B. subtilis 168 trpC2 (Burkholder and Giles, 1947) was grown at 40◦C in Spizizen medium containing per liter 50 mmol CaCl2, 14 g K2HPO4, 6 g KH2PO4, 2 g (NH4)2SO4, 1 g C6H5Na3O7.2H2O, 2 mmol MgSO4, 11 mg Fe III citrate, 10µmol MnCl2, 1µmol FeSO4, 4 mg FeCl3, 2g D-glucose, 100 mg tryptophan, and 1 g casein hydrolysate. 50 ml of culture were shaken at 240 r.p.m. in a 250 ml Ehrlenmeyer flask in a Buhler incubator. After a 1:50 dilution of an overnight culture, bacteria were grown for 3 h to reach exponential phase (in which the mass doubling time was measured at OD600 as 42 min). At the start of this exponential phase, 13C6-D-glucose (99% 13C, ISOTEC, USA) and 15NH4Cl (98% 15N, ISOTEC, USA) were added to give final ratios of 13C-D-glucose:12C-D-glucose of 1 and 15NH4Cl:14NH4Cl of 1. Samples of 10 mL were then taken at 90 min and 120 min. Growth was stopped by adding 10 mL of Spizizen medium at 0◦C. All subsequent manipulations were then performed at or below 4◦C. Samples were centrifuged in a Sigma 3K18C at 6000 r.p.m. (8,700 g) for 10 min and the pellets were resuspended in 0.1 mol/L cacodylate plus 0.04% MgCl2 (to help avoid autolysis), then centrifuged again at 6000 r.p.m. for 15 min.

### Preparation of Cells

Samples were centrifuged in a Sigma 3K18C at 6000 r.p.m. (8,700g) for 10 min. Then one of two fixation methods was used. Cells were fixed either with ethanol but without formaldehyde (in the case of all Figures shown except **Figure 1I** and Supplementary Figure 2G) or with formaldehyde but without ethanol (in the case of **Figure 1I** and Supplementary Figure 2G):

Method 1: The pellets were resuspended in 0.2 mol/L cacodylate plus 0.04% MgCl2 (to help avoid autolysis), then centrifuged again at 6000 r.p.m. for 15 min. The following dehydration step was performed at 4◦C for both B. subtilis and E. coli. The pellet was resuspended in 25% ethanol (in 0.2 mol/L cacodylate) for 20 min, then 50% ethanol (in 0.2 mol/L cacodylate pH 7.4) for 15 min and centrifuged at 6,000 r.p.m. (8700g) for 10 min. The pellet was then suspended in 10µL of 50% ethanol (in 0.2 mol/L cacodylate) and 2µL were deposited on a silicon chip (which had previously been cleaned by sonication during successive immersions in distilled water, 96% ethanol, acetone, 3:1 H2SO4:H2O2, water, and finally 96% ethanol). The dehydrated bacteria on the chip were then baked in a vacuum at 50◦C for B. subtilis and at 40◦C for E. coli.

Method 2: The pellets were resuspended in 1 mL of formaldehyde 3.7% buffered with 0.2 mol/L cacodylate pH 7.4 for 30 min, then centrifuged again at 12,000g for 10 min at 4◦C. The pellet was rinsed for 10 min in 0.2 mol/L cacodylate pH 7.4; this rinsing was then repeated once. The pellet was then suspended in 20µL of 0.2 mol/L cacodylate and 2µL were deposited on a silicon chip. The dehydrated bacteria on the chip were then baked in a vacuum at 50◦C.

### SIMS Analysis

The bacteria on the chip were analyzed using a Zeiss Scanning Electron Microscope (which showed the bacterial morphology expected with no evidence for lysis, Supplementary Figure 3) and a CAMECA NanoSIMS 50. The latter entails bombardment of the sample with a high energy primary beam of ions that

arrowheads show cells with mass doubling times estimated using isotope incorporation as 61 min and 82 min in (G), and as 70 min and 80 min in (H). Scale bars 2µm.

is focused on the surface and rastered across it; this fragments molecules contained in a nanovolume and sputters them out; a proportion of these fragments is collected as secondary ions and counted with a mass spectrometer (the useful yield), hence the name Secondary Ion Mass Spectrometry. These ions counted in each nanovolume can then be displayed as an image. In our conditions, the chip was bombarded with a primary beam of Cs+ ions (the impact energy for the primary ions was set to 16 kV and the primary current to 1.3 pA or to 0.8 pA for the B. subtilis and E. coli experiments, respectively), the mass spectrometer set at a mass resolution of 5000 = M/1M (where 1M is the difference between the masses of two secondary ions with the same mass number, M); five out of the six following secondary ions were counted simultaneously from the same place on the surface of the sample: 12C, 13C, 12C 14N, 12C 15N or 13C 14N, and 13C 15N.

In a "sputtered section," all material containing carbon and nitrogen is sputtered out and a proportion then detected (the useful yield); the labeled nitrogen is sufficient to allow detection of membrane and cytoplasm which both contain this element. Three successive sections of the cells were sputtered out in our conditions and all the counts in the middle section were used for imaging in all the figures shown except for Supplementary Figures 4, 5. Assuming that the density of the sample is 5.1022 atoms/cm3 and that the sputtering yield (the number of sputtered particles per incident primary ion) was in the range 2 to 10, about 5 nm of material was sputtered out per section of B. subtilis (1.3 pA, raster 8µm, pixels 128, dwell time 40 ms/px) and about 5 nm per section of E. coli (0.8 pA, raster 10µm, pixels 256, dwell time 10 ms/px). ImageJ was then used for the extraction, treatment and assembly of the data (Mutterer and Zinck, 2013). Excel and OpenOffice were also used for data collection and statistical analysis. The naturally occurring ratios of 15N:14N and 13C:12C are 0.0036 and 0.011, respectively. Isotope fractions such 12C 15N/(12C 14N+12C 15N) were used to both estimate mass doubling times (as equaling 15N/(14N+15N), see below) and avoid artifacts due to matrix effects (Lhuissier et al., 2000; Peteranderl and Lechene, 2004); the latter arise because of different useful yields in different matrices and can be allowed for by using the fraction 12C 15N/(12C 14N+12C 15N) because both ions are affected in the same way.

The mean values of 13C 15N/(13C 15N+12C 15N) obtained for E. coli, 0.0112 ± 0.0003 after 128 min, were the same as the natural fraction of 13C, as expected. To ensure the validity of our results, we excluded those images in the E. coli experiments that had a bacterium with a mean value of 13C 15N/(13C 15N+12C 15N) >0.016 or <0.006, that is, those bacteria with fractions out of the range 0.011 ± 0.005. This choice was because 0.005 is the standard deviation typically obtained from the pixels of a single bacterium.

Instrumental isotope fractionation in SIMS might be considered a potential problem. This fractionation means that, globally, two different numbers of secondary ions are counted for the same number of atoms of two different isotopes of the same element present in the sample (so that, for example, the isotope ratio of 12C 15N/(12C 14N+12C 15N) counted by the machine differs from the actual ratio). The reason this might be a problem is because the fractionating effects during ionization and analysis are likely to contribute to the isotope ratios on which our conclusions are based. However, a biological variability in a reference sample (a single cell with the normal isotope ratio) has been reported as being <2% (Kopf et al., 2015). Hence the differences between the real and the counted numbers of isotopes are small compared with (1) the differences due to the Poisson nature of the counting process [the percentage error = 100/(square root of the mean)] which in the case of the 32 min and 128 min labeling, for example, is typically 7 and 3%, respectively, on the same image), (2) the differences between the bacteria (in the text we report that the %CV ranges from 5 to 20) and (3) the level of diversity apparent even within an image, which is acquired with the same settings within a short period of time; admittedly, the counting artifact could become more significant when comparing acquisitions obtained at different times and on different samples. However, it should be noted that the diversity we describe below is apparent even within a single image. And even if there were a significant error in our fractionation, this would apply to all the bacteria and simply shift the mean of the distribution of isotope ratios without affecting its variance (i.e., the growth rate diversity we observe would still be true); this is because there is a linear relationship between the isotope ratio of 12C 15N/(12C 14N+12C 15N) measured by SIMS and bulk isotope ratios measured by other techniques (see Table 1 and Figure C8 in Kopf et al., 2015). Finally, our estimation of the mass doubling time of the population based on the isotope ratios matches the estimation of the mass doubling time based on optical density, which it would not do if there were a significant error in our measurements.

### Mass Spectrometry

Electrospray ionization-mass spectrometry (ESI-MS) experiments were performed using a Bruker HCT Ultra ETD II quadrupole ion-trap (QIT) equiped with an ESI source, Esquire control 6.2 and Data Analysis 4.0 software package (Bruker, Bremen, Germany). The ESI parameters were capillary and end plate voltages set respectively to +3.5 and +3.0 kV in negative ion mode. The skimmer and the capillary exit voltages were set respectively to −30 and −80V and the injection low mass cut-off (LMCO, corresponding to the "trap drive" parameter) value was m/z 40. Sample solutions were either infused into the ESI source at a flow-rate of 3µL.min−1 by means of a syringe pump (Cole-Palmer, Vernon Hills, Illinois, USA) for full scan data or injected into the ESI source at a flow-rate of 100µL.min−1 by means of a liquid chromatographic (LC) system (Agilent 1200 chromatographic system equipped with a G1379B degasser, a G1312A high pressure binary pump and a G1329A autosampler, Agilent technologies, Waghaeusel-Wiesental, Germany) for full scan data as well for the more specific precursor-to-product ion transition termed selected reaction monitoring (SRM). For infusion, the nebulizer gas (N2) pressure, drying gas (N2) flow rate and drying gas temperature were 10 psi, 7.0 L/min and 300◦C, respectively. For LC-injection, the nebulizer gas (N2) pressure, drying gas (N2) flow rate and drying gas temperature were 30 psi, 10.0 L/min and 300◦C, respectively. Helium pressure in the ion trap was 1.1 × 10−5 mbar. Full scan spectra were acquired in the m/z 50–1,000 range, using a scan rate of 8100 m/z per second ("Standard-Enhanced" mode). The number of ions entering the ion trap was automatically adjusted by controlling the accumulation time with the ion charge control (ICC) mode (target 50,000) with a maximum accumulation time of 50 ms. The values of spectra averages and rolling average were 5 and 2 respectively. The selected reaction monitoring (SRM) experiments using glutamate specific precursor-to-product ion transition "m/z 146 → m/z 128" were carried out by collision induced dissociation (CID) using a resonant excitation frequency with ramping amplitude from 0.06 to 0.42 Vp−p, helium as the collision gas, isolation width of 1 m/z unit for the precursor ion (m/z 146) and 1 m/z unit for the product ion (m/z 128).

### Estimation of Population Growth Rate

We took growth to be the increase in cell mass (irrespective of its molecular nature) as determined by nitrogen uptake thereby disregarding the question of whether cells have significantly different densities (Poole, 1977; Martinez-Salas et al., 1981; Makinoshima et al., 2003). After addition of 15N to the E. coli culture, we measured the isotope fraction of nitrogen, x(t), using 12C 15N/(12C 14N+12C 15N) after different periods, t. Two independent series of experiments, A and B, were performed and from 2 to 11 images were acquired and analyzed after each period of labeling. The values at 4 min for one of the B series were discarded (see SIMS analysis above). A total of 2174 bacteria in 83 images were analyzed. Given that an unsynchronized population of bacteria grows exponentially:

$$\mathbf{x}(\mathbf{t}) = \langle \mathbf{x}\_{\text{ext}} \times (2^{\mathbf{t}/\mathbf{g}} - 1) + \mathbf{x}\_{\text{0}} \rangle / 2^{\mathbf{t}/\mathbf{g}} \tag{1}$$

Where t is the period of labeling, g is the mass doubling time, x0 is the natural level of 15N (0.00366), and xext is the fraction of 15N to total N in the final medium after addition of the 15N. xext is obtained from adding 9 volumes of 98% 15N (as estimated by the supplier) to 1 volume of natural nitrogen hence xext = (0.9 × 0.98) + (0.1 × 0.00366) = 0.882

To be coherent with the model of exponential growth of a population of bacteria where it is the population mass which is considered, all the values of the x(t) measured at each pixel of all the bacteria visible in a given image j were pooled and the mean of the isotope fractions of the pixels in an image, xj(t), and standard deviations were calculated. In exponential growth where γ = (ln2)/g, the growth of the population can be expressed as:

$$\mathbf{Y}(\mathbf{t}) = (\mathbf{x}\_{\text{ext}} - \mathbf{x}\_0) / (\mathbf{x}\_{\text{ext}} - \mathbf{x}(\mathbf{t})) = \exp(\mathbf{y}\mathbf{t}) \, \mathbf{z} \tag{2}$$

Hence, this can be expressed in terms of the image means:

$$\overline{\mathbf{Y\_j(t)}} = \langle \mathbf{x\_{ext}} - \mathbf{x\_0} \rangle / \langle \mathbf{x\_{ext}} - \overline{\mathbf{x\_j(t)}} \rangle \ = \exp(\gamma \mathbf{t}) \tag{3}$$