The recently developed concept of cancer stem cells (CSC) sheds new light on various aspects of tumor growth and progression. Here, we present a mathematical model of malignancies to investigate how a hierarchical organized cancer cell population affects the fundamental properties of solid malignancies. We establish that tumors modeled in a CSC context more faithfully resemble human malignancies and show invasive behavior, whereas tumors without a CSC hierarchy do not. These findings are corroborated by in vitro studies. In addition, we provide evidence that the CSC model is accompanied by highly altered evolutionary dynamics compared with the ones predicted to exist in a stochastic, nonhierarchical tumor model. Our main findings indicate that the CSC model allows for significantly higher tumor heterogeneity, which may affect therapy resistance. Moreover, we show that therapy which fails to target the CSC population is not only unsuccessful in curing the patient, but also promotes malignant features in the recurring tumor. These include rapid expansion, increased invasion, and enhanced heterogeneity. Cancer Res; 70(1); 46–56

Major Findings

This research shows that the hierarchical organization of malignant clones, as advocated in the CSC concept, has major implications for tumor biology. CSC-driven tumor growth intrinsically orchestrates tumor invasion, influences clonal selection, and has crucial consequences for the development of successful cancer treatments.

Quick Guide to Equations and Assumptions
Proliferation

The stem cellular automaton (SCA) model is a hybrid cellular automaton (4) in which a biological cell is a point (10 × 10 μm) in a lattice. Each point can be a normal cell, a cancer cell, or a necrotic cell, and has the attributes in Table 1. After a cell division, an empty place is created by shifting the surrounding cells outward.

Table 1.

Cellular and microenvironmental variables in the SCA model

VariablesSymbolValueReference/justification
Proliferation speed (average cell cycle duration) T 20 h (44)
O2 diffusion coefficient Dc 10−5 cm−2 s−1 (16)
O2 concentration, healthy tissue ù 10−4 g cm−3 (16, 45)
O2 consumption, proliferative cells êp 10−6 g cm−3 s−1 (46)
O2 consumption, senescent cells ês 5 × 10−7 g cm−3 s−1 (47)
O2 concentration resulting in necrosis θ 5 × 10−7 g cm−3 s−1 Senescent consumption as minimum concentration for cell survival
Random mobility Dn 10−10 cm−2 s−1 (48)
Maximum number of generations generated by DCCs H We choose to fix H to 5 and vary PS to generate various CSC fractions (3)
VariablesSymbolValueReference/justification
Proliferation speed (average cell cycle duration) T 20 h (44)
O2 diffusion coefficient Dc 10−5 cm−2 s−1 (16)
O2 concentration, healthy tissue ù 10−4 g cm−3 (16, 45)
O2 consumption, proliferative cells êp 10−6 g cm−3 s−1 (46)
O2 consumption, senescent cells ês 5 × 10−7 g cm−3 s−1 (47)
O2 concentration resulting in necrosis θ 5 × 10−7 g cm−3 s−1 Senescent consumption as minimum concentration for cell survival
Random mobility Dn 10−10 cm−2 s−1 (48)
Maximum number of generations generated by DCCs H We choose to fix H to 5 and vary PS to generate various CSC fractions (3)

NOTE: Variables used in the SCA model. See the Quick Guide and Materials and Methods for details.

Major Assumptions

Proliferation concentrates in the proximity of the tumor borders where oxygen levels are higher, the pH is more normal, and the pressure is lower (1214). To model this, we assume the probability of a cell to divide to be linear, from 1 at the tumor edges to 0 at a distance λ = 600 μm (60 cells) from the tumor margins.

Metabolism
$∂c∂t=Dc∇2c−κn$
(A)

Oxygen is critical for cell survival, it diffuses into the tumor and it is consumed by cancer cells. Upon discretization of Eq. (A), cancer cells n consume oxygen at rate κ.

Major Assumptions

Oxygen around the tumor is kept constant by the vascular system. Cell quiescence occurs below an oxygen threshold κs and necrosis below an oxygen threshold θ.

Migration
$∂c∂t=Dn∇2n$
(B)

During tumor progression, cancer cells lose cell-to-cell attachment and invade the surrounding tissues (15). This process is modeled using the Hybrid Discrete-Continuum Technique (16), in which cancer cells n disperse with coefficient Dn. To simulate adhesion, we allow cell movement if the number of neighboring tumor cells gα, where α is the adhesion coefficient of the cell in the range 0 to 4.

Major Assumptions

Cells move according to a random motion coefficient Dn and an adhesion coefficient α.

CSCs

CSCs divide symmetrically with probability PS and asymmetrically with probability 1 − PS. Two new CSCs result from the former, a differentiated cancer cell (DCC) and a CSC from the latter. CSCs possess unlimited replicative potential and self-renewal whereas DCCs can divide up to H times. A CSC growth model has small PS values, whereas for PS = 1, we simulate the classical model of malignancies.

Major Assumptions

We fix H = 5 and vary PS to simulate different CSC frequencies. As suggested experimentally for both CSCs and normal stem cells, we have restricted migration to CSCs (1720).

Tumor Evolution

We introduced tumor phenotypical heterogeneity in some experiments by assuming that every CSC self-renewal division has a chance PMut for the new CSC to acquire a different phenotype (see Supplementary Materials and Methods).

Major Assumptions

Only CSCs contribute to tumor evolution in the long run.

