 Research article
 Open Access
 Published:
Stochastic modeling of aging cells reveals how damage accumulation, repair, and celldivision asymmetry affect clonal senescence and population fitness
BMC Bioinformatics volume 20, Article number: 391 (2019)
Abstract
Background
Asymmetry during cellular division, both in the uneven partitioning of damaged cellular components and of cell volume, is a cell biological phenomenon experienced by many unicellular organisms. Previous work based on a deterministic model claimed that such asymmetry in the partitioning of cell volume and of agingassociated damage confers a fitness benefit in avoiding clonal senescence, primarily by diversifying the cellular population. However, clonal populations of unicellular organisms are already naturally diversified due to the inherent stochasticity of biological processes.
Results
Applying a model of aging cells that accounts for natural celltocell variations across a broad range of parameter values, here we show that the parameters directly controlling the accumulation and repair of damage are the most important factors affecting fitness and clonal senescence, while the effects of both segregation of damaged components and division asymmetry are frequently minimal and generally contextdependent.
Conclusions
We conclude that damage segregation and division asymmetry, perhaps counterintuitively, are not necessarily beneficial from an evolutionary perspective.
Background
Even though the somatic cells of multicellular organisms accumulate significant amounts of agingrelated damage throughout the lifetime of the organism, their young progeny generally start with low levels of protein damage. A number of mechanisms have been proposed to explain this phenomenon, generally involving some special way of eliminating damage in germ cells, or elimination of damaged germ cells [1,2,3,4,5]. A major hallmark of aging in higher organisms is the depletion or dysfunction of stem cells, which accumulate various forms of molecular damage during the aging process [6,7,8]. For unicellular organisms such as the budding yeast Saccharomyces cerevisiae or the fission yeast Schizosaccharomyces pombe undergoing mitosis, there is no somatic/germ cell distinction. Yet both S. cerevisiae and S. pombe exhibit lineagespecific aging [9,10,11,12,13,14]. In the budding yeast, for instance, the mother cell is known to accumulate agingrelated damage markers such as extrachromosomal rDNA circles (ERCs) as it ages and eventually enters replicative senescence and dies, while the daughter cells that bud off from the mothers are mostly protected from the accumulated ERCs and generally enjoy a full replicative lifespan even if born from old mother cells [15]. Similar segregation of damaged proteins has been observed during the binary fission of S. pombe, where oxidatively damaged carbonylated proteins are concentrated in one of the two daughter cells – in this case, the one carrying the previous birth scar [9]. Lineagespecific aging has also been observed in the bacteria Caulobacter crescentus and Escherichia coli [10, 16,17,18,19], leading some to analogize the lineagespecific aging behavior in unicellular organisms to the somatic/germ cell distinction in higher ones.
The observation of damagepartitioning behavior even in unicellular species naturally raises the question of whether there is any selective advantage resulting from it. Using a deterministic model based on ordinary differential equations (ODEs) with fixed parameter values, Erjavec and colleagues examined two forms of asymmetry during the celldivision process, which we will denote as division asymmetry and damage segregation, respectively. The former refers to the asymmetric partitioning of cell volume (as naturally seen in S. cerevisiae) between the two cells after division, while the latter refers to the asymmetric partition of damaged cellular components. They concluded that such behavior indeed confers a fitness advantage: under their model, both damage segregation and division asymmetry allowed the population to survive a higher level of protein damage without entering clonal senescence than it otherwise would have [9]. The researchers attributed this effect to the ability of these mechanisms to diversify individuals within a population; in the absence of these mechanisms, the population of cells in their simulations are entirely homogeneous [9].
The fact that real cells are not homogeneous at all raises questions about the reliability of these predictions made based on such a fully deterministic model. It is well known that the expression level of genes can fluctuate substantially, even among cells that are genetically identical or indeed in the same cell over time [20,21,22,23,24,25,26]. This kind of fluctuations, commonly known as noise, can come from a variety of sources: celltocell variations in the abundance of transcription and translation machinery (such as RNA polymerase, general transcription factors, and ribosomes), for instance, or the stochastic nature of transcription events that take place in any single cell [20, 27]. Indeed, it has been shown that stochastic noise can cause drastic differences between reality and what a deterministic model predicts [28].
Bringing a more realistic approach to the study of aging cells in terms of damage accumulation, segregation behavior, and their effects on clonal senescence and fitness, here we investigate whether, and under what circumstances, damage segregation and division asymmetry confer fitness advantages in freely dividing unicellular organism populations when noise is taken into account. We focus on two forms of fitness advantages: resistance against clonal senescence, and increased rates of population growth. We find that damage mitigation and the rate of damage accumulation play major roles in determining the fitness of the cells.
Materials and methods
Modeling of cell growth and division in aging cells
We consider a cell that grows exponentially in volume during the cell cycle [29] and accumulates damage as it grows (Fig. 1a). Cells are assumed to accumulate damage (D) at a constant rate r_{dmg}, and reduce damage via two sources, actively by repair and passively by dilution due to cell division and volume growth:
Here, D is an abstract value representing the concentration of damaged cellular components and other harmful artifacts of aging. For simplicity, we assume that the forms of damage represented by D are freely diffusible, segregable, and repairable, and that there are no other sources of damage affecting cell growth. The rate of damage repair r_{repair} as a function of D is assumed to follow MichaelisMenten kinetics with parameters v_{max} and k_{m}:
The cell volume (V) grows exponentially at a rate that is slowed by accumulated damage:
where r_{growth} is the maximum growth rate constant and α is a nonlinearity coefficient.
In the initial population of cells, each cell starts at an initial volume V_{i}. A cell is assumed to divide when it reaches a generationdependent critical volume V_{crit} (Fig. 1b). This critical volume increases linearly with replicative age (Additional file 1: Table S1), consistent with the observations on single budding yeast cells [26]. The parameters of volume growth during cell cycle were selected to roughly correspond with the microscopic growth dynamics measured in budding yeast cells (Additional file 1: Table S1) [30], with an approximate expected damagefree doubling time of 100 min for symmetrically dividing cells. We separated global noise into two categories: noise in cell volume control (n_{v}), and noise in damage and its repair (n_{d}). In each case, global noise was simulated as a random perturbation applied to each corresponding parameter: the initial parameter value of each individual cell was sampled from a normal distribution N(μ = p, σ = np), where p is the selected mean parameter value from Additional file 1: Tables S1S2 and n is the applicable noise level.
During cell division, we consider the original cell (“mother”) to retain its identity and produce a new daughter cell, for ease of reference. The accumulated damage is distributed between mother and daughter cells as follows. Let D be the damage level of the mother cell before division, then the damage level of the newly produced daughter cell is equal to (1 − s)D, where s in the range [0, 1] is the parameter quantifying the extent of damage segregation, and the damage level of the mother cell after division is \(\frac{D\left(R+s\right)}{R}\), where R is the ratio of the volume between the mother and daughter components after division. In other words, the total amount of the damage (equal to the damage level times volume) is constant across the specific cell division event.
At each cell division, we assumed that each daughter cell partially inherits its mother’s parameter values for all parameters listed in Additional file 1: Tables S1S2. For each parameter p, the new value p_{new} follows the relationship p_{new} = c p_{mother} + (1 − c) p_{fresh}, where p_{fresh} is a parameter value freshly sampled from the normal distribution applicable to that parameter, p_{mother} is the parameter value for the mother cell, and c is a constant in the range [0, 1] characterizing the level of inheritance: when c = 0, the daughter gets a new parameter value from the same distribution used to generate the parameter values used for the initial population of cells, while when c = 1, the daughter perfectly inherits its mother’s parameter value.
Since the simulated cells, just like real ones grown without nutrient limitations, exhibit exponential growth and can easily overwhelm the computing capacity if left unchecked (Fig. 1d), we kept the population size low by means of periodic resampling as a mimicry of using a turbidostat [24]. Because cells slow to divide due to damage are expected to be rapidly overtaken by fasterdividing cells, we did not include a separate procedure for killing cells due to accumulated damage.
To keep the generality of the model intact, we determined the parameter values to use in our model by combinatorially selecting from a grid spanning a wide range (Additional file 1: Table S2), with a total of 7.185 million sets of parameter values tested. The parameter bounds are handselected to cover the arguably biologically plausible range. The noise parameters were capped at 10% because each parameter is perturbed independently, and so the noise in each parameter is expected to combine to produce higher noise in the overall phenotype. We chose the range of r_{dmg} so that the maximum will virtually immediately block cell growth in the absence of repair, and then chose the range of the repair parameters to match the range of r_{dmg}. These parameters are also sampled on a logarithmic scale so as to capture a wide variety of damage strengths.
For each parameter set, we recorded its population size trajectory over the course of the simulation. From these numbers we calculated the average doubling time of the population. If the calculated average doubling time was greater than 1500 min, it was treated as 1500 min for the purpose of analysis. Each simulation was repeated three times and the average doubling time resulting from the three repeats was calculated. For parameter sets producing reasonable fitness levels (< 400 min doubling time), we do not observe significant changes in the computed doubling time if the initial 600 min of the trajectory is omitted. This is expected since these populations relax rapidly to the steady state and the initial conditions have very limited impact when averaged over the long course of the simulation.
Simulations of the stochastic model
The model as described above was implemented in C++ using customwritten code, utilizing the SUNDIALS library [31]. The complete set of model parameters are summarized in Additional file 1: Tables S1S2.
Simulations for each parameter set chosen according to the tables above were performed from an initial set of 2000 cells. Every 40 min, the number of cells in the population was recorded and the foldchange in population growth from the previous time point was calculated, then the population was randomly resampled down to 2000 cells.
For the illustration of exponential growth as shown in Fig. 1e, the simulation was performed as described above, except that

