Toward Decoding the Principles of Cancer Metastasis Circuits

Understanding cell fate decisions during embryonic development and tumorigenesis remain a major research challenge in modern developmental and cancer biology (1). In recent years, we have witnessed rapid progress in mapping the genetic regulatory networks that determine the fate of different cells. Examples include networks that govern the transition between epithelial and mesenchymal phenotypes, networks that control the differentiation of pluripotent stem cells into different lineages and/or progenitor cells, networks that are involved in the transition into and from cancer stem–like cells (CSC), as well as networks that are involved in cellular dedifferentiation during the formation of induced pluripotent stem cells (iPSC; refs. 2–5). In all these examples, a cell's fate is orchestrated by changes in the expression of transcription factors (TF) and microRNAs (miRNA; miR) that in turn govern downstream regulatory networks, ultimately leading to the genome-wide gene expression patterns and corresponding protein levels specific to a particular cell lineage (fate).

An important representative example of cell fate decision is seen during metazoan embryonic development, when sessile epithelial cells reversibly attain mesenchymal-like characteristics that allow them to migrate and invade adjacent tissues. Successive rounds of these forward and backward transitions, namely epithelial–mesenchymal transition (EMT) and mesenchymal–epithelial transition (MET), play a crucial role in forming many internal organs. Cancer cells in various carcinomas (including lung, breast, prostate, colon, pancreas, and ovaries) adopt this embryonic process of EMT and MET during invasion and metastasis. This aberrant activation of EMT by carcinomas is considered to be an important hallmark of cancer metastasis, an outcome that is responsible for more than 90% of cancer deaths (6, 7). The core decision network that regulates carcinoma EMT and MET also allows for transition into a hybrid epithelial–mesenchymal (E/M) phenotype, which has combined epithelial (cell–cell adhesion) and mesenchymal (motility) traits. These mixed characteristics of the hybrid E/M phenotype enable cells to migrate collectively, as compared with the single-cell migration of purely mesenchymal cells (8). The hybrid E/M state, also referred to as partial EMT (pEMT), can revert to the epithelial phenotype when required, as seen during the final stages of branching morphogenesis and wound healing (9–12). More recently, this hybrid state has been observed at different stages of epiboly (13), organ morphogenesis (14), and metastasis (15–17).

Carcinoma metastases typically begin when epithelial cells from the primary tumor lose their apical-basal polarity and cell–cell adhesion and acquire migratory and invasive mesenchymal traits. The newly transformed motile cells navigate through the extracellular matrix (ECM) toward blood vessels. It has been shown that cancer cells can migrate as dispersed individuals with mesenchymal character or in a coordinated collective motion of cells with a hybrid character (10, 11, 17). Collective migration obviates the need for all cells to be able to detect extrinsic signals for migration (12), enabling the cohort of cells to adeptly adapt to different microenvironments. Thus, hybrid E/M cells have an important functional role for successful migration through the ECM. After reaching blood vessels, some of the cells succeed in penetrating into the bloodstream (intravasation) and stay as circulating tumor cells (CTC) until they exit (extravasation) to reach an appropriate colonization site in a distant organ. Subsequently, the cells undergo the MET and regain their epithelial characteristics, growing later into macrometastases (6, 16).

Notably, while hybrid E/M cells have been observed during forward EMT, they have not been observed during the reverse, MET (18, 19). This implies an innate asymmetry in the underlying mechanisms regulating the epithelial and mesenchymal cell fate. Also, many CTCs coexpress both epithelial and mesenchymal markers (17), hence demonstrating the significance of cellular plasticity in metastasis and colonization. The tumor cells that undergo EMT are usually resistant to both chemotherapy and radiotherapy (20, 21), and have properties similar to that of CSCs (22–24). Thus, understanding the backward and forward transitions between these three cell phenotypes [epithelial (E), mesenchymal (M), and E/M] could inform strategies against metastases.

