# Computational Modeling of the Crosstalk Between Macrophage Polarization and Tumor Cell Plasticity in the Tumor Microenvironment

^{1}Center for Theoretical Biological Physics, Rice University, Houston, TX, United States^{2}Department of Bioengineering, Rice University, Houston, TX, United States^{3}Medical Scientist Training Program, Baylor College of Medicine, Houston, TX, United States^{4}The James Buchanan Brady Urological Institute, Johns Hopkins University School of Medicine, Baltimore, MD, United States^{5}Department of Physics and Astronomy, Rice University, Houston, TX, United States^{6}Department of Physics, Northeastern University, Boston, MA, United States

Tumor microenvironments contain multiple cell types interacting among one another via different signaling pathways. Furthermore, both cancer cells and different immune cells can display phenotypic plasticity in response to these communicating signals, thereby leading to complex spatiotemporal patterns that can impact therapeutic response. Here, we investigate the crosstalk between cancer cells and macrophages in a tumor microenvironment through *in silico* (computational) co-culture models. In particular, we investigate how macrophages of different polarization (M_{1} vs. M_{2}) can interact with epithelial-mesenchymal plasticity of cancer cells, and conversely, how cancer cells exhibiting different phenotypes (epithelial vs. mesenchymal) can influence the polarization of macrophages. Based on interactions documented in the literature, an interaction network of cancer cells and macrophages is constructed. The steady states of the network are then analyzed. Various interactions were removed or added into the constructed-network to test the functions of those interactions. Also, parameters in the mathematical models were varied to explore their effects on the steady states of the network. In general, the interactions between cancer cells and macrophages can give rise to multiple stable steady-states for a given set of parameters and each steady state is stable against perturbations. Importantly, we show that the system can often reach one type of stable steady states where cancer cells go extinct. Our results may help inform efficient therapeutic strategies.

## Introduction

Cancer has been largely considered as a cell-autonomous disease, but recent investigations have highlighted the crucial role of the tumor microenvironment in determining cancer progression (1). Cancer cells can communicate bi-directionally through various mechanical and/or chemical ways with their neighboring cancer cells (2, 3), and/or with other components of the tumor microenvironment such as macrophages and fibroblasts, driving aggressive malignancy (4–6). The interconnected feedback loops formed by these interactions can often generate many emergent outcomes for the tumor. Interestingly, many of the latest therapeutic innovations such as immunotherapy are aimed at targeting aspects of the tumor microenvironment instead of the cancer cells (7).

Tumor-associated macrophages (TAMs) are one of the most abundant immune cell populations in the microenvironment (8, 9). They have been shown to promote cancer progression in many ways, such as promoting angiogenesis, suppressing function of cytotoxic T lymphocytes, and assisting extravasation of cancer cells (8–12). Generally, the secretome and functions of TAMs have been shown to be close to that of the so-called alternatively activated macrophages (M_{2}) (13). In the case of pathogen infections, macrophages that can engulf the pathogen and present processed antigens to adaptive immune cells, are generally characterized as the classically activated ones (M_{1}) (14). M_{1} and M_{2} macrophages have different roles during wound healing: while M_{1} macrophages initiate inflammatory responses, M_{2} macrophages contribute to tissue restoration (13). In the context of cancer, M_{1} macrophages have been generally considered anti-tumor (15–17), whereas M_{2} macrophages have been considered as pro-tumor (10).

However, macrophage polarization is not as rigid as the differentiation of T lymphocytes (18); instead, M_{1}, M_{2}, and any intermediate state(s) of macrophage polarization are quite plastic (13, 19, 20). Thus, the idea that reverting TAMs in the cancer microenvironment to its cancer-suppressing counterpart is tempting, the proof of principle of which has been demonstrated in mice models (21–25).

Not only TAMs, but also cancer cells themselves can be extremely plastic, a canonical example of which is epithelial-mesenchymal plasticity, i.e., cancer cells can undergo varying degrees of Epithelial-Mesenchymal Transition (EMT) and its reverse Mesenchymal-Epithelial Transition (MET) (26, 27). EMT/MET has been associated with metastasis (28), chemoresistance (29), tumor-initiation potential (30), resistance against cell death (31), and evading the immune system (32).

Importantly, macrophages and cancer cells can interact with and influence the behavior of one another, as shown by many *in vitro* experiments. Specifically, some epithelial cancer cells are capable of polarizing monocytes into M_{1}-like macrophages (33). Forming a negative feedback loop, these M_{1}-like macrophages can decrease the confluency of the cancer cells that polarized them (33). Moreover, pre-polarized M_{1} macrophages can induce senescence and apoptosis in human cancer cell lines A549 (34) and MCF-7 (35). Intriguingly, factors released by apoptosis of cancer cells can convert M_{1} macrophages into M_{2}-like macrophages (35), thus switching macrophage population from being tumor-suppressive to being tumor-promoting. On the other hand, mesenchymal cancer cells can polarize monocytes into M_{2}-like macrophages (33, 36, 37), that can in turn assist EMT (37, 38). Thus, the interaction network among macrophages and cancer cells is formidably complex, and the emergent dynamics of these interactions can be non-intuitive (39) yet are often crucial in deciding the success of therapeutic strategies targeting cancer and/or immune cells. For example, even if TAMs at some time can be converted exogenously to M_{1}-like macrophages, if most cancer cells still tend to polarize monocytes to TAMs, other coordinated perturbations may be needed at different time-points to restrict the aggressiveness of the disease.

Here, we develop mathematical models to capture the abovementioned set of interactions among cancer cells in varying phenotypes (epithelial and mesenchymal) and macrophages of different polarizations (M_{1}-like and M_{2}-like). We characterize the multiple steady states of the network that can be obtained as a function of different initial conditions and key parameters, and thus analyze various potential compositions of cellular populations in the tumor microenvironment. This *in silico* co-culture system can not only help explain *in vitro* multiple experimental observations and clinical data, but also help acquire novel insights into designing effective therapeutic strategies aimed at cancer cells and/or macrophages.

## Results

### Crosstalk Among Cancer Cells and Macrophages Can Lead to Two Distinct Categories of Steady States

We first considered the following interactions in setting up our mathematical model: (a) proliferation of epithelial and mesenchymal cells (but not that of monocytes M_{0}, or macrophages M_{1} and M_{2}), (b) EMT promoted by M_{2}-like macrophages and MET promoted by M_{1}-like macrophages, (c) polarization of monocytes (M_{0}) to M_{1}-like cells aided by epithelial cells, and that to M_{2}-like cells aided by mesenchymal cells, (d) induction of senescence in epithelial cells by M_{1}-like macrophages (Figure 1A). No inter-conversion among M_{1}-like and M_{2}-like macrophages or cell-death of macrophages is considered here in this model (hereafter referred to as “Model I”; see section Methods).

**Figure 1**. Cancer-immune interplay can give rise to the co-existence of two types of steady states. **(A)** Interaction network for Model I. Conversions between cells are indicated by solid lines. Cell proliferation is indicated by dashed lines. Inhibition (in black) and activation (in red) is indicated by dotted lines. **(B)** Steady states of the epithelial population are plotted as a function of M_{2}-like macrophage population. Stable steady states are plotted in solid blue lines and unstable steady states are plotted in dotted red line. The key parameters are as indicated. **(C)** As the cooperativity of epithelial cancer cells in their resistance of M_{2}-promoted EMT is reduced, i.e., *k* = 4, 3, 2 (blue, cyan, and light-green line), the overlapped region between state I and II shrinks and then disappears. Increasing the cooperativity of M_{2}-promoted EMT or M_{1}-promoted MET, i.e., m and n increase from 1 to 2, can expand the overlapped region (between state I and II) of the system (dark-green lines). Stable steady states are plotted in solid lines and unstable steady states are plotted in dotted lines. **(D)** As the total number of macrophages (M_{c}) increases, the overlapped region (between state I and II) of the system shrinks and then disappears.

Furthermore, this model also considers that mesenchymal cells can secrete soluble factors, such as TGFβ, that can induce or maintain the mesenchymal state in autocrine or paracrine manners (40, 41). Therefore, we hypothesized that mesenchymal cells can resist M_{1}-promoted MET. However, this resistance might not fully suppress MET, as MET still happens in the presence of M_{1} macrophages (38). Therefore, in Equation (1), we assume a Hill-like function to represent the resistance to M_{1}-promoted MET and the resistance term saturates to a finite value as a function of mesenchymal population M. And the corresponding half-saturation constant is ${\text{K}}_{\text{M}}^{0}$. Similarly, epithelial cells adhere to each other via E-cadherin, sequestering β-catenin on the membrane, thus interfering with the induction of EMT (42). Thus, we hypothesized that the M_{2}-promoted EMT can be inhibited by epithelial cells in a cooperative manner, because this inhibition of EMT requires direct physical cell-cell contact and hence involves multiple epithelial cells. Therefore, in Equation (1), we assume a Hill-like function to represent this resistance to M_{2}-dependent EMT and the resistance term saturates to a finite value as a function of epithelial population E. And the corresponding half-saturation constant is ${\text{K}}_{\text{E}}^{0}$. Furthermore, the epithelial-cell-dependent term has a cooperativity parameter k to represent the effects of E-cadherin-β-catenin interaction between multiple epithelial cells.

Note that in all of our calculations, we assumed that there is a carrying-capacity of cells (N_{max}, including both cancer cells and macrophages) in the co-culture system *in vitro*. We vary the initial number of different cells while keeping the total number of macrophages (M_{1}+M_{2}+M_{0}) to be a constant (M_{c}). Therefore, the maximum number of cancer cells will be N_{max} -M_{c}.

In Model I, the final populations of E, M, M_{1}, and M_{2} cells are simply determined by the following equations:

where η_{em} and η_{me} are the EMT and MET rate constants, respectively. The above equations give two categories of stable steady states: (a) state I, dominated by epithelial cancer cells, and (b) state II, dominated by mesenchymal cancer cells. The final steady state of the system depends on the number of M_{2} macrophages in that state; since there are only three equations for four unknowns, M_{2} can be used as a parameter specifying (possibly discrete set of) solutions. As soon as M_{2}/N_{max} crosses a threshold, the system switches from state I to state II (Figure 1B). This prediction is largely robust to parameter variation (Figure S1).

Next, we explored what factors could change the qualitative behavior of the model, i.e., enable a smoother and continuous transition of epithelial and mesenchymal percentages as a function of the M_{2}-macrophage population. We identified that reducing the cooperativity of epithelial cancer cells in their resistance of M_{2}-promoted EMT can lead to a loss of the feature with two-types of steady states (Figure 1C, blue, cyan, and light-green lines). Conversely, increasing the cooperativity of M_{2}-promoted EMT or M_{1}-promoted MET can expand the region of overlapping between state I and II of the system (Figure 1C, green lines). Another factor that can alter the behavior of the model is the initial number of monocytes in the system. At high enough number of monocytes, the absolute number of cancer cells will be small, thus the effect of the cooperativity between epithelial cancer cells will be reduced, and consequently, as discussed above, the overlap of the two types of steady states of the system will disappear (Figure 1D). Together, this increased propensity of multi-stability in the system upon including cooperative effects in the interactions among different species (or variables) is reminiscent of observations in models of biochemical networks (43).

Note that using the M_{2}/N_{max} as the “control variable” is specific for Model I, because there is no interconversion between M_{1} and M_{2} macrophages. For a given set of parameters, the steady state level of M_{2} can vary continuously from 0 to M_{c} (Figures 1B–D), and in practice is determined by the initial conditions (initial number of E, M, M_{1}, M_{2}, and M_{0}). For models in the following sections, M_{2}/N_{max} will be shown to be fixed to discrete allowed values (instead of continuously varying) for a given set of parameters. For example, if we add a very small conversion rate between M_{1} and M_{2}, the steady state of the system will collapse to only one steady state (Figure S2): on the shorter time scale, the trajectory of the system (on the phase plane of E and M_{2}) will first converge to one steady state in Model I (with the same M_{2}/N_{max}); on the longer time scale, as determined by the value of inter-conversion rate between M_{1} and M_{2}, the system will slowly evolve to the steady state (blue dot in Figure S2), following along the now-approximate steady states in Model I.

### Cancer-Cell Enhanced Interconversion Between M_{1}- and M_{2}-Like Macrophages Lead to Bi-Stability

In the next iteration of our model (hereafter referred to as “Model II”), we included the possibility of interconversion between M_{1}-like and M_{2}-like macrophages, as reported in the literature (13, 35, 44, 45) (Figure 2A, see section Methods). We first assume constant interconversion between M_{1} and M_{2} (with rates denoted as ${\eta}_{21}^{0}$ and ${\eta}_{12}^{0}$) and investigate the effects of varying the rate of conversion of mesenchymal cells to epithelial cells (η_{me}). At small values of η_{me}, the system has small number of epithelial cells, which is equivalent to state II in Model I; with increasing η_{me}, a threshold is crossed, and the system can switch to states with larger number of epithelial cells, which is equivalent to state I in Model I (solid black lines, Figure 2B). Thus, we can observe bi-stability of cancer cell population in this system. However, the populations of macrophages stay constant as a function of η_{me}, since the M_{1}/M_{2} ratio is simply determined by ${\eta}_{21}^{0}$/${\eta}_{12}^{0}$ according to Equation (2) when η_{12} = η_{21} = 0). In this model, we focused on increased cooperativity of M_{2}-assisted EMT and M_{1}-assisted MET (i.e., m = n = 2), because the interconversion between M_{1} and M_{2}, in absence of such cooperativity (as considered in model I with m = n = 1; Figure 1B), gives rise to a narrower bi-stable region (Figure S3). Note again that in our calculations, we assumed that there is a carrying-capacity of cells (N_{max}, including both cancer cells and macrophages) in the co-culture system *in vitro*. Again, the total number of macrophages is constant, since no cell death or cell division of macrophages is considered.

**Figure 2**. Effects of interconversion between M_{1} and M_{2} cells, mediated by cancer cells. **(A)** Network of Model II. **(B,C)** Steady states of the epithelial **(B)** or M_{1} **(C)** populations are plotted as a function of η_{me}. Stable steady states are plotted in solid lines and unstable steady states are plotted in dotted lines. The key parameters in **(B,C)** are: η_{12} = η_{21} = 0 (for black lines) or 1/72 h^{−1} (for blue lines), *m* = 2, *n* = 2, *k* = 4, ${\eta}_{21}^{0}$ = ${\eta}_{12}^{0}=$1/72 h^{−1}. **(D)** In this plot, the key parameters are: η_{21} = 1/72 h^{−1} (for blue lines) or 1/36 h^{−1} (for green lines), *m* = 2, *n* = 2, *k* = 4, ${\eta}_{21}^{0}$ = ${\eta}_{12}^{0}$ = 1/72 h^{−1}, and η_{12} = 1/72 h^{−1}. In **(B,D)**, the gray line is the maximum fraction of cancer cells (=0.7) as M_{c} is set as 0.3.

We chose η_{me} as a control parameter, because it can act as a bottleneck for the transition between the two types of stable steady states of the system. Lowering the transition threshold of η_{me} can be helpful in a sense that the system can potentially stay at state I (high epithelial state, supposed to be less aggressive) only. Due to the inherent symmetry in the network, the effect of lowering the threshold of η_{me} can be recapitulated via other perturbations, such as a smaller rate of M_{1}-assisted MET (η_{em}), lower M_{1} to M_{2} conversion rate (η_{12}), or higher M_{2} to M_{1} conversion rate (η_{21}) (Figure S4).

We next consider the case where interconversion between M_{1}-like and M_{2}-like macrophages is enhanced by cancer cells, i.e., mesenchymal cells enhance the M_{1}-to-M_{2} transition, while epithelial cells enhance M_{2}-to-M_{1} transition. In this case, we observe the existence of a bi-stable region, and the range of the epithelial-cell-low solution (equivalent to state II) is wider than that for the previous case without any effects of epithelial (E) and mesenchymal (M) cancer cells on M_{1}-M_{2} interconversion (Figures 2B,C, solid blue lines). The reason for the expanded epithelial-cell-low region is that the positive feedback loop between mesenchymal and M_{2} cells makes the mesenchymal- and M_{2}-dominated state more stable. Therefore, a higher η_{me} is required to compensate for the effects of low M_{1} population. For therapeutic purposes, state I is favored, for it is believed that epithelial cells are typically less aggressive than mesenchymal ones. Therefore, symmetrical mutual enhancement might not be helpful for the therapy because of the expanded epithelial-low (mesenchymal-high) region. For the same reason, the lower threshold of η_{me} to switch from state II to state I also shift to a higher value, which means that for a small η_{me}, the system will only have one stable steady state, i.e., state II with high number of mesenchymal cells.

In order to drive the system to be bi-stable at a smaller η_{me}, we tested the case with asymmetrical interconversion rates between M_{1} and M_{2}, i.e., higher conversion rate from M_{2} to M_{1} enhanced by epithelial cells (η_{21}). In this case, the lower threshold of η_{me} can indeed shift to a smaller value (Figure 2D), which makes it theoretically possible to switch the system from state II (low epithelial cells) to state I (high epithelial cells) at a fairly small η_{me} value.

In summary, our results for Model II suggest that one can switch the system from state II (low epithelial cells) to state I (high epithelial cells) by increasing both η_{me} and the effective conversion from M_{2} to M_{1}. Note that the interactions in Model 2 are rather symmetric: the interconversion rates between M_{1} and M_{2} are either constants or enhanced by the corresponding cancer-cells. In the following section, we will consider the asymmetric case where only the M_{1} to M_{2} conversion is enhanced by factors released by apoptotic epithelial cancer cells.

### Cell Apoptosis-Induced M1-to-M2 Conversion Leads to Symmetry Breaking in the Cancer-Immune Interaction Network

Finally, we incorporated one other set of experimentally documented interactions, i.e., M_{1} macrophages can drive the apoptosis of epithelial cells, and factors released during cancer cell apoptosis can drive M_{1}-to-M_{2} conversion (referred as Model III, see section Methods) (35). The newly introduced interactions are highlighted in thick lines in Figure 3A. This interaction induces “symmetry breaking” (46) in the system, as previously most of the interactions considered were of a “symmetric” nature, i.e., M_{1} cells driving MET and M_{2} cells driving EMT, and epithelial cells driving M_{1} maturation while mesenchymal cells driving M_{2} maturation. With this new interaction, the system is now biased against epithelial cells, because (a) epithelial but not mesenchymal cells, are killed by M_{1}-macrophages, and (b) their dead counterparts may convert M_{1} to M_{2} cells that can, in turn, convert some epithelial cells to mesenchymal ones (EMT). In addition, the conversion from M_{2} to M_{1} is essential to bring the system back from the mesenchymal-cancer-cell biased state and it can be enhanced by introducing IL-12 (45) or Type 1 T helper cells (44).

**Figure 3**. Steady state solutions for a model including the M_{1} to M_{2} conversion assisted by apoptotic epithelial cancer cells**. (A)** Network of Model III. In this network, M_{1} induces the apoptosis of epithelial cancer cells. And the factors released by apoptotic cancer cells enhance the M_{1} to M_{2} conversion. These interactions are highlighted by thicker lines. **(B–E)** Steady states of the epithelial **(B)**, mesenchymal **(C,D)** and M_{1} **(E)** populations are plotted as a function of η_{me}. **(D)** is a zoomed-in version of **(C)**. Stable steady states are plotted in solid lines and unstable steady states are plotted in dotted lines. The model-specific parameters are: β = 1/36 h^{−1}, β_{c} = 1/1200 h^{−1}, η_{12} = η_{21} = 1/72 h^{−1}. When η_{me} is higher than a critical value ${\eta}_{\text{me}}^{\text{a}}$, the steady state E = M = 0 becomes stable **(D)**. When η_{me} is higher than a critical value ${\eta}_{\text{me}}^{\text{c}}$, the steady state E = M = 0 is the only stable steady state of the system. **(F)** For a given λ_{M}, the critical value ${\eta}_{\text{me}}^{\text{c}}$ is plotted.

Thus, in the parameter region investigated, there are 3 types of stable steady states, which are represented by solid blue, red and black lines, respectively. At small values of the control parameter η_{me} (rate of M_{1} macrophage assisted MET), there is only one type of stable steady state: the system is biased toward the mesenchymal-dominated state (solid blue curves in Figures 3B–E), whereas the cancer-extinction state (E = M = 0, dashed black curves in Figures 3B–E) is unstable. As η_{me} increases and goes across a critical value ${\eta}_{\text{me}}^{\text{a}}=$0.0626 h^{−1}), the extinction state (E = M = 0) becomes stable as well as a new set of steady states emerges (red lines in Figures 3B–E). Between η_{me} = 0.0663 h^{−1} and η_{me} = 0.0747 h^{−1}, the steady states depicted as a solid red line (Figures 3B–E) are stable; here both populations of epithelial and mesenchymal cancer cells are at a low level. This phenomenon can be understood in the following sense: at higher η_{me}, proliferating mesenchymal cells are continuously being converted to epithelial cells, which will be killed by M_{1} macrophages. This effect brings down the overall fraction of cancer cells. Note that in this region, three types of stable steady states co-exist in the system for a given set of parameters. As η_{me} increases and goes across ${\eta}_{\text{me}}^{\text{b}}=$0.0747 h^{−1}, the stable steady states in solid red lines disappears and the other two types of stable steady states co-existed (solid blue and black lines in Figures 3B–E). As η_{me} further increases and goes across ${\eta}_{\text{me}}^{\text{c}}=$0.1532 h^{−1}, there is only one stable steady state: cancer cells are necessarily eliminated from the system. We note in passing that the instability that drives the red solution unstable as η_{me} is lowered past 0.0663 h^{−1} is a Hopf bifurcation to an unstable limit cycle.

Furthermore, for the breast cancer cell line MDA-MB-231 used in Yang et al. (38), η_{me} is estimated to be around 1/120 h^{−1} (~0.0083 h^{−1}). In order to explore the conditions for the absolute extinction of cancer cells around the estimated η_{me}, we varied different parameters (such as λ_{M}, η_{em}, and η_{21}) to get the corresponding ${\eta}_{\text{me}}^{\text{c}}$ where the extinction state is the only stable steady state of the system. We found that lowering the growth rate of mesenchymal cells (λ_{M}) can reduce ${\eta}_{\text{me}}^{\text{c}}$ to the nominal experimental value 1/120 h^{−1} (Figure 3F) whereas lowering down η_{em} or increasing η_{21} might not (Figure S5). Therefore, for therapeutic purposes, increasing η_{me} and decreasing λ_{M} can be a promising combination to help to eliminate cancer cell populations. In addition, increasing the number of macrophages in the co-culture system can also shift ${\eta}_{\text{me}}^{\text{c}}$ to a lower value (Figure S6) because of a potentially higher level of cancer-killing M_{1}-like macrophages.

### EMT Scores Correlate With Multiple Genes Upregulated in M2 Macrophages

As a proof of principle for the predictions of our model, we investigated multiple TCGA datasets (see section Methods for the source of datasets), using our previously developed EMT scoring metric (47). This metric quantifies the extent of EMT in a particular sample, and correspondingly assigns a score between 0 (fully epithelial) and 2 (fully mesenchymal). We calculated the correlation coefficients for EMT scores with various genes reported to be differently regulated in M_{1} or M_{2} macrophages relative to M_{0} macrophages. The list of those genes is taken from Jablonski et al. (48). Out of 52 genes investigated, we observed that many genes upregulated in M2 and downregulated in M_{1} macrophages–ACTN1, FLRT2, MRC1, PTGS1, RHOJ, TMEM158–correlated positively with the EMT scores (Table 1, first 6 rows) across multiple cancer types. The higher the EMT scores, the higher the levels of those genes, including the canonical M_{2} macrophage marker CD206 (MRC1). On the other hand, out of the 52 genes investigated, many genes upregulated in M_{1} macrophages but downregulated in M_{2} macrophages–FIIR, STAT1, RSAD2, TUBA4A, and XAF1–showed either a negative or an overall weak positive correlation with EMT scores (Table 1, last 5 rows). The only exception observed in this trend was that for ARHGP24 which correlated positively with EMT scores across cancer types. Together, these correlation results in multiple TCGA datasets offer a promising initial validation of our model predictions that a state dominated by epithelial cells typically has higher number of M_{1} macrophages, while the other state dominated by mesenchymal cells typically has higher number of M_{2} macrophages.

**Table 1**. Correlation coefficients of expression levels of specific genes with EMT scores, across many TCGA datasets.

Furthermore, to go beyond the correlation analysis between a single gene and the EMT-score, we investigated whether genes that are expressed higher in M_{2}-like macrophages will be enriched in EMT-score-high tumors and whether genes that are expressed higher in M_{1}-like macrophages will be enriched in EMT-score-low tumors. Again, the lists of genes that are expressed higher in M_{1}- or M_{2}-like macrophages are from Jablonski et al. (48) and the EMT-scores were calculated for each sample of the four TCGA dataset investigated in this work. We define a tumor to be EMT-score-high if its EMT-score is higher than the median score of the specific dataset and vice versa. Then, the expressions of a gene in EMT-score-high and low tumors in each dataset are used to determine whether this gene has a significantly high/low expression in EMT-score-high/low tumors. As expected, the genes that are expressed higher in M_{2}-like macrophages (but lower in M_{1}-like macrophages) are enriched in EMT-score-high tumors across the 4 types of cancers investigated (Table 2). This result is consistent with our model predictions: M_{2}-like macrophages and mesenchymal cancer cells tend to stably co-exist. However, somewhat to our surprise, the genes that expressed higher in M_{1}-like macrophages (but lower in M_{2}-like macrophages) are also enriched in EMT-score high tumors across the 4 types of cancers investigated (Table S1). This result is not in line with our model prediction where M_{1}-like macrophages are expected to stably co-exist with epithelial cancer cells. Possible reasons for this inconsistency will be discussed later. Nevertheless, the gene enrichment analysis in multiple TCGA datasets offer a promising initial validation of (most of) our model predictions.

**Table 2**. M_{2}-associated gene set enrichment: Reported are all genes (out of 31 total) having statistically significant differential expression (samples segregated based on median EMT score) across TCGA colorectal adenocarcinoma (Adeno COL; *n* = 286), prostate adenocarcinoma (Adeno PCA; *n* = 497), invasive breast cancer (BCA invasive; *n* = 1097), and lung adenocarcinoma (Lung Adeno; *n* = 515).

## Discussion

The tumor microenvironment involves multiple cell types that interact among each other in diverse ways, thus giving rise to a complex adaptive ecological system (49–51). For such a system, mathematical models can help reveal the mechanisms underlying the emergent behavior and eventually aid in designing effective therapeutic strategies to modulate those aspects of the microenvironment that fuel disease aggressiveness (52–58).

Here, we focused on the interactions among macrophages of different polarizations (M_{1} and M_{2}) and cancer cells with different phenotypes (epithelial and mesenchymal). Based on the literature, we focused on two types of models: with (Model II and III) or without (Model I) interconversion between M_{1} and M_{2} macrophages. All three models share a common feature: with a given set of parameters, multiple types of stable steady states can co-exist. More specifically, in Model I (without M_{1}-M_{2} interconversion), with a given set of parameters, the system can converge to continuous range of states depending on the initial condition. However, these steady states belong to two categories: state I with a higher epithelial population and state II with a lower epithelial population. After the system reaches a steady state, perturbations applied only to cancer cell populations cannot drive the system out of the original steady state whereas perturbations on macrophage populations will drive the system out of the original steady state. However, the perturbation on macrophage populations might not drive the system out of the original category of steady states unless the perturbation is sufficiently strong. In Model II and III, with a given set of parameters, the system can again converge to two types of steady states depending on the initial condition: state I with a higher epithelial population and state II with a lower (or even zero for Model III) epithelial population. Now, however, the states are discrete. After the system reaches a steady state, perturbations of any single cell population might not drive the system out of the original steady state. Therefore, in general, it is not easy to switch the cancer-macrophage system from a mesenchymal- and M_{2}-dominated state to an epithelial- and M_{1}-dominated state.

Mathematical approaches similar to ours may be useful in explaining, and even predicting, the efficacy of different therapeutic intervention(s) and their combinations. For example, various efforts have been made to switch populations of macrophages from M_{2}- into M_{1}-like (25), such as depletion of TAMs, reprogramming of TAMs toward M_{1}-Like macrophages, inhibition of circulating monocyte recruitment into tumor, etc. The polarization of M_{2} macrophages can also be altered via inhibiting the endothelial-mesenchymal transition (EndMT), a process similar to EMT (59). However, it is unclear that how effective these types of strategies would be. Our modeling studies suggest that when we consider the interactions in Model III, the efforts on manipulating the M_{1}-M_{2} interconversion might not be effective to eliminate cancer cells; whereas in Model II increasing the M_{2}-to-M_{1} conversion rate can help the system to switch to the epithelial cancer cell and M_{1} macrophage dominated state, which is believed to be less aggressive. Therefore, the effective therapeutic strategies strongly depend on the type of interactions present between cancer cells and macrophages.

Furthermore, we attempted to investigate whether the model-predicted stable co-existence of M_{1}-like macrophages and epithelial cancer cells (and likewise M_{2}-like macrophages and mesenchymal cancer cells) can be seen in the TCGA dataset. Indeed, we can observe the enrichment of many genes that are expressed higher in M_{2}-like macrophages in more mesenchymal tumors. This observation is consistent with our mathematical modeling analysis. However, we also detected the enrichment of genes that are expressed higher in M_{1}-like macrophages in more mesenchymal tumors. There can be several reasons for the enrichment of M_{1} genes in mesenchymal tumors: (i) the list of genes used here to identify the ones preferentially expressed in M_{1}-like macrophages is from a study of mouse macrophages (48). For human macrophages, not as many markers have been validated as those in mice; future experiments can help characterize these better; (ii) for the TCGA sample used, the gene expression data that are used to calculate the EMT-score could be not exclusively that of cancer cells but may incorporate other types of cells, such as fibroblasts. Future experiments that can separate and measure gene expression of cancer cells separately, through techniques such as laser microdissection, will be more precise in quantifying the EMT status of cancer cells.

In addition, it is also important to recognize that our model suffers from limitations. For instance: (a) it does not consider spatial aspects of these interactions, for instance, mesenchymal cells may migrate and invade through the matrix, thus changing the interactions among the cells considered in our framework; (b) it does not consider the effects of senescence on epithelial cell growth (60); and (c) it considers EMT and macrophages polarization as a binary process, whereas emerging reports support the notion that in both cases there is likely to exist a spectrum of states/ phenotypes (61, 62). A more realistic model that can overcome the abovementioned assumptions and thus reflect the dynamics of tumor microenvironment more accurately can be used to help designing a more effective way to switch and stably maintain the system into a less aggressive state.

Despite the abovementioned limitations, our model can contribute to identifying key parameters of the system. For example, it suggests that to design an effective therapy to maintain the system in a M_{1}-dominated and cancer-free steady state, not only the conversion rate from mesenchymal to epithelial cells should be significant, but also the growth rate of mesenchymal cells should be low enough. In other words, MET-inducing and cell-growth-suppressing mechanisms can together synergistically restrict disease aggressiveness.

In summary, our results show that the interaction network between tumor cells and macrophages may lead to multi-stability in the network: one state dominated by epithelial and M_{1}-like cells, the other dominated by mesenchymal and M_{2}-like cells. We also identify that inducing MET and inhibiting cancer-cell growth might be much more efficient in “normalizing” (1) the tumor microenvironment that can otherwise be engineered by cancer cells to support their growth (63).

## Methods

### Computational Models

According to *in vitro* experiments in literature, we construct three mathematical models. In Model I, factors secreted by epithelial (mesenchymal) cancer cells will polarize monocytes into M_{1}-like (M_{2}-like) macrophages. M_{1}-like macrophages will induce senescence of epithelial cancer cells and convert mesenchymal cancer cells to epithelial ones. On the other hand, M_{2}-like macrophages will convert epithelial cancer cells to mesenchymal ones. In addition, a mesenchymal cancer cell can secrete soluble factors to help maintain the mesenchymal state of itself and its neighbors, whereas an epithelial cancer cell can maintain the epithelial state through being in contact with its neighboring epithelial cancer cells. There is assumed to be no cell growth or death for macrophages and there is a maximum “carrying capacity” of cells in the system. The figure that illustrate this model is in Figure 1A. The equations to describe this system is as follows:

The description and the value of each parameter are given in Table 3.

In Model II, interconversions between M_{1} and M_{2} macrophages are included. Furthermore, the interconversions can be enhanced by corresponding cancer cells. The figure that illustrate this model is in Figure 2A. The equation to describe this system is as follows:

In the third model, additional interactions were introduced as illustrated in Figure 3A. M_{1}-like macrophages can induce apoptosis of epithelial cancer cells and factors released by apoptotic cancer cells can convert M_{1}-like macrophages into M_{2}-like macrophages. In order to restore the symmetry of the system, we further consider the therapeutic interaction: M_{2}-like macrophages can be re-polarized back to M_{1}-like macrophage by Type 1 T helper cells or IL-12. The equation to describe this system is as follows:

The description and the corresponding values of additional parameters are given in Table 3.

### Correlation Analysis

For a fixed dataset, linear regression was performed for each gene of interest against the predicted EMT score (47). In each case samples were discretized into EMT-high or EMT-low based on median EMT score for hypothesis testing. The linear correlation coefficient was recorded in each case. We performed statistical analysis under the null hypothesis of zero correlation between gene expression fold-change and EMT score and recorded the corresponding *p*-values at significance level α = 0.05 using unpaired *t*-test. All datasets (colon adenocarcinoma, *n* = 286; lung adenocarcinoma, *n* = 525; prostate adenocarcinoma, *n* = 497; breast invasive carcinoma, *n* = 1,097) were obtained from the R2: Genomics Analysis and Visualization Platform (http://r2.amc.nl).

## Author Contributions

XL, MKJ, KJP, and HL designed research. XL, MKJ, and JTG performed research. MKJ and JTG analyzed data. XL, MKJ, KJP, and HL wrote the paper.

## Funding

This work was sponsored by the National Science Foundation NSF grant PHY-1427654 (Center for Theoretical Biological Physics). XL was supported by Stand Up to Cancer and The V Foundation. MKJ was also supported by a training fellowship from Gulf Coast Consortium as computational cancer biology training grant (CPRIT RP1705593). JTG was supported by National Cancer Institute of NIH (F30CA213878). KJP is supported by National Cancer Institute NCI grant CA093900. HL is also a Cancer Prevention and Research Institute of Texas Scholar in Cancer Research of the State of Texas at Rice University.

## Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2019.00010/full#supplementary-material

## References

1. Jain RK. Normalizing tumor microenvironment to treat cancer: bench to bedside to biomarkers. *J Clin Oncol*. (2013) 31:2205–18. doi: 10.1200/JCO.2012.46.3653

2. Bocci F, Jolly MK, Tripathi SC, Aguilar M, Hanash SM, Levine H, et al. Numb prevents a complete epithelial-mesenchymal transition by modulating Notch signaling. *J R Soc Interface* (2017) 14:20170512. doi: 10.1098/rsif.2017.0512

3. Peng DH, Ungewiss C, Tong P, Byers LA, Wang J, Canales JR, et al. ZEB1 induces LOXL2-mediated collagen stabilization and deposition in the extracellular matrix to drive lung cancer invasion and metastasis. *Oncogene* (2016) 36:1925–38. doi: 10.1038/onc.2016.358

4. Wang T, Liu G, Wang R. The intercellular metabolic interplay between tumor and immune cells. *Front Immunol*. (2014) 5:358. doi: 10.3389/fimmu.2014.00358

5. del Pozo Martin Y, Park D, Ramachandran A, Ombrato L, Calvo F, Chakravarty P, et al. Mesenchymal cancer cell-stroma crosstalk promotes niche activation, epithelial reversion, and metastatic colonization. *Cell Rep.* (2015) 13:2456–69. doi: 10.1016/j.celrep.2015.11.025

6. Nagarsheth N, Wicha MS, Zou W. Chemokines in the cancer microenvironment and their relevance in cancer immunotherapy. *Nat Rev Immunol*. (2017) 17:559–72. doi: 10.1038/nri.2017.49

7. Pitt JM, Marabelle A, Eggermont A, Soria JC, Kroemer G, Zitvogel L. Targeting the tumor microenvironment: removing obstruction to anticancer immune responses and immunotherapy. *Ann Oncol.* (2016) 27:1482–92. doi: 10.1093/annonc/mdw168

8. Mantovani A, Bottazzi B, Colotta F, Sozzani S, Ruco L. The origin and function of tumor-associated macrophages. *Immunol Today* (1992) 13:265–70. doi: 10.1016/0167-5699(92)90008-U

9. Noy R, Pollard JW. Tumor-associated macrophages: from mechanisms to therapy. *Immunity* (2014) 41:49–61. doi: 10.1016/j.immuni.2014.06.010

10. Mantovani A, Schioppa T, Porta C, Allavena P, Sica A. Role of tumor-associated macrophages in tumor progression and invasion. *Cancer Metastasis Rev*. (2006) 25:315–22. doi: 10.1007/s10555-006-9001-7

11. Allavena P, Sica A, Solinas G, Porta C, Mantovani A. The inflammatory micro-environment in tumor progression: The role of tumor-associated macrophages. *Crit Rev Oncol Hematol*. (2008) 66:1–9. doi: 10.1016/j.critrevonc.2007.07.004

12. Solinas G, Germano G, Mantovani A, Allavena P. Tumor-associated macrophages (TAM) as major players of the cancer-related inflammation. *J Leukoc Biol*. (2009) 86:1065–73. doi: 10.1189/jlb.0609385

13. Sica A, Mantovani A. Macrophage plasticity and polarization: *in vivo* veritas. *J Clin Invest*. (2012) 122:787–95. doi: 10.1172/JCI59643

14. Mantovani A, Sozzani S, Locati M, Allavena P, Sica A. Macrophage polarization: tumor-associated macrophages as a paradigm for polarized M2 mononuclear phagocytes. *Trends Immunol*. (2002) 23:549–55. doi: 10.1016/S1471-4906(02)02302-5

15. Edin S, Wikberg ML, Dahlin AM, Rutegård J, Öberg Å, Oldenborg PA, et al. The Distribution of macrophages with a M1 or M2 phenotype in relation to prognosis and the molecular characteristics of colorectal cancer. *PLoS ONE* (2012) 7:e47045. doi: 10.1371/journal.pone.0047045

16. Ohri CM, Shikotra A, Green RH, Waller DA, Bradding P. Macrophages within NSCLC tumour islets are predominantly of a cytotoxic M1 phenotype associated with extended survival. *Eur Respir J*. (2009) 33:118–26. doi: 10.1183/09031936.00065708

17. Ma J, Liu L, Che G, Yu N, Dai F, You Z. The M1 form of tumor-associated macrophages in non-small cell lung cancer is positively associated with survival time. *BMC Cancer* (2010) 10:112. doi: 10.1186/1471-2407-10-112

18. Murphy KM, Stockinger B. Effector T cell plasticity: flexibility in the face of changing circumstances. *Nat Immunol*. (2010) 11:674–80. doi: 10.1038/ni.1899

19. Das A, Sinha M, Datta S, Abas M, Chaffee S, Sen CK, et al. Monocyte and macrophage plasticity in tissue repair and regeneration. *Am J Pathol*. (2015) 185:2596–606. doi: 10.1016/j.ajpath.2015.06.001

20. Moghaddam AS, Mohammadian S, Vazini H, Taghadosi M, Esmaeili S-A, Mardani F, et al. Macrophage plasticity, polarization and function in health and disease. *J Cell Physiol*. (2018) 233:6425–40. doi: 10.1002/jcp.26429

21. Guiducci C, Vicari AP, Sangaletti S, Trinchieri G, Colombo MP. Redirecting *in vivo* elicited tumor infiltrating macrophages and dendritic cells towards tumor rejection. *Cancer Res*. (2005) 65:3437–46. doi: 10.1158/0008-5472.CAN-04-4262

22. Rolny C, Mazzone M, Tugues S, Laoui D, Johansson I, Coulon C, et al. HRG inhibits tumor growth and metastasis by inducing macrophage polarization and vessel normalization through downregulation of PlGF. *Cancer Cell* (2011) 19:31–44. doi: 10.1016/j.ccr.2010.11.009

23. Xu M, Liu M, Du X, Li S, Li H, Li X, et al. Intratumoral delivery of IL-21 overcomes anti-Her2/Neu resistance through shifting tumor-associated macrophages from M2 to M1 phenotype. *J Immunol*. (2015) 194:4997–5006. doi: 10.4049/jimmunol.1402603

24. Geeraerts X, Bolli E, Fendt SM, Van Ginderachter JA. Macrophage metabolism as therapeutic target for cancer, atherosclerosis, and obesity. *Front Immunol*. (2017) 8:289. doi: 10.3389/fimmu.2017.00289

25. Genard G, Lucas S, Michiels C. Reprogramming of tumor-associated macrophages with anticancer therapies: radiotherapy versus chemo- and immunotherapies. *Front Immunol*. (2017) 8:828. doi: 10.3389/fimmu.2017.00828

26. Nieto MA, Huang RY, Jackson RA, Thiery JP. EMT: 2016. *Cell* (2016) 166:21–45. doi: 10.1016/j.cell.2016.06.028

27. Jolly MK, Ware KE, Gilja S, Somarelli JA, Levine H. EMT and MET: necessary or permissive for metastasis? *Mol Oncol*. (2017) 11:755–69. doi: 10.1002/1878-0261.12083

28. Krebs AM, Mitschke J, Losada ML, Schmalhofer O, Boerries M, Busch H, et al. The EMT-activator Zeb1 is a key factor for cell plasticity and promotes metastasis in pancreatic cancer. *Nat Cell Biol*. (2017) 19:518–29. doi: 10.1038/ncb3513

29. Siebzehnrubl FA, Silver DJ, Tugertimur B, Deleyrolle LP, Siebzehnrubl D, Sarkisian MR, et al. The ZEB1 pathway links glioblastoma initiation, invasion and chemoresistance. *EMBO Mol Med.* (2013) 5:1196–212. doi: 10.1002/emmm.201302827

30. Jolly MK, Jia D, Boareto M, Mani SA, Pienta KJ, Ben-Jacob E, et al. Coupling the modules of EMT and stemness: A tunable “stemness window” model. *Oncotarget* (2015) 6:25161–74. doi: 10.18632/oncotarget.4629

31. Huang RYJ, Wong MK, Tan TZ, Kuay KT, Ng AHC, Chung VY, et al. An EMT spectrum defines an anoikis-resistant and spheroidogenic intermediate mesenchymal state that is sensitive to e-cadherin restoration by a src-kinase inhibitor, saracatinib (AZD0530). *Cell Death Dis*. (2013) 4:e915. doi: 10.1038/cddis.2013.442

32. Tripathi SC, Peters HL, Taguchi A, Katayama H, Wang H, Momin A, et al. Immunoproteasome deficiency is a feature of non-small cell lung cancer with a mesenchymal phenotype and is associated with a poor outcome. *Proc Natl Acad Sci USA*. (2016) 113:E1555–64. doi: 10.1073/pnas.1521812113

33. Hollmén M, Roudnicky F, Karaman S, Detmar M. Characterization of macrophage–cancer cell crosstalk in estrogen receptor positive and triple-negative breast cancer. *Sci Rep*. (2015) 5:9188. doi: 10.1038/srep09188

34. Yuan A, Hsiao YJ, Chen HY, Chen HW, Ho CC, Chen YY, et al. Opposite effects of M1 and M2 macrophage subtypes on lung cancer progression. *Sci Rep*. (2015) 5:14273. doi: 10.1038/srep14273

35. Weigert A, Tzieply N, von Knethen A, Johann AM, Schmidt H, Geisslinger G, et al. Tumor cell apoptosis polarizes macrophages role of sphingosine-1-phosphate. *Mol Biol Cell* (2007) 18:3810–9. doi: 10.1091/mbc.E06-12-1096

36. Sousa S, Brion R, Lintunen M, Kronqvist P, Sandholm J, Mönkkönen J, et al. Human breast cancer cells educate macrophages toward the M2 activation status. *Breast Cancer Res*. (2015) 17:101. doi: 10.1186/s13058-015-0621-0

37. Su S, Liu Q, Chen J, Chen J, Chen F, He C, et al. A Positive feedback loop between mesenchymal-like cancer cells and macrophages is essential to breast cancer metastasis. *Cancer Cell* (2014) 25:605–20. doi: 10.1016/j.ccr.2014.03.021

38. Yang M, Ma B, Shao H, Clark AM, Wells A. Macrophage phenotypic subtypes diametrically regulate epithelial-mesenchymal plasticity in breast cancer cells. *BMC Cancer* (2016) 16:419. doi: 10.1186/s12885-016-2411-1

39. Mobius W, Laan L. Physical and mathematical modeling in experimental papers. *Cell* (2015) 163:1577–83. doi: 10.1016/j.cell.2015.12.006

40. Gregory PA, Bracken CP, Smith E, Bert AG, Wright JA, Roslan S, et al. An autocrine TGF-beta/ZEB/miR-200 signaling network regulates establishment and maintenance of epithelial-mesenchymal transition. *Mol Biol Cell* (2011) 22:1686–98. doi: 10.1091/mbc.E11-02-0103

41. Scheel C, Eaton EN, Li SHJ, Chaffer CL, Reinhardt F, Kah KJ, et al. Paracrine and autocrine signals induce and maintain mesenchymal and stem cell states in the breast. *Cell* (2011) 145:926–40. doi: 10.1016/j.cell.2011.04.029

42. Schmalhofer O, Brabletz S, Brabletz T. E-cadherin, beta-catenin, and ZEB1 in malignant progression of cancer. *Cancer Metastasis Rev*. (2009) 28:151–66. doi: 10.1007/s10555-008-9179-y

43. Jia D, Jolly MK, Harrison W, Boareto M, Ben-Jacob E, Levine H. Operating principles of tristable circuits regulating cellular differentiation. *Phys Biol*. (2017) 14:035007. doi: 10.1088/1478-3975/aa6f90

44. Heusinkveld M, de Vos van Steenwijk PJ, Goedemans R, Ramwadhdoebe TH, Gorter A, Welters MJ, et al. M2 macrophages induced by prostaglandin E2 and IL-6 from cervical carcinoma are switched to activated M1 macrophages by CD4+ Th1 cells. *J Immunol*. (2011) 187:1157–65. doi: 10.4049/jimmunol.1100889

45. Wang Y, Lin Y-X, Qiao S-L, An H-W, Ma Y, Qiao Z-Y, et al. Polymeric nanoparticles promote macrophage reversal from M2 to M1 phenotypes in the tumor microenvironment. *Biomaterials* (2017) 112:153–63. doi: 10.1016/J.BIOMATERIALS.2016.09.034

46. Frost JJ, Pienta KJ, Coffey DS. Symmetry and symmetry breaking in cancer: a foundational approach to the cancer problem. *Oncotarget* (2017) 9:11429–40doi: 10.18632/oncotarget.22939

47. George JT, Jolly MK, Xu S, Somarelli JA, Levine H. Survival outcomes in cancer patients predicted by a partial EMT gene expression scoring metric. *Cancer Res*. (2017) 77:6415–28. doi: 10.1158/0008-5472.CAN-16-3521

48. Jablonski KA, Amici SA, Webb LM, Ruiz-Rosado JDD, Popovich PG, Partida-Sanchez S, et al. Novel markers to delineate murine M1 and M2 macrophages. *PLoS ONE* (2015) 10:e0145342. doi: 10.1371/journal.pone.0145342

49. Whiteside TL. The tumor microenvironment and its role in promoting tumor growth. *Oncogene* (2008) 27:5904–12. doi: 10.1038/onc.2008.271

50. Gajewski TF, Schreiber H, Fu YX. Innate and adaptive immune cells in the tumor microenvironment. *Nat Immunol*. (2013) 14:1014–22. doi: 10.1038/ni.2703

51. Augsten M. Cancer-associated fibroblasts as another polarized cell type of the tumor microenvironment. *Front Oncol*. (2014) 4:62. doi: 10.3389/fonc.2014.00062

52. Kuznetsov VA, Knott GD. Modeling tumor regrowth and immunotherapy. *Math Comput Model* (2001) 33:1275–87. doi: 10.1016/S0895-7177(00)00314-9

53. Knútsdóttir H, Pálsson E, Edelstein-Keshet L. Mathematical model of macrophage-facilitated breast cancer cells invasion. *J Theor Biol*. (2014) 357:184–99. doi: 10.1016/j.jtbi.2014.04.031

54. Knutsdottir H, Condeelis JS, Palsson E. 3-D individual cell based computational modeling of tumor cell–macrophage paracrine signaling mediated by EGF and CSF-1 gradients. *Integr Biol.* (2016) 8:104–19. doi: 10.1039/C5IB00201J

55. Serre R, Benzekry S, Padovani L, Meille C, Andre N, Ciccolini J, et al. Mathematical modeling of cancer immunotherapy and its synergy with radiotherapy. *Cancer Res.* (2016) 76:4931–40. doi: 10.1158/0008-5472.CAN-15-3567

56. Guo Y, Nie Q, MacLean AL, Li Y, Lei J, Li S. Multiscale modeling of inflammation-induced tumorigenesis reveals competing oncogenic and oncoprotective roles for inflammation. *Cancer Res*. (2017) 77:6429–41. doi: 10.1158/0008-5472.CAN-17-1662

57. Mahlbacher G, Curtis LT, Lowengrub J, Frieboes HB. Mathematical modeling of tumor-associated macrophage interactions with the cancer microenvironment. *J Immunother Cancer* (2018) 6:10. doi: 10.1186/s40425-017-0313-7

58. Michor F, Beal K. Improving cancer treatment via mathematical modeling: surmounting the challenges is worth the effort. *Cell* (2015) 163:1059–63. doi: 10.1016/j.cell.2015.11.002

59. Choi SH, Kim AR, Nam JK, Kim JM, Kim J-Y, Seo HR, et al. Tumour-vasculature development via endothelial-to-mesenchymal transition after radiotherapy controls CD44v6+ cancer cell and macrophage polarization. *Nat Commun*. (2018) 9:5108. doi: 10.1038/s41467-018-07470-w

60. Storer M, Mas A, Robert-Moreno A, Pecoraro M, Ortells MC, Di Giacomo V, et al. Senescence is a developmental mechanism that contributes to embryonic growth and patterning. *Cell* (2013) 155:1119–30. doi: 10.1016/j.cell.2013.10.041

61. Jolly MK, Ward C, Eapen MS, Myers S, Hallgren O, Levine H, et al. Epithelial–mesenchymal transition, a spectrum of states: Role in lung development, homeostasis, and disease. *Dev Dyn* (2018) 247:346–58. doi: 10.1002/dvdy.24541

62. Aras S, Raza Zaidi M. TAMeless traitors: macrophages in cancer progression and metastasis. *Br J Cancer* (2017) 117:1583–91. doi: 10.1038/bjc.2017.356

Keywords: MET–mesenchymal-to-epithelial transition, EMT–epithelial-to-mesenchymal transition, M1-/M2-polarized macrophages, interaction network, multi-stability

Citation: Li X, Jolly MK, George JT, Pienta KJ and Levine H (2019) Computational Modeling of the Crosstalk Between Macrophage Polarization and Tumor Cell Plasticity in the Tumor Microenvironment. *Front. Oncol.* 9:10. doi: 10.3389/fonc.2019.00010

Received: 03 October 2018; Accepted: 03 January 2019;

Published: 23 January 2019.

Edited by:

Ubaldo Emilio Martinez-Outschoorn, Thomas Jefferson University, United StatesReviewed by:

Rafael Coelho Lopes De Sa, University of Massachusetts Amherst, United StatesVirendra K. Chaudhri, Cincinnati Children's Hospital Medical Center, United States

Copyright © 2019 Li, Jolly, George, Pienta and Levine. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Herbert Levine, herbert.levine@rice.edu

^{†}These authors have contributed equally to this work

^{‡}Present Address: Mohit Kumar Jolly, Centre for BioSystems Science and Engineering, Indian Institute of Science, Bangalore, India