Cancer in a host induces responses that increase the ability of the microenvironment to sustain the growing mass, for example, angiogenesis, but cancer cells can have varying sensitivities to these sustainability signals. Here, we show that these sensitivities are significant determinants of ultimate tumor fate, especially in response to treatments and immune interactions. We present a mathematical model of cancer–immune interactions that modifies generalized logistic growth with both immune-predation and immune-recruitment. The role of a growing environmental carrying capacity is discussed as a possible regulatory mechanism for tumor growth, and this regulation is shown to modify cancer–immune interactions and the possibility of achieving immune-induced tumor dormancy. This mathematical model qualitatively matches experimental observations of immune-induced tumor dormancy as it predicts dormancy as a transient period of growth that necessarily ends in either tumor elimination or tumor escape. As dormant tumors may exist asymptomatically and may be easier to treat with conventional therapy, an understanding of the mechanisms behind tumor dormancy may lead to new treatments aimed at prolonging the dormant state or converting an aggressive cancer to the dormant state. Cancer Res; 73(12); 3534–44. ©2013 AACR.

Major Findings

We show, using a mathematical model, how the sensitivity of tumor cells to immune-mediated environmental signals can significantly alter tumor dynamics and, thus, treatment outcomes. Moreover, immune-induced tumor dormancy is predicted to be a transient period of tumor growth that must necessarily end in either tumor elimination or tumor escape, in agreement with several experimental observations.

Quick Guide to Equations and Assumptions

As we focus on cancer dormancy induced by immune predation, the cancer–immune interactions considered are immune predation of, and immune recruitment by, cancer cells. The cancer, C(t), and immune, I(t), populations are assumed to exhibit generalized logistic growth. Without immune predation, C(t) will grow until it reaches the carrying capacity KC. We consider the carrying capacity to be dictated by the growth regulatory signal emitted by the environment and sensed by the cells. The immune compartment I(t) impedes tumor growth through cytotoxic actions targeted at the cancer cells. These actions are contained in the predation term Ψ(I,C), which acts to modulate the growth rate of C(t). The equation governing cancer growth is thus,

formula

Here, μ and α are parameters that capture the growth rate and sensitivity of cancer cells to environmental regulatory signals.

The immune population is assumed to consist of cytotoxic cells, for example, CD8+ T cells, NK cells, and macrophages that are known to directly target cancer cells. The immune response is assumed to grow toward a theoretical limit, or carrying capacity KI through proliferation and recruitment from the blood, spleen, and bone marrow. With parameter r, we allow for both immune cell– and cancer cell–mediated recruitment; the latter occurs through, for example, the production of danger-associated molecules and necrosis (1). The equation governing immune growth is thus,

formula

where λ is the immune growth-rate parameter.

Immune predation of cancer cells is assumed to contain both fast and slow dynamics, expressed as by definition. The fast saturation behavior observed in some cancer–immune lysis assays (2, 3) is described by the first, ratio-dependent term (4). Although this form has been suggested to better capture the cytotoxic effects of T cells, this term alone neglects innate immunity, which should not exhibit saturation (compare, e.g., the percent lysis curves for cytotoxic T cells vs. NK cells in refs. 2, 3, 5). The second term accounts for innate immunity, and allows for a slow increase in the saturation limit with increasing immune presence, as seen when the lysis assays are conducted at larger ratios (3, 5, 6). The logarithm transforms the cytotoxic actions of the innate immunity into a range appropriate for Ψ. For small immune populations, we assume that the innate and adaptive effects can be combined into the ratio-dependent term, but for large populations, innate immunity should still have a small cytotoxic effect. Thus, we assume .

Tumor growth is defined by interactions with the local environment. Immune cells can both promote tumor growth through inflammatory mechanisms and inhibit growth through direct cytotoxic effects (7). Other stromal cells contribute to a tumor-promoting environment by increasing the local vasculature and by forming a supportive niche for cancer cell proliferation and invasion (8). The intercellular signaling between cancer and stromal cells can trigger prosurvival pathways in both cancer cells and host stroma, including endothelial cells, fibroblasts, and immune cells (8).

Inherently, there is variability in the intercellular signaling, which may modulate macroscopic measurements such as tumor diameter. This variability may arise from many factors, including cell sensitivities to intercellular signals and the rates at which cells respond to those signals. One example of environmental regulation is shown by dormant and fast-growing tumor clones that exhibit similar growth kinetics in vitro but display strikingly different growth rates in vivo (9). Here, we show that intercellular signaling may act to regulate tumor growth. We focus on the sensitivity of cancer cells to growth regulatory signals from the tumor microenvironment and the altered tumor dynamics resulting from treatment-induced disruption of these signals. We show that the variable sensitivities of cancer cells to stromal intercellular signaling may fundamentally control tumor dynamics, even to the point of inducing an immune-induced dormant state, with clear implications for therapeutic efficacies.

Immune-induced cancer dormancy is a state of cancer progression in which the cancer is maintained in a viable, but nonexpanding, state (10), often described as an “equilibrium” phase in immunoediting nomenclature (11). Although this state may persist for days to decades, its immunologic realization is one of transience, i.e., the eventual elimination of the disease or the development of immune-resistance followed, some time later, by tumor escape.

In contrast, mathematical models typically describe the dormant state by a stable equilibrium point (or limit cycle) with a basin of attraction (12–17). This implies that the dormant state can attract tumor trajectories and maintain itself for long times. Such analyses neglect the transient nature of the dormant state, however, and require external perturbations to the system to explain the eventual escape from dormancy. Recent mathematical explorations of possible escape mechanisms include random fluctuations in immune presence (18), intercellular communication of learned cancer cell resistance (19), and immunoediting or evolution in cancer cell phenotypes (20).

Here, without considering specific escape mechanisms, we present a formalism that contains one long-term dormancy-associated equilibrium and predicts tumor dormancy that will generally end in either tumor elimination or tumor escape. The equilibrium point is a saddle node with a separatrix that divides two attractor regions of ultimate tumor fate. We show that the duration of dormancy is determined by tumor–immune dynamics and the proximity of the tumor trajectory to the separatrix.