Cell fate determination between these three phenotypes is regulated by many internal and external signals such as hypoxia-inducible factor 1 (HIF1), p53, TGFβ, HGF, FGF, EGF, Notch, and Wnt (25). These signals converge on a core regulatory unit composed of four components, two families of transcription factors that induce EMT, SNAIL and ZEB, and two families miRNAs that inhibit EMT, miR34, and miR200 (5). The epithelial phenotype corresponds to high levels of miR34 and miR200, whereas the mesenchymal phenotype corresponds to high levels of SNAIL and ZEB. SNAIL and ZEB repress epithelial-specific gene expression including E-cadherin, the hallmark of the epithelial phenotype, and promote the expression of mesenchymal markers such as N-cadherin and vimentin (26). The EMT-inducing signals (HIF1, TGFβ, HGF, FGF, EGF, Notch, and Wnt) activate SNAIL and ZEB (24, 26, 27, 28), whereas EMT-repressing signals (p53) activate miR34 and miR200 (29). The input signals have elaborate mutual interactions such as the degradation of p53 via TGFβ-mediated activation of MDM2 (Fig. 1; ref. 30).

The core decision unit is also known to play an important role as “a motor of cellular plasticity” in various human carcinomas, as it is coupled to a variety of other key cellular features, including stemness, cell-cycle arrest, apoptosis and senescence, resistance to chemotherapy, and cell–cell communication, as discussed in the last section (31, 32).

From a clinical perspective, SNAIL and ZEB are overexpressed at the tumor–stroma interface and correlate significantly with many interrelated features including overall survival, poor prognosis, as well as tumor subtypes and grades associated with worse outcomes (6, 26). These observations reinforce the importance of EMT in carcinoma metastasis, as SNAIL and ZEB are the major drivers of EMT.

The four components of the core regulatory unit form two highly interconnected modules—the miR34/SNAIL and the miR200/ZEB modules (Fig. 1). Each of these modules is a double-negative or a mutually inhibiting miRNA–TF chimeric feedback loop (33–35), also referred to as a chimera switch (36). In the miR34/SNAIL module, miR34 (μ₃₄) binds to two conserved sites on the 3′-untranslated region (UTR) of the SNAIL (S) mRNA, whereas SNAIL represses miR34 transcriptionally through binding at one site in its promoter region (35). In the miR200/ZEB module, the miR200 family (μ₂₀₀)—mir-141, miR200a/b/c, and miR429—has a total of eight binding sites in the ZEB1 3′-UTR and nine in the ZEB2 3′-UTR, whereas ZEB1 and ZEB2 (considered together as the ZEB family) bind to the miR200 family promoter region (33, 34). Because stable expression of miR200c alone is sufficient to reverse EMT and restore E-cadherin expression (37), in Fig. 1, we considered six binding sites of the miR200 (μ₂₀₀) family to ZEB (Z) mRNA, the number of miR200c binding sites on the ZEB family, and three binding sites for ZEB on the miR200 promoter, the number of ZEB binding sites on miR200c. SNAIL activates ZEB (38) and is transcriptionally self-inhibiting (39). ZEB has two potential binding sites in its promoter (40) and activates itself indirectly by stabilizing the SMAD complexes (41), thus we represent ZEB as a transcriptional self-activating gene of rank 2 (two sites in its promoter region). It may be noted that these interactions are not specific for a particular cancer, but rather hold true for most human carcinomas undergoing EMT (see Supplementary Data for a list of these interactions as identified in various carcinoma cell lines).

Aimed for readers from diverse disciplines, we include in this section more detailed description of modeling miRNA-based circuits. Understanding the mathematical details is not essential for understanding the conclusions reached and their implications about the crucial role of epithelial plasticity, a conserved phenomenon among carcinomas.

Because miRNAs are key components of the EMT regulatory circuits, we developed a theoretical framework for modeling miRNA-based chimera (MBC) circuits (36). This framework captures the essential features of miRNA-mediated regulation by considering the binding/unbinding chemical reactions for both TF–promoter and miRNA–mRNA complexes (Fig. 2). In the model, when one or multiple miRNAs bind to mRNA, the miRNAs are able to silence mRNA translation by inhibiting the translation process (translation inhibition) and/or actively degrading the mRNA-containing complex (Fig. 2A). Consider the case in which miRNA (μ) molecules target the mRNA (m) of a protein (B), and there are n distinct miRNA binding sites on the mRNA. In this case, the deterministic equations are given by