The population size was recorded every 20 min;

Resampling was not performed until the population size reached 1000 times the initial sample size (i.e., 2 million cells);

The population size for the resampling was 50 times the initial sample size (100,000 cells).
Analysis of simulation results
We quantified the effect of changing the value of a parameter P on fitness (resistance against clonal senescence or increased rate of population growth) as follows. For simplicity, we denote the nine parameters of the model P_{1}, … , P_{9}, which can take N_{1}, … , N_{9} possible values, respectively (Additional file 1: Table S2). Without loss of generality, let P_{1} be the parameter P we want to examine. Then, we partition the 7.185 million combinations into \({G}_1=\prod \limits_{i=2}^9{N}_i\) groups of N_{1} combinations each, where the combinations in each group only differ in the value of P_{1} (and have the same values of P_{2}, … , P_{9}). For each group, we then computed a minimum and a maximum doubling time, from which we determined whether changing the value of P_{1} for this particular set of parameter value combination could cause a significant change in fitness (for clonal senescence, a difference in outcome; for growth rate, defined as more than 5% difference between minimum and maximum doubling time). Repeating this for all G_{1} groups, we found that, in M_{1} of them, changing the value of P_{1} resulted in a significant change in fitness. Then, the frequency at which a change in the value of P_{1} could significantly alter fitness was M_{1}/G_{1}. The total number of groups for all parameters combined is G = ∑ G_{i} = 12.56 million. For clonal senescence, we found changes in M = ∑ M_{i} = 1.41 million groups. For fitness, we found significant changes in M = ∑ M_{i} = 2.38 million groups.
Quantification of relative abundance
A key aspect of our analysis is the quantification of the relative abundance of each possible value of a parameter Q among the parameter combinations where changes in the value of a different parameter P has a significant effect (> 5%) on fitness. We denote the possible values of Q as Q_{1}, … , Q_{n}, and also denote the number of combinations where Q = Q_{i} and changes in P can cause a change (> 5%) in fitness as \({Z}_i^P\). We further partition \({Z}_i^P\) into three groups (colored blue, green and red in the figure panels) based on the value of P at which maximum fitness is attained (i.e., \({Z}_i^P={Z}_i^P(blue)+{Z}_i^P(green)+{Z}_i^P(red)\)), blue if P is at its largest possible value, red if P is at its smallest possible value, and green if P is at an intermediate value.
Since the possible values of Q are arbitrarily selected from a large grid, the values are not necessarily equally responsive to fitness changes. Thus, as a normalization measure, we normalized the value of \({Z}_i^P\) (and the partitioned) by the total number of combinations where Q = Q_{i} and changes in the value of any other parameter can cause a change (> 5%) in fitness. In other words, the normalized value is \({S}_i^P=\frac{Z_i^P}{\sum \limits_P{Z}_i^P}\). Similarly, for each color C the normalized value is \({S}_i^P(C)=\frac{Z_i^P(C)}{\sum \limits_P{Z}_i^P}\).
The value of \({S}_i^P\) can vary significantly depending on the identity of the parameter P. To make the abundance value more uniform across panels, we further multiplied \({S}_i^P\) and \({S}_i^P(C)\) for each color by a scaling factor equal to \(\frac{100}{\sum \limits_i{S}_i^P}\) to produce the normalized abundance of Q_{i} plotted in each panel. Thus, the normalized abundance points within in each panel add up to a constant value of 100.
Results
Dissecting the key parameters affecting the clonal senescence outcome
In our model, a cell reaches senescence when its growth rate is slow enough as to virtually stop dividing. A clonal population of cells exhibits clonal senescence if every single cell in the population reach senescence. For the purposes of our analysis, we classified a cell population as clonally senescent if it exhibits an average doubling time greater than 1000 min over the course of the simulation, which is more than 10 times the expected doubling time in damagefree conditions.
To determine the degree of importance of a model parameter for clonal senescence, we examined how likely it is for changes in the value of one parameter to alter the senescence outcome. More formally, for each parameter P, we partitioned the 7.185 million parameter value combinations into disjoint groups such that the combinations in each group only differ in the value of P, and calculated the fraction of groups whose combinations diverge in the clonal senescence outcome. While the absolute value of this fraction is necessarily dependent on the values of the other parameters we picked in the study, the relative value is still indicative of whether the parameter is relevant generally, or only when the other parameters are in a relatively narrow region of the parameter space.
Changes in the damage rate r_{dmg} caused a different senescence outcome in 93% of parameter combinations, and changes in the repair rate v_{max} caused a different outcome in 60% of cases. These observations were intuitive and reaffirmed the validity of the model setup, as the strongest effect was exerted on the amount of accumulated damage, with the possible values of the parameters spanning a wide range. As expected, we find the mostfit combination to be when the damage rate is lowest and the repair rate is highest (Fig. 2a, as indicated by red and blue coloring of their respective bars).
Changes in damagerelated noise (n_{d}), the Michaelis constant (k_{m}) for repair, and the nonlinearity of damage’s effect on growth (α) are less likely to affect the clonal senescence outcome. In only 4.6% of the parameter combinations did a change in α affect the clonal senescence outcome; for k_{m}, the number is slightly higher at 6.1%, while for n_{d}, it is lower at 2.3%. These parameters also affect the amount of accumulated damage or its effect on the cell, but the effects are weaker and less direct. When changes in these parameters did affect the senescence outcome, the direction is essentially uniform: in almost all of the cases, clonal senescence is avoided by having high noise, low k_{m}, or low α (Fig. 2a, color). This again makes sense: a lower k_{m} means a higher repair rate, while a lower α means a higher growth rate (when D > 1, which is necessary for clonal senescence to even come into play because if D < 1 then the volume growth rate can’t fall below half of the maximum growth rate). In the borderline cases where noise matters, moreover, higher noise means a better chance to come across good parameter values that could sustain the population.
Damage segregation and division asymmetry only affected the clonal senescence outcome in a very small fraction of parameter combinations – 1.6 and 0.4% respectively (Fig. 2a). We did find that damage segregation is overwhelmingly beneficial in the few cases where it did matter: in 99% of the cases in which segregation made a difference on the senescence outcome, some damage segregation (represented as the blue and green portions of the bar) was needed to avoid clonal senescence (Fig. 2a). On the other hand, division asymmetry is more likely to be detrimental, if not overwhelmingly so: in 60% of the cases where asymmetry made a difference in the senescence outcome (represented by the red portion of the bar), lack of asymmetry is necessary to avoid clonal senescence, while in the remaining 40% some level of asymmetry is necessary.
We therefore conclude that damage segregation and division asymmetry are not the main effectors of the senescence outcome. Interestingly, neither mechanism is capable of altering the total amount of damage accumulated in the population, which appears to be the key determinant. Thus, changing the damage accumulation rate and the maximum repair rate are most likely to cause (or avoid) clonal senescence.
Characterizing the effect of ageassociated damage on population fitness
Clonal senescence, which implies a complete loss of fitness, is a drastic outcome, and it is certainly conceivable that a parameter might affect population fitness incrementally without causing the entire population to become senescent. We therefore examined the ability of each model parameter to affect the growth rate (or fitness) of the population (Fig. 2b). For this analysis, we calculated the doubling time output of our model using the parameter sets determined as described in the previous section, and examined the cases where the change in the value of one parameter could lead to a significant change (> 5%) in doubling time. For each parameter P, we again partitioned the 7.185 million parameter value combinations into disjoint groups such that the combinations in each group only differ in the value of P, and calculated the fraction of groups whose maximum and minimum doubling time are different by more than 5%.
The parameters most likely to affect the senescence outcome are also most likely to have strong fitness effects
As in the output of clonal senescence, we found that r_{dmg} and v_{max} were the two parameters most likely to cause a fitness differential (> 5% in terms of doubling time). Changing the damage accumulation rate r_{dmg} is capable of significantly altering fitness in more than 99% of all parameter combinations used, while changes in the maximum repair rate v_{max} significantly altered fitness in 78% of the parameter combinations (Fig. 2b). Other damagerelated parameters are also more likely to affect fitness: changes in k_{m} and α are each capable of affecting fitness in about onefifth (23 and 22%, respectively) of the parameter combinations tested, compared to 3% for the volumemodule noise, the parameter that turned out to be the least likely to cause fitness differences (Fig. 2b). For three of the four parameters, moreover, the effect of a parameter value change on fitness is monotonic (as indicated by the prevalence of one color in the figures): a higher v_{max} (Additional file 1: Figure S2) virtually always improved fitness, while a higher r_{dmg} (Additional file 1: Figure S4) or k_{m} (Additional file 1: Figure S3) decreased fitness. This is expected given the functional forms linking these parameters to the model. α, on the other hand, turned out to be a parameter with a doubleedged impact on fitness, though unsurprisingly (Figs. 2b, 3, blue and red color). Introducing ultrasensitivity means that some damage levels will have less impact on fitness while others will have more. Depending on the other parameter values, then, the impact of α on fitness could and did go in both directions (Fig. 3).
Effect of damage segregation on fitness
Somewhat surprisingly, we found that changes in damage segregation affected fitness in only 7% of the parameter value combinations used, compared to 6% for division asymmetry, and 7% for inheritance level and damagerelated noise (Fig. 2b). To gain additional insights into what other parameters may interact with damage segregation to produce a fitness effect, we next examined those cases in which damage segregation could cause a significant change in population fitness. We found that most of the cases where damage segregation produced a fitness effect were seen when the damage rate was low and the repair rate was even lower (Fig. 4e, g), meaning that dilution is the primary method of damage reduction instead of active damage repair, giving significantly more prominence to the ability to sequester damage in the mother compartment. Another interesting observation related to these cases of parameter values is that higher values of α caused damage segregation to behave more like a doubleedged sword: when α = 4, there were as many cases when damage segregation reduced fitness as when it improved fitness (Fig. 4f; compare red vs blue). This can be explained by the fact that the ultrasensitivity of fitness to damage level caused by the high nonlinearity diminishes the impact of damage segregation on fitness when the damage level is low but amplifies the effect when a threshold is crossed – in either direction. Once the damage rate becomes higher, however, damage segregation becomes more of a doubleedged sword, causing fitness decreases about as often as it causes fitness increases. As a modular validation of the modeling approach we took in this study, consistent with previous reports [32], we also found that the highest damage rate (1 min^{− 1}) means that segregation is more likely to be beneficial compared to lower damage rates (0.1 or 0.2 min^{− 1}) (Fig. 4e).
Effect of division asymmetry on fitness
We next examined the effect of introducing division asymmetry on population doubling time or fitness. Introduction of division asymmetry caused significant changes (> 5%) in doubling time in only 6% of the parameter value combinations used, making it the parameter the second least likely to alter the fitness outcome (Fig. 2b). Unlike damage segregation, the effects of division asymmetry on fitness are more likely to be doubleedged: in 40% of the cases where a change in the division asymmetry parameter significantly altered fitness, the highest fitness was seen when there was no division asymmetry (Fig. 2b, red color). Like damage segregation, this behavior was mostly seen when the damage rate was low and the repair rate was even lower (Fig. 5e, g), corresponding to situations where dilution is the primary means of reducing damage levels. Moreover, such detrimental effects from division asymmetry were only seen when some damage segregation was present (Fig. 5c, red line). We interpret this result as follows. The smaller daughter cell usually takes longer to reach the volume threshold before it is ready to divide for the first time, which drags down the volume growth (and therefore dilution) rate at the population level. This effect is more pronounced when the damage segregation mechanism enriched the larger mother compartment with damage, slowing the growth of the mother compartment and dilution of damage. At higher rates of damage accumulation and repair, division asymmetry also becomes predominantly beneficial (Fig. 5e, g, blue line).
Effect of noise on fitness
Just as in the case for clonal senescence, we found that in most cases, higher noise levels had a beneficial effect on population fitness (Figs. 6 and 7, blue); this was particularly pronounced when the inheritance level was high (Figs. 6d, 7d). We interpret this as due to the high inheritance level permitting the propagation of “good” sets of parameter values selected by chance (and is more likely to be chosen if the noise value is high).
Summary
Overall, we found that model parameters directly affecting the accumulation of damage and its effect on the growth rate are the most likely to affect fitness. The introduction of damage segregation or division asymmetry, on the other hand, had a significant effect on population fitness in only a small fraction of the cases and in a contextdependent manner; however, when there was an effect, it was far more likely to be beneficial than detrimental. The fitness impact of asymmetric partitioning of cell volume was similarly contextdependent: when dilution was the predominant mechanism for damage removal, it was more likely to be detrimental, whereas if active damage repair was predominant, it was more likely to be beneficial.
Discussion
The expansive explorations above demonstrate that damage segregation and division asymmetry impacts clonal senescence only in a small fraction of borderline cases once we account for the natural diversity of the cell population. We note that our model differs from the model used in the previously published work [9] in more ways than just the use of randomized coefficients. For instance, our model keeps track of singlecell volume explicitly, while the previous model used the protein count as a proxy for volume. Using protein count for cell volume inevitably led to some questionable assumptions where the partitioning of damaged proteins necessitated the inverse partitioning of undamaged proteins to maintain the protein count (and so the volume) of the daughter cell. Similarly, the sole mechanism for damaged protein to exert a detrimental effect in the previous model was by negative feedback on the production of new proteins, which required the amount of damaged proteins to be roughly comparable to that of intact proteins to have a meaningful effect, in turn requiring likely unrealistic amounts of damage. Our model avoids these problems by using a more abstract “damage level” concept.
Uneven distribution of aging factors between daughter cells following cell division is a wellknown phenomenon that has been observed in a variety of unicellular organisms, and asymmetric partitioning of volume has similarly also been observed in many unicellular organisms. In the present work, we comprehensively examine the fitness impact of these asymmetric damage and volume partitioning schemes and find that they, perhaps counterintuitively, have minimal fitness impact most of the time, as long as the natural diversity of the population is taken into account. When the repair rate was low and dilution was the predominant form of damage elimination, we found damage segregation to be more likely to be beneficial for fitness but division asymmetry to be generally detrimental. On the other hand, when active damage repair was the predominant damage elimination mechanism operating with a high damage repair rate, division asymmetry was found to be generally beneficial for fitness, while damage segregation had a doubleedged impact, becoming more beneficial when the damage accumulation rate was very high.
Overall, our results here indicate that parameters governing the accumulation and elimination of cellular damage are the most important determinants of population fitness. Even though asymmetric partitioning of either cell volume or ageassociated damage might seem beneficial at first glance, neither mechanism actually eliminates any damage on the population level, and, as we show here, they are far from being consistently beneficial evolutionarily. Why, then, are these mechanisms seen in some real organisms? To start with, the fitness impact of both mechanisms is contextdependent; depending on the values of other parameters, representing conditions both intrinsic and extrinsic to the cell, introduction of damage segregation and/or division asymmetry may be either beneficial or detrimental. Thus, it is possible that the organisms exhibiting damage segregation and/or division asymmetry fall within the section of the parameter space where such mechanisms are beneficial rather than detrimental, at least for some portion of their lifespan. And even if this region of the parameter space is not a common occurrence, the importance of avoiding irreversible senescence when it is encountered may be sufficient to preserve the mechanism evolutionarily, similar to how obscure nutrient pathways are preserved due to their critical importance under certain growth conditions. Second, in many cases, introduction of these mechanisms resulted in minimal fitness impact, but even minimal levels of fitness impact can accumulate over time and drive evolutionary selection. Moreover, several natural forms of damage are resistant to active repair, and thus may fall within the region of the parameter space where damage segregation is beneficial. For instance, carbonylated proteins can form aggregates that are resistant to proteasome digestion [33, 34], and ERCs are selfreplicating, suggesting that their effective repair rate – accounting for such selfreplication – is probably also low [15, 35]. Finally, for asymmetric partitioning of damage in particular, it has been reported that some organisms like S. pombe and E. coli only exhibit this behavior during high levels of external stress [14, 17, 32, 36], suggesting that they may actually have evolved mechanisms to activate or inactivate damage partitioning depending on the region of parameter space they are in, just like other stress response mechanisms that are only activated in the presence of stress and can be epigenetically inherited by daughter cells [37, 38].
Conclusions
In this study, we comprehensively examine the effects of ageassociated damage accumulation and removal on the phenotypes of clonal senescence and population fitness. Contrary to the results from a previous study that was based on a fully deterministic model with fixed parameter values [9], we found that neither damage segregation nor division asymmetry played a major role in the avoidance of clonal senescence once the natural diversity of the population is taken into account. Only in a small fraction of borderline cases did introduction of damage segregation eliminated clonal senescence or improve population fitness. As for division asymmetry, not only did it have an effect on clonal senescence or population fitness in an even smaller fraction of cases, but that effect also is frequently detrimental. While we acknowledge that the exact fraction will depend on the set of values and range chosen for the parameters, we believe that the relative value is still a good and useful indicator of the approximate importance and effect size of the parameters. We conclude that damage segregation and division asymmetry, perhaps counterintuitively, are not necessarily beneficial from an evolutionary perspective.
Availability of data and materials
The data and materials related to this work are available upon request from the corresponding author.
References
 1.
