Epigenetic feedback and stochastic partitioning during cell division can drive resistance to EMT

Epithelial-mesenchymal transition (EMT) and its reverse process mesenchymal-epithelial transition (MET) are central to metastatic aggressiveness and therapy resistance in solid tumors. While molecular determinants of both processes have been extensively characterized, the heterogeneity in the response of tumor cells to EMT and MET inducers has come into focus recently, and has been implicated in the failure of anti-cancer therapies. Recent experimental studies have shown that some cells can undergo an irreversible EMT depending on the duration of exposure to EMT-inducing signals. While the irreversibility of MET, or equivalently, resistance to EMT, has not been studied in as much detail, evidence supporting such behavior is slowly emerging. Here, we identify two possible mechanisms that can underlie resistance of cells to undergo EMT: epigenetic feedback in ZEB1/GRHL2 feedback loop and stochastic partitioning of biomolecules during cell division. Identifying the ZEB1/GRHL2 axis as a key determinant of epithelial-mesenchymal plasticity across many cancer types, we use mechanistic mathematical models to show how GRHL2 can be involved in both the abovementioned processes, thus driving an irreversible MET. Our study highlights how an isogenic population may contain subpopulation with varying degrees of susceptibility or resistance to EMT, and proposes a next set of questions for detailed experimental studies characterizing the irreversibility of MET/resistance to EMT.


INTRODUCTION
Epithelial-Mesenchymal Transition (EMT) is a cell biological process involved in driving cancer metastasis and therapy resistance-the two grand clinically unsolved challenges. EMT and its reverse process Mesenchymal-Epithelial Transition (MET) are believed to enable cancer cell dissemination from the primary tumor, facilitate survival in the bloodstream, and are implicated in extravasation and the formation of macrometastases at multiple distant organs [1]. Thus, understanding the dynamics of EMT and MET is essential to develop novel therapeutic interventions.
Recent studies have highlighted that EMT and MET are not binary processes as thought earlier [1]. Instead, besides the epithelial and mesenchymal phenotypes, cells can acquire and stably maintain one or more hybrid epithelial/mesenchymal (E/M) phenotypes. These hybrid E/M phenotypes may drive collective cell migration as clusters of tumor cells (CTCs) and can be more aggressive