where g_m and g_μ are the synthesis rates of m and μ, respectively, and these could depend on the concentration of some external signals. k_m, k_B, and k_μ are the innate degradation rates of m, B, and μ, respectively, and g_B is the innate mRNA's translation rate for protein B. Equation 1 also contains three μ-dependent functions to quantify the miRNA–mRNA coupling, that is, the translational inhibition term (L), the mRNA active degradation term (Y_m), and the miRNA active degradation term (Y_μ). The formulae and their derivation can be found in the Supplementary Data of ref. 36. As shown in Fig. 2, as μ increases, L is markedly repressed (Fig. 2B), while Y_μ increases slightly and Y_m increases markedly (Fig. 2C). Compared with the previously derived models for miRNA–TF chimeric circuits (42–50), the new approach is more consistent with the biologic mechanism. Because miRNAs are typically more stable than mRNAs (49), we sometimes can reduce Equation 1 to two coupled equations for B and μ by assuming that m is always at steady states (36).

Applying the theoretical framework for miRNA-based circuits, we studied the dynamics of the two miRNA–TF chimera toggle switches in the core regulatory unit. Typically, the operating characteristics of circuit modules depend on the nonlinearity of the inhibitory regulations by TFs/miRNAs and the nonlinearity of the auto-regulations of the TFs (Fig. 1). To better understand the relationship between nonlinearity and multistability, we examined the behavior of several similar circuit modules.

We started with the miR200/ZEB module, which is a typical miRNA–TF chimera switch (Fig. 3A). To simplify the problem, we first omitted the ZEB self-activation in the initial stage. The EMT transcription factor SNAIL transcriptionally regulates both miR200 (inhibition) and ZEB (activation). Considering SNAIL as an external control input signal S, we obtained the following deterministic equations for miR200 (μ₂₀₀), ZEB mRNA (m_Z), and ZEB protein (Z) as

Here, the shifted Hill function, defined as |$H^S (X,\lambda) = H^ - (X) \,+\, \lambda H^ + (X)$|⁠, where |$H^ - (X) = 1/[1 \, +\, (X/X_0)^{n_X} ]$|⁠, |$H^ + (X) = 1 - H^ - (X)$|⁠, and n_X is the Hill rank for X (which is usually associated with the number of binding sites of TF X on the promoter); λ is a positive number, which quantifies the “fold change” of the synthesis rate caused by X. According to this definition (51), the activation corresponds to λ > 1, while the inhibition corresponds to λ < 1.

From Fig. 3A, the miR200/ZEB switch is seen to be bistable (two-way switch) for a certain set of biologically relevant parameters (see also Definitions box). We also tested the multistability of the system by adjusting a wide range of parameters. The system can be monostable for some parameters, but not tristable (three-way switch).

Given the experimental knowledge about the miR200/ZEB module, ZEB mRNA has six miRNA binding sites for the miR200 family (including miR141, miR200a/b/c, and miR429; see Fig. 1), so this translation inhibition is highly nonlinear. Similarly, ZEB transcriptionally inhibits miR200 with high nonlinearity (Hill rank is 3). Such high nonlinearity renders the circuit to have multiple steady states.

We further examined the effect of the nonlinearity (the value of the Hill rank) on a generic chimera toggle switch. We found that (see Supplementary Data in ref. 36), in general, when there is one miRNA binding site, the system is monostable, even for a very high rank of transcriptional inhibition. For the case of two miRNA binding sites, the system is monostable when the Hill rank is less than or equals to three for the TF transcription inhibition and is bistable when the Hill rank is 4 or more. The system is easily bistable (even for low Hill rank of transcriptional inhibition) when there are more miRNA binding sites.

Next, we considered the miR200/ZEB module with the ZEB self-activation (Fig. 3B). The deterministic equations for miR200 (μ₂₀₀), ZEB mRNA (m_Z), and ZEB protein (Z) are

As mentioned earlier, the module can be bistable when the ZEB self-activation is omitted (Fig. 3A). Typically, self-activation on both sides of the toggle switch can render the circuit tristable (ternary or three-way) switch (52–55). However, self-activation on one side can also render the circuit tristable for a wide range of parameters, as seen in both TF–TF toggle switches and miRNA–TF chimera toggle switches (36). This is consistent with our modeling on the miR200/ZEB module (Fig. 3B).