Kirkwood TB. Evolution of ageing. Nature. 1977;270:301–4. https://doi.org/10.1038/270301a0.
 2.
Kirkwood TBL, Holliday R. The evolution of ageing and longevity. Proc R Soc B Biol Sci. 1979;205:531.
 3.
Hernebring M, Brolén G, Aguilaniu H, Semb H, Nyström T. Elimination of damaged proteins during differentiation of embryonic stem cells. Proc Natl Acad Sci U S A. 2006;103:7700–5. https://doi.org/10.1073/pnas.0510944103.
 4.
Song R, Sarnoski EA, Acar M. The Systems Biology of SingleCell Aging. iScience. 2018;7:154–69. https://doi.org/10.1016/j.isci.2018.08.023.
 5.
Holliday R. Growth and death of diploid and transformed human fibroblasts. In: Biology of aging and development. Boston: Springer US; 1975. p. 97–106. https://doi.org/10.1007/9781468426311_11.
 6.
Espada L, Ermolaeva MA. DNA damage as a critical factor of stem cell aging and organ homeostasis. Curr Stem Cell Rep. 2016;2:290–8. https://doi.org/10.1007/s4077801600526.
 7.
Moskalev AA, Shaposhnikov MV, Plyusnina EN, Zhavoronkov A, Budovsky A, Yanai H, et al. The role of DNA damage and repair in aging through the prism of Kochlike criteria. Ageing Res Rev. 2013;12:661–84.
 8.
