The role of backward cell migration in two-hit mutants’ production in the stem cell niche

It has been discovered that there are two stem cell groups in the intestinal crypts: central stem cells (CeSCs), which are at the very bottom of the crypt, and border stem cells (BSCs), which are located between CeSCs and transit amplifying cells (TAs). Moreover, backward cell migration from BSCs to CeSCs has been observed. Recently, a bi-compartmental stochastic model, which includes CeSCs and BSCs, has been developed to investigate the probability of two-hit mutant production in the stem cell niche. In this project, we improve this stochastic model by adding the probability of backward cell migration to the model. The model suggests that the probability of two-hit mutant production increases when the frequency of backward cell migration increases. Furthermore, a small non-zero probability of backward cell migration leads to the largest range of optimal values for the frequency of symmetric divisions and the portion of divisions at each stem cell compartment in terms of delaying 2-hit mutant production. Moreover, the probability of two-hit mutant production is more sensitive to the probability of symmetric divisions than to the rate of backward cell migrations. The highest probability of two-hit mutant production corresponds to the case when all stem cell’s divisions are asymmetric.


Introduction
Studying stem cell dynamics is important for determining the origin of many diseases including cancer, and it may also suggest ways to obtain optimal treatments for these diseases. Stem cell therapy has been used for treating several diseases such as cancer [1]. Recently, scientists are trying to use umbilical cord stem cells (USCs), which are a source of mesenchymal stem cells (HUCMSCs) that promote tissue repair and modulate immune responses, to treat solid tumors. There is evidence that co-culture of rUSCs with Lewis lung carcinoma cells causes cancer cells to remain at the G0/G1 phase [2]. Saliently, in an in vivo study, the injection of rat umbilical cord SCs (rUSCs) could completely abolish rat mammary carcinomas [3].
Knowledge of stem cell division patterns such as their division and death rates, and the rate at which they divide symmetrically or asymmetrically can suggest ways to alter the stem cell niche in order to minimize the number of mutant cells in a tissue. Moran models, which assume a constant number of cells at each updating time step, are commonly used to study cell a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 dynamics [4][5][6][7][8][9], because the number of cells in normal adult tissues stays approximately constant. For instance, it has been observed that the total number of cells in the normal intestinal and colon crypts stays approximately constant [10,11], and because of the fairly simple structure of colon and intestinal crypts, many computational models have been developed to investigate cell dynamics in the crypts [12][13][14][15][16][17][18][19]. Additionally, several mathematical models have been designed to study the interplay between mutants and normal cells [20][21][22][23][24][25][26][27][28][29].
Tissue cells are categorized into two general groups, stem cells and non-stem cells. Stem cells are characterized by their ability to divide both symmetrically and asymmetrically. There are two types of stem cell symmetric divisions: proliferation (two newborn cells are SCs) and differentiation (two newborn cells are TAs). It has been suggested that stem cells in many tissues, including hair, blood, intestine, and brain [30], follow a bi-compartmental structure, which includes 'border stem cells' (BSCs) and 'central stem cells' (CeSCs). Lately, Ristma et al. [31] provided more details about how the two SC compartments, where each consists of approximately 7 SCs, work together to maintain a constant cell population in the mouse intestinal crypt. They observed that the BSCs, which are located between the transit amplifying cells (TAs) and the CeSCs, mostly differentiate in order to control the number of non-stem cells. Additionally, the CeSCs, which are located at the base of the crypt, mostly proliferate to control the total number of SCs. They also found that central stem cells can divide and migrate to the BSC compartment to replace cells in that region. Moreover, a small number of migrations of BSCs to CeCS was observed.
There are several mathematical models suggesting that stem cell symmetric division delays the production of two-hit mutants [6,28,32]. Two-hit mutant production is important because inactivation of tumor-suppressor genes resulting from double-hit mutations is one of the most common causes of carcinogenesis [33]. Recently, computational models have been designed to investigate the role of the bi-compartmental structure of the stem cell niche in the production and spread of mutants [7,8]. We follow a model developed by Shahriyari and Komarova [7], which provides optimal division patterns in the SC niche in terms of minimizing the rate of double-hit mutations. This model does not consider the possibility of the migration of stem cells from the BSC compartment to the CeSC compartment. In the present work, we improve this model by incorporating this backward cell migration between the two stem cell compartments. We use both stochastic numerical simulations and analytical methods to investigate the role of backward cell migration in the probability of two-hit mutant production.