A typical phase space diagram is shown in Fig. 3B. The circuit is seen to have three coexisting stable states (three filled circles at points of intersection of the curves), which correspond to (i) high miR200 and low ZEB [denoted as (1, 0) and furthest to the right], (ii) medium miR200 and medium ZEB [denoted as (½, ½) and in the middle], and (iii) low miR200 and high ZEB [denoted as (0, 1) and furthest to the left]. The states (1, 0) and (0, 1) correspond to the E and M phenotypes, respectively (33, 34). We suggest that the intermediate (½, ½) state should be associated with the hybrid E/M phenotype. Thus, the miR200/ZEB module serves as the three-way decision switch that enables the cells to adopt a hybrid E/M state or a completely epithelial or mesenchymal state (51). The Hill rank for ZEB self-activation was chosen to be 2 for the reasons described earlier. It was shown that both roles of miRNA—translational inhibition and active mRNA degradation (see Fig. 1)—are required for tristability and for the ZEB/miR200 module to operate as a three-way switch (51).

The miR34/SNAIL module (Fig. 3C) represents a chimera toggle switch with TF self-inhibition (in this case, SNAIL). In the model, the module is driven by an external input signal I to capture the effect of various signals that target the gene (e.g., TGFβ and HIF1). The deterministic equations for miR34 (μ₃₄), SNAIL mRNA (m_S), and SNAIL protein (S) are

A typical phase space diagram is shown in Fig. 3C. As seen, there exists only one stable steady state, which corresponds to fixed levels of SNAIL and miR34. When the signal I increases, the SNAIL level increases smoothly and the miR34 level decreases smoothly.

There are only two miR34 binding sites on SNAIL mRNA and ZEB inhibits miR34 transcriptionally with weak nonlinearity (Hill rank is 1; ref. 35). Thus, the modeling results are consistent with the findings mentioned above, where we considered the chimera toggle switch without self-inhibition. Further investigations showed that self-inhibition usually makes the system even less likely to have multistability. Moreover, the self-inhibition reduces the effect of external noise in various EMT-inducing signals and determines the sensitivity threshold to those signals (56). The miR34/SNAIL module serves as a noise-buffer signal integrator (51). Such noise-buffering feature and relatively weak dependency of the stable sate enables the module to be more reliable and prevents erroneous activation of EMT by some transient signals. The proposed role of miR34/SNAIL also explains why the epithelial phenotype is stable (57).

Toward characterizing the combined system of the two chimera modules here, we analyzed the complete core regulatory unit (Fig. 3D), which includes the decision module (miR200/ZEB) and the integrator module (miR34/SNAIL). The combined circuit is also driven by an external signal I that activates SNAIL transcriptionally. The deterministic equations for miR200 (μ₂₀₀), ZEB mRNA (m_Z), ZEB protein (Z), miR34 (μ₃₄), SNAIL mRNA (m_S), and SNAIL protein (S) are

As shown in Fig. 3D, the combined core decision unit acts as a three-way switch (has three stable states). In Fig. 4, bifurcation curves were plotted for the cases with and without feedback from ZEB to miR34. When the ZEB mRNA levels of the steady states were plotted with respect to the SNAIL protein levels, as shown in Fig. 4A, the feedback did not change the bifurcation curve. When they were plotted with respect to the external signal levels I, as shown in Fig. 4B, the feedback (inhibition of miR34 by ZEB) shifted the bifurcation curve to the left, especially for high levels of mRNA. Thus, the inhibitory feedback plays important roles in amplifying the symmetry breaking between the forward (EMT) and backward (MET) transitions. During EMT, cells are more likely to go through the hybrid (E/M) phenotype, while during MET, cells are likely to go directly from M to E, without attaining the hybrid state en route. This result is consistent with the experimental observations that the hybrid phenotype has not yet been observed during MET, for example, in reprogramming to iPSCs (18, 19), but frequently observed during EMT, for example, tumor invasion and metastasis (15–17).

A gene regulatory network, although complicated, usually consists of many circuit modules that interact with each other. In the current study, we found that the core regulatory circuit of EMT contains the miR34/SNAIL module, which serves as a noise-buffer signal integrator, and the miR200/ZEB module, which acts as a three-way decision circuit or switch. In the previous sections, we have shown how we first analyzed the behavior of each individual circuit module and then characterized the whole circuit by adding the coupling between the two modules. There are a total of three regulatory links between miR34/SNAIL and miR200/ZEB, that is, inhibition of miR200 by SNAIL, activation of ZEB by SNAIL, and inhibition of miR34 by ZEB (Fig. 1). To simplify the problem, we included the regulation of miR200 and ZEB by SNAIL into the circuit module of miR200/ZEB, and regarded SNAIL as an external signal of this module. By doing this, we reduced the number of links between two modules to one without sacrificing the accuracy in describing the full circuit. From the analysis mentioned in the previous section, we found that the feedback from ZEB to miR34 does not dramatically affect the functions and stand-alone dynamics of the two modules. Instead, the feedback adds more asymmetry between forward (EMT) and backward (MET) transitions. We also found that the regulation of the miR200/ZEB module on SNAIL is unaffected by the feedback, which validates the decomposition approach we adopted.