Siudeja K, Nassari S, Gervais L, Skorski P, Lameiras S, Stolfa D, et al. Frequent somatic mutation in adult intestinal stem cells drives neoplasia and genetic mosaicism during aging. Cell Stem Cell. 2015;17:663–74.
 9.
Erjavec N, Cvijovic M, Klipp E, Nyström T. Selective benefits of damage partitioning in unicellular systems and its effects on aging. Proc Natl Acad Sci U S A. 2008;105:18764–9. https://doi.org/10.1073/pnas.0804550105.
 10.
Florea M. Aging and immortality in unicellular species. Mech Ageing Dev. 2017;167:5–15. https://doi.org/10.1016/J.MAD.2017.08.006.
 11.
Liu P, Acar M. The generational scalability of singlecell replicative aging. Sci Adv. 2018;4:eaao4666. https://doi.org/10.1126/sciadv.aao4666.
 12.
Liu P, Young TZ, Acar M. Yeast Replicator: A HighThroughput Multiplexed Microfluidics Platform for Automated Measurements of SingleCell Aging. Cell Rep. 2015;13:63444.
 13.
Barker MG, Walmsley RM. Replicative ageing in the fission yeast Schizosaccharomyces pombe. Yeast. 1999;15:1511–8. https://doi.org/10.1002/(SICI)10970061(199910)15:14<1511::AIDYEA482>3.0.CO;2Y.
 14.
