Robustness in phenotypic plasticity and heterogeneity patterns enabled by EMT networks

Abstract

Epithelial-mesenchymal plasticity (EMP) is a key arm of cancer metastasis and is observed across many contexts. Cells undergoing EMP can reversibly switch between three classes of phenotypes: epithelial (E), mesenchymal (M), and hybrid E/M. While a large number of multistable regulatory networks have been identified to be driving EMP in various contexts, the exact mechanisms and design principles that enable robustness in driving EMP across contexts are not yet fully understood. Here, we investigated dynamic and structural robustness in EMP networks with regard to phenotypic heterogeneity and plasticity. We use two different approaches to simulate these networks: a computationally inexpensive, parameter-independent continuous state space Boolean model, and an ODE-based parameter-agnostic framework (RACIPE), both of which yielded similar phenotypic distributions. While the latter approach is useful for measurements of plasticity, the former model enabled us to extensively investigate robustness in phenotypic heterogeneity. Using perturbations to network topology and by varying network parameters, we show that multistable EMP networks are structurally and dynamically more robust compared with their randomized counterparts, thereby highlighting their topological hallmarks. These features of robustness are governed by a balance of positive and negative feedback loops embedded in these networks. Using a combination of the number of negative and positive feedback loops weighted by their lengths, we identified a metric that can explain the structural and dynamical robustness of these networks. This metric enabled us to compare networks across multiple sizes, and the network principles thus obtained can be used to identify fragilities in large networks without simulating their dynamics. Our analysis highlights a network topology-based approach to quantify robustness in the phenotypic heterogeneity and plasticity emergent from EMP networks.

Figure 1

Measurement of robustness in multistable biological networks. (a) Depiction of plasticity as the ability of a system (defined by a set of parameters) to achieve multiple phenotypes depending on the initial state. (b) EMP networks analyzed, namely GRHL2 (4 nodes, 7 edges), OCT4 (5 nodes, 10 edges), and NRF2 (8 nodes, 16 edges). (c) Demonstration of how JSD between two frequency distributions varies; JSD ranges from 0 to 1. (d) Different ways to perturb the topology/structure of a network: deletion of an edge, addition of an edge, edge nature change from activation to inhibition and vice versa. (e) Depiction of dynamic robustness where different parameter spaces converge to largely overlapping solution space.

Figure 2

High dynamical robustness of EMP networks. (a) Histogram of average JSD for random networks of size 4, with corresponding EMP networks marked with arrows with different shapes at their tail. Percentiles for WT networks reported as mean ± standard deviation. (b) Similar to (a), but for average fold change in plasticity. (c and d) Similar to (a) and (b) but for random networks of size 5, with corresponding EMP network marked with an arrow. (e and f) Similar to (a) and (b), but for random networks of size 8, with corresponding EMP network marked with an arrow.

Figure 3

Continuous state space Boolean (CSB) formalism captures the dynamic robustness of EMP networks. (a) Visualization of the CSB framework. Networks are simulated in a continuous state space and the resultant final states are discretized. (b) Comparison of JSD of the phenotypic distributions obtained by RACIPE from the Boolean framework (dotted bars) and CSB framework (cross-patterned bars) for various EMP networks. (c) Comparison of JSD of the phenotypic distributions obtained from RACIPE against those obtained from CSB formalism (x axis) and Boolean formalism (y axis). The line represents $x=y$. (d) Histogram of average JSD for random networks of size 4, with corresponding EMP networks marked with arrows. Percentiles for WT networks reported as mean $\pm$ standard deviation. (e) Same as (d), but for networks of size 5. (f) Same as (e), but for networks of size 8.

Figure 4

Structural robustness and feedback loops. (a) Comparison of perturbation JSD obtained via RACIPE (y axis) and CSB (x axis) for OCT4 network. Each dot represents one perturbed network. The correlation coefficient and the residual for linear fit has been reported, with the p value for correlation to be interpreted as: ^∗∗∗p < 0.001. (b) Correlation coefficient between perturbation JSD from RACIPE and perturbation JSD from Boolean (dotted bar) and perturbation JSD from CSB (criss-cross patterned bar). (c) Histogram of average perturbation JSD for networks of size 4, with corresponding EMP networks marked with arrows. Percentiles for WT networks reported as mean $\pm$ standard deviation. Inset shows the distribution of perturbation JSD for GRHL2 network. (d) Same as (c), but for random networks of size 5. (e) Same as (c), but for average fold change in plasticity. (f) Same as (e), but for networks of size 5.