We propose the following multistep approach to investigate more elaborated networks: (i) identify and study the “stand-alone” dynamics of the basic modules; (ii) formulate/devise “solvability” conditions of the mutual feedbacks (constraints) between the modules when functioning as a combined unit. More specifically, first, the combined network is decomposed into different circuit modules in such a way that each module only has a few elements (TFs or miRNAs) linked with the other modules. The most straightforward way (although not necessarily the most efficient) is to simply have two elements in each module. Second, each module is analyzed by the so-called two-signal bifurcation technique, which provides information about the multistability with respect to various possible combinations of the two signals. Third, the coupling among different modules is included by solving proper solvability conditions, which ensure consistency in the dynamics of the various modules. By consistency, we mean that the resultant levels of expression of the elements (TFs and miRNAs) are consistent with the levels of expression obtained from the stand-alone dynamics of the modules. Here, we elaborate below on the two-signal bifurcation for the specific example of the miR200/ZEB module.

The idea is illustrated in Fig. 5 in which we show the phenotypic phase diagram for the miR200/ZEB module when driven by two input signals representing the action of SNAIL in the combined circuit. Therefore, miR200 is driven by an inhibitory signal S₁ and ZEB is driven by an activator signal S₂. The resultant phase diagram (shown in Fig. 5A) shows the existence of seven different phases. From a dynamical system perspective, each of the phases corresponds to a different nullcline describing monostability or multistability (coexistence) of different states (nullclines of three different phases are shown in Fig. 5B–D).

Using our framework for miRNA-based chimeric circuits, we decipher the operation principles of the EMT regulatory network (Fig. 3). First, the miR34/SNAIL acts as a noise-buffer signal integrator that prevents aberrant activation of EMT due to transient signals and explains the stability of epithelial phenotype (57). Second, the miR200/ZEB acts as a three-way (ternary) switch, which explains the existence of three different phenotypes—the canonical E and M ones, and the more recently discovered hybrid phenotype (E/M).

Our theoretical predictions are consistent with experimental observations. For example, the activation of SNAIL can initiate EMT by repressing CDH1 (the gene for E-cadherin), but ZEB1 is required for complete inhibition of E-cadherin (56), and hence the completion of EMT. Similarly, it has been observed that a complete reversal to an epithelial phenotype requires a strong inhibition of ZEB1 (58), as the knockdown of SNAIL is not sufficient (59). Furthermore, the cells that attain high ZEB levels, for example, by being continuously treated with the signal TGFβ, do not immediately revert to an epithelial phenotype when the signal is removed. However, the cells with remarkably low levels of ZEB do revert, indicating that miR200/ZEB is the module that acts as the commitment point for cells undergoing EMT (27). Finally, ZEB1 transcriptionally inhibits most of the genes that are downregulated during EMT (60), suggesting that ZEB is the master regulator for cell fate decision-making during EMT.

Recent studies of gastrulation in Drosophila embryos (61) show that collectively migrating cells, the hallmark of the hybrid E/M state, coexpress ZEB1 and E-cadherin, thus validating that the hybrid state has intermediate levels of both miR200 and ZEB. Also, the (½, ½) state for the hybrid phenotype allows for cell–cell communication via the Jag1–Notch–Delta system (32), as observed during collective migration in wound healing (62). Because Jag1 mRNA has five binding sites for miR200 and is thus strongly repressed by miR200 (32), its expression in the hybrid state is possible only when the hybrid state has either low or intermediate levels of miR200, instead of the high levels that are characteristic of epithelial cells. This implies that the level of expression of miR200 in the hybrid phenotype of collectively migrating cells that maintain the Notch–Jag1 signaling has to be lower than its high level of expression in the epithelial phenotype. Thus, the results validate our hypothesis that the observed hybrid phenotype corresponds to the theoretical predicted (½, ½) state.