This model is simple enough to analytically investigate, yet complex enough to capture all qualitative behaviors of tumor growth, including tumor dormancy. Using parameter sets estimated by a Markov Chain Monte Carlo (MCMC) algorithm, we show that the sensitivity of cancer cells to environmental signals is a prominent factor in determining tumor fate. We generate four parameter sets: one assuming a constant environmental signal, producing a static carrying capacity, and three assuming a variable environmental signal, giving rise to a dynamic carrying capacity. These four sets fit the experimental tumor growth data equally well, but when the variable signals are disrupted, the differing cancer cell sensitivities predict different tumor growth fates. Interestingly, these four sets predict drastically different results for the same treatment, ranging from rapid tumor elimination to tumor escape.

These results may explain both the high variability of treatment success among patients as well as the observation of accelerated repopulation after treatment (21, 22), a phenomenon arising here from disruption of growth regulatory signals in the microenvironment. These same signals are shown to ultimately determine whether dormancy will be possible, with higher tumor sensitivities to these signals corresponding to decreased likelihood of dormancy when these signals are disrupted, as may happen from therapy or chronic microenvironmental perturbations such as low-dose radiation exposure (23).

Experimental data

Basic tumor growth dynamics in the absence of an anticancer immune response are estimated from experimental measurements of a subcutaneous fibrosarcoma induced by 3-methylcholanthrene. Tanooka and colleagues (24) induced sarcomas in wild-type WBB6F1/J and C3H/He mice. Volume measurements were taken when the resulting tumors reached a palpable size. We note that WBB6F1/J mice are mast cell deficient, which inhibits the host's inflammatory response.

Estimation of model parameters

Parameters are estimated using a MCMC algorithm (25, 26). From an initial estimation, a Markov chain of permitted parameter sets is created by randomly perturbing the previous parameter set and accepting this perturbed set with a probability determined by a measure of goodness of fit. Here, we use the sum of squared deviations between model prediction and experimental measurements. Each parameter is perturbed and tested for acceptance independently, except parameters μ and α, which are handled together because they are inherently related in Eq. A. The algorithm is repeated 10 times, with each run having 20,000 iterations. Final parameters are those that provide the best fits.

Parameters are estimated under the assumption of no immune predation |$(\Psi = 0)$|⁠. The results assuming a constant or dynamic carrying capacity are listed in Tables 1 and 2, respectively. Immune size and growth rate values are selected to approximately match those of the tumor as methods to measure immune system checkpoint blockades and other size-determining factors remain elusive.

Table 1.

Parameter values for model simulations

ParameterValueInterpretation
μ 0.16 days−1 Cancer growth 
α 0.72 Cancer growth 
λ 0.22 days−1 Immune growth 
r 1·10−3 Immune recruitment 
θ 2.5 Predation strength 
β 0.5 Predation saturation 
φ 50 Predation saturation 
ϵ 0.01 Innate immunity predation strength 
K_C 3.92 × 1010 Cell num Cancer-carrying capacity 
K_I 3.92 × 1010 Cell num Immune-carrying capacity 
ParameterValueInterpretation
μ 0.16 days−1 Cancer growth 
α 0.72 Cancer growth 
λ 0.22 days−1 Immune growth 
r 1·10−3 Immune recruitment 
θ 2.5 Predation strength 
β 0.5 Predation saturation 
φ 50 Predation saturation 
ϵ 0.01 Innate immunity predation strength 
K_C 3.92 × 1010 Cell num Cancer-carrying capacity 
K_I 3.92 × 1010 Cell num Immune-carrying capacity 

NOTE: C(t) and I(t) measure cell number, but their values are plotted in terms of cellular mass diameter in millimeters according to the formula |${\rm diameter(X)} = 2\left({\frac{X}{{10^6}}}{\frac{3}{{4\pi}}} \right)^{1/3}$| mm.

Table 2.

Tumor growth parameter sets estimated by the MCMC parameter fitting method assuming a dynamic carrying capacity

Growth rate and sensitivityDynamic carrying capacityConstant carrying capacity
 μ α p q K_C(0) |$K_C = \left({\frac{p}{q}} \right)^{3/2}$| 
Set 1 0.28 1.08 0.90 7.97 × 10−8 1.12 × 108 3.79 × 1010 
Set 2 0.23 1.03 34.5 3.14 × 10−6 8.47 × 105 3.64 × 1010 
Set 3 0.19 1.44 × 10−5 0.73 5.98 × 10−8 2.90 × 107 4.27 × 1010 
Growth rate and sensitivityDynamic carrying capacityConstant carrying capacity
 μ α p q K_C(0) |$K_C = \left({\frac{p}{q}} \right)^{3/2}$| 
Set 1 0.28 1.08 0.90 7.97 × 10−8 1.12 × 108 3.79 × 1010 
Set 2 0.23 1.03 34.5 3.14 × 10−6 8.47 × 105 3.64 × 1010 
Set 3 0.19 1.44 × 10−5 0.73 5.98 × 10−8 2.90 × 107 4.27 × 1010 

NOTE: The last column gives the maximum value of K_C⁠, a value used when a constant carrying capacity is forced on the system.

Mathematical analysis

The nonlinear system given by Eqs. A, B, and C describes a 2-compartment model for cancer–immune interactions. The general form |$\Ω = 1 + \Psi$| allows for immune modulation of tumor growth, which may be stimulatory or inhibitory, depending on Ω. It should be noted that this form fundamentally differs from classic predator–prey-type models where immune effectors only inhibit tumor growth and are stimulated directly by tumor presence (4, 13, 17, 27, 28). In contrast to decaying oscillations or limit cycle behaviors, our model predicts tumor elimination or escape, with only transient periods of dormancy. This is partially due to the assumed form of the immune response, which, once initiated, will progress to a maximal response that is independent of the cancer population.