Malignancies arise after sequential accumulations of mutations in oncogenes and tumor suppressor genes (1). During the process of malignant transition, fundamental regulatory mechanisms are lost. The cancerous cell population has unlimited growth potential, evades the immune system and apoptosis, and often acquires the ability to breach tissue boundaries and expand into foreign environments (1). The cancer stem cell (CSC) concept sheds new light on all of these features (2, 3). In this study, we investigate the influence of a hierarchical organization of malignant clones on tumor growth, evolution, invasion, and morphology using a multi-scale cellular automaton–based computer model (4).

### CSCs

Malignancies are highly heterogeneous tissues containing largely diverse cancer cell populations as well as other cells such as fibroblasts and macrophages (1). The dominant genetic view of malignancies explains the heterogeneity in cancer cells with the presence of genetically diverse clones emerging from the continuous acquiring of additional genetic lesions by cancer cells. Such clones compete for resources, resulting in a highly dynamic process known as Tumor Darwinism (5). In this article, we refer to this view of malignancies as the classical model. Although this model greatly contributes to our understanding of malignancies, recent experimental evidence suggests an additional layer of complexity. The heterogeneity present in tumors could, in part, be the result of the diversity in differentiation grade of genetically identical cells (3, 6). For example, in glioblastoma multiforme, cells with an immature phenotype expressing the cell surface marker AC133 are the cells that fuel tumor growth and have the exclusive capacities to self-renew, differentiate, and transplant the malignancy into severe combined immunodeficiency mice (7). Self-renewal and spin-off of differentiated cells are features shared with normal stem cells, and therefore, cells with such features in malignancies are defined as CSCs (6, 8).

### Modeling Tumor Growth

Tumor growth is generally accepted to be the result of several highly complex interacting processes. Fundamental cellular characteristics such as genetic and epigenetic features influence signal transduction route activities that in turn control cellular functions such as mitosis, apoptosis, and cell migration. In addition, environmental factors including nutrients and growth factor concentrations interplay with these processes. To study the emergent properties of such systems regarding proliferation speed, infiltrative growth, and phenotypical evolution of cancer, advanced mathematical models have been developed (911). Here, we apply computational modeling techniques to investigate the consequences of hierarchically organized cancer cell populations on solid tumor growth dynamics and progression. We describe that implementing the developing concept of CSCs in a mathematical tumor growth model directly results in an invasive morphology. Moreover, we found that hierarchical organized malignant clones have highly altered evolutionary dynamics. Most strikingly, the CSC organization promotes phenotypical heterogeneity, a feature that could have immediate consequences for therapy resistance.

### SCA Model

We developed a hybrid tumor growth model based on cellular automata (4) and partial differential equations. We refer to this model as the SCA model. In the SCA model, the individual cancer cell is the fundamental unit of the tumor, we simulate its proliferation, metabolism, migration, stemness, and differentiation (see Quick Guide and Supplementary Materials and Methods).

Implementing the CSC model of malignancies means to simulate cancer cells with different replicative potential within the tumor. For simplicity, we assume that in our model there are only two types of cells: CSCs and DCCs. CSCs possess unlimited replicative potential and could either generate new CSCs (with a probability PS) or DCCs. DCCs can divide for a maximum of H generations before stopping to proliferate irreversibly. This method yields a hierarchy with CSCs at the top and DCCs at the bottom. We simulate the classical model of malignancies, in which all cells possess tumor growth–promoting capacities, by simply setting PS = 1. In such a situation, all cells possess stem cell characteristics. Using this method, we have an intuitive way to compare the flat, classical tumor model with the hierarchically organized CSC model.

### Emergent Invasive Morphology

Computational modeling allows the exploration of highly complex systems, such as tumor growth. In this study, we have used a computational tumor growth model to test the consequences of hierarchical organized clones on a vast range of areas of tumor biology. These include tumor growth dynamics and morphology, but also tumor evolution and therapeutic effects. The applied computational modeling technique allows us to get more insight into the underlying dynamics of these aspects of cancer growth and progression that would be impossible in a conventional experimental biological setting.

First, we investigated how tumor growth dynamics are altered upon varying the CSC fraction. In the SCA model, this corresponds to changing the variable PS. For high values of PS, we expect to model the classical interpretation of tumors because all cell divisions are symmetrical and all cells are therefore clonogenic. In contrast, low values of PS represent the CSC model. We simulated the growth of tumors with PS = 1, PS = 0.1, and PS = 0.03 with the variables described in Table 1.

Figure 1A shows the fraction of CSCs on the total amount of tumor cells for different values of PS. The selected PS values correspond to CSCs populations comprising roughly 100%, 1%, or 0.1% of the total tumor volume and therefore cover mainly the CSC fractions observed in a variety of solid malignancies (2). In the CSC model (PS = 0.1 and PS = 0.03), small (early) lesions have relatively high fractions of CSCs, although this number decreases and tends to stabilize when the tumor progresses. This observation is supported by in vivo studies that find increased numbers of CSC marker–bearing cells in micrometastases compared with larger tumors (21). This indicates that even with fixed self-renewal rates (PS), this phenomenon is intrinsic to lesions initiated by a single CSC although environmental factors influencing self-renewal frequencies are also likely to contribute. Tumor growth curves for various PS values all display the classical Gompertzian-like growth kinetics. However, as expected with equal cell cycle durations, the self-renewal rates of the stem cell fraction influences proliferation rate greatly, hence, low self-renewal rates (small PS) correspond to slow tumor growth (Fig. 1B). Interestingly, the decrease in accumulation of tumor volume is accompanied by a stabilization of the fraction of CSCs, suggesting an intimate relationship between these two processes.