Of note, after the first version of this manuscript was submitted, a phenotypic study of 43 well-characterized ovarian carcinoma cell lines identified an EMT spectrum among them. Twenty-six of the 43 cell lines had a hybrid phenotype, were highly aggressive, resistant to anoikis, had enhanced sphere-forming ability, and expressed both ZEB1 and E-cadherin, thus indicating collective cell migration (63). This suggests that most carcinoma cells undergo partial or incomplete EMT, and a complete EMT may not be necessary for cell survival during metastasis. Thus, characterizing the hybrid state is critical for developing therapeutic targets to reverse EMT in a selective group of patients.

We use the notation of “0,” “1,” and “½” to denote the expression levels of miR200 and ZEB in the three different phenotypes. It may be noted that “0,” “1,” or “½” expression levels of a given element in one cell line or context may be different from those in another cell line or context, due to phenotypic heterogeneity and nonheritable variability pertaining to different cell lines. Thus, the absolute expression levels of many elements in the EMT regulatory network need to be measured quantitatively. Our modeling approach, which incorporates a detailed analysis based on biologically realistic parameters (see Supplementary Data in ref. 51 for details on range of parameters), might be well suited to capture this phenotypic variability, as opposed to a simplistic Boolean framework.

We also discuss the possibility of coexistence of different phenotypes among cells, through phase plane analysis (Fig. 5). From a biological perspective, each phase corresponds to (the coexistence of) different phenotypes. In three out of seven phases, only one of the three phenotypes (epithelial, mesenchymal, or hybrid) can exist (denoted as {E}, {M}, and {E/M}). In three other phases, two of the phenotypes can coexist (denoted as {E, M}, {E, E/M}, and {E/M, M}) and in one phase, we see the coexistence of all three phenotypes (denoted as {E, E/M, M}).

The abundant information from these identified phases can also help to elucidate the experimental data on the diversity of behaviors associated with SNAIL in triggering EMT in different contexts. SNAIL can induce complete EMT (64) or partial or transient EMT (28, 65) and is also required for re-epithelialization of keratinocytes after wound healing (66). Depending on the rate of induction of SNAIL by different signals (e.g., by TFGβ or hypoxia, both of which also regulate ZEB), the cells would undergo a different trajectory in the phase diagram and therefore exhibit different phenotypic transitions. Also, the largest area in the phase diagram is the monostable phase with only epithelial phenotype, which agrees with the observation that epithelial type is the default phenotype (57).

In Fig. 6, we show examples of how the EMT decision unit is coupled to other key cellular processes including programmability, cellular motility, genome plasticity, and metabolism.

miR200 inhibits the oncogenic driver LIN28 (67), which forms a double-negative feedback loop with let-7 (68), which is another family of miRNAs that inhibit EMT (69). The LIN28/let7 loop has some indirect auto-regulations (70–72) that target the pluripotency factor OCT4 (73), which forms a complex with SOX2 and activates itself indirectly, at least in human embryonic stem cells (hESC; ref. 3). Also, OCT4 and SOX2 activate miR200 during MET while reprogramming of fibroblasts to iPSCs (74). In addition, let-7 also mediates the double-negative feedback loop between metastasis repressor RKIP and metastasis-promoting self-inhibiting gene BACH1. RKIP activates let-7 indirectly, which inhibits BACH1, and BACH1 inhibits RKIP directly. A similar mutual repression also exists between RKIP and SNAI, again via let-7 (75, 76).

miR200 and miR34 are coupled to the regulatory network of Rac1 and RhoA (77, 78), the two auto-regulatory GTPases that determine the mode of cancer cell invasion (79). High levels of active Rac1 lead to the formation of cell protrusions and a mesenchymal phenotype, which degrade the ECM and are “path generators,” while high levels of active RhoA lead to actomyosin contractility and amoeboid mode of invasion, which squeeze through the gaps in ECM and are referred to as “path finders” (79, 80). miR200 family members differentially regulate these two major modes of cell invasion by reducing both the formation of cell protrusions and actomyosin contractility (78).

miR200 forms a chimera toggle switch with SIRT1 (81), a key player in linking cell metabolism to stress response. SIRT1 is inhibited by itself and miR34 as well (82), and regulates hypoxia response through HIF1 and HIF2 (83, 84). Besides, SIRT1 maintains proper chromatin structure, thus being crucial for genome integrity and DNA damage response (85).

