Introduction

Asymmetric cell division refers to the mechanism by which a mother cell splits into two daughter cells with distinct cellular fates. As first shown by Edwin Conklin [1], this process is usually achieved by asymmetric inheritance of intrinsic determinants of cell fate, such as specific proteins or RNA, and plays a crucial role during the development of organisms, favoring cell differentiation and self-renewal [2, 3].

Increasing evidence is vouching for a more ubiquitous presence of fluctuations in the partitioning of cellular elements. In fact, asymmetric segregation has been observed in non differentiating cell populations as different as bacteria [4], yeast [5], and tumor cells [6, 7]. In bacterial populations, differential segregation of cellular components has been linked to antibiotic resistance [8]. Similarly, Katajisto et al. [9], studying mitochondrial partition in stem-like mammalian epithelial cells, found that asymmetric segregation of older mitochondria enables cells to protect themselves from aging by preventing the accumulation of misfolded proteins. In addition, asymmetric centrosome partition has been observed in stem cells from human neuroblastoma and colorectal cancer, controlled by polarity factors and influenced by the segregation of subcellular vesicles during cell division [10]. To quantify the strength of fluctuations at division, the common route goes through fluorescence microscopy measurements of either fixed or live cells, where fluorescent dyes are used as a proxy for the cellular element of interest. The partition statistics are then estimated by looking at the fraction of fluorescence intensity inherited by the two daughter cells [11, 12]. This approach is, in principle, very informative and a lot of work has been done in recent years to write theoretical models of stochastic gene expression [1316] that can account for the coupling of noise sources and can even be extended to time-lapse microscopy data of single cell’s growing dynamics [17]; however, relying on microscopy data implies the need of the (often manual) identification of hundreds of division events and in some cases even complex experimental set-ups and techniques where cells (usually bacteria) have to be followed for a great number of generations [17]. These are highly time-consuming procedures that dampen the possibility of a wide characterization of the partition statistics of different cell types.

To cope with these limitations, we propose using high-throughput flow cytometry measurements to quantify fluctuations in the partitioning of cellular components in adherent cells. In particular, flow cytometry has already been applied to study asymmetric division in specific cellular systems: Yang et al. [5], sorted budding yeast cells to analyze the differential segregation of proteins between daughters cells, while Peruzzi and coworkers first applied multi-color flow cytometry to show that leukemia cells divide mitochondria and membrane proteins asymmetrically [18]. Here, we generalized the protocol to adherent cell types and probed its accuracy by comparing the measured fluctuations to those observed with extensive time-lapse microscopy experiments on a panel of normal and cancerous epithelial cells.

Our work shows that (i) there is an exact analytical expression linking the inherited cellular component distribution dynamic to the specifics of the partition process. Based on this result, we (ii) proposed an experimental protocol to tag with fluorescent dyes cellular components and measure the dynamics of fluorescent distribution, whose analysis allows accurate estimates of the degree of asymmetry in the partition process. Applying our method to a panel of normal and cancerous human colon cells, we found that (iii) different lines of colon adenocarcinoma display very distinct extent of fluctuations in cytoplasm partition, reflecting in (iv) different degrees of asymmetric division of their size.

Results

Modeling cell population dynamic in presence of partitioning noise

In this section, we propose a theoretical framework that is able to describe the evolution of the fluorescence intensity of a cell population stained with live fluorescent marker, by uncovering its dependence on the underlying partitioning process to which single cells are subjected. The rationale is that, if cells follow a common division rule, knowing the shape of the mother distribution, it should be able to determine the daughters’ one. By iterating this computation, one can trace the whole proliferation. The specificity of our approach lies in its specific design for live fluorescent markers. The selected component is tagged through a staining procedure with fluorophores that persist in the cell and its daughter cells, as long as the targeted protein is not degraded.

Monitoring the time evolution of the abundance of specific cellular elements over multiple divisions in singlecell experiments, it is observed [19] that the counts of a component have a qualitatively similar dynamic to those shown in Figure 1(b). Indeed, variations in the components number are determined by: (i) the production and degradation processes along the duration of the cell growth phase; and (ii) non identical apportioning of the components between the two daughter cells. Here, we will assume that the number of components varies only due to partition events and look for a quantitative description of the evolution of the element distribution. By neglecting variability in division times and the intercurrent production and degradation of components, we are reducing the dynamics to a progressive dilution. The stability of our model with respect to these assumptions will be studied in the following section and more in detail in the Supplementary Information. Note that it is possible to model the complete growth and division process (see for instance refs. [14, 16, 2022]), but it is more difficult to separate the different sources of noise. Here, we start by considering a population of cells characterized by a certain initial distribution, P(m)0, of element m. Each cell i of the population divides into two daughter cells that we can label unambiguously as 2i and 2i + 1, which inherit respectively m2i and m2i+1 components. The process, exactly at division, must conserve the total number of elements, i.e. mi = m2i + 1 and we call partitioning fraction the quantity . To model noise in partitioning, we assume that at each division a random value for f is extracted from a probability distribution function, ∏(f). Given these assumptions, we can write the probability distribution of the number of components of a daughter cell m after g divisions from the mother M as:

where Pk(M) is the components distribution at generation k and it is given by the sum of the distributions of all siblings cells of the k-th generation. After g generations, the cell lineage generated from a single cell forms a lineage tree composed of ng = 2g cells. Without losing generality, the division probability can be expressed, as:

Therefore, we can write the probability distribution, P(m), of the daughter cells inheriting a fraction f, of the mother elements:

As is often the case, it is convenient to characterize the distributions in terms of their moments [21]. We identify with m2i (m2i+1) the subpopulation of cells inheriting a fraction f(1—f). In particular, the mean of the daughter cell component distribution can be expressed, for 2i, as :

Where E[·] stands for the average. Note that every computation can be symmetrically done for the siblings subpopulation 2i + 1. Thus, we have that after a single generation, the mean number of inherited elements of the daughter cells depends on the asymmetry of the ∏(f) and on the mean of the P(mi).

Modeling partitioning noise at cell division.

a) Schematic representation of the growth and division process of a mother cell. The cell undergoes an initial phase of duplication, where its internal elements are multiplied, followed by a division phase, where these elements are partitioned between the two daughter cells. One daughter inherits a fraction p of the mother elements, while the other receives the remaining fraction q = 1 − p. b) (blue) Idealized behavior of the components number of a cell element over time, considering three different noise terms: fluctuations in compound counts during growth, uncertainty in the timing of division, and noise in the compound’s partition fraction. (red) Same description but considering only production and degradation processes. c) Microscopy images of a single cell division event for an HCT116 cell, whose cytoplasm is stained with Celltrace Yellow. d) Time evolution of the distribution of the components number of a cellular element for a population of cells subjected only to partitioning noise. The distribution used for the simulation is a Gaussian distribution with mean and standard deviation σ = 0.07. e) Examples of partition distributions II, with increasing coefficient of variation. f) Mean (μg) and variance (σg) of the components number as a function of the generations, g, for the proliferation of a population subject to different partition noise distributions. Different colors correspond to distributions with varying coefficients of variation, as represented in the top panels. The dashed line represents the theoretical behavior obtained from the model. Dots are colored according to the distributions shown in panel e). g) Behavior of the mean (left) and standard deviation (right) of the simulated dynamics compared to the expected theoretical behavior for a proliferating population assuming a sizer division strategy (see SI).

In Figure 1(d), we show the evolution of the components number distribution at different generations. The ∏(f) chosen for the simulation is a Gaussian distribution with mean and standard deviation σ = 0.7.

The total distribution is the sum of 2g sub-populations which derives from the division process. For example, if we assume , after one division we would have a population distributed around and one at (1 − f)μ. After two division the populations would grow to four sub-populations centered in f2μ, f(1 − f)μ, (1 − f) and (1 − f)2μ and so on.

To compute the mean relative to an entire generation, we need to sum over all the different sub-populations:

where is the mean of k-th sub-population at generation g and 2g is the total number of sub-populations for that generation. This equation [5] can be simplified into the form (see Supplementary Note for the full computation):

where μ0 is the mean of the initial population. It is reievant to notice, that independently of the characteristic of ∏(f), the mean always halves itself at each generation. This result is simply understood considering that for each cell that inherits a fraction f of the compound, its sibling inherits a fraction 1 − f. Hence, at every division, the compound counts is diluted by an average factor of 2 in the population, and partition distributions with different properties lead to the same behavior of μg (see Figure 1(f)). Therefore, no relevant insights on the ∏(f) can be obtained from the first moment. We thus moved to consider the second moment, i.e. the variance. Similarly to what is done for the mean, the variance of the inherited fraction of elements can be expressed as:

and in turn, the variance of a generation, starting from the consideration of the variance of a mixture of distributions can be written as (see Supplementary materials for details):

with

It is possible to obtain a recursive equation for Ag (see Supplementary Notes), Ag = Ag−1E[f2] = A0E[f2]g which leads to a compact form for equation [8]:

where .

Finally, with equation 10 we can link the variance of the components distribution for the entire cell population at generation g, to the properties of ∏(f). It is interesting to highlight that there is no dependence on the bias asymmetry of the process. The distribution ∏(f) is necessarily symmetric, as daughter cells are indistinguishable. Therefore, by bias we refer to a situation where the most probable outcome of the division process is asymmetric. To this order, the only relevant feature shaping the dynamics is the extent of the fluctuations, regardless of their origin. Indeed, as shown in Figure 1(e), the same ∏(f), but with a different variance shows a different behavior for .

Our analytical results show that the dynamical behavior of the population variance depends non-trivially on the second moment of the underlining partitioning distribution. The next section will be devoted to demonstrating how this condition can be exploited in experiments to measure the partition noise component and to evaluate the stability of our assumptions with respect to concurrence of the other noise sources (see also Figure 1g).

Measuring fluctuations via flow cytometry

As shown by Peruzzi et al. [18] for leukemia cells growing in suspension, it is possible to use flow cytometry to follow the evolution of properly marked cellular elements in time. Here, we generalized the protocol considering a panel of colon cell lines, i.e. Caco2, a human epithelial colorectal adenocarcinoma cell line; HCT116, a human colorectal carcinoma cell line; and CCD-I8C0, a normal human colon fibroblast cell line. To measure the degree of asymmetric division, we marked cells’ cytoplasm via a fluorescent dye and followed its dilution through successive generations. As explained in detail in the Materials and Methods section, the experimental procedure we used is based on the following main steps: staining of the selected compound, sorting of the initial population, plating and acquisition of the distinct populations in the following days. The sorting allows us to obtain an initial population with a narrow peak at a specific fluorescence intensity, and it plays a dual role (a full explanation can be found in the Supplementary Information): it indirectly induces cell synchronization and enables a clearer reconstruction of the coexisting generations.

To achieve this, we used fluorescence-activated cell sorting (FACS) to isolate a specific cell population based on cell morphology and fluorescence intensity. The sorted population is divided into distinct wells, one per acquisition, in equal numbers and kept under identical growth conditions to mitigate inoculum-dependent variability among samples [23]. Figure 2(a) shows an example of the time course dynamics obtained with HCT116 cells (see Methods for details). Each distribution represents a different acquisition of a whole well at increasing times from plating (0 to 84 hours from the bottom to the top). It is possible to neatly observe the succession of peaks associated with different generations that succeed one another and the average decrease of the fluorescence intensity. Note that different generations can coexist at the same time point. One can spot the flow of cells towards higher generations and the behavior of the peak heights, which grow over time, reach dominance and then decay. The shift on the y-axis allows to follow the evolutionary dynamics while being able to see the overlapping of the same generation at different time points. As in each measurement, the distribution of fluorescent intensity I, is a mixture of distributions corresponding to different generations:

Quantification of partitioning noise via population-level measurements.

a) Time evolution of CellTrace-Violet™ fluorescence intensity distribution measured in a flow cytometry time course experiment for a population of HCT116 cells. Time progresses from the darker shade of blue to the lightest, spanning from [0,84] hours (bottom line to top). b) Snapshots of the evolution of the distribution of CellTrace-Violet fluorescence intensity measured in a flow cytometry time course experiment for a population of HCT116 cells. Experimental data are represented by the light blue histogram, while the best-fit Gaussian Mixture Model is displayed as lines, with different colors representing different generations. c) Mean (left) and variance (right) of the intensity of the fluorescent markers as a function of generations, normalized to the initial population values. Each replica of the experiment is identified by a different color and point marker’s shape. The experiment showed are conducted on Caco2 cells and the points correspond to the mean values for each generation. Dashed lines represent the best fit according to Equations 10 and 6, respectively. d) Division asymmetry of the ∏(f) obtained by fitting Equation 10 to data for all experiments and cell lines. The division asymmetry is measured via the percentage of the coefficient of variation.

To obtain the necessary information on the properties of each generations we performed a fitting of the data via a Gaussian Mixture Model (GMM) combined with an Expectation Maximization algorithm (details are explained in Materials and Methods). Indeed, assuming that each generation is log-normally distributed with mean μg and variance of the corresponding Gaussian distribution, equation 11 can be rewritten as:

An example of the outcome of the GMM algorithm is shown in Figure 2(b) for two acquisitions (48h and 60 h from plating). Blue-shaded distributions represent experimental data, while the best solution of the Gaussian mixture is shown with a solid line. Each generation’s fit is identified with a different color. Upon an overnight, we can identify generations 0, 1, and 2 in both the acquisitions. One can see that their relative fractions changed: generation 3 became dominant overnight and generation 4 appeared.