Research Paper
Oncotarget 2612 www.oncotarget.com than cells in pure epithelial or mesenchymal phenotypes [2]. Importantly, tumor cells may switch among different phenotypes: E, M, and hybrid E/M [3][4][5]. Such dynamic and reversible switching can help tumor cells to overcome various challenges during disease progression such as anoikis [6], and assaults by the immune system [7]. Thus, epithelial-mesenchymal plasticity (EMP)-a combination of EMT and MET-needs to be utilized with spatiotemporal precision to drive metastasis. For example, a failure to undergo MET at a metastatic site may compromise colonization [8,9]. The reversibility of EMT and MET, mediated by multiple interconnected feedback loops regulating a balance between epithelial and mesenchymal traits, thus forms the backbone of metastasis [10].
Are EMT and MET always reversible? Recent experiments decoding the dynamics of EMT/MET using live-cell imaging and/or induction and withdrawal of various EMT-inducing external signals such as TGFβ or tuning the levels of EMT-specific transcription factors (EMT-TFs) have provided important insights into the reversibility of EMT and MET. Cells induced to undergo EMT for shorter durations (~2-6 days) may revert to an epithelial state after withdrawal of the signal/stimulus. However, some cells exposed to EMT-inducing signals for longer durations (~10 days or more) may get "locked" in a mesenchymal state, making EMT largely irreversible, at least for the timescale observed experimentally [11][12][13][14]. The possibility of an irreversible EMT is also supported by multiple phenomenological observations [15][16][17]. Multiple mechanisms have been proposed to explain the existence of a "tipping point"-a time point beyond which continued treatment with EMT inducing signals can drive an irreversible EMT. These include self-stabilizing feedback loops [18][19][20][21] in the regulatory circuits for EMT/ MET and/or epigenetic alterations [13,22,23]. However, similar investigations about the irreversibility of MET, or in other words, the resistance of epithelial cells to undergo EMT in response to EMT-inducing signals, remain to be done. Some sporadic observations about the resistance of epithelial cells to undergo EMT have been reported [14,24], but a causative mechanistic understanding still remains elusive.
Here, we propose two independent mechanism that may explain the resistance of epithelial tumor cells to undergo EMT: 1) epigenetic feedback mediated via GRHL2-an MET-inducing transcription factor (MET-TF) [25][26][27]; and 2) stochastic partitioning of parent cell biomolecules among the daughter cells at the time of cell division [28][29][30]. GRHL2 and miR-200 both form mutually inhibitory feedback loops with ZEB1-a key EMT-inducing transcription factor (EMT-TF) [31]. Previously, we have shown that incorporating an epigenetic feedback term acting on the inhibition of miR-200 by ZEB1 could drive an irreversible EMT [13]. This epigenetic feedback term was incorporated at a phenomenological level to represent the idea that the longer a gene is turned on, the easier it becomes for it to stay transcriptionally active; thus, epigenetic feedback can modulate the thresholds for the influence of a transcription factor on its downstream target [32,33]. Conversely, here, we show that incorporating this epigenetic feedback loop acting on the inhibition of ZEB1 by GRHL2 can cause an irreversible MET. Cells undergoing irreversible MET may exhibit resistance in undergoing EMT when exposed to EMT-inducing signal. Thus, our results offer a conceptual framework to decode the impact of various epigenetic mechanisms in terms of modulating the reversibility of EMT and MET in a cancer cell population. Further, our previous analysis illustrated how epithelial-mesenchymal heterogeneity can be generated from a phenotypically homogeneous population by stochastic partitioning of molecules during cell division [3]. Here, we demonstrate how this stochasticity in the presence of GRHL2 can lead to an irreversible MET or, in other words, a resistance to undergo EMT. Together, our results describe how an isogenic cellular population may contain these different subpopulations with varying degrees of susceptibility and resistance to undergo EMT in response to an EMTinducing signal.