The miR200/ZEB loop links hypoxia to EMT as HIF1 and HIF2 induce SNAIL and ZEB, respectively (28, 86). The most intriguing coupling between hypoxia, metabolism, and EMT is through reactive oxygen species (ROS), which not only mediates the mutual inhibition between HIF1 and HIF2 (87), but also activates both miR200 and ZEB (88, 89). Because HIF1 and HIF2 are believed to govern the response to acute and chronic hypoxia, respectively (90, 91), this coupling can unravel how epithelial–mesenchymal plasticity is linked to changes in tumor metabolism and tumor angiogenesis, as induced by hypoxia.

In addition, the role of hypoxic niches in maintaining CSCs in glioblastoma through the Notch signaling pathway (92, 93) underline the interplay between hypoxia, cell–cell communication, and stemness. miR200 inhibits EMT and metastasis by altering the tumor microenvironment (94) and blocking tumor angiogenesis (95). Furthermore, ZEB1 forms feedback loops with transcription factors GRHL2 (96), OVOL1, and OVOL2 (97) to regulate EMT/MET during metastasis. We hypothesize that the hybrid phenotype is associated with transition into CSC. The transition can start during migration toward the blood vessels, during circulation, and also after seeding in the new niche.

miR203 forms a chimera toggle switch with SNAIL and has similar interconnections with the miR200/ZEB loop as that for the miR34/SNAIL loop (98), thus indicating that both miR203 and miR34 form noise-buffering integrators with SNAIL.

Epithelial–mesenchymal plasticity is crucial during embryonic development and cancer metastasis (6). The core regulatory network for these transitions is the miR200/ZEB chimera toggle switch coupled with the miR34/SNAIL loop. Also, the miR200/ZEB decision circuit is coupled to many key cellular properties such as stemness, cell motility, cell–cell communication, metabolism, and resistance to apoptosis (31, 32, 99). Despite its widespread influence, this decision circuit has had limited theoretical attention.

Here, we studied the involvement of miRNAs in EMT decision-making circuit, using the recent theoretical framework devised by Lu and colleagues (36). Our approach incorporates in detail both modes of translational silencing by miRNAs and includes the effects of number of binding sites of miRNA on mRNA, which have not been considered in previous studies (42–50). Using this framework, we unraveled the modular design principles of the core regulatory network for EMT/MET. We found that the miR34/SNAIL module is monostable and acts as a noise-buffering integrator, whereas miR200/ZEB module is tristable and functions as a three-way switch. The three stable states of miR200/ZEB correspond to the epithelial [high miR200/low ZEB, denoted as (1, 0)], mesenchymal [low miR200/high ZEB, denoted as (0, 1)], and hybrid E/M [medium miR200/medium ZEB, denoted as (½, ½)] phenotypes.

We elaborated on the above to clarify that our hypothesis for the (½, ½) hybrid state (medium miR200/medium ZEB) is different from that proposed by another recent study on EMT circuitry modeling, which neither included the self-regulations of SNAIL and ZEB nor distinguished between the transcription and translation inhibition processes (they were both modeled by inhibitory Hill functions of rank 2; ref. 45). Under these assumptions, the two modules have similar bistable dynamics and act as binary switches, and the authors proposed high miR200/low miR34 to be the hybrid state. The experimental data mentioned above is more consistent with our medium miR200, medium ZEB hypothesis. However, it might be possible that different cell lines exhibit different hybrid phenotypes, as multiple hybrid states have been proposed to be en route during EMT in certain contexts.

In Fig. 6, we show examples of how the EMT decision unit is coupled to other key cellular processes. Future theoretical investigations into these circuits hold the key to valuable new insights about how these cellular characteristics are modified during backward and forward transitions between epithelial and mesenchymal cells.

A better understanding of the transitions involving the hybrid phenotype is essential for a better comprehension of cancer progression. Also, it can help answer the long-standing fundamental question—the difference between EMT that happens during embryonic development and wound healing and EMT that happens during metastasis and organ fibrosis. From the perspective of developing better therapeutics, we require normal wound healing to continue during antimetastatic therapies and also at the same time, restrict wound healing augmentation from promoting metastasis in case a malignancy is present (100).

No potential conflicts of interest were disclosed.

The authors have benefited from useful discussions with Herbert Levine, Ilan Tsarfaty, and M. Cindy Farach-Carson.

This work was supported by the NSF Center for Theoretical Biological Physics (NSF grant no. PHY-1308264) and by the Cancer Prevention and Research Institute of Texas (CPRIT). E. Ben-Jacob was also supported by a grant from the Tauber Family Funds and the Maguy-Glass Chair in Physics of Complex Systems.