Comparative Study of Transcriptomics-Based Scoring Metrics for the Epithelial-Hybrid-Mesenchymal Spectrum

Abstract

The Epithelial-mesenchymal transition (EMT) is a cellular process implicated in embryonic development, wound healing, and pathological conditions such as cancer metastasis and fibrosis. Cancer cells undergoing EMT exhibit enhanced aggressive behavior characterized by drug resistance, tumor-initiation potential, and the ability to evade the immune system. Recent in silico, in vitro, and in vivo evidence indicates that EMT is not an all-or-none process; instead, cells can stably acquire one or more hybrid epithelial/mesenchymal (E/M) phenotypes which often can be more aggressive than purely E or M cell populations. Thus, the EMT status of cancer cells can prove to be a critical estimate of patient prognosis. Recent attempts have employed different transcriptomics signatures to quantify EMT status in cell lines and patient tumors. However, a comprehensive comparison of these methods, including their accuracy in identifying cells in the hybrid E/M phenotype(s), is lacking. Here, we compare three distinct metrics that score EMT on a continuum, based on the transcriptomics signature of individual samples. Our results demonstrate that these methods exhibit good concordance among themselves in quantifying the extent of EMT in a given sample. Moreover, scoring EMT using any of the three methods discerned that cells can undergo varying extents of EMT across tumor types. Separately, our analysis also identified tumor types with maximum variability in terms of EMT and associated an enrichment of hybrid E/M signatures in these samples. Moreover, we also found that the multinomial logistic regression (MLR)-based metric was capable of distinguishing between "pure" individual hybrid E/M vs. mixtures of E and M cells. Our results, thus, suggest that while any of the three methods can indicate a generic trend in the EMT status of a given cell, the MLR method has two additional advantages: (a) it uses a small number of predictors to calculate the EMT score and (b) it can predict from the transcriptomic signature of a population whether it is comprised of "pure" hybrid E/M cells at the single-cell level or is instead an ensemble of E and M cell subpopulations.

Keywords: EMT; EMT quantification; EMT score; MET; hybrid epithelial/mesenchymal; tumor heterogeneity.

FIGURE 1

General outline of all three EMT scoring methods. (A) 76GS score is calculated by weighted sum of 76 genes, where EMT score_i is score for ith sample, w_j is correlation of jth gene (G_j) with CDH1 gene in that dataset to which the sample i belongs, G_ij is the jth gene’s normalized expression in ith sample. (B) MLR utilizes log₂(VIM)/log₂(CDH1) and log₂(CLDN7) space to predict categorization of a sample into E, E/M, or M category. Where P_E, P_H, and P_M are the probabilities of a sample falling into each phenotype. EMT score_i is the score for ith sample, which is defined in relation to P_E, P_H, and P_M. Figure adapted from George et al. (2017) with permission. (C) KS score is estimated by the empirical cumulative distributions of epithelial and mesenchymal gene set, denoted by ecdf (GS_mes) and ecdf (GS_epi), respectively. EMT score_i is the maximum vertical distance between the ecdf (GS_mes) and ecdf (GS_epi) (given by Eq. 1 in the section “Materials and Methods”) for a given sample i.

FIGURE 2

Scatter plot depicting the correlation between the EMT scores of cancer cell line samples calculated via three methods. Each pairwise relation is estimated by a linear regression line (red), Pearson’s correlation coefficient (R), and p-value (p) reported in each plot. (A) NCI60 dataset and (B) CCLE dataset.

FIGURE 3

Concordance across all three EMT scoring methods in quantification of EMT and survival prediction of tumor patients. Each pairwise relation is estimated by linear regression (red), Pearson’s correlation coefficient (R), and p-value (p), reported in each plot. (A) TCGA ovarian cancer dataset, (B) TCGA sarcoma dataset, (C) TCGA breast cancer dataset. (D) Correlation between EMT score (high vs. low) and overall survival (OS) in breast cancer patients. Kaplan–Meier survival analysis is performed to estimate differences in survival of 76GS^high, MLR^high, KS^high and 76GS^low, MLR^low, KS^low groups, respectively, in GSE1456. p-values (p) reported are based on the log rank test. HR (hazard ratio) and confidence interval (95% CI) reported are estimated using cox regression.

FIGURE 4

Plots depicting pairwise comparisons of all three EMT scores. (A) Volcano plots showing the correlation of different EMT scoring methods across 85 different GEO microarray datasets along with the p-values for the respective correlation coefficient values. In each case, –log10(p-value) is plotted as a function of Pearson’s correlation coefficient. Thresholds for correlation (R < –0.3 or R > 0.3; vertical blue lines) and p-values (p < 0.05; horizontal red line) are denoted. (B) Bar plots for different categories based on the correlation sign and statistical significance of all three pairwise comparisons across 85 datasets. p < 0.05 and R < −0.3 or R > 0.3. (C) Venn diagram showing the common GEO datasets across all pairwise comparisons that are significantly correlated in the expected direction.

FIGURE 5

Bar plots showing EMT scores of different cell lines calculated using the three EMT scoring methods. (A) EMT induction is shown in three cell lines – A549, HCC87, and NCIH358 (GSE49664). (B) EMT induction in MCF7 cell line (GSE58252). (C) EMT induction in EpRas cells (GSE59922). (D) EMT induction by different EMT-inducing transcription factors. “a” denotes statistical significant difference (p < 0.05, n = 3, two-tailed Student’s t-test) for pairwise comparison of a given set with untreated (first column), “b” denotes the same when compared with empty vector (EV; second column) (GSE43495). (E) MET induction by GRHL2 in MDA-MB-231 cell line (GSE36081). (F) Two cell lines of hepatocellular carcinoma with varying EMT status (GSE26391). Each control case has been compared to EMT/MET induced case (*p < 0.05, n = 3, two-tailed Student’s t-test; error bars represent standard deviation).

FIGURE 6

Variance and mean of EMT scores in CCLE samples grouped by tumor subtype, highlighting the most variable tumor types (circled). (A) 76GS EMT scores, (B) MLR EMT scores, and (C) KS EMT scores. (D) Venn diagram showing the overlap between each tumor type based on the abundance of hybrid samples as defined by the MLR method, where #EM > 10 denote the cases where the absolute number of hybrid E/M samples in a tumor subtype is >10; %EM > 20 denote the cases where the percentage of cell lines identified as hybrid E/M in a given tumor subtype is >20%.

FIGURE 7

Distinguishing between hybrid E/M cells vs. mixtures of E and M cells. (A) Scatter plot showing CCLE cell lines that display a hybrid E/M phenotype (red) on the mixture curve (dotted curve) determined by the mean of 35 pure E (orange) and pure M (blue) reference samples in CCLE dataset. (B) Scatter plot showing the 100 farthest (purple) and 100 closest (green) samples based on the distance from the mixture curve. (C) Bar plots showing EMT scores of N (10, 20, 50, and 100) closest and farthest hybrid E/M samples from mixture curve. (D) Mesenchymal proportion (%M) distribution of the 100 closest and farthest hybrid samples from mixture curve. *p < 0.05, N = 10, 20, 50 and 100, two-tailed Student’s t-test; error bars represent standard deviation.