RESULTS
GRHL2 correlates with a more epithelial phenotype across many cancer types GRHL2 has been identified as an MET inducer in breast cancer [25][26][27], where it forms a mutually inhibitory feedback loop with ZEB1, an EMT-TF. Overexpression of GRHL2 suppresses EMT induced by TGF-β or Twist by directly binding to the ZEB1 promoter, and inhibits various other properties associated with a partial or complete EMT such as higher mammosphere-forming efficiency and anoikis resistance [6,25,34,35]. In general, the GRHL2/ ZEB1 feedback loop was identified as a key regulator of EMP and associated traits in breast cancer, lung cancer [36], colorectal cancer [37] and ovarian cancer [38,39]. Consistently, GRHL2 was shown to inhibit EMT in gastric cancer [40], oral cancer [41] and pancreatic cancer [42]. These reports drove us to investigate the correlation of levels of GRHL2 with EMT/MET across many cancer types both in the Cancer Cell Line Encyclopedia (CCLE) cohort [43] and in many TCGA datasets.
GRHL2 levels correlated positively with CDH1 (E-cadherin) levels and negatively with ZEB1 in the CCLE dataset and TCGA datasets from breast cancer, ovarian cancer and colorectal cancer ( Figure 1). Given that GRHL2 is one of the top transcriptional activators of CDH1 and ZEB1 is one of its strongest transcriptional repressors, ZEB1 and CDH1 correlated negatively (Supplementary Figure 1). We also investigated the correlation of GRHL2 levels with three transcriptomics-based EMT scoring www.oncotarget.com algorithms: MLR, KS, and 76GS. While MLR and KS methods assign higher scores to mesenchymal samples, the 76GS method assigns higher scores to epithelial samples [44]. As expected, GRHL2 levels correlated positively with EMT scores calculated via the 76GS method and negatively with EMT scores calculated by MLR or KS methods across these TCGA datasets (Supplementary Figure 1). Consistently, a constitutive expression of GRHL2 in an inducible EMT model system-HMLE cells that contain a Twist-ER fusionled to reduced ZEB1 levels and corresponding changes in EMT scores as identified by all three abovementioned EMT scoring metrics (Supplementary Figure 2).
Next, in the CCLE dataset, we calculated the pairwise correlations of various canonical epithelial and mesenchymal markers and regulators. We found that GRHL2 correlates positively with its family members GRHL1 and GRHL3, its downstream target OVOL2 and corresponding family member OVOL1, and with CDH1, while it correlates negatively with TWIST1/2, SNAI1, VIM, and ZEB1/2 ( Figure 2A, Supplementary Figure  3A). On the other hand, ZEB1 correlates negatively with GRHL1/2/3, OVOL1/2 and positively with SNAI1/2, TWIST1/2 and VIM ( Figure 2A). With these consistent observations regarding the antagonistic roles of ZEB1 and GRHL2 in regulating EMT/MET, we next calculated the correlation of all genes in CCLE with GRHL2 and with ZEB1, and observed that most genes that showed significant correlation with both of them were either positively correlated with GRHL2 and negatively with ZEB1 or vice-versa. The relatively smaller set of genes that correlated either positively or negatively with both ZEB1 and GRHL2 showed relatively weak correlations (R < 0.3) ( Figure 2B, Supplementary Figure 3B). Together, these observations indicate that GRHL2 associates with epithelial traits across cancer types.