Figure 1.

Emergent invasive behavior in the CSC model. A, different PS values result in different CSC fractions. B, growth curves for different PS values. C, quantitative measure for invasiveness shows increasing invasive behavior with declining PS. A to C, bars, SD (n = 16). D, hierarchical organization in the SCA model affects tumor morphology. Tumors for different values of self-renewal probability (PS) and different volumes. Dark blue, cells which have divided within the last 48 h (depicted larger); light blue, nondividing cells; brown, necrotic center. Right, localization of CSCs in the tumor mass. Gray, tumor mass; red, CSCs (depicted larger). In all figures, 6 × 6 mm of tissue are represented. See also Movies S1 and S2.

Figure 1.

Emergent invasive behavior in the CSC model. A, different PS values result in different CSC fractions. B, growth curves for different PS values. C, quantitative measure for invasiveness shows increasing invasive behavior with declining PS. A to C, bars, SD (n = 16). D, hierarchical organization in the SCA model affects tumor morphology. Tumors for different values of self-renewal probability (PS) and different volumes. Dark blue, cells which have divided within the last 48 h (depicted larger); light blue, nondividing cells; brown, necrotic center. Right, localization of CSCs in the tumor mass. Gray, tumor mass; red, CSCs (depicted larger). In all figures, 6 × 6 mm of tissue are represented. See also Movies S1 and S2.

Close modal

Spatially, all experiments display a three-layered structure consisting of an external proliferative area, an inner senescent layer, and a necrotic core. However, tumor morphologies for different PS values are remarkably dissimilar (Fig. 1D). For PS = 1, in which there is no hierarchy, a symmetrical, sphere-like tumor morphology is generated that closely resembles early tumor growth models (22, 23). In contrast, the shape of the tumors generated with low PS values is highly irregular (especially PS = 0.03). CSC-driven tumors yield highly invasive morphologies with fingering fronts and clusters of cancer cells beyond the tumor margin, driven by the mobility and the exclusive proliferative properties of CSCs (Movies S1 and S2).

From Fig. 1C, it is evident how hierarchically organized tumors (PS = 0.03 and PS = 0.1) generate a higher degree of invasiveness (see Supplementary Materials and Methods), compared with tumors in the classical model (PS = 1). It is important to note that the intrinsic properties of the cells in the classical tumor and the stem cells in the CSC-driven tumor are completely identical.

Magnification of a tumor border in a CSC-fueled tumor growth model (PS = 0.03) shows how CSCs migrate beyond the margins of the tumor mass. CSCs colonize the surrounding tissue and expand locally, forming small satellites that grow back into, and are engulfed by, the main tumor mass (Fig. 2A). This is accompanied by the observation that the invasive behavior of individual experiments that make up Fig. 1C occurs in a wave-like fashion (Supplementary Fig. S1). These results are paralleled by recent findings in a different model system in which high migration levels in a small subset of cells give rise to small proliferating extratumoral lesions and therefore tumors are “conglomerates of self-metastasis” as the authors propose (24).

Figure 2.

Invasive behavior in silico, in vivo, and in vitro. A, close-up of tumor border showing invasive behavior for PS = 0.03. CSCs (red) infiltrate surrounding tissue and spin-off DCCs that proliferate (dark blue). Small satellites are formed in the surrounding normal tissue (white) and grow back to and are engulfed by the main tumor mass. Nondividing cells (light blue). B, representative figures of various malignancies. All images reveal island-like formations at the rim of the main tumor mass (T). These findings are in line with tumor expansion as predicted by the SCA model that implements a CSC hierarchy. C, cell lines (n = 24) have been plated at clonal density in Matrigel. Simultaneously, the clonogenic fraction of the lines has been determined by limiting dilution analysis. A significant (P = 0.02) inverse relationship exists between clonogenicity and invasion as quantified by the measure of invasiveness we defined. Two examples of cell lines are shown and invasiveness (IM) and clonogenicity are indicated (right). See Supplementary Figs. S2 to S5 and Supplementary Table S1 for details.

Figure 2.

Invasive behavior in silico, in vivo, and in vitro. A, close-up of tumor border showing invasive behavior for PS = 0.03. CSCs (red) infiltrate surrounding tissue and spin-off DCCs that proliferate (dark blue). Small satellites are formed in the surrounding normal tissue (white) and grow back to and are engulfed by the main tumor mass. Nondividing cells (light blue). B, representative figures of various malignancies. All images reveal island-like formations at the rim of the main tumor mass (T). These findings are in line with tumor expansion as predicted by the SCA model that implements a CSC hierarchy. C, cell lines (n = 24) have been plated at clonal density in Matrigel. Simultaneously, the clonogenic fraction of the lines has been determined by limiting dilution analysis. A significant (P = 0.02) inverse relationship exists between clonogenicity and invasion as quantified by the measure of invasiveness we defined. Two examples of cell lines are shown and invasiveness (IM) and clonogenicity are indicated (right). See Supplementary Figs. S2 to S5 and Supplementary Table S1 for details.

Close modal