Figure 2(c) shows the behaviors of μg and σg as functions of the generation number g for three independent replicas of the time course experiments with CaCo2 cells. Dashed lines represent the best fit of the data against Equation 5 for the mean and Equation 10, which instead yields the variance σf of the partition distribution probability, ∏(f). The outcomes of the fitting procedure on the time course dynamic for all the replicas of three studied cell lines are reported in Figure 2(d). To evaluate the degree of asymmetry we have defined the division asymmetry as the percentage of the coefficient of variation of the obtained partition distribution (CV: ). In particular, Caco2 and CCDI8C0 cells show a higher level of asymmetry with respect to HCT116. As already stated in the description of the theoretical model, we consider only partitioning noise and neglect other potential sources: (i) production and degradation processes, and (ii) variability in cell cycle length. Since we use non-endogenous live fluorescent dyes, production processes are excluded by definition, while degradation results in a negligible time-dependent decrease in mean fluorescence intensity as used dyes are optimized to maintain stable fluorescence levels for the duration of the experiment. Indeed, by looking at the data in Figure 2 (c) the values of the measured μg seem to only fluctuate around the theoretical line, but with the same slope, indicating that degradation can be neglected without major consequences. The effects of correlations between the inherited fractions of cellular components and the fluctuations in the duration of the cell cycle are in principle more difficult to disentangle. In this respect, it is widely known that cells adjust their growth/division strategy to mitigate size fluctuations at division [24, 25]. According to the ‘sizer’ model, they divide upon reaching a certain size so that, smaller cells, i.e. those that inherit a smaller fraction of their mother’s volume, will have longer division times and vice-versa. This results in fluctuations in the cell division times that may be coupled with partitioning noise.

To assess the robustness of our model under such conditions, we performed simulations assuming that cells follow a sizer-like division strategy, while the marked component partition perfectly matches cell size partition (see Supplementary Information for more details). This scenario maximizes the coupling among the two noise sources. As can be observed in Figure 1(g), the behaviors of the mean value μg and the variance σg deviate from the the case in which division time is uncoupled with partition noise. However, the values obtained are still within the estimated experimental error (see Supplementary Figure 3). It is worth noticing that the presence/absence of the increase in the mean, μg, for newer generations may highlight the presence of the coupling between the two noise sources.

Validation of the method via live fluorescence microscopy

The flow cytometry-based protocol combined with the developed theoretical model provides insights into the properties of the cell partition function (see Figure 2). To validate our approach, we performed extensive live time-lapse microscopy experiments to measure the partition fraction of a labelled compound during proliferation, providing a direct measure of the ∏(f). Specifically, we aimed to follow a growing cell colony over time, detect divisions, evaluate the fluorescence intensity of mother and daughter cells, and calculate the inherited fraction. The detailed experimental procedure is described in the Materials and Methods. Briefly, cells are stained for the cytoplasm and are plated at low density (approximately 103 cells/well) on IBIDI cell imaging chambers (μ-Slide 4 and 8 wells). After 24 hours (sufficient time for complete cell adhesion to the slide) the wells are washed, the media is renewed and live imaging acquisition time-lapse is started. Multiple fields are followed in each well, and for each field, we acquire a z-stack bright-field and fluorescence image (a single plane is shown in Figure 3a) at constant intervals (20 min) for 3 days. The fields are manually selected at the beginning of the acquisition based on cell density and fluorescence intensity. The recorded time lapses have been then manually analyzed to identify divisions. A sample outcome of the analysis procedure is shown in Figure 3(b). The total fluorescence of the cells is displayed as a function of time, before and after the division. The mother cell fluorescence decays exponentially over time and at division it halves due to components partitioning between the two daughters. In this specific case, we can observe that the fraction of cytoplasm inherited by the two cells is not equal, indicating an asymmetrical division. To compute the partition fraction f, we fit the fluorescence intensity of the daughters’ cells vs time with:

where I(t) is the fluorescence intensity at time t, m accounts for the decaying process and 2i and 2i + 1 respectively refer to the two daughter cells. We constrain the slope to be the same between the two cells. The fraction of inherited fluorescence, f , is obtained as:

To ensure the highest possible accuracy, we (i) measured the total fluorescence of the mother cell and those of the two daughter cells for at least 2 hours before and after the division event to increase the determination of the fluorescence splitting for each detected mitosis. Moreover, (ii) we removed all dynamics that showed an absolute Pearson correlation value between luminosity and pixel size higher than 0.9 as a function of time (see Supplementary Information). In fact, single cell’s pixel size, which corresponds to the cell size projection on the focal plane, is found to shrink before mitosis and enlarges after division. Fluorescent intensity is instead proportional to the number of stained cytoplasmic proteins, which are expected to remain constant (except for a constant decay of the fluorescence). Therefore, a high correlation between pixel size and fluorescence intensity indicates a high noise-to-signal ratio. Figure 3(c) (bottom) displays the strip plot of the obtained partition fraction, f, for Caco2, HCT116, and CCDI8C0 cells, randomly spread on the y-axis, and their corresponding fitted ∏(f) distributions (top). Note that via microscopy imaging, it is possible to measure the whole partition distribution and not just its moments; however, hundreds of events are required to get a reliable estimate. Here, we characterize it by fitting the data with a double Gaussian distribution of the form:

that allows for the measure of both 〈f〉 and f2f − 1/4.

Quantification of partitioning noise via single-cell measurements.

a) Example of a recorded cell colony of HCT116 cells in brightfield (left) and on CTFR fluorescence (right). b) Cell cytoplasm fluorescence intensity as a function of time for a cell before and after division. Dark green circles correspond to the fluorescence intensity of the mother cell up to the division frame, and then to the sum of the daughters’ fluorescence. Lighter green triangles and squares represent the fluorescence intensity of the daughter cells. Solid Unes are the Unear fit of the points. The intercepts of the linear fit are used to compute the fraction of tagged components inherited by the daughter cells. The time is counted from the start of the experiment. c) (bottom) Strip plot of the distribution of the inherited fraction of cytoplasm for the different cell lines. The points are randomly spread on the y-axis to avoid overlay. (top) Fit of the inherited fraction distribution with the sum of two Gaussians with mean symmetric to 1/2. d) Comparison of division asymmetry obtained with time-lapse fluorescent microscopy measures (striped bars) with the ones obtained from flow cytometry experiments (plain bars). The flow cytometry bars are obtained as the mean over the multiple conducted experiments.

The fitted distributions clearly show how HCT116 cells are the most symmetric with only one central and narrow peak, while for CCDI8C0 and CaCO2 cells two asymmetric peaks are visible, associated to higher fluctuations. In Figure 3(d), we compare the outcomes of the flow cytometry vs microscopy measurements. For all three cell lines, the found degrees of asymmetry are statistically consistent between the two approaches.

Size division bias accounts for cytoplasmic fluctuations

In the above sections, we have proposed and validated a method to evaluate fluctuations in the division of adherent cells. The model makes no assumptions about the shape of the partitioning distribution and, in its general form, accounts only for the effect of its variance in shaping the population dynamics. Indeed, we cannot yet distinguish between divisions which are biased, defined by a partitioning distribution with peak values different from 1/2 (as shown in figure 3(C) for CaCO2 and CCDI8C0 cells) or symmetric ones (as shown in figure 3(c) for HCT-116 cells) which though have higher fluctuations. Also, apart from measuring it, we want to propose a qualitative interpretation of the obtained level of asymmetry. We start by recalling that our experimental protocol uses dyes that bind a-specifically to cytoplasmic amines. Thus, (i) it is expected to be uniformly distributed in the cellular cytoplasmic space and (ii) the number of labeled cytoplasmic components can be considered large. With this hypothesis, the least complex, model one can assume is a binomial one with parameter p, measuring the bias in the process. Via equation 10 we got a direct link between the measurable variance of the population and the second moment of the underlining partition probability distribution, 〉(f). Henceforth we sought for a relationship between the second moment of n〉(f) and the parameters of the binomial distribution. To begin with, we can write the partition distribution for the fraction of inherited component, assuming a level of asymmetry p (q = 1 − p), as:

As we have already observed in the general model, due to symmetry, 〈f〉 = 1/2 for any form of the ∏(f). No information on the system can be obtained by looking at the first moment.

For the second moment, we note that the variance of a single branch of the ∏(f) is:

where P(Ni) is the probability distribution of the mother’s components number at the moment of division, while , is the variance of the ∏(f)p given the number of internal component Ni in the mother cell.

Recalling that m2i is the number of inherited intracellular components of one of the two daughter cells, the partitioning fraction f can be expressed as f = m2i/Ni which leads to:

Therefore, by substituting it into equation 16 we obtain:

Which explicitly depends on mother’s component probability distribution, P(Ni), and coherently returns for P(Ni) = δ(NiN).