Epigenetic feedback on self-activation of GRHL2 does not largely affect EMT/MET dynamics
We have previously analyzed the dynamics of the EMT/MET regulatory network that incorporates the connection of GRHL2 with the two double negative loops that are central to EMT/MET dynamics: miR-200/ ZEB1 and miR-34/SNAIL [36] ( Figure 3A) miR-34 and miR-200 are EMT-inhibiting microRNAs that can inhibit the translation of EMT-TFs SNAIL and ZEB1, thus safeguarding an epithelial phenotype. ZEB and SNAIL can repress E-cadherin and other epithelial genes, and/or drive the expression of mesenchymal genes [31,45]. The knockdown of GRHL2 drives EMT and impairs collective cell migration, while its overexpression may drive an MET [46], the reverse is true for ZEB1 [31,47]. Finally, ZEB1 and GRHL2 can both promote their own expression, albeit indirectly [31,48].
The bifurcation diagram ( Figure 3B) shows that this network can enable cells to exhibit three distinct phenotypes: epithelial (E), hybrid epithelial/mesenchymal (H) and mesenchymal (M). At low levels of an external EMT-inducing signal I, cells are in an epithelial phenotype; they switch to a hybrid E/M phenotype and finally a mesenchymal phenotype as the levels of I increase. For Oncotarget 2614 www.oncotarget.com certain range of values of I, more than one phenotype may co-exist, enabling stochastic cell-fate transitions. Under these conditions, the phenotype exhibited by a cell may switch spontaneously under the influence of biological noise ( Figure 3B) [49], leading to non-trivial phenotypic distributions.
Recent experiments showed how cells in an isogenic population can exhibit varying degrees of susceptibility or resistance to EMT in response to EMT induction [14]. A large subset of single-cell clones (SCCs) (72%) derived from this population successfully underwent EMT and displayed irreversible EMT upon withdrawal of the signal. The remaining 28% SCCs, in response to EMT inducing signal, did not exhibit a reduction in epithelial traits and only underwent a reversible and partial EMT. To capture this behavior in our simulations, we started from all cells in an epithelial state, but with each of them with different random values of threshold of the Hills function denoting the transcriptional inhibition of ZEB by GRHL2. Thus, in this simulation, each cell in this population has different extents of inhibition of ZEB by GRHL2. For a fixed value of EMT-inducing signal I = 75 × 10 3 , molecules in 100 heterogenous cells, only 58% of them underwent an EMT ( Figure 3C). This percentage was higher for larger values of I, but a small subset of cells remained in H or E state, confirming that some cells can intrinsically resist EMTinducing signals ( Figure 3D).
To depict these transitions, we plotted the stochastic dynamics of a population of 1000 cells. The mean value of the signal was fixed at 75 × 10 3 molecules and all cells were initially in an epithelial state. We observed a stable phenotypic distribution with 58% E, 37% hybrid E/M, and 5% M cells ( Figure 4B). Next, we added a strong epigenetic feedback governing the threshold value of GRHL2 that governs its self-activation. The threshold value, instead of being a static variable, is now governed by a dynamical equation with a lower steady-state value which is proportional to the feedback factor α. This feedback has a minimal effect on the bifurcation diagram of the circuit (compare the black and blue curves in Figure  4A); this feature is true for all biologically relevant values of α. We next investigated the effect of this feedback on the equilibrium phenotypic distribution. We started with a scenario where the entire population (n = 1200) exhibits an epithelial phenotype, and tracked the dynamics when a noise term was added to induce spontaneous transitions among the different states (SI sections 1-3). We noticed that the phenotypic distribution seen for the case with epigenetic feedback (55% E, 37% H and 8% M) is largely similar to the scenario without any such feedback (58% E, 37% H, and 5% M) ( Figure 4B).
These results strengthen the observations made from the bifurcation diagram. Moreover, a comparative analysis of the dynamics of two cases show minimal differences; in both cases, instances of partial EMT/MET and complete EMT/MET can be observed, and cells can continue to transition among all three phenotypes ( Figure 4C and 4D). Together, these results suggest that epigenetic feedback acting only on self-activation of GRHL2 has only a weak effect on the dynamics of EMT/MET. Epigenetic feedback on the inhibition of ZEB by GRHL2 can stabilize an epithelial state Next, we examined the effect of adding an epigenetic feedback on the inhibition of ZEB by GRHL2. Unlike the scenario of epigenetic feedback on GRHL2 self-activation, incorporating epigenetic feedback on the inhibition of ZEB by GRHL2 significantly alters the bifurcation diagram.
Compared to the case without any epigenetic feedback (black curve), the bifurcation curve for a case of strong epigenetic feedback case (blue lines) shifts to the right ( Figure 5A), suggesting that a higher external stimuli is required to induce EMT. In other words, this epigenetic feedback can stabilize an epithelial state, or in other words, offer resistance to undergo an EMT. Consistently, large changes were noted in the population distribution as well-the equilibrium population distribution for the case including epigenetic feedback on inhibition of ZEB by GRHL2 was 79% E, 20% H, and 1%M-thus, compared to the control case, the epithelial population increased by 22%, while H and M populations both decreased ( Figure  5B). Moreover, corresponding dynamical trajectories show that in the presence of a strong epigenetic feedback on inhibition of ZEB by GRHL2, it becomes exceedingly difficult for cells to reach a mesenchymal state. Thus, most cells stay much more robustly in an epithelial state ( Figure  5C and 5D). Hence, a sufficiently strong epigenetic feedback on the inhibition of ZEB by GRHL2 can make cells resistant to undergoing a full EMT.
The abovementioned analysis was also conducted on the EMT circuit without miR-34 (because miR-34/ SNAIL has been proposed to be a noise-buffering integrator for EMT). In this simplified model, SNAIL serves as inducing signal. The simulation showed similar results that epigenetic feedback of GRHL2 on the inhibition of ZEB can stabilize an epithelial state, while the epigenetic feedback on GRHL2's self-activation does not largely change EMT/MET dynamics (SI section 5; Supplementary  Figures 4-6).
Next, we studied whether adding an epigenetic feedback on the inhibition of ZEB by GRHL2 can affect the reversibility of EMT. To mimic the experiment where cells are treated with an EMT-inducing signal (say TGF-β; represented by I in simulations) for varying time durations, here, we increased the value of I = 125 × 10 3 molecules and then decreased it back to original value of I = 71 × 10 3 molecules. In the case without any epigenetic feedback case, for a short duration increase in the levels of I, cell can quickly revert to epithelial state upon the removal of the signal ( Figure 6A). For somewhat longer treatment, the cell can stay in hybrid E/M phenotype for a long time and eventually revert to being epithelial ( Figure 6B-6C). A further increase in the duration of cells to EMTinduction can render the induced EMT irreversible; i.e., cells can stay in mesenchymal state and not revert to being epithelial even after the EMT-inducing signal is effectively withdrawn. However, the presence of epigenetic feedback on the inhibition of ZEB by GRHL2 can alter the abovementioned dynamics. In presence of such feedback, a cell can revert to being epithelial rapidly soon after the external signal is withdrawn, irrespective of the duration for which cells were exposed to this signal ( Figure 6D-6F). All these results reveal that epigenetic feedback on the inhibition of ZEB by GRHL2 can significantly stabilize an epithelial state and thus increase the reversibility of EMT.