In an endeavor to validate these observations from our computational model, we investigated human tumor specimens. Close examination of the histology of different highly diverse human malignancies displayed a relatively confined large tumor mass with clearly detached tumor cells forming small lesions in the surrounding normal tissue (Fig. 2B). This exemplifies that human tumor histology contains indications of a stepwise infiltration and colonization of the surrounding tissue as our model predicts would follow from a hierarchical organized malignancy. Next, we attempted to determine the relationship between CSC fraction and invasive properties, as predicted by the model, using in vitro cell culture (Fig. 2C; Supplementary Figs. S2–S5 and Table S1). We determined the clonogenicity of a set of various cell lines (n = 24) as a surrogate for their CSC fraction. Additionally, we quantified the invasive properties of clonally derived structures with our measure of irregular morphology both on adherent plates and in Matrigel for all these lines. A significant inverse relationship between clonogenic fraction and the invasive properties of these lines exists. This implicates that tumor structures driven by a small fraction of clonogenic cells tend to generate a more irregular and invasive morphology, a finding that corroborates the predictions of our model. Combined, we take this as evidence that the SCA model based on the CSC concept closely resembles tumor growth patterns and morphology.

### Expansion of the Model in Three Dimensions Retains Invasive Morphology

To show that our results are not limited to two dimensions only, we expand the SCA model to three dimensions (see Supplementary Materials and Methods). The three-dimensional SCA model reproduces the invasive tumor morphology induced by the hierarchical organization, yielding three-dimensional fingering tumor fronts and clusters of cancer cells beyond the borders of the main tumor mass (Fig. 3A; Movie S3). On the contrary, again a classical model with equal volume does not display any apparent invasive pattern and instead exhibits a spherical morphology (Fig. 3B; Movie S4). We therefore propose that the CSC model, in contrast to the flat model, more faithfully reproduces the human tumor morphology for both the two-dimensional and the three-dimensional implementations.

Figure 3.

Expansion of the model in three-dimensions retains infiltrative morphology. Three-dimensional representation of CSC-driven tumor growth (A, PS = 0.03) and of the classical tumor model (B, PS = 1). A and B, 7 × 106 cells are represented. Inset, region that is enlarged in the right image. Normal tissue (black) and tumor cells (white). Cube represents 4003 cell lattice. See also Supplementary Movies S3 and S4.

Figure 3.

Expansion of the model in three-dimensions retains infiltrative morphology. Three-dimensional representation of CSC-driven tumor growth (A, PS = 0.03) and of the classical tumor model (B, PS = 1). A and B, 7 × 106 cells are represented. Inset, region that is enlarged in the right image. Normal tissue (black) and tumor cells (white). Cube represents 4003 cell lattice. See also Supplementary Movies S3 and S4.

Close modal

### CSC Organization Stimulates Tumor Heterogeneity

The dominant model evoked to explain the advancement of malignancies is clonal evolution, as was first proposed by Nowell (5). The term effective population size is used in population genetics to indicate the fraction of total individuals in a population that effectively contributes to the next generation and are therefore evolutionarily relevant. Hence, the effective population size in a CSC-driven malignancy is smaller than in the classical model because only mutations in the CSC compartment contribute to the evolutionary process (2). To implement clonal evolution in our model, we assume that, at each symmetrical division a CSC has a probability PMut to acquire a genetic hit and generate a daughter cell with a different phenotype selected from a randomly generated pool of 30 phenotypes (Table S2 and S3; see Supplementary Materials and Methods). Under an equal mutation rate (PMut = 0.1), the CSC model exhibits a slower acquisition of new phenotypes due to its smaller effective population size (Fig. 4A), but also shows a radically different selection process. Strikingly, although the rate of emergence of new phenotypes is much lower than observed in the classical model, a wide range of newly acquired phenotypes expands and contributes to the malignancy (Supplementary Fig. S6A and B).

Figure 4.

Tumor evolution and phenotypical selection in a cancer stem cell context. Phenotypes are randomly generated and possess traits listed in Supplementary Table S3. A, with equal mutation rates (PMut = 0.1), the cumulative amount of phenotypes that emerge is higher in the classical model (PS = 1) compared with the CSC model (PS = 0.03). B, we adapted PMut to obtain an equal rate of emergence of phenotypes. C, under equal conditions, the phenotypes in the CSC model are more diverse compared with the classical model. A representative example is shown, see also Supplementary Movies S5 and S6. A and B, bars, SD (n = 8). C, fraction of cells for each newly generated phenotype. The original phenotype “1” is ignored. For an average of eight experiments, see Supplementary Figs. S6B and S6C.

Figure 4.

Tumor evolution and phenotypical selection in a cancer stem cell context. Phenotypes are randomly generated and possess traits listed in Supplementary Table S3. A, with equal mutation rates (PMut = 0.1), the cumulative amount of phenotypes that emerge is higher in the classical model (PS = 1) compared with the CSC model (PS = 0.03). B, we adapted PMut to obtain an equal rate of emergence of phenotypes. C, under equal conditions, the phenotypes in the CSC model are more diverse compared with the classical model. A representative example is shown, see also Supplementary Movies S5 and S6. A and B, bars, SD (n = 8). C, fraction of cells for each newly generated phenotype. The original phenotype “1” is ignored. For an average of eight experiments, see Supplementary Figs. S6B and S6C.

Close modal

To compare the phenotypical selection more closely, we synchronize the pace at which new traits emerge. Evolutionarily, a CSC model with PS = 0.03 and PMut = 0.5 corresponds to a classical model with PS = 1 and PMut = 0.001 (Fig. 4B). Intriguingly, despite this adaptation of the mutation pace, the two models differ substantially in the way they exercise clonal selection. From a representative example in Fig. 4C, it is evident that the CSC model allows for a much higher phenotypical heterogeneity whereas the classical model seems to select for a small number of aggressive clones (average result in Supplementary Fig. S6C). Quantification of the heterogeneity clearly underscores this (Supplementary Fig. S7A and Supplementary Materials and Methods). It is important to realize that the selective pressure from the environment is equal in both models. This, combined with the equal rate of new phenotype occurrence, suggests that the geometric properties of the CSC model promote heterogeneity. We argue that the intrinsic properties of the CSC model might propel an alternative process to natural selection, referred to as genetic drift (25). In populations with small effective population sizes, sampling errors are frequent and might allow for the expansion of clones with no clear survival benefit. Interestingly, the invasive properties of the CSC model might fuel such a sampling error promoting mechanism. The phenomenon we observe at the tumor margins during invasion: a CSC migrating out of the tumor initiating a small satellite lesion, resembles the Founder Effect (26). The Founder Effect occurs when a new population is established by a low number of individuals from a larger population, this process is often accompanied by a loss of phenotypical variation in the new population due to sampling errors. In such a scenario, the new generation can differ substantially from the previous and is not in direct competition for nutrients and space. This implies that the pattern of tumor growth in the CSC model stimulates genetic drift at the infiltrative edges of the tumor and therefore promotes phenotypical heterogeneity.

### Dynamics of Therapeutic Interventions in Hierarchical Organized Malignancies

After having established the effects of CSC-driven tumor growth on the invasion and evolution of malignancies, we now explore the crucial topic of therapeutic intervention and tumor relapse. CSCs have been suggested to be more resistant to therapeutic interventions such as chemotherapy or irradiation compared with their differentiated counterparts (27, 28). Also tumors that relapse after seemingly successful therapy are believed to regrow from the CSC that survived the therapeutic regimen (27). Here, we investigate the dynamics associated with therapeutic interventions that are either selective for non-CSCs or equally efficient against both cell types. We find that the morphology and growth kinetics of relapses for both types of therapeutic interventions are very much different (Fig. 5; see Supplementary Materials and Methods). Regrowth after therapy that specifically targets DCCs is accompanied by enhanced invasive growth patterns whereas relapsing tumors after stochastic tumor cell killing are similar to the malignancy before treatment (Fig. 5A). Simultaneously, in case CSCs are resistant to therapy, the pace at which the malignancy relapses is greatly enhanced due to the relatively high fraction of CSCs directly following therapy (Fig. 5B). Also, the invasiveness of the recurrent tumors is markedly increased following intervention which is not effective against CSCs (Fig. 5C). These findings are in line with a range of clinical observations describing increased growth speed and enhanced invasion in the relapsing malignancy that are mostly attributed to the selection of more aggressive clones by the drug (2931). However, our data now indicate that these observations could be partially explained by the failure of conventional therapies to eradicate the CSC compartment and the subsequent relapse dynamics in CSC-driven tumors. Evaluation of evolutionary dynamics during relapse after both types of intervention revealed significant differences as well. Following therapy which is ineffective against CSCs, relapsing tumors display a marked increase in heterogeneity, whereas therapy that does target CSCs results in a dramatic decrease of heterogeneity (Supplementary Fig. S8). This latter scenario is related to the fact that relapses are very much different compared with the primary malignancy with respect to the clonal lineages that contribute to the relapse of the tumor (Fig. 5D). Combined, this implicates that applying therapy that is ineffectively targeting the CSC population is not only unsuccessful in curing the patient but also promotes malignant features including rapid expansion, increased invasion, and stimulates heterogeneity directly after therapy.

Figure 5.

Dynamics of therapeutic intervention in hierarchical organized malignancies. A, effects of therapy in simulated tumors considering CSCs to be resistant to treatment (top) or sensitive to treatment (bottom). Left, pretreated tumor at VT = 25,000 cells. Right, after treatment and relapsed tumor for VRelapse = VT. B, growth curves of tumor relapse. Whereas therapy-sensitive CSCs recapitulate the initial growth pace, resistant CSCs yield a faster growth at relapse. C, a marked increase of invasive morphology just after treatment occurs particularly in the case of resistant CSCs. B and C, arrows, therapy (n = 8). D, phenotype distribution in the two cases of therapy response assuming no further mutation to occur after therapy (see Supplementary Fig. S8 for an average of eight experiments). Although a treatment-resistant CSC-driven tumor is able to recapitulate the clonal features of the initial malignancy by conserving the initially developed heterogeneity, for therapy-sensitive CSCs, the recurrent tumor is highly altered in terms of phenotypes.

Figure 5.

Dynamics of therapeutic intervention in hierarchical organized malignancies. A, effects of therapy in simulated tumors considering CSCs to be resistant to treatment (top) or sensitive to treatment (bottom). Left, pretreated tumor at VT = 25,000 cells. Right, after treatment and relapsed tumor for VRelapse = VT. B, growth curves of tumor relapse. Whereas therapy-sensitive CSCs recapitulate the initial growth pace, resistant CSCs yield a faster growth at relapse. C, a marked increase of invasive morphology just after treatment occurs particularly in the case of resistant CSCs. B and C, arrows, therapy (n = 8). D, phenotype distribution in the two cases of therapy response assuming no further mutation to occur after therapy (see Supplementary Fig. S8 for an average of eight experiments). Although a treatment-resistant CSC-driven tumor is able to recapitulate the clonal features of the initial malignancy by conserving the initially developed heterogeneity, for therapy-sensitive CSCs, the recurrent tumor is highly altered in terms of phenotypes.

Close modal

We have applied a computational tumor growth model to investigate the effects of hierarchical organization on a range of areas in tumor biology. The obtained results provide novel insights into the influence of CSC-driven tumor growth on tumor morphology and invasion. Importantly, we have been able to partially validate these findings in a biological experimental setting. We have also established a significant effect of hierarchical organized malignant clones on tumor evolution, and consequently, on tumor progression and relapse after therapy. Our results stress the need to develop therapeutic interventions that efficiently target the CSC compartment because drugs that fail to do so are not only unable to cure the patient—but even enhance the malignant properties of tumors. Therefore, we conclude that computational oncology has the potential to greatly enhance our understanding of fundamental cancer biology and improve future therapeutics.

### CSC-Driven Tumor Growth

The capacity to display infiltrative growth is a pathognomic feature of malignant cells (1). Crossing tissue boundaries is the first step in the process of metastasis and is therefore of special interest. Surprisingly, this aspect of tumor growth has rarely been investigated in mathematical models. However, recent important work of Anderson and colleagues (11) and Bearer and colleagues (32) illuminate how selective pressure, applied by the microenvironment, in combination with tumor cell heterogeneity, could result in an invasive morphology. Here, we report how, in the SCA model, invasion emerges from the hierarchical organization of malignant clones and need not be driven by external interaction with a heterogeneous microenvironment. Hence, invasive morphology is an intrinsic feature of CSC-driven malignancies and has to be considered when studying tumor invasion.

Simultaneously, high levels of clonal heterogeneity within malignancies are associated with rapid progression of the disease (33), poorer survival (34), and occurrence of therapy resistance (35). Understanding this phenomenon in more detail might improve the design of more effective treatments and chemopreventives. Our results show how CSC-driven tumor growth stimulates heterogeneity in malignancies. This is related to the reduced effective population size in these malignancies and the typical growth pattern that results in the segregation of clones by formation of small invasive lesions. These satellite structures are not in direct competition with the dominant clone, allowing suboptimal clones to expand. We speculate that interfering with these dynamics, for example by inhibiting migration, and thereby segregation, might be a potent means by which reduced heterogeneity accompanied by the slower progression of the disease and better response to therapy can be achieved either in a chemopreventive or therapeutic setting.

Additionally, our results underlie the need to develop therapeutic modalities that successfully eradicate the CSC compartment. We show how therapeutic interventions that do not target CSCs efficiently are not only ineffective as a cure, but might even contribute to an increase in malignant properties of the relapsing tumor. We find, for example, enhanced heterogeneity and invasion in the relapsing tumors that are treated with interventions that are only directed against the differentiated cells. This is not due to the selection of the most malignant clones present in the cancer tissue but are intrinsic to the dynamics of regrowth in CSC-driven malignancies.

### Experimental Validation and Future Directions

The model presented provides a range of predictions and implications regarding CSC-fuelled tumor growth that can be tested and exploited to investigate the properties of hierarchical organized malignancies. Here, we report that the clonogenicity of cell lines is connected with the invasive properties of these lines in vitro as suggested by our model. In future research, it would be interesting to expand this finding to human tumor specimens or in established CSC lines in which the CSC fraction can be manipulated. Additionally, our model illuminates how tumor features, such as clonal heterogeneity, are related to the CSC fraction (Supplementary Fig. S7B). Therefore, we speculated that careful analysis of these tumor properties might lead to optimized techniques that establish the CSC fraction in malignancies without the need for transplantation studies which are currently highly disputed (36, 37).

The current version of the SCA model clearly shows the dynamics of CSC-driven tumor growth and its consequences for tumor morphology and evolution. However, further research efforts will undoubtedly lead to increased insight into the nature of the hierarchical organization of tumor cells. If so, the SCA model can be easily adapted to implement potential new information and subsequently come to even more accurate description of tumor growth. For instance, we have assumed several variables such as CSC self-renewal frequency, mutation rate, and migration speed to be constant. However, in in vivo tumors, these variables might be more dynamic both throughout tumor progression and as a consequence of microenvironmental interactions. For example, CSCs have been suggested to reside in a so-called “CSC niche” that is formed by secreted factors and direct cell-cell interactions (38, 39). This niche is believed to protect CSC from differentiation and to promote CSC self-renewal and therefore introduces dynamic PS values. Moreover, physical circumstances such as oxygen concentration might have differential effects on CSCs and DCCs, and therefore, might play a role in dynamic CSC functions (40). In addition, mutation rates are potentially subject to changes during tumor development (41) and intriguingly, might even be influenced by environmental conditions, such as that reported for bacteria (42). Currently, it remains elusive as to what extent this latter phenomenon, which is referred to as “adaptive mutation,” is involved in human malignancies (43). All these considerations justify further investigations, but will be dependent on thorough experimental examination of the variables involved.

Importantly, from the current formulation of the SCA model, it is indisputable that hierarchical organization of malignancies significantly contributes to the invasive morphology and increased heterogeneity of tumors and is therefore a crucial issue for better understanding tumor biology and to improve current anticancer treatments.

No potential conflicts of interest were disclosed.

We thank M.R. Sprick, D.J. Richel, and F. de Sousa Mello for useful comments and R. Belleman and J.H. de Jong for technical assistance.

Grant Support: Academic Medical Center (J.J.C. Verhoeff, J.P. Medema, and L. Vermeulen), ZonMW VICI program (J.P. Medema), and the Faculty of Science of the University of Amsterdam (A. Sottoriva, L. Naumov, and P.M.A. Sloot).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.

1
Hanahan
D
,
Weinberg
RA
.
The hallmarks of cancer
.
Cell
2000
;
100
:
57
70
.
2
Vermeulen
L
,
Sprick
MR
,
Kemper
K
,
Stassi
G
,
Medema
JP
.
Cancer stem cells—old concepts, new insights
.
Cell Death Differ
2008
;
15
:
947
58
.
3
Reya
T
,
Morrison
SJ
,
Clarke
MF
,
Weissman
IL
.
Stem cells, cancer, and cancer stem cells
.
Nature
2001
;
414
:
105
11
.
4
Sloot
PMA
,
Hoekstra
AG
.
Modeling dynamic systems with cellular automata
. In:
Fishwick
PA
, editor.
Handbook of dynamic system modeling
.
Chapman & Hall/CRC
;
2007
, p.
1
20
.
5
Nowell
PC
.
The clonal evolution of tumor cell populations
.
Science
1976
;
194
:
23
8
.
6
Vermeulen
L
,
Todaro
M
,
de Sousa Mello
F
, et al
.
Single-cell cloning of colon cancer stem cells reveals a multi-lineage differentiation capacity
.
Proc Natl Acad Sci U S A
2008
;
105
:
13427
32
.
7
Singh
SK
,
Hawkins
C
,
Clarke
ID
, et al
.
Identification of human brain tumour initiating cells
.
Nature
2004
;
432
:
396
401
.
8
Clarke
MF
,
Dick
JE
,
Dirks
PB
, et al
.
Cancer stem cells—perspectives on current status and future directions: AACR Workshop on Cancer Stem Cells
.
Cancer Res
2006
;
66
:
9339
44
.
9
Araujo
RP
,
McElwain
DL
.
A history of the study of solid tumour growth: the contribution of mathematical modelling
.
Bull Math Biol
2004
;
66
:
1039
91
.
10
Anderson
AR
,
Quaranta
V
.
Integrative mathematical oncology
.
Nat Rev Cancer
2008
;
8
:
227
34
.
11
Anderson
AR
,
Weaver
AM
,
Cummings
PT
,
Quaranta
V
.
Tumor morphology and phenotypic evolution driven by selective pressure from the microenvironment
.
Cell
2006
;
127
:
905
15
.
12
Sutherland
RM
.
Cell and environment interactions in tumor microregions: the multicell spheroid model
.
Science
1988
;
240
:
177
84
.
13
Folkman
J
,
Hochberg
M
.
Self-regulation of growth in three dimensions
.
J Exp Med
1973
;
138
:
745
53
.
14
Vaupel
P
,
Hockel
M
.
Blood supply, oxygenation status and metabolic micromilieu of breast cancers: characterization and therapeutic relevance
.
Int J Oncol
2000
;
17
:
869
79
.
15
Zaman
MH
,
Trapani
LM
,
Sieminski
AL
, et al
.
Migration of tumor cells in 3D matrices is governed by matrix stiffness along with cell-matrix adhesion and proteolysis
.
Proc Natl Acad Sci U S A
2006
;
103
:
10889
94
.
16
Anderson
AR
.
A hybrid mathematical model of solid tumour invasion: the importance of cell adhesion
.
Math Med Biol
2005
;
22
:
163
86
.
17
Brabletz
T
,
Jung
A
,
S
,
Hlubek
F
,
Kirchner
T
.
Opinion: migrating cancer stem cells—an integrated concept of malignant tumour progression
.
Nat Rev Cancer
2005
;
5
:
744
9
.
18
Hermann
PC
,
Huber
SL
,
Herrler
T
, et al
.
Distinct populations of cancer stem cells determine tumor growth and metastatic activity in human pancreatic cancer
.
Cell Stem Cell
2007
;
1
:
313
23
.
19
Laird
DJ
,
von Andrian
UH
,
Wagers
AJ
.
Stem cell trafficking in tissue development, growth, and disease
.
Cell
2008
;
132
:
612
30
.
20
Sheridan
C
,
Kishimoto
H
,
Fuchs
RK
, et al
.
CD44+/CD24− breast cancer cells exhibit enhanced invasive properties: an early step necessary for metastasis
.
Breast Cancer Res
2006
;
8
:
R59
.
21
Balic
M
,
Lin
H
,
Young
L
, et al
.
Most early disseminated cancer cells detected in bone marrow of breast cancer patients have a putative breast cancer stem cell phenotype
.
Clin Cancer Res
2006
;
12
:
5615
21
.
22
Dormann
S
,
Deutsch
A
.
Modeling of self-organized avascular tumor growth with a hybrid cellular automaton
.
In Silico Biol
2002
;
2
:
393
406
.
23
Jiang
Y
,
Pjesivac-Grbovic
J
,
Cantrell
C
,
Freyer
JP
.
A multiscale model for avascular tumor growth
.
Biophys J
2005
;
89
:
3884
94
.
24
Enderling
H
,
Hlatky
L
,
Hahnfeldt
P
.
Migration rules: tumours are conglomerates of self-metastases
.
Br J Cancer
2009
;
100
:
1917
25
.
25
Merlo
LM
,
Pepper
JW
,
Reid
BJ
,
Maley
CC
.
Cancer as an evolutionary and ecological process
.
Nat Rev Cancer
2006
;
6
:
924
35
.
26
Mayr
E
.
Change of genetic environment and evolution
. In:
Huxley
J
,
Hardy
AC
,
Ford
EB
, editors.
Evolution as a process
.
Princeton
:
Princeton University Press
;
1954
, p.
157
80
.
27
Jordan
CT
,
Guzman
ML
,
Noble
M
.
Cancer stem cells
.
N Engl J Med
2006
;
355
:
1253
61
.
28
Rich
JN
.
Cancer stem cells in radiation resistance
.
Cancer Res
2007
;
67
:
8980
4
.
29
Huff
CA
,
Matsui
W
,
Smith
BD
,
Jones
RJ
.
The paradox of response and survival in cancer therapeutics
.
Blood
2006
;
107
:
431
4
.
30
Spiegl-Kreinecker
S
,
Pirker
C
,
Marosi
C
, et al
.
Dynamics of chemosensitivity and chromosomal instability in recurrent glioblastoma
.
Br J Cancer
2007
;
96
:
960
9
.
31
El Sharouni
SY
,
Kal
HB
,
Battermann
JJ
.
Accelerated regrowth of non-small-cell lung tumours after induction chemotherapy
.
Br J Cancer
2003
;
89
:
2184
9
.
32
Bearer
EL
,
Lowengrub
JS
,
Frieboes
HB
, et al
.
Multiparameter computational modeling of tumor invasion
.
Cancer Res
2009
;
69
:
4493
501
.
33
Maley
CC
,
Galipeau
PC
,
Finley
JC
, et al
.
Genetic clonal diversity predicts progression to esophageal adenocarcinoma
.
Nat Genet
2006
;
38
:
468
73
.
34
Flyger
HL
,
Larsen
JK
,
Nielsen
HJ
,
Christensen
IJ
.
DNA ploidy in colorectal cancer, heterogeneity within and between tumors and relation to survival
.
Cytometry
1999
;
38
:
293
300
.
35
Iwasa
Y
,
Nowak
MA
,
Michor
F
.
Evolution of resistance during clonal expansion
.
Genetics
2006
;
172
:
2557
66
.
36
Hill
RP
.
Identifying cancer stem cells in solid tumors: case not proven [discussion 0]
.
Cancer Res
2006
;
66
:
1891
5
.
37
Quintana
E
,
Shackleton
M
,
Sabel
MS
,
Fullen
DR
,
Johnson
TM
,
Morrison
SJ
.
Efficient tumour formation by single human melanoma cells
.
Nature
2008
;
456
:
593
8
.
38
Gilbertson
RJ
,
Rich
JN
.
Making a tumour's bed: glioblastoma stem cells and the vascular niche
.
Nat Rev Cancer
2007
;
7
:
733
6
.
39
Borovski
T
,
Verhoeff
JJ
,
ten Cate
R
, et al
.
Tumor microvasculature supports proliferation and expansion of glioma-propagating cells
.
Int J Cancer
2009
;
125
:
1222
30
.
40
Li
Z
,
Bao
S
,
Wu
Q
, et al
.
Hypoxia-inducible factors regulate tumorigenic capacity of glioma stem cells
.
Cancer Cell
2009
;
15
:
501
13
.
41
Bielas
JH
,
Loeb
KR
,
Rubin
BP
,
True
LD
,
Loeb
LA
.
Human cancers express a mutator phenotype
.
Proc Natl Acad Sci U S A
2006
;
103
:
18238
42
.
42
Hastings
PJ
,
Bull
HJ
,
Klump
JR
,
Rosenberg
SM
.
Adaptive amplification: an inducible chromosomal instability mechanism
.
Cell
2000
;
103
:
723
31
.
43
Rosenberg
SM
.
.
Nat Rev Genet
2001
;
2
:
504
15
.
44
Tomita
K
,
Plager
JE
.
In vivo cell cycle synchronization of the murine sarcoma 180 tumor following alternating periods of hydroxyurea blockade and release
.
Cancer Res
1979
;
39
:
4407
11
.
45
Sherwood
L
.
Human physiology: from cells to systems
.
Belmont
:
;
1997
.
46
Casciari
JJ
,
Sotirchos
SV
,
Sutherland
RM
.
Variations in tumor cell growth rates and metabolism with oxygen concentration, glucose concentration, and extracellular pH
.
J Cell Physiol
1992
;
151
:
386
94
.
47
Freyer
JP
,
Tustanoff
E
,
Franko
AJ
,
Sutherland
RM
.
In situ oxygen consumption rates of cells in V-79 multicellular spheroids during growth
.
J Cell Physiol
1984
;
118
:
53
61
.
48
Bray
D
.
Cell movements
.
New York
:
Garland Pub.
;
1992
.