An equivalent computation can be done for the symmetric branch and therefore the variance of the full distribution is given by:

Note that this relationship depends on the integral:

which requires the knowledge of the mother element distribution P(Ni) to be computed. To get an idea of its behavior, we observe that assuming a population of cells with a fixed number of elements before division, i.e. P(Ni) = δ(N − Ni), the integral simply goes as 1/N; in a more realistic scenario, the elements can be considered log-normally distributed, since they are the product of multiple division processes, hence:

In this case, one gets the following expression for :

where μ = 〈log(Ni)〉.

The quantity can be computed from flow cytometry or microscopy measurements. Indeed, if we consider the relation I = NiF, where I is the total fluorescence intensity, Ni the number of components and F the fluorescence of each component, then = . Since is accessible from the data, so is . On the other hand, the values of μ are not directly accessible, limiting the possibility of an exact analytical computation. However, by computing the value of ∑(N) as a function of for fixed values of (Figure 4(a)), one can see that already for N = 100 the role of ∑(N) is negligible. In general, for the partitioning of components such as cytoplasm which are present in much larger number in the cell, this factor can be neglected. In this regime, the expression for becomes:

By inverting this relationship, we obtained a general mapping between the variance and the level of biased asymmetry of the binomial distribution.

In Figure 4(c), we validated the analytical expression to the values obtained via microscopy measurements, showing again good comaptibility of the results. CaCO2 cells confirm to be the most asymmetric among the cell lines and since they are known for their heterogeneity in cell morphology [26], we explored the hypothesis of shape heterogeneity and cytoplasmic asymmetry to be linked. Figure 4(d) displays three sample cases of division examined through time-lapse fluorescent microscopy following cytoplasm. Different time frames are overlaid to depict the entire dynamic at once, with time progressing from the bottom left to the top right of each image.

Cytoplasm partition fluctuations vs cells size.

a) Behavior of the integral term ∑(N) (Eq.19) as a function of the number of dividing elements, N, at fixed = 0.8. Vertical lines mark typical values of cellular elements, like mitochondria, b) Theoretical value of the variance of ∏(f) in the binomial assumption for different levels of asymmetries and increasing values of N. c) Asymmetry of the partitioning distribution ∏(f) in the binomial limit as measured by the binomial bias, p. d) Sample cases of volume asymmetric division for different cell lines. Images show the overlay of consecutive times during the division dynamic. Time flows from the bottom left to the top right.

We computed the pixel size ratio between the daughter cells in the last analyzed frame and compared it to the fraction of cytoplasm inherited from the mother cell. Although pixel size is only a rough proxy for cell volume, we observe that the cell inheriting a larger fraction of cytoplasm also appears larger in size. This finding suggests that asymmetry in cytoplasmic partitioning is linked to size fluctuations at birth, where cell size reflects the bias observed in the segregation process.

Discussions

Cell division is orchestrated by hundreds of molecular interactions, which are intrinsically stochastic processes [27] . The presence of such stochasticity results in the in-surgence of variability among the cells phenotypes sharing the same genome, as it is the case for cancer cells. In this respect, gene expression noise has been extensively characterized [28] incorporating cell cycle effects [2931], as well as the different strategies that cells evolved to regulate the extent of the noise in these channels [32, 33].

Beside such noise sources, proliferating cells are subject to the fluctuations originated by the partition of cellular elements at division. Analysis of this partition noise highlighted that under some circumstances, its contribution exceeds the others in creating heterogeneity [19, 34]. In fact, asymmetrically dividing cells can produce daughters that differ in size, cellular components, and, in turn, fate [35].

In this framework, the origin of such asymmetry may be ascribed to the basal stochasticity of molecular processes taking place at different levels during the division and/or to specific, evolutionary-conserved mechanisms used by cells to control cell fate and generate diversity [3, 36]. In fact, asymmetric partitioning of components alters the initial counts of molecules in the following cell cycle, which can lead the system to a different phenotypic state [37, 38]. Here, we show that an accurate determination of partition fluctuations can be obtained via standard flow cytometry measurements of properly marked cell populations. Our approach is easier and faster with respect to the use of microscopy imaging, which requires the acquisition and analysis of extensive recording of the cells, ultimately limiting the reachable statistics. In addition, it is specifically designed to isolate partitioning noise from other sources via the usage of non-endogenous dyes to target cell components.

Although, the need to use live fluorescent dyes to track proliferation surely limits the extent of applicability of our method to markable components, we anticipate that there are numerous contexts in which our approach could give a determining insight. For instance, organelles such as mitochondria, the endoplasmic reticulum, lysosomes, peroxisomes, and centrosomes can be easily marked [18] and their unequal distribution between daughter cells has been observed to provoke functional differences that influence their fate. In this respect, asymmetric division of mitochondria in stem cells is associated with the retention of sternness traits in one daughter cell and differentiation in the other [9] and, although the functional consequences are not yet fully determined, to therapeutic resistance in Cancer Stem Cells [39, 40]. Asymmetric inheritance of lysosomes, together with mitochondria, appears to play a crucial role in determining the fate of human blood stem cells [41]. Furthermore our approach could be extended to investigate cell growth dynamics and the origins of cell size homeostasis in adherent cells [21, 24, 4244],

Here, we measured the segregation statistics of the cytoplasm of different colon cell types, finding that they are neatly distinct. CaCo2 cells exhibits the highest degree of cytoplasmic partition fluctuations, while HCT116 division is the one with lowest asymmetry. It is interesting to note that normal colon fibroblasts sit in between, suggesting that the cancer-associated dis-regulation of the normal cell activities leads to the reactivation of asymmetric cell division in a cancer-specific manner. Indeed, CaCo2 cells are known to display a markedly heterogeneous distribution of phenotypes [26]. Notably, comparing the degree of asymmetry in the division of the cytoplasm in the three different analysed cell types, we found that microscopy-derived fluctuations are systematically higher than those measured by flow cytometry. While this could be linked to the different used protocols, it seems more probable that the difference is due to a systematic overestimation of the fluctuations measured by microscopy due to its lower statistics.

We have assessed the stability of our model with respect to perturbations arising from the co-occurrence of multiple noise sources. In particular, we explored the case where fluctuations at division and variability in cell cycle length are dependent. In particular, we investigated a strongly coupled scenario in which division is triggered upon reaching a specific size (sizer strategy) and partitioning fluctuations correspond to size inheritance. We observed that -even under this extreme condition- our approach yields reliable parameter estimates that remain consistent with the experimental error. Moreover, we observed coupling-dependent increase in μg for newer generations. This observation opens new questions and potential applications for our method.

Finally, to get insights on the possible mechanisms behind the observed asymmetries, we looked for specific shape of the partition distribution. Following a maximum entropy approach [27, 45, 46], we started by considering the distribution that employs the fewest assumptions about the system, i.e., we consider independent segregation, which directly maps into binomial partitioning. Indeed, this is a common assumption in segregation models that traces back to the pioneering works of Berk [47] and Rigney [48]. Notably, our results indicated biased segregations with colon cancer cells displaying different levels of biases. To explain their origin, we looked at relative sizes between the daughter cells right after division. The found correlation between daughter sizes and the fractions of inherited fluorescence vouches size as the origin of the observed biased segregation.

In conclusion, by combining experimental data with statistical modeling, we show how it is possible to use flow cytometry data to extract reliable estimates of the strength of fluctuations during cell division. Our approach has the potential to be applied to different cell types where a quantification of the level of division asymmetries daughter cells can experience may provide new insights into the mechanisms of asymmetric cell division and its role in cancer heterogeneity and plasticity.

Materials and Methods

Cell Culture

Caco-2, colorectal adenocarcinoma cell line, was purchased from ATCC (Manassas, VA, USA, HBT-37) and maintained in complete culture media DMEM (D6046) containing 20% FBS, penicillin/streptomycin plus glutamine, nonessential amino acids (NEAA) and sodium pyruvate, all 1/100 dilution. HCT116 VIM RFP, colorectal carcinoma cell line, was purchased from ATCC (Manassas, VA, USA, CCL-247EMT) and maintained in complete culture media McCoy’s 5A (M9309) containing 10% FBS, penicillin/streptomycin plus glutamin.

CCD-I8C0, colon fibroblast cell line, were purchased from ATCC (Manassas, VA, USA, CRL-1459) and maintained in complete culture media DMEM containing 10% FBS, penicillin/streptomycin plus glutamin, nonessential amino acids (NEAA) and sodium pyruvate, all 1/100 dilution.

All cells were kept in culture at 37° C in 5% CO2 and passaged according to the experimental protocol. For all experiments cells were harvested, counted, and washed twice in serum-free solutions and resuspended in room temperature (RT) PBS w/o salts for further staining.

In detail, to detach cells from culture flasks we used Trypsin EDTA 0.25%. Culture media was removed and the flasks were washed once with PBS, then Trypsin EDTA was added and allowed to work for 1 min in an incubator at 37°C. Once the cells appeared detached the flasks were mechanically agitated to facilitate cells’ loss of adhesion. Trypsin was then inactivated with one volumes of complete media and cells were collected, pelleted, and washed a second time in PBS. To determine cell viability, prior to dye staining, the collected cells were counted with the hemocytometer using Trypan Blue, an impermeable dye not taken up by viable cells.

Flow Cytometry and Cell Sorting

To track cell proliferation by dye dilution for establishing the progeny of a mother cell, cells were stained with CellTrace Violet™ (CTV). The CTV dye staining (C34557, Life Technologies, Paisley, UK), used to monitor multiple cell generations, was performed according to the manufacturer’s instruction, diluting the CTV 1/1000 in 0.5-1 ml of PBS for 20 min in a water bath at 37°C’, mixing every 10 min. Afterward, 5x complete media was added to the cell suspension for an additional 5 min incubation before the final washing in PBS.

Labeled cells were sorted using a FACSArialll (Becton Dickinson, BD Biosciences, USA) equipped with Near UV 375nm, 488nm, 561nm, and 633nm lasers and FACSDiva software (BD Biosciences version 6.1.3). Data were analyzed using FlowJo software (Tree Star, version 10.7.1). Briefly, cells were first gated on single cells, by doublets exclusion with morphology parameters area versus width (A versus W), both side and forward scatter. The unstained sample was used to set the background fluorescence. The sorting gate was set around the max peak of fluorescence of the dye distribution. In this way, the collected cells were enriched for the highest fluorescence intensity. Following isolation, an aliquot of the sorted cells was analyzed with the same instrument to determine the post-sorting purity and population width, resulting in an enrichment 〉 97 % for each sample.

To monitor multiple cell division, the sorted cell population was seeded in 12 wells plates (Corning, Kennebunk, ME, USA) at cell density bewteen [30—70] · 103 cells/well, according to the experiment, and kept in culture for up to 84 h. Each well corresponds to a time point of the acquisition and cells in culture were analyzed every 24, 36, 48, 60, 72, 84 h by the LSRFortessa flow cytometer. In order to set the time zero of the kinetic, prior culturing, a tiny aliquot of the collected cells was analyzed immediately after sorting at the flow cytometer. The unstained sample was used to set the background fluorescence as described above.

Time-lapse Microscopy

To better investigate the cell proliferation dynamics, we performed time-lapse experiments, for up to 3 days. In our previous paper [18] we verified that cell growth is not affected by different dyes combination. Therefore, in the present work, we used both CellTrace Yellow™ (C34573 A, Life Technologies, Paisley, UK) and Cell-Trace Far Red™ (C34572, Life Technologies, Paisley, UK), to stain the cytoplasm for microscopy analysis. Low passage cells at around 60 % confluency were counted and stained with CellTrace dyes 1/500 and plated on IBIDI cell imaging chambers (μ-Slide 4 and 8 wells) at low cell density of 103 cells/well. After overnight incubation, the chamber is transferred to the inverted microscope adapted with an incubator to keep cells in appropriate growing conditions. Brightfield and fluorescent confocal image stacks were acquired with a 20x air objective (Olympus, Shinjuku, Japan) and Zen Microscopy Software (Zeiss, Oberkochen, Germany), every 20 min. Time-lapse images were analyzed using ImageJ and inhouse Python programs.

Data Availability

The data that support the findings of this study are available from the corresponding author upon reasonable request.

Code Availability

All codes used to produce the findings of this study are available from the corresponding author upon request. The code for the Gaussian Mixture algorithm is available at https://github.com/ggosti/fcGMM.

Acknowledgements

This research was partially funded by grants from ERC-2019-Synergy Grant (ASTRA, n. 855923); EIC-2022-PathfinderOpen (ivBM-4PAP, n. 101098989); Project ‘National Center for Gene Therapy and Drugs based on RNA Technology’ (CN00000041) financed by NextGeneration EU PNRR MUR—M4C2—Action 1.4—Call ‘Potenziamento strutture di ricerca e creazione di campioni nazionali di R&S’ (CUP J33C22001130001).

Additional information

Author contributions statement

M.M. and G.G. conceived research; D.C., G.P., and S.S. performed cell biology experiment and flow cytometry measurements. D.C., C.G., V.d.T., and M.M. performed timelapse microscopy experiments; F.G. and G.R. contributed additional ideas and helped interpret results; D.C. analyzed data; M.M. performed analytical calculations; D.C. and G.G. performed numerical simulations and statistical analysis; all authors analyzed results, and wrote and revised the paper.