Coelho M, Dereli A, Haese A, Kühn S, Malinovska L, DeSantis ME, et al. Fission yeast does not age under favorable conditions, but does so after stress. Curr Biol. 2013;23:1844–52.
 15.
Sinclair DA, Guarente L. Extrachromosomal rDNA circlesa cause of aging in yeast. Cell. 1997;91:1033–42. https://doi.org/10.1016/s00928674(00)804936.
 16.
Stewart EJ, Madden R, Paul G, Taddei F. Aging and death in an organism that reproduces by morphologically symmetric division. PLoS Biol. 2005;3:e45. https://doi.org/10.1371/journal.pbio.0030045.
 17.
Rang CU, Peng AY, Chao L. Temporal Dynamics of Bacterial Aging and Rejuvenation; 2011.
 18.
Ackermann M. Senescence in a bacterium with asymmetric division. Science. 2003;300:1920. https://doi.org/10.1126/science.1083532.
 19.
Ackermann M, Schauerte A, Stearns SC, Jenal U. Experimental evolution of aging in a bacterium. BMC Evol Biol. 2007;7:126. https://doi.org/10.1186/147121487126.
 20.
Elowitz MB, Levine AJ, Siggia ED, Swain PS. Stochastic gene expression in a single cell. Science. 2002;297:1183–6. https://doi.org/10.1126/science.1070919.
 21.