Noise in the partitioning of parent cell biomolecules among the daughter cells can cause a seemingly irreversible MET at the population level
So far, we have described an epigenetic-based mechanisms which may underlie irreversibility of MET. We next investigated the effect of another stochastic behavior in population of cancer cells-partitioning of molecules during cell division [28][29][30]-independent of epigenetic feedback. Thus, we investigated the dynamics of EMT/MET at a population level, where we incorporated cell division with an average doubling time of 38 hours, typical of cancer cells [50]. At every instance of cell division, we incorporated some noise in the partitioning of parent cell biomolecules among the daughter cells.
Our previous analysis has shown that such noise in the distribution of biomolecules during cancer cell division can generate epithelial-mesenchymal heterogeneity in an initially homogeneous population of cancer cells [3]. Thus, a daughter cell may or not have the same phenotype (E, hybrid E/M, or M) as the parent cell. Further, this noisy partitioning can lead to some cells in the population undergoing a seemingly irreversible EMT if the cells are treated with an EMT-inducing signal followed by the withdrawal of the EMT-inducing signal. We investigated if the noisy partitioning model described previously can also lead to a seemingly irreversible MET at the population level.
Starting with a population of all mesenchymal cells on day 0 (all cells had high concentration of the EMT-inducing signal initially), fixed dosages of the EMT-inducing signal were withdrawn each day for a period of 10 days (to simulate MET). Day 11 onwards, fixed dosages of the EMT-inducing signal were added for another 10 days. We carried out the simulation in the absence and presence of GRHL2, and in both cases, ~40% of the cells had undergone MET by day 10. In the absence of GRHL2, upon treatment with the EMT-inducer starting on day 11, almost all the cells in the population that had undergone MET returned to a mesenchymal state by day 10. Thus, in this scenario, MET was reversible. However, in the scenario when GRHL2 was present, it was observed that > 10% of cells still exhibited an epithelial phenotype on day 20 (Figure 7). These cells thus represented a subpopulation that had under-gone an irreversible MET; i.e., a subpopulation which would exhibit resistance to EMT upon exposure to EMT-inducing signal.

DISCUSSION
With an increasing appreciation of the nonlinear dynamics of EMT/MET at the single-cell level [51][52][53] and its implications for metastatic aggressiveness [9], questions regarding the degree of reversibility/ irreversibility of EMT/MET in different contexts and the molecular determinants of these processes have gained importance. Recent studies have illustrated that the trajectories taken by individual cells in the highdimensional molecular and/or morphological landscape of EMP en route to EMT and MET may be different [51,54]; thus, MET cannot be simply thought of as a mirror image of EMT. There may be molecular and/or morphological changes happening at different stages of EMT/MET to varying degrees, hence making it difficult to identify the molecular mechanisms that may render the dynamics of EMT/MET as reversible or irreversible.
The dynamics of EMT has been studied much more in detail as compared to that of MET [51][52][53][54][55][56]; therefore, it is not surprising that irreversible EMT has been reported more frequently. The degree of reversibility of EMT has been proposed to be largely a function of the timescale Oncotarget 2618 www.oncotarget.com of EMT induction and corresponding epigenetic changes. However, a causative role of epigenetic changes in regulating the irreversibility of EMT remains to be firmly established [11][12][13][14]. Here, we propose two independent mechanisms that may enable an irreversible MET: 1) epigenetic feedback mediated by the inhibitory action of GRHL2 on the promoter of ZEB1; and 2) noise in the partitioning of biomolecules during cell division.
GRHL2 can inhibit and reverse EMT and associated molecular and/or morphological traits. Specifically, overexpression of GRHL2 has been shown to induce epithelial gene expression, inhibit mesenchymal gene expression, restore metabolic reprogramming caused by EMT, and suppress tumor cell migration/invasion, at least in breast and ovarian cancer cell lines [46]. It can also activate directly or indirectly other drivers of an epithelial phenotype such as p63 [57,58], OVOL2 [59,60], and miR-200 [38,61]. While epigenetic changes induced by EMT-TFs have been reported extensively [62], recent studies have pointed out the possibility of GRHL2 contributing to epigenetic control of genes involved in EMT/MET [39]. GRHL2 can employ both DNA methylation and histone modification to inhibit and/or reverse EMT, and also act as a pioneer transcription factor that can regulate chromatin accessibility at epithelial enhancers [39,63,64].
The experimental observations described above offer a possible underlying mechanism by which GRHL2 overexpression can resist EMT, as demonstrated by our model simulations. Indeed, a global epigenetic program that limits the action of ZEB1 was found to underlie the retention of epithelial traits in cells exposed to persistent Twist1 activation for 21 days [14]. In single-cell clones established from an inducible HMLE-Twist1 population, Oncotarget 2619 www.oncotarget.com two subsets responded quite differently (E-SCCs, M-SCCs); while both E-SCCs and M-SCCs exhibited an upregulation of mesenchymal markers to a similar degree upon Twist1 induction, but E-SCCs did not display any reduction in epithelial genes even after 28 days of Twist1 induction. Upon inactivation of Twist1, E-SCCs reverted to an epithelial phenotype, while M-SCCs did not. Changes in chromatin accessibility seen in E-SCCs also reverted upon Twist1 deactivation, but not in M-SCCs, suggesting a strong correlation of transcriptional changes induced in E-SCCs/M-SCCs and dynamic epigenetic status. Intriguingly, this study highlighted that isogenic/ clonal cells can also exhibit variability in terms of their susceptibility or resistance to EMT/MET, indicating non-genetic mechanisms at play [65], as indeed captured via our simulations. Our simulations also offer a causal connection between the two axes phenomenologically associated with one another in this study: epigenetic changes driven by GRHL2 and consequent resistance to EMT.
More generally, resistance to EMT has also been observed across many cell lines, where their treatment with 5 ng/ml TGFβ for 48 hours led to varying degrees of changes in molecular and morphological axes of EMT: spindle-shape acquisition, loss of junctional E-cadherin and/or ZO-1, and actin stress fiber formation [24]. Another recent study in multiple breast cancer cell lines showed that 100pM TGFβ treatment for 9 days need not be sufficient for at least a subset of cells within a cell line to lose their E-cadherin expression [52]. While longertime measurements would probably be more appropriate to fully assess the degree of resistance to EMT (in other words, the degree of irreversibility of MET), both these studies emphasize the possible implications of non-genetic heterogeneity prevalent in multistable regulatory networks [66] as seen in EMT/MET and its associated traits such as stemness [67]. Different degrees of couplings between these multistable networks may enable co-occurrence of partial or full EMT with these associated traits. However, a mechanism-based mathematical modeling of these networks suggests that this association is not likely to be universal [68], hence offering a reconciliatory framework to integrate various contradictory results associating the epithelial, mesenchymal and hybrid E/M phenotypes with degrees of stemness [69][70][71][72][73][74][75].
Aside from epigenetic feedback, our results indicate one additional mechanism possibly contributing to irreversible MET (or resistance to EMT), namely stochasticity in the partitioning of molecules during cell division. Such noise in the distribution of molecules may affect cell-fate and drive non-genetic heterogeneity [28][29][30], leading to different phenotypic distributions in terms of EMT [3]. While EMT is believed to repress the cell cycle [76,77], this association remains controversial [78]. Thus, noise in the partitioning of parent cell biomolecules among the daughter cells can further alter the subpopulation structure, and may underlie different bimodal distributions of surface CDH1 expression seen in breast cancer cell lines [52]. Another factor contributing to these non-genetic differences explored here could be the different EMT vs. MET trajectories. When a cell undergoes MET, it need not return to the precise epithelial (E) state coordinates that it started form [51], possibly due to hysteresis in EMT dynamics [79]. Even if it returns to the same state, the minimum action path (MAP) from E to M need not be the same as the MAP from the M to the E state [80].
There are various limitations to our analysis. First, we do not consider the detailed molecular mechanisms underlying epigenetic changes; instead, our treatment is phenomenological; we assume the net effect of Oncotarget 2620 www.oncotarget.com these changes is the feature that if gene i can regulate the expression of gene j, when gene i is expressed, the activation of gene j is more likely to become stronger [32]. A more detailed molecular mechanism based epigenetic model would be needed in order to study in more detail how epigenetics regulates EMT, including for example, the PRC2-mediated histone methylation epigenetic framework [81]. This will become possible as more details emerge regarding epigenetic regulation of EMT [62]. Also, we have not considered any spatial effects in our model; these may become more relevant in crowded regions in terms of enabling access to signaling molecules. Spatially extended frameworks of EMT [79] can be integrated with our epigenetic framework to investigate the effect of crowding and nutrient/oxygen access.
Together, our results offer mechanistic insights into two possible mechanisms that may drive varying degrees of susceptibility and resistance to undergoing EMT in response to an EMT-inducing signal in a given isogenic population. Future efforts should decode the molecular mechanisms of any such epigenetic feedback of GRHL2 on ZEB1 expression as well as track the distribution of molecules during cell divisions happening while cells are being induced to undergo EMT/MET. Moreover, our study calls for concerted efforts to map the single-cell dynamics of MET induction.

MATERIALS AND METHODS
The analysis shown in Figures 1 and 2 was carried out using publicly available datasets. TCGA data was obtained from https://xenabrowser.net/datapages/. Different EMT scoring metrics were calculated as described previously [44].
The set of differential equations used to simulate the dynamics of the EMT/MET regulatory circuit are enumerated in the SI. The SI also includes all the model parameters and a description of how epigenetic feedback was incorporated into the mathematical model of EMT/ MET regulation (Supplementary Tables 1-3). Finally, the analysis shown in Figure 7 was carried out using the population-level model described previously (3). The computer code used to generate the data shown in Figure 7 is available online (https://www.github.com/ st35/cancer-EMT-heterogeneity-noise/tree/master/ ExternalIntervention; see files "DailyIntervention_MET. cpp" and "DailyIntervention_GRHL2_MET.cpp").

Author contributions
MKJ designed research; MKJ and HL supervised research; WJ, ST, PC, AC conducted research; AR, MKJ, HL analyzed data; all authors participated in writing and revision of the manuscript.