Set-up
Recently, Shahriyari and Komarova [6] developed a bi-compartmental stochastic model, which includes two stem cell groups BSCs and CeSCs, for the stem cell niche, and they obtained the optimal parameters' values, which delay two-hit mutant's production. Here, we modify that model by adding the probability of backward cell migration from BSCs to CeSCs. In this Moran model, at each updating time step, two TA cells die and two stem cells divide to keep the total number of cells constant. We choose two cells in order to accommodate the symmetric division of stem cells. Note that two dead TAs could be replaced by divisions of TAs, however, here we only investigate cell dynamics in the stem cell niche. Therefore, we only focus on divisions in the SC niche. At each updating time step, with probability σ, two SCs will divide symmetrically (a differentiation is coupled with a proliferation), or with probability 1 − σ two SCs will divide asymmetrically. When SCs divide asymmetrically, one of newborn cells is a stem cell and the other one is a TA cell. In each division, we assume with a very small probability that a mutation occurs in one of the newborn children. In the asymmetric division, we assume with probability ν the mutation occurs in the stem daughter cell because of the immortal DNA strand hypothesis [34]. The hypothesis suggests that, upon asymmetric division, the DNA of a stem cell does not segregate randomly, but instead the daughter stem cell retains the parental strand. As a result, SCs pass mutations arising from errors in DNA replication on to their non-stem daughter cells. Although several groups support the immortal DNA strand hypothesis [35,36], some scientists have a concern that there is no convincing experimental confirmation [37].
In this model, asymmetric divisions and differentiations can only happen in the BSC compartment because BSCs are close to the TA compartment. However, proliferations will happen in the CeSC compartment with probability γ, or in the BSC compartment with probability 1 − γ. When a proliferation occurs in the CeSC compartment, we assume a random cell migrates to the BSC compartment in order to keep the number of cells constant at each compartment. Since Ritsma et al. [31] observed backward cell migration from BSCs to CeSCs, we assume with probability α a random cell will migrate from BSCs to CeSCs when a BSC proliferates. When a random cell migrates from BSCs to CeSCs, then another migration happens from CeSCs to BSCs to keep the number of cells constant at each compartment.
We assume that the fitness of all wild-type stem cells is one, while the relative fitness of onehit mutant cells is r. That means that when a division occurs in the CeSC (or BSC) compartment, with probability re Ã re Ã þe (or rb Ã rb Ã þb ) a one-hit mutant cell is chosen to divide, or with probability e re Ã þe (or b rb Ã þb ) a wild-type cell is chosen to divide. Symbols e Ã and e are respectively the number of one-hit CeSC mutants and wild-type CeSCs, and b Ã and b are respectively the number of one-hit BSC mutants and wild-type BSCs. We stop the algorithm as soon as the first two-hit mutant appears in the system. In summary, at each updating time step, we follow the below algorithm, which is presented in Fig 1, and the parameters of the model are given in Table 1.
The algorithm. Although we only keep track of divisions in the stem cell niche, the process begins with the death of two TA cells. These lost TA cells are replaced by either symmetric or asymmetric stem cell divisions in the BSC compartment; stem cells are randomly chosen to divide based on their fitness. With probability 1 − σ two asymmetric divisions will occur, and with probability σ a symmetric differentiation will occur. When there is a symmetric differentiation, further changes in the stem cell compartments must follow to maintain the constant cell population. With probability 1 − γ a second stem cell in the BSC compartment proliferates. And with probability α a migration between the two stem cell compartments will follow with one random cell from the BSC compartment migrating to the CeSC compartment, and one random CeSC migrating to the BSC compartment. The stem cell proliferation will occur in the CeSC compartment with probability γ, along with the migration of one CeSC to the BSC compartment.
Each time a wild type stem cell proliferates, there is a probability 1 − u 1 that another wild type cell is created, and probability u 1 that one of the daughter cells is a single-hit mutant. When a single-hit mutant cell proliferates, there is a probability 1 − u 2 that another single-hit mutant is created and probability u 2 that one of the daughter cells is a double-hit mutant. We chose similar values to those used in [38] for u 1 and u 2 to be able to compare the results. Note that there are some mathematical restrictions for these values in order to be in the tunneling regime, which aligns with reality. The assumption of the tunneling rate is that second-hit mutants are generated before one-hit mutants take over the entire tissue. In Moran models, in order to be in the tunneling regime the fixation rate must be smaller than the tunneling rate. For more details please see the Methods section.
If a mutation occurs when a stem cell divides asymmetrically, it is possible for the mutation to occur in the SC or in the differentiated daughter cell. We include the parameter ν to represent the probability of the mutation occurring in the stem cell.