Raser JM, O’Shea EK. Control of stochasticity in eukaryotic gene expression. Science. 2004;304:1811–4. https://doi.org/10.1126/science.1098641.
 22.
Raser JM, O’Shea EK. Noise in Gene Expression: Origins, Consequences, and Control. Science. 2005;309:2010–3. https://doi.org/10.1126/science.1105891.
 23.
Liu P, Song R, Elison GL, Peng W, Acar M. Noise reduction as an emergent property of singlecell aging. Nat Commun. 2017;8:680.
 24.
Acar M, Mettetal JT, van Oudenaarden A. Stochastic switching as a survival strategy in fluctuating environments. Nat Genet. 2008;40:471–5.
 25.
Acar M, Becskei A, van Oudenaarden A. Enhancement of cellular memory by reducing stochastic transitions. Nature. 2005;435:228–32. https://doi.org/10.1038/nature03524.
 26.
Sarnoski EA, Song R, Ertekin E, Koonce N, Acar M. Fundamental Characteristics of SingleCell Aging in Diploid Yeast. iScience. 2018;7:96–109. https://doi.org/10.1016/j.isci.2018.08.011.
 27.
Swain PS, Elowitz MB, Siggia ED. Intrinsic and extrinsic contributions to stochasticity in gene expression. Proc Natl Acad Sci. 2002;99:12795–800.
 28.