Stability analysis of the system equilibrium points and possible cancer–immune dynamics.

The dynamical system described by Eqs. A and B is simple enough to allow for analytic investigation of system stability. The system has four equilibrium points: |$(I,C) = (0,0)$|⁠, |$(K_I,0)$|⁠, |$(K_I,K_C)$|⁠, and |$(K_I,\bar C)$| where |$\bar C$| satisfying |$1 + \Psi (K_I,\bar C) = 0$| is given by

formula

and represents a cancer equilibrium state that, once obtained, will persist for a long time. We note that Ψ(I, C) is not differentiable at (0, 0), but that we can classify this point based on physical arguments and a topologically equivalent system. The equilibria stability, based on the Jacobian, are as follows: (0, 0) is an unstable point, (KI, 0) is a stable point if |$1 + \Psi (K_I,\,0) < 0$|⁠, otherwise it is a saddle point, and |$(K_I,\,K_C)$| is a stable point if |$1 + \Psi (K_I,K_C)\,> \,\,0$|⁠, otherwise it is a saddle point. The point |$(K_I,\bar C)$| with |$0 < \bar C < K_C$| is a saddle point if (KI,0) and (KI, KC) are both stable. The two restrictions guaranteeing the attractor states bound immune predation strength as follows, ensuring a dormancy-associated equilibrium point:

formula

Under this condition, a large immune response is strong enough to eliminate small tumors but not large tumors, resulting in a possible dormant tumor of size |$\bar C$|⁠. It guarantees the existence of all possible tumor outcomes: tumor elimination, tumor escape, and long-term tumor dormancy followed by elimination or escape. Three sample cancer–immune phase portraits are shown in Fig. 1. The long-term dormant state increases as immune efficacy increases from θmin to θmax, and trajectories that pass near the separatrix have transient dormant states.

Figure 1.

The progression of the immune–cancer phase plane as immune efficacy θ increases. A, if θ < θmin, immune predation is weak and all cancers reach the carrying capacity K_C⁠. B, if |$\theta _{{\rm min}} < \theta < \theta _{{\rm max}}$|⁠, a dormant cancer state exists as a saddle point between two stable equilibria creating two basins of attraction, one representing cancer elimination (bottom, blue region) and one representing cancer escape (top, red region). The separatrix is the black curve separating the two basins of attraction, and leads to the dormancy-associated equilibrium point |$(K_I,\bar C)$|⁠. C, finally, if |$\theta \,>\, \theta _{{\rm max}}$|⁠, all cancers are eventually eliminated. In D, a bifurcation diagram shows the stability of the cancer equilibrium points as immune efficacy θ increases. In E, the sensitivity of the size of the dormant state, |$\bar C$|⁠, to changes of 10% in the immune model parameters K_I⁠, θ, φ, and β, is shown. Analyses were completed with parameter values from Table 1 unless otherwise stated.

Figure 1.

The progression of the immune–cancer phase plane as immune efficacy θ increases. A, if θ < θmin, immune predation is weak and all cancers reach the carrying capacity K_C⁠. B, if |$\theta _{{\rm min}} < \theta < \theta _{{\rm max}}$|⁠, a dormant cancer state exists as a saddle point between two stable equilibria creating two basins of attraction, one representing cancer elimination (bottom, blue region) and one representing cancer escape (top, red region). The separatrix is the black curve separating the two basins of attraction, and leads to the dormancy-associated equilibrium point |$(K_I,\bar C)$|⁠. C, finally, if |$\theta \,>\, \theta _{{\rm max}}$|⁠, all cancers are eventually eliminated. In D, a bifurcation diagram shows the stability of the cancer equilibrium points as immune efficacy θ increases. In E, the sensitivity of the size of the dormant state, |$\bar C$|⁠, to changes of 10% in the immune model parameters K_I⁠, θ, φ, and β, is shown. Analyses were completed with parameter values from Table 1 unless otherwise stated.

Close modal

Bifurcation and sensitivity analysis of tumor dormancy.

If |$\theta < \theta _{{\rm min}}$| (or |$1 + \Psi (K_I,0) \,> \,0$|⁠), then the immune response is weak and |$(K_I,0)$| is a saddle point with all tumors escaping immune control (Fig. 1A). In this case, |$\bar C$| is negative and, thus, unphysical. If |$\theta \,> \,\theta _{{\rm max}}$| (or |$1 + \Psi (K_I, K_C) < 0$|⁠), then the immune response is strong and |$(K_I,K_C)$| is a saddle point, implying that all tumors are eliminated (Fig. 1C). In this case, |$\bar C \,> \,K_C$|⁠, and because we only consider tumor growth in the positive quadrant bounded by the carrying capacity, we neglect this dormant equilibrium point. Furthermore, if |$\theta = \theta _{{\rm min}}$| (or |$1 + \Psi (K_I,0) = 0$|⁠), then |$\bar C = 0$|⁠, and if |$\theta = \theta _{{\rm max}}$| (or |$1 + \Psi (K_I,K_C) = 0$|⁠), then |$\bar C = K_C$|⁠. Thus, θ is a bifurcation parameter with two transcritical bifurcations. A bifurcation diagram showing the transitions and stability of the equilibrium points is shown in Fig. 1D. Because the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| guarantees the existence of all possible tumor outcomes, we focus on it for the remainder of the discussion.

Ranges of immune efficacy guaranteeing tumor dormancy.

Figure 2 shows the dependence of the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| on parameters KI, KC, β, and φ. As the immune-carrying capacity increases, the upper bound on immune efficacy decreases because the same amount of cancer cell predation can be mediated by either a large and weakly effective immune presence or a small but strongly effective immune presence. Similarly, as the cancer-carrying capacity increases, the upper bound on immune efficacy increases, as stronger immune predation is required to eliminate the larger tumor. The effects of the immune predation shape parameters, (β,φ), on the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| depends on the relative size of the two carrying capacities. Notice that if |$K_C = K_I$|⁠, then β has no effect on |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$|⁠. If this is not the case, however, then if |$K_I < K_C$|⁠, increasing β increases θmax, and if |$K_I \,> \,K_C$|⁠, increasing β decreases θmax. This behavior results from the shape of the fast saturation kinetic in the cell-kill term, which is a function of the ratio |${K_I}\slash{K_C}$| [see Eq. C]. Parameter φ also controls the shape of the fast saturation kinetic. If |$\phi = 0$|⁠, then |$\theta _{{\rm min}} = \theta _{{\rm max}}$|⁠, and as φ increases, the upper bound also increases. Notice that for all these cases, the lower bound, θmin, is essentially constant, suggesting that there exists a strictly positive threshold for immune predation efficacy, below which long-term tumor dormancy cannot exist.

Figure 2.

The dependence of the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| on model parameters K_I⁠, K_C⁠, β, and φ. In all figures, the 3 shaded areas represent regions describing distinct tumor outcomes: all tumors are eliminated in the (light blue) upper regions; tumor elimination, tumor escape, and long-term tumor dormancy are all possible in the (green) middle regions; and all tumors escape in the (pink) lower regions. The range of θ permitting a long-term dormant cancer state decreases as immune carrying capacity increases (A) and increases as cancer carrying capacity increases (B), for fixed (β, θ) in the predation term, C–H, the effects of changes in β and φ on the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| depend on both K_I and K_C⁠. Unless otherwise stated, parameter values are from Table 1.

Figure 2.

The dependence of the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| on model parameters K_I⁠, K_C⁠, β, and φ. In all figures, the 3 shaded areas represent regions describing distinct tumor outcomes: all tumors are eliminated in the (light blue) upper regions; tumor elimination, tumor escape, and long-term tumor dormancy are all possible in the (green) middle regions; and all tumors escape in the (pink) lower regions. The range of θ permitting a long-term dormant cancer state decreases as immune carrying capacity increases (A) and increases as cancer carrying capacity increases (B), for fixed (β, θ) in the predation term, C–H, the effects of changes in β and φ on the range |$[\theta _{{\rm min}},\theta _{{\rm max}} ]$| depend on both K_I and K_C⁠. Unless otherwise stated, parameter values are from Table 1.

Close modal

Numerical simulations of tumor growth are presented with and without immune predation. Simulations are computed numerically in MAPLE (www.maplesoft.com) using the parameter values in Tables 1 and 2 as appropriate.

Tumor growth is regulated by signals from the microenvironment

We investigate the effect of environmental signaling on tumor dormancy. When a tumor grows, it modifies the local environment, through angiogenesis, to increase the nutrient supply to that area. In this way, the carrying capacity increases with tumor mass. This effect can be modeled by a differential equation that describes the growth of the cancer-carrying capacity, KC(t), in response to the growing tumor, C(t). The functional form, as first proposed by Hahnfeldt and colleagues (29), is

formula

where p and q are growth stimulation and inhibition constants, respectively.

Replacing the constant carrying capacity in Eq. A with a dynamic capacity as described earlier, creates an extension of the mathematical model given by Eqs. A, B, and C. By comparing these two models, we show that disturbance of the support signals from the tumor environment can significantly alter tumor growth dynamics. This formulation attempts to capture some of the effects stromal cells (including immune cells) can have on the proliferation rate, viability, differentiation state, and invasiveness of cancer cells (30).

Parameter estimation for tumor growth in the absence of immune predation is conducted with the experimental data of Tanooka and colleagues (24) using an MCMC parameter estimation algorithm. Parameter sets are estimated (i) for the tumor growth model assuming a constant carrying capacity given by Eq. A; and (ii) the tumor growth model assuming a dynamic carrying capacity given by Eqs. A and F. Table 1 lists the estimated parameter values used for the constant carrying capacity model, referred to as set 0 (μ, α, and KC from Table 1), and Table 2 lists 3 possible parameter sets for the dynamic carrying capacity model, referred to as sets 1, 2, and 3.

As can be seen from Fig. 3, both models and all parameter sets fit the experimental data equally well. When microenvironmental regulation of tumor growth is disturbed, however, sets 1 and 3 no longer fit the data. Microenvironmental signaling is disrupted in the model by forcing the dynamic carrying capacity to be constant and equal to the maximum possible value for each parameter set (these values are listed as KC in Table 2). Doing so removes early damping of tumor growth, providing a strongly protumor environment. Its effect is different for the three parameter sets. Parameter set 2 still fits the data well, but sets 1 and 3 predict much faster tumor growth. More specifically, the rate at which the environment increases either its vascularization, the sensitivity of the stromal cells to the proangiogenic or other growth factors, or the sensitivity of the tumor cells to the signaled levels of environmental support, all contribute to tumor-specific inherent variabilities that determine basic tumor growth dynamics.

Figure 3.

Tumor microenvironmental growth signals modulate tumor growth. A, tumor growth is predicted by fitting growth parameters to the experimental data of Tanooka and colleagues (24). Set 0 (Table 1) assumes a constant carrying capacity, whereas sets 1, 2, and 3 (Table 2) allow for a dynamic carrying capacity. All curves fit the experimental data equally well. B, the same parameter sets are used, but a constant carrying capacity is forced on sets 1, 2, and 3 after the parameter fitting using |$K_C = \left({{\textstyle{p \over q}}} \right)^{{\textstyle{3 \over 2}}}$| (values from Table 2). Note that set 2 still fits the data well, whereas sets 1 and 3 do not. This indicates that the tumor described by set 2 is not sensitive to the local environmental growth cues, whereas sets 1 and 3 are very sensitive. Cancer–immune phase portraits for sets 2 and 3 show the same general behavior when the dynamic carrying capacity is allowed (C and D) but show significantly different behavior when a constant carrying capacity is enforced (E and F) simulating a disrupted environmental signal. After disruption, set 2 still allows immune-induced transient periods of tumor dormancy, whereas set 3 no longer predicts any significant dormant periods.

Figure 3.

Tumor microenvironmental growth signals modulate tumor growth. A, tumor growth is predicted by fitting growth parameters to the experimental data of Tanooka and colleagues (24). Set 0 (Table 1) assumes a constant carrying capacity, whereas sets 1, 2, and 3 (Table 2) allow for a dynamic carrying capacity. All curves fit the experimental data equally well. B, the same parameter sets are used, but a constant carrying capacity is forced on sets 1, 2, and 3 after the parameter fitting using |$K_C = \left({{\textstyle{p \over q}}} \right)^{{\textstyle{3 \over 2}}}$| (values from Table 2). Note that set 2 still fits the data well, whereas sets 1 and 3 do not. This indicates that the tumor described by set 2 is not sensitive to the local environmental growth cues, whereas sets 1 and 3 are very sensitive. Cancer–immune phase portraits for sets 2 and 3 show the same general behavior when the dynamic carrying capacity is allowed (C and D) but show significantly different behavior when a constant carrying capacity is enforced (E and F) simulating a disrupted environmental signal. After disruption, set 2 still allows immune-induced transient periods of tumor dormancy, whereas set 3 no longer predicts any significant dormant periods.

Close modal

These results highlight the important fact that the same tumor growth curve can be observed for three different tumors, each with different cellular proliferation rates and different sensitivities to the regulatory signals originating from their different environments. At the same time, when these signals are altered, the resulting changes in growth dynamics can be drastic.

Altering the balance of the tumor and microenvironmental compartments, for example, by selective reduction of tumor burden by therapy, leaving the microenvironmental support nearly intact (or enhanced; ref. 23), may accelerate tumor repopulation, a problematic phenomenon familiar to clinicians that remains poorly understood (21, 22, 31). In this regard, Fig. 3 shows how the rate of tumor regrowth depends on the signaling sensitivities of the tumor cells.

Another possible interpretation of our disrupted regulatory signals comes from Hu and colleagues (32). They found that when tumor-associated and arthritis-associated (inflammatory) fibroblasts were coinjected with breast cancer cells, tumor weight was increased. This contribution may explain our simulations for how enhanced capacity leads to accelerated growth in sets 1 and 3. Typically, these stromal cells would develop their protumor associations over time as the tumor grows (modeled here by the dynamic capacity), but in these in vivo experiments, the stromal cells are already protumorigenic and are thus capable of supporting a tumor faster (modeled here by the forced-constant capacity).

In addition, the ability of the immune response to induce tumor dormancy may be lost when environmental regulatory signals are disrupted. Using parameter set 2, simulations of tumor growth with disrupted signals result in a growth curve similar to the data, and the resulting tumor–immune phase portraits, shown in Fig. 3C and E, are similar for both the dynamic and the (disrupted) forced-constant carrying capacity models. This suggests that the tumor represented by parameter set 2 is not very sensitive to regulation signals produced by the microenvironment and that the cancer cells proliferate at a rate consistent with the growth curves. This allows the immune response to achieve transient dormant states in both the dynamic and (disrupted) forced-constant capacity models.

In contrast, the cancer–immune phase portraits for set 3 in Fig. 3D and F show a lost ability of the immune response to induce tumor dormancy, due to the tumor being very sensitive to growth regulation signals emanating from the local microenvironment. In addition, disrupting the regulatory signals for set 3 causes the tumor to grow much faster than the experimental data. These observations suggest that immune-mediated tumor dormancy is dependent on the sensitivity of cancer cells to environmental regulatory signals.

Success of immunotherapy depends on tumor sensitivity to environmental regulatory signals

To elucidate the role of variability in signaling sensitivity to treatment success, we now consider immunotherapies that increase the immune effectiveness of the cytotoxic T cells, achieved conceivably, for instance, by checkpoint blockade-targeted antibodies to cytotoxic T-lymphocyte antigen 4 (CTLA-4) and programmed cell death protein 1 (PD-1). These therapies release the immunologic block that limits the immune response, allowing the development of a large specifically targeted T-cell response to the cancer. We simulate the effects of these treatments by allowing a large immune response to develop with an increased recruitment ability (r) and an increased cytotoxic efficacy (θ). Treatment occurs after the tumor is initially detected (at |$t = 0$|⁠) with no appreciable initial immune presence [|$I(0) = 0$|]. In Fig. 4A and B, the experimental data for tumor growth without predation is included as a reference to the initial growth trend.

Figure 4.

Immunotherapy treatment is simulated by an increased immune recruitment rate, r, and an increased predation efficacy, θ. Treatment is given at the time of tumor detection (here, |$t = 0$|⁠) and is assumed to be immediately effective with no decay. Any immune presence at detection is ignored [|$I(0) = 0$|]. A and B, the data points show the growth trend with no immune interaction (as in Fig. 3). The dynamic carrying capacity model with parameter sets 1, 2, and 3 are shown in (A) along with the constant capacity model with parameter set 0. In B, the dynamic capacity models are forced to have constant capacities, disrupting the regulatory signal to the tumor. All curves are modified to different extents by immune predation when compared with the data (no predation). C, tumor fate is indicated for increasing initial immune presence, I(0)⁠. From left to right, the pink regions represent tumor escape with no dormancy, the red regions represent tumor escape after dormancy, the blue regions represent tumor elimination after dormancy, and the light blue regions represent tumor elimination with no dormancy. Dormancy is defined as a period of non-growth lasting at least 30 days. The width of each dormancy range (red + blue region) is labeled to the right of each range. Parameter sets 1, 2, and 3 are used with the dynamic capacity (D) and constant capacity (C) models. Parameter values are |$r = 0.47$| and |$\theta = 4$|⁠, together with the values from Tables 1 and 2.

Figure 4.

Immunotherapy treatment is simulated by an increased immune recruitment rate, r, and an increased predation efficacy, θ. Treatment is given at the time of tumor detection (here, |$t = 0$|⁠) and is assumed to be immediately effective with no decay. Any immune presence at detection is ignored [|$I(0) = 0$|]. A and B, the data points show the growth trend with no immune interaction (as in Fig. 3). The dynamic carrying capacity model with parameter sets 1, 2, and 3 are shown in (A) along with the constant capacity model with parameter set 0. In B, the dynamic capacity models are forced to have constant capacities, disrupting the regulatory signal to the tumor. All curves are modified to different extents by immune predation when compared with the data (no predation). C, tumor fate is indicated for increasing initial immune presence, I(0)⁠. From left to right, the pink regions represent tumor escape with no dormancy, the red regions represent tumor escape after dormancy, the blue regions represent tumor elimination after dormancy, and the light blue regions represent tumor elimination with no dormancy. Dormancy is defined as a period of non-growth lasting at least 30 days. The width of each dormancy range (red + blue region) is labeled to the right of each range. Parameter sets 1, 2, and 3 are used with the dynamic capacity (D) and constant capacity (C) models. Parameter values are |$r = 0.47$| and |$\theta = 4$|⁠, together with the values from Tables 1 and 2.

Close modal

Including strong immune predation in the tumor growth simulations of Fig. 4 results in the parameter sets all predicting different outcomes for the same treatment, mimicking real-world interpatient variability. A striking comparison can be seen between the case with no immune predation, where all curves fit the data equally well (Fig. 3A), and the case with strong immune predation, where each curve predicts a different outcome (Fig. 4A). Note also that while set 2 indicates elimination happens faster than in the constant capacity model (set 0), set 1 indicates a short dormant period before escape, and set 3 indicates escape at only a slightly slower rate than the experimental reference. Thus, with environmental growth regulation signals intact, these simulations suggest that the tumor may escape, be eliminated, or be transiently dormant before elimination or escape. Such results may explain why patient responses to immunotherapies are quite varied, and shed light on the connection between regulatory sensitivities and treatment response.

With the loss of environmental growth regulation, however, all tumors are predicted to escape (Fig. 4B). Disrupting the regulatory signals by forcing the dynamic capacity to be constant results in tumor escape: (i) with a short dormant period (set 2); (ii) at a slightly slower rate than the reference (set 1); or (iii) at a faster rate than the reference (set 3). Again, forcing the dynamic carrying capacity to be constant simulates a reduction of tumor burden with little modification of environmental support. When regulation is disrupted, only tumor escape is predicted, implying that a stronger dose or combinations of therapies may be required to eliminate the tumor upon recurrence.

Dormancy depends on tumor cell-inherent sensitivities and parameters of the immune response

As discussed above, immune-mediated tumor dormancy is dependent not only on immune effectiveness but also on the sensitivity of cancer cells to environmental regulatory signals. To show this, we examine the initial immune presence required to induce a period of dormancy in our simulations. We define dormancy as a period of at least 30 days within the first 150 days of growth, wherein the effective tumor growth rate is small, i.e., |$\left({\frac{1}{C}}{\frac{{{\rm d}C}}{{{\rm d}t}}} \right) < 0.04$|⁠, and the tumor size is moderate, i.e., between 1 and 30 mm in diameter. This definition excludes dormant regions where the tumor has effectively been eliminated or where the tumor has grown to the maximal size.

In Fig. 4C, predicted tumor fate for an initial immune presence is indicated by color for each model (i.e., for “constant” and “dynamic” carrying capacities) for parameter sets 0, 1, 2, and 3 as appropriate. The dormancy windows, i.e., ranges of initial immune presence guaranteeing tumor dormancy, are indicated by the red (escape afterwards) and blue (elimination afterwards) regions. Pink regions correspond to tumor escape and light blue regions to elimination. Parameter values are |$r = 0.47$| and |$\theta = 4$| as used in the rest of Fig. 4 to simulate immunotherapy. As can be seen, the constant carrying capacity model with parameter set 0 predicts only elimination, with dormancy occurring for small initial immune presence. When a dynamic carrying capacity is used, however, the three parameter sets (sets 1D, 2D, and 3D) predict different outcomes with different dormancy windows. After disrupting the growth regulation signals (sets 1C, 2C, and 3C), these dormancy windows tend to shift to the right, indicating that a higher initial immune presence is required to obtain dormancy. For comparison, the left-most color of Fig. 4C corresponds to outcomes predicted in Fig. 4A and B. Because parameter set 3 requires a much higher initial immune presence to obtain dormancy, and the size of the window is small compared with those predicted by sets 1 and 2, we conclude that dormancy is harder to achieve in this instance due to high sensitivity to environmental regulation.

Dormancy is explored further through the parameter sensitivity analysis shown in Fig. 5 using the constant carrying capacity model and parameters from Table 1. Each immune parameter is varied, and the resulting tumor fate is indicated by color for ranges of initial immune presence. Under the assumption that large dormancy windows at low levels of initial immune presence correspond to more easily achieved dormancy, these sensitivities may suggest treatment strategies targeted at obtaining and maintaining the dormant state.

Figure 5.

Parameter sensitivity for immune system parameters [λ (A), r (B), K_I (C), θ (D), φ (E), and β (F)] using the model with constant carrying capacities and parameter values from Table 1 (set 0). The lower (pink) regions indicate tumor escape, the upper (light blue) regions indicate tumor elimination, and the mid-regions indicate that the tumor either escaped (red) or was eliminated (blue) after a period of at least 30 days of dormancy. Vertical axes indicate the initial immune presence relative to the initial cancer presence. Theoretically, tumor dormancy should be easier for the immune response to achieve when the required initial presence I0 is small and the width of the dormancy region (red and blue, middle regions) is large. The width of the dormancy regions (measured in initial immune presence) is indicated at the left and right end points of the curves as well as on the inset line graphs.

Figure 5.

Parameter sensitivity for immune system parameters [λ (A), r (B), K_I (C), θ (D), φ (E), and β (F)] using the model with constant carrying capacities and parameter values from Table 1 (set 0). The lower (pink) regions indicate tumor escape, the upper (light blue) regions indicate tumor elimination, and the mid-regions indicate that the tumor either escaped (red) or was eliminated (blue) after a period of at least 30 days of dormancy. Vertical axes indicate the initial immune presence relative to the initial cancer presence. Theoretically, tumor dormancy should be easier for the immune response to achieve when the required initial presence I0 is small and the width of the dormancy region (red and blue, middle regions) is large. The width of the dormancy regions (measured in initial immune presence) is indicated at the left and right end points of the curves as well as on the inset line graphs.

Close modal

The analysis presented here shows that, in most cases, immune-induced tumor dormancy is a transient state that necessarily ends in either elimination or escape. This result contrasts with behavior predicted by many predator–prey-type models wherein dormancy is described by decaying oscillations converging to a limit cycle or steady-state, that, once obtained, is maintained for a long time. The transient nature of the dormant state described here complements the immunoediting hypothesis because the equilibrium stage must end in either tumor elimination or escape, and qualitatively matches real-world observations (33).

Cancer and immune cell behaviors are constantly evolving through autocrine and paracrine cytokine feedback loops in the microenvironment. Tumor escape from immune-induced dormancy can involve selective editing of cancer cell immunogenicity (34) and accumulated resistance to apoptotic signals from immune cells or chemotherapeutics (6, 35). In addition, the local cytokine milieu can determine the antitumor or protumor polarization of the immune response (5, 36, 37). These phenomena are now being investigated through mathematical models of both tumor escape mechanisms (18–20) and immune polarization. The model presented here does not explicitly address these phenomena, yet captures the transient nature of dormancy in a comparatively simple framework.

With this model, we investigated the range of immune efficacy required to guarantee the existence of a dormant state. Obtaining this state, however, is not guaranteed as we also showed how growth regulatory signals emanating from the tumor microenvironment play a significant role in determining tumor fate. The regulatory mechanism of a carrying capacity that grows with the tumor may be disturbed after debulking treatments that leave the microenvironment unaltered or possibly enhanced. We found that, for tumors that are highly sensitive to these regulatory signals, tumor regrowth may be accelerated, offering one possible rationale for the phenomenon of accelerated repopulation frequently observed in the clinic following treatment.

Variability in tumor sensitivities, in turn, may account for unpredictable outcomes following therapy and thus may serve as a proxy for patient variability. The mathematical model discussed here is the first, to our knowledge, to incorporate such variability and predict outcomes that qualitatively match those seen in the clinic. Tracking variability in sensitivity to environmental regulatory signals may well improve our understanding of why treatments work for some patients but not others.

No potential conflicts of interest were disclosed.

Conception and design: K.P. Wilkie

Development of methodology: K.P. Wilkie

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): P. Hahnfeldt

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): K.P. Wilkie, P. Hahnfeldt

Writing, review, and/or revision of the manuscript: K.P. Wilkie, P. Hahnfeldt

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): K.P. Wilkie

Study supervision: K.P. Wilkie, P. Hahnfeldt

The authors thank Dr. M. La Croix for his help with computation and illustrations.

This work was financially supported by the National Cancer Institute under Award Number U54CA149233 (L. Hlatky) and by the Office of Science (BER), U.S. Department of Energy, under Award Number DE-SC0001434 (P. Hahnfeldt).

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.
Grivennikov
SI
,
Greten
FR
,
Karin
M
. 
Immunity, inflammation, and cancer
.
Cell
2010
;
140
:
883
99
.
2.
Dudley
ME
,
Wunderlich
JR
,
Robbins
PF
,
Yang
JC
,
Hwu
P
,
Schwartzentruber
DJ
, et al
Cancer regression and autoimmunity in patients after clonal repopulation with antitumor lymphocytes
.
Science
2002
;
298
:
850
4
.
3.
Diefenbach
A
,
Jensen
ER
,
Jamieson
AM
,
Raulet
DH
. 
Rae1 and H60 ligands of the NKG2D receptor stimulate tumour immunity
.
Nature
2001
;
413
:
165
71
.
4.
de Pillis
LG
,
Radunskaya
AE
,
Wiseman
CL
. 
A validated mathematical model of cell-mediated immune response to tumor growth
.
Cancer Res
2005
;
65
:
7950
8
.
5.
Saudemont
A
,
Jouy
N
,
Hetuin
D
,
Quesnel
B
. 
NK cells that are activated by CXCL10 can kill dormant tumor cells that resist CTL-mediated lysis and can express B7-H1 that stimulates T cells
.
Blood
2005
;
105
:
2428
35
.
6.
Saudemont
A
,
Quesnel
B
. 
In a model of tumor dormancy, long-term persistent leukemic cells have increased B7-H1 and B7.1 expression and resist CTL-mediated lysis
.
Blood
2004
;
104
:
2124
33
.
7.
de Visser
KE
,
Eichten
A
,
Coussens
LM
. 
Paradoxical roles of the immune system during cancer development
.
Nat Rev Cancer
2006
;
6
:
24
37
.
8.
Liotta
L
,
Kohn
E
. 
The microenvironment of the tumour–host interface
.
Nature
2001
;
411
:
375
9
.
9.
Naumov
GN
,
Bender
E
,
Zurakowski
D
,
Kang
SY
,
Sampson
D
,
Flynn
E
, et al
A model of human tumor dormancy: an angiogenic switch from the nonangiogenic phenotype
.
J Natl Cancer Inst
2006
;
98
:
316
25
.
10.
Teng
MWL
,
Swann
JB
,
Koebel
CM
,
Schreiber
RD
,
Smyth
MJ
. 
Immune-mediated dormancy: an equilibrium with cancer
.
J Leukoc Biol
2008
;
84
:
988
93
.
11.
Dunn
GP
,
Old
LJ
,
Schreiber
RD
. 
The three Es of cancer immunoediting
.
Annu Rev Immunol
2004
;
22
:
329
360
.
12.
Page
K
,
Uhr
JW
. 
Mathematical models of cancer dormancy
.
Leuk Lymphoma
2005
;
46
:
313
27
.
13.
Kuznetsov
VA
. 
Mathematical modeling of the development of dormant tumors and immune stimulation of their growth
.
Cybern Syst Anal
1988
;
23
:
556
64
.
14.
Kirschner
D
,
Panetta
JC
. 
Modeling immunotherapy of the tumor–immune interaction
.
J Math Biol
1998
;
37
:
235
52
.
15.
Kuznetsov
VA
,
Makalkin
IA
,
Taylor
M
,
Perelson
AS
. 
Nonlinear dynamics of immunogenic tumors: parameter estimation and global bifurcation analysis
.
Bull Math Biol
1994
;
56
:
295
321
.
16.
Arciero
J
,
Jackson
TL
,
Kirschner
DE
. 
A mathematical model of tumor–immune evasion and siRNA treatment
.
Discret Contin Dyn S
2004
;
4
:
39
58
.
17.
Wilkie
KP
. 
A review of mathematical models of cancer–immune interactions in the context of tumor dormancy
.
Adv Exp Med Biol
2013
;
734
:
201
34
.
18.
d'Onofrio
A
. 
Bounded-noise-induced transitions in a tumor–immune system interplay
.
Phys Rev E
2010
;
81
:
021923
:
1
7
.
19.
d'Onofrio
A
,
Ciancio
A
. 
Simple biophysical model of tumor evasion from immune system control
.
Phys Rev E
2011
;
84
:
031910
:
1
7
.
20.
Al-Tameemi
MM
,
Chaplain
MA
,
d'Onofrio
A
. 
Evasion of tumours from the control of the immune system: consequences of brief encounters
.
Biol Direct
2012
;
7
:
31
.
21.
Withers
H
. 
Treatment-induced accelerated human tumor growth
.
Semin Radiat Oncol
1993
;
3
:
135
43
.
22.
Sharouni El
SY
,
Kal
HB
,
Battermann
JJ
. 
Accelerated regrowth of non-small-cell lung tumours after induction chemotherapy
.
Br J Cancer
2003
;
89
:
2184
9
.
23.
Nguyen
DH
,
Oketch-Rabah
HA
,
Illa-Bochaca
I
,
Geyer
FC
,
Reis-Filho
JS
,
Mao
J-H
, et al
Radiation acts on the microenvironment to affect breast carcinogenesis by distinct mechanisms that decrease cancer latency and affect tumor type
.
Cancer Cell
2011
;
19
:
640
51
.
24.
Tanooka
H
,
Tanaka
K
,
Arimoto
H
. 
Dose response and growth rates of subcutaneous tumors induced with 3-methylcholanthrene in mice and timing of tumor origin
.
Cancer Res
1982
;
42
:
4740
3
.
25.
Robert
CP
,
Casella
G
. 
Metropolis–Hastings algorithms
. 
Introducing Monte Carlo methods with R
.
Springer Science+Business Media, New York, NY
; 
2010
.
p.
167
97
.
26.
Cirit
M
,
Haugh
JM
. 
Data-driven modelling of receptor tyrosine kinase signalling networks quantifies receptor-specific potencies of PI3K- and Ras-dependent ERK activation
.
Biochem J
2012
;
441
:
77
85
.
27.
Sotolongo-Costa
O
,
Morales
Molina L
,
Rodriguez
Perez D
,
Antoranz
J
,
Chacon
Reyes M
. 
Behavior of tumors under nonstationary therapy
.
Physica D
2003
;
178
:
242
53
.
28.
d'Onofrio
A
. 
A general framework for modeling tumor–immune system competition and immunotherapy: mathematical analysis and biomedical inferences
.
Physica D
2005
;
208
:
220
35
.
29.
Hahnfeldt
P
,
Panigrahy
D
,
Folkman
J
,
Hlatky
LR
. 
Tumor development under angiogenic signaling: a dynamical theory of tumor growth, treatment response, and postvascular dormancy
.
Cancer Res
1999
;
59
:
4770
5
.
30.
Polyak
K
,
Kalluri
R
. 
The role of the microenvironment in mammary gland development and cancer
.
Cold Spring Harb Perspect Biol
2010
;
2
:
a003244
.
31.
Kim
JJ
,
Tannock
IF
. 
Repopulation of cancer cells during therapy: an important cause of treatment failure
.
Nat Rev Cancer
2005
;
5
:
516
25
.
32.
Hu
M
,
Yao
J
,
Carroll
DK
,
Weremowicz
S
,
Chen
H
,
Carrasco
D
, et al
Regulation of in situ to invasive breast carcinoma transition
.
Cancer Cell
2008
;
13
:
394
406
.
33.
Matsushita
H
,
Vesely
MD
,
Koboldt
DC
,
Rickert
CG
,
Uppaluri
R
,
Magrini
VJ
, et al
Cancer exome analysis reveals a T-cell-dependent mechanism of cancer immunoediting
.
Nature
2012
;
482
:
400
4
.
34.
Koebel
CM
,
Vermi
W
,
Swann
JB
,
Zerafa
N
,
Rodig
SJ
,
Old
LJ
, et al
Adaptive immunity maintains occult cancer in an equilibrium state
.
Nature
2007
;
450
:
903
7
.
35.
Saudemont
A
,
Hamrouni
A
,
Marchetti
P
,
Liu
J
,
Jouy
N
,
Hetuin
D
, et al
Dormant tumor cells develop cross-resistance to apoptosis induced by CTLs or imatinib mesylate via methylation of suppressor of cytokine signaling 1
.
Cancer Res
2007
;
67
:
4491
8
.
36.
Maldonado
RA
,
Irvine
DJ
,
Schreiber
R
,
Glimcher
LH
. 
A role for the immunological synapse in lineage commitment of CD4 lymphocytes
.
Nature
2004
;
431
:
527
32
.
37.
DeNardo
DG
,
Andreu
P
,
Coussens
LM
. 
Interactions between lymphocytes and myeloid cells regulate pro- versus anti-tumor immunity
.
Cancer Metast Rev
2010
;
29
:
309
16
.