Tunneling rate
In order to calculate the probability of two-hit mutant production, we obtained tunneling rate -the rate at which the system of a given size produces two-hit mutants (assuming that one-hit mutants drift at relatively low levels). In other words, this is the rate at which the first double-  The role of backward cell migration in two-hit mutants' production in the stem cell niche hit mutant is produced before all cells become one-hit mutants. We found that the tunneling rate R 0 ! 2 satisfies the following formula; the detail of calculations is provided in the Methods section.
Where x and y are the solutions of the following equations.
We also ran stochastic simulations and compared the results of the simulations with the results of the formulas. In the numerical simulations, at each time step, two stem cells were chosen to divide following the above algorithm. We ran the algorithm until producing a twohit mutant or reaching the maximum number of time-steps, T. Then, we checked if a two-hit mutant has been generated in T number of time steps. In order to obtain the probability of two-hit mutant production, we repeated this process for 100 times and obtained the ratio of the number of runs that lead to the production of a two-hit mutant to the total number of runs, which was 100. We repeated this procedure for 10 times to find the mean and standard deviation of the probabilities. In other words, we obtained the probability of two-hit mutant production for 10 times and then we calculated the mean and standard deviation of these 10 obtained probabilities. Note that we chose 100 and 10 because we get a reasonable standard deviation in a fairly fast computing time period using our limited computational power.

Backward cell migration increases the probability of 2-hit mutant production
Although the probability of 2-hit mutant production is an increasing function of α when γ > 0, the probability of double-hit mutant production is not very sensitive to α, when α > 0 (Fig 2). When γ = 1, no cell divides symmetrically in the BSC compartment. Therefore, when γ = 1, the probability of 2-hit mutant production is independent of α, because backward migration occurs when border SCs divide symmetrically. When γ = 0, all stem cell divisions occur in the BSC compartment, and no division occurs in the CeSC compartment. In this case the probability of 2-hit mutant production is a decreasing function of α, because when a mutant migrates to the CeSC it will never have a chance to divide. In other words, when γ = 0 and α is a high number, there is a high probability that one-hit BSC mutants migrate to the CeSC compartment and never divide to generate a double-hit mutant. However, when 0 < γ < 1, the probability of two-hit mutant production is an increasing function of α. This means when stem cells are able to divide symmetrically in both SC compartments, the optimal value for α is zero in terms of delaying 2-hit mutant production. Note, in general the probability of two-hit mutant production is not very sensitive to α as long as α is not zero (Fig 2).
The variation of α changes the optimal values for σ and γ Fig 3 shows a trend of changes in the optimal values for γ and σ when α varies. This figure reveals that a small non-zero percentage of backward cell migration leads to a higher range (or a higher interval) for the optimal values of σ and γ in terms of delaying two-hit mutant production. In other words, the regimes in which probability of two-hit mutant production is small, i.e. the dark areas, are much larger for α = 0.001 than for any other values of α in Fig 3. However, α = 1 corresponds to the worst case, which gives the largest intervals for the values of σ and γ causing a high probability of double-hit mutant production. Fig 3 also shows that the optimal case corresponds to the parameter values, γ % 0.3, σ = 1, and α = 0. In other words, the probability of two-hit mutant production is minimized when SCs only divide symmetrically while around 70% of SC proliferations happen in the BSC compartment, and no migration occurs from BSCs to CeSCs.
Symmetric division delays two-hit mutant production when proliferations mostly occur in the BSC compartment and α > 0.2 Stochastic simulations and analytic calculations show that the probability of two-hit mutant production is crucially affected by the frequency of symmetric divisions (parameter σ) . Fig 2(a) and 2(c) show that when the probability of backward cell migration, α, is more than 0.2 and γ = 0.5, then the minimum of the probability of 2-hit mutant production corresponds to σ = 1. Moreover, when α > 0.2, the probability of two-hit mutant production decreases when σ increases (Fig 4). However, when γ = 0.5 and there are no backward cell migration in the niche (α = 0), then the optimal value for σ is around 0.3 in terms of delaying two-hit mutant production (Fig 4). Importantly, the probability of two-hit mutant production maximized when stem cells divide only asymmetrically (σ = 0) for all values of α (Fig 4). More precisely, if at least 50% The role of backward cell migration in two-hit mutants' production in the stem cell niche The role of backward cell migration in two-hit mutants' production in the stem cell niche PLOS ONE | https://doi.org/10.1371/journal.pone.0184651 September 20, 2017 of times mutations occur in the stem cell offspring in the asymmetric divisions (ν ! 0.5) and stem cells only divide asymmetrically, then the probability of two-hit mutant production is very high (Fig 4).
When all SC proliferations occur in the CeSC compartment, i.e. γ = 1, then SCs should divide both symmetrically and asymmetrically, but mostly asymmetrically to delay two-hit mutant production  4(d) show that the minimum value of the probability of double-hit mutant production corresponds to σ % 0.2, when α = 0.5 and γ = 1. When γ = 1, then no stem cell proliferation occurs in the BSC compartment, and all stem cell proliferations happen in the CeSC compartment. Since mutants in the CeSC compartment have a higher chance of survivability, the probability of production of two-hit mutants minimizes when the rate of proliferation in the CeSC compartment decreases. Thus, when γ = 1, in order to minimize the number of divisions in the CeSC compartment, σ, probability of symmetric division, should be small. Note that the maximum of the probability of two-hit mutant production always corresponds to σ = 0, because when stem cells divide purely asymmetrically, mutant cells will never have a chance to differentiate and be removed from the stem cell niche. However, when stem cells divide The role of backward cell migration in two-hit mutants' production in the stem cell niche symmetrically, there is always a chance that mutant cells differentiate into two mutant nonstem cells and become removed from the stem cell niche.
Probability of two-hit mutant production is not sensitive to the percentage of cells at each compartment indicates that the probability of two-hit mutant production is not very sensitive to the percentage of cells at each compartment. Additionally, like other figures, it clearly shows that when σ = γ = 0.5, the probability of two-hit mutant production is an increasing function of α. Moreover, when α = γ = 0.5, the probability of two-hit mutant production is a decreasing function of σ. In general, Fig 5 indicates that the probability of two-hit mutant production is not very sensitive to the percentage of cells at each compartment, and it is more sensitive to σ than to α.
Note that in the analytical calculations, we assume the number of cells is large in order to be able to derive the formulas. More precisely, we simplify the transition probabilities by assuming the number of mutants is much less than the total number of cells. This simplification leads to a noticeable difference between the result of the formulas and numerical simulations when the number of cells is small, particularly when the number of CeSCs is small (See Fig 5). However, the trend stays the same, and therefore the results do not change.

The immortal DNA strand hypothesis
It has been hypothesized that, in asymmetric divisions, the DNA of a SC does not segregate randomly, but instead the stem daughter cell retains a distinct template set of DNA strands (called the parental strand). Therefore, SCs pass mutations arising from errors in DNA replication onto their TA daughters, which soon terminally differentiate to minimize the probability of mutants' production in the stem cells.
We incorporate this mechanism into our model by introducing a parameter ν, which quantifies the probability of a mutation happing in the stem daughter cell of an asymmetrically dividing stem cell rather than in its TA offspring. Figs 2-5 show the results for ν = 0.5, which The role of backward cell migration in two-hit mutants' production in the stem cell niche corresponds to a complete symmetry between stem and TA offspring (i.e. relaxing the assumption of immortal DNA strand hypothesis). Fig 6 presents the probability of double-hit mutant generation as a function of ν. In the case of ν = 0, the asymmetric division is optimal, because no mutations can occur in SC daughter cells when SCs divide asymmetrically, while CSs have a chance to acquire mutations in symmetric divisions. On the other hand, if ν > 0.5, then symmetric division is optimal. However, if 0 < ν < 0.5, then a mixture of symmetric and asymmetric divisions comprises the optimal strategy. The counter plots in Fig 6 show how the optimal regimes move toward σ = 1 when ν is increasing. The role of backward cell migration in two-hit mutants' production in the stem cell niche

Discussion
One of the most common patterns of carcinogenesis is the inactivation of tumor suppressor genes. Knudson suggested that two hits (mutation in both alleles) are required to inactivate a tumor suppressor gene [33,39]. Many cell dynamics models have been developed to investigate the process of two-hit mutants' production to answer some of the following questions. What is the optimal division pattern that minimizes the probability of two-hit mutant productions? Are healthy tissues following the optimal patterns? Are the optimal patterns violated in abnormal tissues like tumors? Is there anyway that we can alter the tissue's structure to minimize the size of tumors? Lately, mathematical models show that stem cell symmetric divisions delay 2-hit mutant productions [6,28,32], and importantly it has been observed that adult stem cells in normal tissues divide mostly symmetrically [40]. Here, the model indicates that a small number of backward SC migrations, which have been experimentally observed, enlarge the interval of optimal values for other parameters. More precisely, the probability of two-hit mutant production is not very sensitive to a small perturbation in the optimal values of other parameters like the rate of symmetric divisions or CeSC's proliferation rates when there are a small number of backward migrations in the SC niche.
The model also predicts that if γ < 1 and α > 0, then symmetric division delays the production of two-hit mutants. In other words, if the probability of CeSCs's proliferation is less than one, and there is a migration from CeSCs to BSCs, then SCs should divide mostly symmetrically to delay two-hit mutant production. Ritsma et al. [31] observed backward cell migration from BSCs to CeSCs, i.e. α > 0. They also estimated that the probability of proliferation of CeSCs to be 0.3, i.e. 0 < γ < 1. Importantly, it has been observed that the intestinal stem cells divide mostly symmetrically, and asymmetric stem cell divisions rarely occur [40].
The model suggests that backward cell migration very slightly increases the chance of twohit mutant productions when other parameter values are fixed. Importantly, a small non-zero value for α, the probability of backward cell migration, leads to the largest interval of optimal values for σ, the frequency of symmetric divisions, and γ, the probability of CeSC's proliferations, in terms of delaying 2-hit mutant production. Although a small non-zero probability of backward cell migration expands the range of optimal values for σ and γ, a high probability of backward cell migration enlarge the intervals for values of σ and γ that lead to a high probability of two-hit mutant production. Furthermore, the frequency of SC symmetric and asymmetric divisions has a much higher effect on the probability of two-hit mutant production than the probability of backward cell migration. Moreover, the rate of double-hit mutant production is maximized when stem cells divide only asymmetrically.