To TL, Maheshri N. Noise can induce bimodality in positive transcriptional feedback loops without Bistability. Science. 2010;327:1142–5. https://doi.org/10.1126/science.1178962.
 29.
Di Talia S, Skotheim JM, Bean JM, Siggia ED, Cross FR. The effects of molecular noise and size control on variability in the budding yeast cell cycle. Nature. 2007;448:947–51.
 30.
Ferrezuelo F, Colomina N, Palmisano A, Garí E, Gallego C, CsikászNagy A, et al. The critical size is set at a singlecell level by growth rate to attain homeostasis and adaptation. Nat Commun. 2012;3:1012.
 31.
Hindmarsh AC, Brown PN, Grant KE, Lee SL, Serban R, Shumaker DE, et al. SUNDIALS: suite of nonlinear and differential/algebraic equation solvers. ACM Trans Math Softw. 2005;31:363–96.
 32.
Vedel S, Nunns H, Košmrlj A, Semsey S, Trusina A. Asymmetric damage segregation constitutes an emergent populationlevel stress response. Cell Syst. 2016;3:187–98.
 33.
Nyström T. Role of oxidative carbonylation in protein quality control and senescence. EMBO J. 2005;24:1311–7. https://doi.org/10.1038/sj.emboj.7600599.
 34.