Figure 5

Feedback loops can explain the robustness of EMP networks. (a and b) Violin plots depicting the distribution of negative feedback loops (NFLs) and positive feedback loops (PFLs) on either side of the median of the histograms in Fig. 4,c and e, respectively. Unpaired t test comparing the distribution of feedback loops on either side of the median was performed, with significance depicted on top of the plots. (c and d) Scatterplot of the average fold change in plasticity upon structural perturbation (structural robustness in plasticity) (y axis) and number of PFLs (c) and the number of NFLs (d). Each dot represents a random network. The WT EMP networks are highlighted using points of different shapes. Spearman correlation coefficient is reported in the legend of both panels, with the p value to be interpreted as: ^∗∗∗p < 0.001. (e and f) Same as (c) and (d), but for average perturbation JSD (structural robustness in distribution). (g–i) Scatterplot between fraction of positive cycles (FPCs, x axis) against the robustness measures: (g) JSD between the phenotypic distributions obtained from RACIPE and CSB formalisms (dynamical robustness in distribution), (h) average fold change in plasticity upon structural perturbations, and (i) average perturbation JSD. The WT EMP networks are highlighted using points of different shapes. Spearman correlation coefficient is reported in the legend of each panel. p Values in the figure to be interpreted as: ∗ p < 0.05, ∗∗ p < 0.01, ^∗∗∗p < 0.001.

Figure 6

Fraction of weighted positive cycles as a measure of robustness. (a) Scatterplot of the correlation of plasticity with PFLs and weighted PFLs, along with the $x=y$ line. (b) Absolute value of the Spearman correlation b/w average perturbation JSD and PFLs, NFLs, fraction of positive cycles (FPC), and FWPC for networks of different sizes. Absolute values have been used for ease of comparison across robustness measures. (c) Change in the strength of correlation between the robustness measures—structural robustness in JSD (pJSD), plasticity (pPlast), dynamic robustness in JSD (dJSD), and plasticity (dPlast)—and FWPC, with change in the weight of WPFL in FWPC, for networks of size 4. (d) Heatmap of the correlation coefficients of robustness measures—pJSD, pPlast, dJSD, and dPlast —with PFLs, NFLs, FPC, and FWPC. Statistical significance of correlation is reported as follows: ^∗∗∗p < 0.001.

Figure 7

Feedback loop-based metric explains robustness in large-scale networks. (a and b) The 22 node EMP network (EMT_RACIPE) and 26 node EMP network (EMT_RACIPE2). (c and d) Barplot of the fraction of multistable parameters in RACIPE for EMT_RACIPE, EMT_RACIPE2. Error bars represent the standard deviation of the fraction obtained over 3 replicates. (e and f) Plasticity of the perturbed network versus change in PFLs for each perturbed network (EMT_RACIPE, EMT_RACIPE2) was plotted. The size of each point represents the number of negative feedback loops. (g and h) Same as (e) and (f), but for change in FWPC. Color bar denotes number of negative feedback loops in a given perturbed network topology. Significance of correlation is reported: ^∗∗∗p < 0.001.

Figure 8

Perturbations reducing FWPC reduce network robustness and vice versa. (a and b) Representative barplot showing the edge perturbations that increase or decrease FWPC maximally iteratively. The arrows show the point at which perturbations were stopped, when either the maximum (1) or minimum (0) FWPC were reached. (c and d) Change in perturbation JSD for NRF2 and a representative random network of size 8 with each perturbation described in (a) and (b). (e) Change in the FWPC metric with network perturbation for three additional random networks. (f) Change in structural robustness (perturbation JSD) with network perturbation for three additional random networks. In (c), (d), and (f); mean ± standard deviation across three replicates are reported.