Transition probabilities
In order to calculate the probability of a two-hit mutant production, we found all other nonzero transition probabilities, which are the probabilities of moving from e Ã mutants in the CeSC compartment and b Ã mutants in the BSC compartment to another state in the model. We define S b as the total number of stem cells in the BSC compartment and S c as the total number of stem cells in the CeSC compartment. We denote the number of wild-type CeSCs by e = S c − e Ã and the number of wild-type BSCs by b = S b − b Ã . The transition probabilities of the model are as follows.
The probability that the number of mutant S b cells increases by one while the number of mutant S c cells does not change: The probability that the number of mutant S b cells decreases by one while the number of mutant S c cells stays the same: Prolif: of a w:t: S b into two w:t: ð1 À aÞ Migrations: a w:t:=mut: S b to S c and a w:t:=mut: S c to S b þ g e re Ã þ e ð1 À u 1 Þ e S c |fflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of a w:t: S c into two w:t: and mig: of w:t: The role of backward cell migration in two-hit mutants' production in the stem cell niche The probability that the number of mutant S b cells decreases by one while the number of mutant S c cells increases by one: Prolif: of a w:t: S b to one w:t: and one mut: Migrations : a mut: S b to S c and a w:t: S c to S b þ rðb Ã À 1Þ rðb Ã À 1Þ þ b ð1 À u 2 Þ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of a mut: S b into two mut:

ffl}|ffl ffl{
Migrations: a mut: S b to S c and a w:t: Prolif : of a mut: S c into two mut: and mig: of w:t: Prolif : of a w:t: S c into one w:t: and one mut: and mig: of a w:t:

zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl}|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{
Prolif: of a w:t: S b into two w:t: Migrations : a mut: S b to S c and a w:t: The probability that the number of mutant S b cells stays the same while the number of mutant S c cells decreases by one: Prolif : of a w:t: S b into two w:t: Migrations : a w:t: S b to S c and a mut: S c to S b þg e re Ã þ e ð1 À u 1 Þ e Ã S c |fflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflffl} Prolif: of a w:t: S c into two w:t: and mig: of a mut: The role of backward cell migration in two-hit mutants' production in the stem cell niche The probability that the number of mutant S b cells increases by one while the number of mutant S c cells decreases by one: Prolif: of a w:t: S b into one w:t: and one mut:

zffl ffl}|ffl ffl{
Migrations : a w:t: S b to S c and a mut: Prolif: of a mut: S b into two mut: Migrations: a w:t: S b to S c and a mut:

zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl}|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{
Prolif: of a w:t: S b into one w:t: and one mut: The probability that the number of mutant S b cells stays the same while the number of mutant S c cells increases by one: Diff: of a w:t: S b " ð1 À gÞ b À 1 rb Ã þ ðb À 1Þ u 1 |fflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of a w:t: S b into one w:t: and one mut:

ffl}|ffl ffl{
Migrations : a mut: S b to S c and a w:t: S c to S b þ rb Ã rb Ã þ ðb À 1Þ ð1 À u 2 Þ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of a mut: S b into two mut: Migrations: a mut: S b to S c and a w:t: S c to S b þg e re Ã þ e u 1 e S c |fflfflfflfflfflfflffl{zfflfflfflfflfflfflffl} Prolif: of a w:t: S c into one w:t: and one mut: and mig: of a w:t: Prolif: of a mut: S c into two mut: and mig: of a w:t: The probability that the number of mutant S b cells increases by two while the number of mutant S c cells decreases by one: Prolif: of a w:t: S b into one w:t: and one mut: Migrations : a w:t: S b to S c and a mut: ð1 À u 2 Þ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of a mut: S b into two mut: Migrations: a w:t: S b to S c and a mut: The probability that the number of mutant S b cells decreases by two while the number of mutant S c cells increases by one: Prolif: of a w:t: S b into two w:t: Migrations: a mut: S b to S c and a w:t: The probability to create a double-hit mutant in the S b compartment: re Ã re Ã þ e u 2 1 S c zfflfflfflfflfflfflffl}|fflfflfflfflfflfflffl{ Prolif: of mut: S c into one mut: and one 2 À hit mut: and mig: of 2 À hit mut: to S b # þ 2ð1 À sÞ rb Ã rb Ã þ b u 2 n |fflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl} two asym: mut: The role of backward cell migration in two-hit mutants' production in the stem cell niche The probability to create a double-hit mutant in the S c compartment: Diff: of a w:t: S b ð1 À gÞ rb Ã rb Ã þ ðb À 1Þ u 2 |fflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of a mut: S b into one mut: and one 2 À hit mut: Diff: of a w:t: or mut: S b g re Ã re Ã þ e u 2 1 À 1 S c zfflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl}|fflfflfflfflfflfflfflfflfflfflfflfflfflffl ffl{ Prolif: of a mut: S c into one mut: and one 2 À hit mut: and all mig: less mig: of 2 À hit mut: to S b þ s rb Ã rb Ã þ b |fflfflfflffl ffl{zfflfflfflffl ffl} Diff: of a mut: S b ð1 À gÞ rðb Ã À 1Þ rðb Ã À 1Þ þ b u 2 |fflfflfflfflfflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflfflfflfflfflffl ffl} Prolif: of mut: S b into one mut: and one 2 À hit mut: The probability to create a double-hit mutant: We can define φ e Ã , b Ã (t) as the probability to have e Ã mutants in the S c compartment and b Ã mutants in the S b compartment at time t. Using the transition probabilities at time t − 1, the Kolmogorov forward equation for the function is given by: The role of backward cell migration in two-hit mutants' production in the stem cell niche We define the probability generating function by: The probability to be in one of the states (e Ã , b Ã ) is given by ψ(1, 1; t). Therefore, the probability to transition to state E is P 2 (t) = 1 − ψ(1, 1; t). The probability generating function gives the following equations: We have, À gÞð1 À u 2 Þx 2 þ ð1 À gÞ a r ðð1 À u 1 Þy þ u 1 xy À xÞ À xð2 À g À u 1 Þ þ 1 À u 1 À 2ð1 À sÞu 2 nxÞ % r S 1 s ð1 À gÞð1 À u 2 Þx 2 þ ð1 À gÞ a r ðy À xÞ À ð2 À gÞx þ 1

Tunneling conditions
Here, we investigate the applicability of the tunneling rate, for more information please see [41]. The tunneling rate is practical when one-hit mutants exist at low numbers, and the creation of a two-hit mutant occurs without prior fixation of 1-hit mutants. A condition for such a scenario is that the fixation rate for one-hit mutants is smaller than the tunneling rate. Let us denote the probability of fixation of 1-hit mutants in the S b compartment (or in the S c group), starting from the state with b Ã number of BSC mutants (or eÃ number of CeSC mutants) by p S b b Ã (or p S c e Ã ). Then tunneling occurs if the probability of absorbing in the state of all one-hit mutants (when there is no chance of new mutations i.e u 1 = u 2 = 0) is less than the probability of creation of new mutation, i.e. p S b b Ã and p S c e Ã ( u 2 . The quantities p S b b Ã and p S c e Ã satisfy the following systems: For more details about obtaining fixation probability please see [42]. We first calculate the probability of transforming from a one-hit mutant to the all 1-hit mutant state by assuming no new mutation will be created, i.e. u 1 = u 2 = 0. Therefore, the non-zero transition probabilities are We calculate the probability of transformation from a single one-hit mutant in the S b compartment (i.e. (e Ã , b Ã ) = (0, 1)) to all one-hit mutants in the S b group (b Ã = S b ). By substituting the transition probabilities Eq (23) in Eq (21), we get This system implies Thus, the fixation probability of moving from the state that has only a single BSC mutant to the state that all BSCs are mutants is given by Additionally, we obtain the probability of transformation from single one-hit mutant in the S c compartment (i.e. (e Ã , b Ã ) = (1, 0)) to the state that all CeSCs are one-hit mutants (e Ã = S c ).