Andersson V, Hanzén S, Liu B, Molin M, Nyström T. Enhancing protein disaggregation restores proteasome activity in aged cells. Aging (Albany NY). 2013;5:802–12.
 35.
Kaeberlein M, McVey M, Guarente L. The SIR2/3/4 complex and SIR2 alone promote longevity in Saccharomyces cerevisiae by two different mechanisms. Genes Dev. 1999;13:2570–80.
 36.
Wang P, Robert L, Pelletier J, Dang WL, Taddei F, Wright A, et al. Robust growth of Escherichia coli. Curr Biol. 2010;20:1099–103. https://doi.org/10.1016/J.CUB.2010.04.045.
 37.
Chatterjee M, Acar M. Heritable stress response dynamics revealed by singlecell genealogy. Sci Adv. 2018;4:e1701775. https://doi.org/10.1126/sciadv.1701775.
 38.
Xue Y, Acar M. Mechanisms for the epigenetic inheritance of stress response in single cells. Curr Genet. 2018:1–8. https://doi.org/10.1007/s0029401808491.
Acknowledgements
We thank G. Elison and G. Urbonaite for comments and feedback on the manuscript, and K. MillerJensen, C. O’Hern, and Acar Lab members for their comments on the project in its various stages. The authors also acknowledge the use of the Open Science Grid resources.
Funding
MA acknowledges funding from the National Institutes of Health (1DP2AG050461–01 and 1R01GM127870–01). This work was supported in part by the facilities and staff of the Yale University Faculty of Arts and Sciences High Performance Computing Center, and by the National Science Foundation under grant #CNS 08–21132 that partially funded acquisition of the facilities. The funding bodies did not play any roles in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.
Author information
Affiliations
Contributions
RS and MA designed the project, including designing the model, simulations and their analyses. RS implemented the model in C++ and performed the simulations and analyses. RS and MA interpreted the results, and prepared, read, and approved the manuscript.
Corresponding author
Correspondence to Murat Acar.
Ethics declarations
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Competing interests
The authors declare that they have no competing interests.
Additional information
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Additional file
Additional file 1:
Figure S1. AH. Relative representation of the other parameters in the cases where changing the level of inheritance caused a significant (> 5%) fitness difference. Figure S2. AH. Relative representation of the other parameters in the cases where changing the maximum repair rate caused a significant (> 5%) fitness difference. Figure S3. AH. Relative representation of the other parameters in the cases where changing the Michaelis constant k_{m} for repair caused a significant (> 5%) fitness difference. Figure S4. AH. Relative representation of the other parameters in the cases where changing the damage accumulation rate caused a significant (> 5%) fitness difference. Table S1. Parameters with fixed mean values. Table S2. Parameters with values selected from a grid of values. (PDF 1035 kb)
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Song, R., Acar, M. Stochastic modeling of aging cells reveals how damage accumulation, repair, and celldivision asymmetry affect clonal senescence and population fitness. BMC Bioinformatics 20, 391 (2019). https://doi.org/10.1186/s1285901929213
Received:
Accepted:
Published: