Since its discovery almost three decades ago, the Janus kinase (JAK)-signal transducer and activator of transcription (STAT) pathway has paved the road for understanding inflammatory and immunity processes related to a wide range of human pathologies including cancer. Several studies have demonstrated the importance of JAK-STAT pathway components in regulating tumor initiation and metastatic progression, yet, the extent of how genetic alterations influence patient outcome is far from being understood.
Focusing on 133 genes involved in JAK-STAT signaling, we investigated genomic, transcriptomic and clinical profiles of over 18,000 patients representing 21 diverse cancer types. We identified a core set of 28 putative gain- or loss-of-function JAK-STAT genes that correlated with survival outcomes using Cox proportional hazards regression and Kaplan-Meier analyses. Differential expression analyses between high- and low-expressing patient groups were performed to evaluate the consequences of JAK-STAT misexpression.
We found that copy number alterations underpinning transcriptional dysregulation of JAK-STAT pathway genes differ within and between cancer types. Integrated analyses uniting genomic and transcriptomic datasets revealed a core set of JAK-STAT pathway genes that correlated with survival outcomes in brain, renal, lung and endometrial cancers. High JAK-STAT scores were associated with increased mortality rates in brain and renal cancers, but not in lung and endometrial cancers where hyperactive JAK-STAT signaling is a positive prognostic factor. Patients with aberrant JAK-STAT signaling demonstrated pan-cancer molecular features associated with misexpression of genes in other oncogenic pathways (Wnt, MAPK, TGF-β, PPAR and VEGF). Brain and renal tumors with hyperactive JAK-STAT signaling had increased regulatory T cell gene (Treg) expression. A combined model uniting JAK-STAT and Tregs allowed further delineation of risk groups where patients with high JAK-STAT and Treg scores consistently performed the worst.
Providing a pan-cancer perspective of clinically-relevant JAK-STAT alterations, this study could serve as a framework for future research investigating anti-tumor immunity using combination therapy involving JAK-STAT and immune checkpoint inhibitors.
In their quest to survive and prosper, tumor cells are armored with a unique ability to manipulate the host’s immune system and promote pro-inflammatory pathways. Inflammation can both initiate and stimulate cancer progression, and in turn, tumor cells can create an inflammatory microenvironment to sustain their growth further (Mantovani et al., 2008; Trinchieri, 2012). Cytokines are secretable molecules that influence immune and inflammatory processes of nearby and distant cells. Although cytokines are responsible for inflammation in cancer, spontaneous eradication of tumors by endogenous immune processes rarely occurs. Moreover, the dynamic interaction between tumor cells and host immunity shields tumors from immunological ablation, which overall limits the efficacy of immunotherapy in the clinic.
Cytokines can be pro- or anti-inflammatory and are interdependent on each other’s function to maintain immune homeostasis (Zhang & An, 2007). Discovered as a critical regulator of cytokine signaling, the Janus kinase (JAK) - signal transducer and activator of transcription (STAT) pathway allows cytokines to transduce extracellular signals into the nucleus to regulate gene expression implicated in a myriad of developmental processes including cellular growth, differentiation and host defense (Slattery et al., 2013). JAK proteins interact with cytokine receptors to phosphorylate signaling substrates including STATs. Unlike normal cells which phosphorylate STATs temporarily, several STAT proteins were found to be persistently phosphorylated and activated in cancer (Yu et al., 2009). Studies have shown that persistent activation of STAT3 and STAT5 promote inflammation of the microenvironment, tumor proliferation, invasion and suppress anti-tumor immunity (Yu et al., 2009). Particularly in cancers associated with chronic inflammation such as liver and colorectal cancers (Grivennikov et al., 2009; Park et al., 2010), STAT3 activation by growth factors or interleukins suppresses T cell activation and promotes the recruitment of anti-immunity factors such as myeloid-derived suppressor cells and regulatory T cells (Poschke et al., 2010; Tartour et al., 2011).
Given its fundamental roles in interpreting environmental cues to drive a cascade of signaling events that control growth and immune processes, it is essential to dissect cell type-specific roles of the JAK-STAT pathway in a pan-cancer context. This is made possible by advances in high-throughput sequencing initiatives and as many genetic alterations have become targetable, detailed understanding on genetic variations would be essential to identify particular weaknesses in individual tumors in order to boost therapeutic success. We predict that genetic alterations in JAK-STAT pathway genes do not occur equally between and within cancer types. Moreover, detection of rare alterations would require a large sample size to unravel genes that are altered within specific histological subtypes of cancer. With over 18,000 samples representing 21 cancer types, we took the opportunity to systematically characterize genetic alterations within 133 JAK-STAT pathway genes to uncover shared commonalities and differences. To identify functionally relevant alterations, genetic polymorphisms were overlaid with transcript expression profiles and correlated with clinical outcomes. Our study identifies a core set of candidate JAK-STAT drivers that correlated with tumor progression and that predict overall survival outcomes in brain, renal, lung and endometrial cancers converging on similar downstream oncogenic pathways. This work provides a rich source of cancer type-dependent alterations that could serve as novel therapeutic targets to support underexploited treatment initiatives targeting JAK-STAT signaling in cancer.
All plots were generated using R packages (pheatmap [version 1.0.12] and ggplot2 [version 3.1.0]).
Cancer cohorts and JAK-STAT pathway genes
133 JAK-STAT pathway genes were obtained from the Kyoto Encyclopedia of Genes and Genomes (KEGG) database listed in Additional file 1. Genomic, transcriptomic and clinical datasets of 21 cancer types (n = 18,484) were retrieved from The Cancer Genome Atlas (TCGA) (Weinstein et al., 2013). The following is a list of cancer cohorts and corresponding TCGA abbreviations in parentheses: bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC), cholangiocarcinoma (CHOL), colon adenocarcinoma (COAD), esophageal carcinoma (ESCA), glioblastoma multiforme (GBM), glioma (GBMLGG), head and neck squamous cell carcinoma (HNSC), kidney chromophobe (KICH), pan-kidney cohort (KIPAN), kidney renal clear cell carcinoma (KIRC), kidney renal papillary cell carcinoma (KIRP), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), pancreatic adenocarcinoma (PAAD), sarcoma (SARC), stomach adenocarcinoma (STAD), stomach and esophageal carcinoma (STES) and uterine corpus endometrial carcinoma (UCEC).
Copy number variation analyses
Level 4 GISTIC copy number variation datasets were downloaded from the Broad Institute GDAC Firehose (Mermel et al., 2011). Discrete amplification and deletion indicators were obtained from GISTIC gene-level tables. Genes with GISTIC values of 2 were annotated as deep amplification events, while genes with values of − 2 were annotated as deep (homozygous) deletion events. Shallow amplification and deletion events were annotated for genes with values of + 1 and − 1 respectively (Chang & Lai, 2019a).
Calculating JAK-STAT and regulatory T cell scores
A JAK-STAT 28-gene signature was developed from putative gain- or loss-of-function candidates. For each patient, 28-gene scores were calculated from the average log2 expression of signature genes: IL7, IFNG, MPL, IL11, IL2RA, IL21R, OSMR, IL20RA, IFNGR1, CDKN1A, CISH, SOCS1, IL10, IL10RA, STAT2, IL24, IL23A, PIAS3, IFNLR1, EPO, TSLP, BCL2, IL20RB, IL11RA, PTPN6, IL13, IL17D and IL15RA. Regulatory T cell (Treg) scores were determined by taking the mean expression of 31 Treg genes identified from the overlap of four Treg signatures to generate a more representative gene set that is cell type-independent (De Simone et al., 2016; Plitas et al., 2016; Tirosh et al., 2016; Zheng et al., 2017).
Multidimensional scaling, survival and differential expression analyses
We previously published detailed methods on the above analyses (Chang et al., 2019a) and thus the methods will not be repeated here. Briefly, multidimensional scaling analyses based on Euclidean distance in Fig. 2f were performed using the vegan package in R (version 2.5–6). Permutational multivariate analysis of variance (PERMANOVA) was used to determine statistical significance between tumor and non-tumor samples (Chang & Lai, 2019b; Chang & Lai, 2019c). Survival analyses were performed using Cox proportional hazards regression and the Kaplan-Meier method coupled with the log-rank test (Chang et al., 2019b; Chang & Lai, 2019d) using the R survival (version 2.43–3) and survminer (version 0.4.3) packages. Predictive performance of the 28-gene signature was assessed using the receiver operating characteristic analysis. To determine the prognostic significance of a combined model uniting the JAK-STAT signature and IRF8 expression or Treg scores, patients were separated into four survival categories based on median 28-gene scores and IRF8 expression or Treg scores for Kaplan-Meier and Cox regression analyses. To determine the transcriptional effects of aberrant JAK-STAT signaling, differential expression analyses were performed on patients within the 4th versus 1st survival quartiles (stratified using the 28-gene signature).
Functional enrichment and transcription factor analyses
Differentially expressed genes (DEGs) identified above were mapped to KEGG and Gene Ontology (GO) databases using GeneCodis3 (Tabas-Madrid et al., 2012) to determine significantly enriched pathways and biological processes. DEGs were also mapped to ENCODE and ChEA databases using Enrichr (Chen et al., 2013; Kuleshov et al., 2016) to identify transcription factors that were significant regulators of the DEGs.
Copy number and transcriptome analyses reveal conserved driver mutations in JAK-STAT pathway genes
We interrogated genomic and transcriptomic landscape of 133 JAK-STAT pathway genes in 18,484 patients across 21 cancer types (Additional file 1). To investigate the effects of JAK-STAT pathway genomic alterations on transcriptional output, we first analyzed copy number variation (CNV) of all 133 genes. CNVs were classified into four categories: low-level amplifications, deep amplifications, low-level deletions (heterozygous deletions) and deep deletions (homozygous deletions). Lung squamous cell carcinoma (LUSC) and papillary renal cell carcinoma (KIRP) had the highest and lowest fraction of samples with deleted JAK-STAT pathway genes respectively (Additional file 2). In terms of gene amplification, this was also the highest in lung squamous cell carcinoma (LUSC) and the lowest in pancreatic adenocarcinoma (PAAD) (Additional file 2). To identify pan-cancer CNV events, we focused on genes that were deleted or amplified in at least one-third of cancer types (≥ 7 cancers). We identified 71 and 49 genes that were recurrently deleted and amplified respectively in at least 20% of samples within each cancer type and at least 7 cancer types (Additional file 2). Esophageal carcinoma (ESCA) had 70 genes that were recurrently deleted while only four recurrently deleted genes were found in papillary renal cell carcinoma (KIRP). When considering recurrent gene amplifications, lung adenocarcinoma (LUAD) had the highest number of gains (48 genes), while the lowest number of gene gains was observed in glioma (GBMLGG) (3 genes) (Additional file 2). CNV events associated with transcriptional changes could represent candidate driver genes. Loss-of-function genes can be identified from genes that were recurrently deleted and downregulated at the transcript level. Similarly, genes that were concomitantly gained and upregulated could represent a gain-of-function. Differential expression profiles (tumor vs. non-tumor) were intersected with CNV profiles and we identified 40 driver genes representing potential loss- or gain-of-function (Fig. 1). In at least 7 cancer types, 18 genes were recurrently deleted and downregulated (log2 fold-change < − 0.5, P < 0.01), while a non-overlapping set of 22 genes were recurrently amplified and upregulated (log2 fold-change > 0.5, P < 0.01) (Fig. 1).
JAK-STAT driver genes predict overall survival in diverse cancer types
Focusing on the 40 driver genes identified above, we next analyzed whether expression levels of individual genes were associated with overall survival outcomes. Cox proportional hazards regression analyses demonstrated that all 40 genes harbored prognostic information in at least one cancer type. IL11, PTPN6 and CISH were significantly associated with survival outcomes in patients from 9 cancer cohorts (Fig. 2a). In contrast, IL19, CNTFR and JAK2 were some of the least prognostic genes (Fig. 2a). When deciphering the contribution of individual genes across cancer types, we observed that the glioma (GBMLGG) cohort had the highest number of prognostic genes (31/40 genes) (Fig. 2a). On the other hand, none of the 40 driver candidates were prognostic in esophageal carcinoma (ESCA) and cholangiocarcinoma (CHOL), which suggests the minimal contribution of JAK-STAT signaling in driving tumor progression in these cancer types (Fig. 2a). To identify a core set of prognostic genes denoting pan-cancer significance, we performed Spearman’s correlation analyses on hazard ratio (HR) values obtained from Cox regression and identified 28 highly-correlated genes: IL7, IFNG, MPL, IL11, IL2RA, IL21R, OSMR, IL20RA, IFNGR1, CDKN1A, CISH, SOCS1, IL10, IL10RA, STAT2, IL24, IL23A, PIAS3, IFNLR1, EPO, TSLP, BCL2, IL20RB, IL11RA, PTPN6, IL13, IL17D and IL15RA (Fig. 2b). These genes were collectively regarded as a pan-cancer JAK-STAT signature. To determine the extent of JAK-STAT pathway variation across cancer types, we calculated an activity score based on the mean expression of the 28 genes. When cancers were sorted according to their pathway activity scores, chromophobe renal cell cancer (KICH) had the lowest median score while the highest median score was observed in head and neck cancer (HNSC) (Fig. 2c). Hierarchical clustering of the 28 driver genes demonstrated that they exhibited a wide range of expression depending on the cellular context where they could serve as potential candidates for therapy. For example, IL7, IL15RA, IL21R, IL10, OSMR, IFNGR1, IL10RA, IL2RA, IFNG, IL24, SOCS1, IL20RA, IL11 and IL23A were highly expressed in gastrointestinal cancers (PAAD, STAD and STES) (Fig. 2d). When the 28-gene scores were employed for patient stratification, we observed that the JAK-STAT signature conferred prognostic information in five diverse cancer cohorts (Fig. 2e). Intriguingly, the significance of the signature in predicting overall survival was cancer type-dependent. Kaplan-Meier analyses and log-rank tests revealed that patients with high scores (4th quartile) had higher death risks in glioma (P < 0.0001), pan-kidney (consisting of chromophobe renal cell, clear cell renal cell and papillary renal cell cancers; P < 0.0001) and clear cell renal cell (P < 0.0001) cohorts (Fig. 2e). In contrast, high expression of signature genes was linked to improved survival rates in lung (P = 0.025) and endometrial (P = 0.032) cancers (Fig. 2E). These results were independently corroborated using Cox regression analyses: glioma (HR = 6.832, P < 0.0001), pan-kidney (HR = 3.335, P < 0.0001), clear cell renal cell (HR = 4.292, P < 0.0001), lung (HR = 0.624, P = 0.028) and endometrium (HR = 0.504, P = 0.027) (Additional file 3). Since high expression of signature genes was associated with poor survival outcomes in brain and renal cancers, as expected, we observed a significant increase in expression scores according to tumor stage (Fig. 2g). An opposite trend was observed in lung and endometrial cancers where more aggressive tumors had lower expression scores (Fig. 2f). Lastly, multidimensional scaling analyses of signature genes in the five cohorts revealed significant differences between tumor and non-tumor samples, implying that dysregulated JAK-STAT signaling may serve as a diagnostic marker for early detection in pre-cancerous lesions (Fig. 2F).
The JAK-STAT 28-gene signature is an independent prognostic factor
Multivariate Cox regression analyses confirmed that the signature was independent of tumor, node and metastasis (TNM) stage and age (where information is available): glioma (HR = 4.316, P < 0.0001), pan-kidney (HR = 2.379, P < 0.0001), clear cell renal cell (HR = 2.552, P = 0.00047), lung (HR = 0.636, P = 0.031) and endometrium (HR = 0.434, P = 0.033) (Additional file 3). Given that the signature was an independent predictor of overall survival, we reasoned that its predictive performance could be increased when used in conjunction with TNM staging. Employing the receiver operating characteristic (ROC) analysis, we demonstrated that a combined model uniting the signature and TNM staging could outperform the signature (higher area under the curve [AUC] values) when it was considered alone: pan-kidney (0.838 vs. 0.800), clear cell renal cell (0.836 vs. 0.789), lung (0.724 vs. 0.703) and endometrium (0.760 vs. 0.713) (Fig. 3a). Independently, Kaplan-Meier analyses and log-rank tests confirmed that the signature allowed further delineation of risk groups within similarly-staged tumors: pan-kidney (P < 0.0001), clear cell renal cell (P < 0.0001), lung (P < 0.0001) and endometrium (P < 0.0001) (Fig. 3b).
High JAK-STAT scores were associated with decreased survival rates in glioma patients. We confirmed that this was also true for histological subtypes of glioma: astrocytoma (P = 0.015) and oligoastrocytoma (P = 0.037) (Fig. 3b) but not for patients with glioblastoma multiforme (Additional file 4). Independently confirmed using Cox regression, patients within the 4th survival quartile had lower survival rates: astrocytoma (HR = 2.377, P = 0.018) and oligoastrocytoma (HR = 2.730, P = 0.038) (Additional file 3). In terms of the signature’s predictive performance, ROC analyses revealed that it performed the best in oligoastrocytoma (AUC = 0.951), followed by astrocytoma (AUC = 0.878) and glioma patients when considered as a full cohort (AUC = 0.853) (Fig. 3a).
Consequences of dysregulated JAK-STAT signaling and significant crosstalk with tumor immunity
Since dysregulated JAK-STAT signaling was associated with survival outcomes (Fig. 2 and Fig. 3), we reasoned that patients from diverse cancer types might harbor similar transcriptional defects caused by aberrant activation of JAK-STAT. Differential expression analyses performed between the 4th and 1st quartile patients revealed that a striking number of over 200 differentially expressed genes (DEGs) were shared between all five prognostic cohorts (Fig. 4a; Additional file 5). Significant overlaps were observed in DEGs where 555 genes were found in at least four cohorts, 1009 genes in at least three cohorts and 2034 genes in at least two cohorts (Fig. 4a; Additional file 5). The highest number of DEGs was observed in glioma (2847 genes), followed by pan-kidney (2810 genes), lung (1221 genes), clear cell renal cell (1161 genes) and endometrial (782 genes) cancers (Fig. 4a; Additional file 5). To determine their functional roles, the DEGs were mapped to Gene Ontology (GO) and KEGG databases. All five cohorts exhibited remarkably similar patterns of enriched biological processes (Fig. 4b). Pan-cancer enrichments of ontologies related to inflammation and immune function were observed, e.g., cytokine and chemokine signaling, T cell and B cell receptor signaling, natural killer cell-mediated processes, Toll-like receptor signaling and NOD-like receptor signaling (Fig. 4b). Additionally, genes associated with other oncogenic pathways (Wnt, MAPK, TGF-β, PPAR and VEGF) were frequently altered at both transcriptional and genomic levels (Fig. 4b and c) (Chang & Lai, 2019e; Chang & Lai, 2019f). CNV analyses performed on DEGs that were found to be in common in at least three cohorts demonstrated that transcriptional dysregulation of the aforementioned oncogenic pathways was attributed to activating or inactivating CNVs (Fig. 4c). For instance, except for THBS1 and THBS2, a vast majority of TGF-β DEGs exhibited somatic gains (Fig. 4c). To further explore which upstream transcriptional regulators were involved, we mapped the DEGs to ENCODE and ChEA databases (Fig. 4b). Interestingly, we observed significant enrichment of transcription factors (TFs) involved in modulating immune function: IRF8 (Mattei et al., 2012), RUNX1 (Ono et al., 2007), RELA (Bonizzi & Karin, 2004) and EZH2 (Vasanthakumar et al., 2017) (Fig. 4b). Taken together, our analyses revealed that pan-cancer JAK-STAT drivers underpin numerous aspects of tumor oncogenesis and immunity, which play important roles in tumor progression and ultimately patient prognosis.
IRF8 and JAK-STAT pathway synergistically influence survival outcomes in glioma and renal cancer
IRF8 was among one of the most enriched TFs implicated in the regulation of transcriptional outputs of patients with dysregulated JAK-STAT signaling (Fig. 4b). As a member of the interferon regulatory factor family, IRF8 is needed for the development of immune cells and is often regarded as a tumor suppressor gene since a loss of function is associated with increased metastatic potential (Lee et al., 2008; Yang et al., 2007). Thus, we predict that tumors with low expression of IRF8 would be more aggressive. To evaluate the combined relationship between JAK-STAT signaling and IRF8 expression, patients were categorized into four groups based on median IRF8 and JAK-STAT scores: 1) high-high, 2) low-low, 3) high IRF8 and low 28-gene score and 4) low IRF8 and high 28-gene score. The combined model encompassing JAK-STAT and IRF8 offered an additional resolution in patient stratification: glioma (full-cohort, P < 0.0001), astrocytoma (P = 0.0007), pan-kidney (P < 0.0001) and clear cell renal cell (P < 0.0001) (Fig. 5a). Indeed, patients with low IRF8 and high 28-gene scores performed the worst in cancers where hyperactive JAK-STAT signaling was linked to adverse survival outcomes: glioma (full-cohort: HR = 5.826, P < 0.0001), astrocytoma (HR = 3.424, P = 0.0032), pan-kidney (HR = 5.131, P < 0.0001) and clear cell renal cell (HR = 5.389, P < 0.0001) (Fig. 5b). Our results support a model in which IRF8 influences the behavior of tumors with aberrant JAK-STAT signaling.
Hyperactive JAK-STAT signaling attenuates regulatory T cell-mediated tumor immunity
Given the wide-ranging effects of JAK-STAT signaling on immune-related functions, we hypothesized that JAK-STAT activity would correlate with immune cell infiltration. We retrieved genes implicated in regulatory T cell (Treg) function from four studies and isolated 31 genes that were common in all four Treg signatures (De Simone et al., 2016; Plitas et al., 2016; Tirosh et al., 2016; Zheng et al., 2017). Treg scores were calculated for each patient based on the mean expression levels of the 31 genes. Remarkably, we observed strong positive correlations between Treg and JAK-STAT scores, suggesting that tumor tolerance was enhanced in patients with hyperactive JAK-STAT signaling: glioma (rho = 0.85, P < 0.0001), pan-kidney (rho = 0.88, P < 0.0001) and clear cell renal cell (rho = 0.83, P < 0.0001) (Fig. 6a). As in the previous section, patients were stratified into four categories based on median JAK-STAT and Treg scores for survival analyses. Log-rank tests and Cox regression analyses confirmed that elevated Treg activity further exacerbated disease phenotypes in patients where JAK-STAT scores were already high: glioma (full cohort, HR = 6.183, P < 0.0001), astrocytoma (HR = 3.035, P = 0.00042), pan-kidney (HR = 3.133, P < 0.0001) and clear cell renal cell (HR = 2.982, P < 0.0001) (Fig. 6b and c). Taken together, elevated JAK-STAT signaling may increase the ability of tumors to escape immunosurveillance, resulting in more aggressive tumors and increased mortality rates.
Tumor-promoting and tumor-suppressing roles of JAK-STAT signaling is very much cell type-dependent (Fig. 1, Fig. 2, Fig. 3). JAK-STAT activation appears to drive oncogenic progression in liver cancer (Calvisi et al., 2006) and infection with hepatitis viruses could also induce pathway activation (Arbuthnot et al., 2000; Yoshida et al., 2002). Moreover, liver tumors with downregulated SOCS proteins (inhibitors of JAK-STAT signaling) are associated with poor prognosis (Calvisi et al., 2006). JAK and STAT3 promote cell proliferation, invasion and migration in colorectal cancer through the regulation of cell adhesion molecules and growth factors (Xiong et al., 2008). In contrast, phosphorylated STAT5 promotes cellular differentiation and inhibits invasive properties in breast cancer cells (Sultan et al., 2008; Sultan et al., 2005). A decrease in STAT5 is also linked to poorly differentiated morphology and advance histological grades in breast tumors (Nevalainen et al., 2004; Peck et al., 2011). Similarly, in rectal cancers, patients with tumors positive for phosphorylated STAT3 had improved survival outcomes (Monnien et al., 2010). In contrast, high levels of phosphorylated STAT3 is associated with reduced survival rates in glioblastoma (Birner et al., 2010) and renal cancer (Horiguchi et al., 2002), which independently corroborates our findings on the tumor-promoting effects of JAK-STAT signaling in these cancer types (Figs. 2 and 3). Given its ambiguous role, understanding the function of JAK-STAT signaling in a pan-cancer context would be necessary to increase the success of therapy in tumors with abnormal pathway activity. In an integrated approach employing genomic, transcriptomic and clinical datasets, we elucidated pan-cancer patterns of JAK-STAT signaling converging on a core set of candidate driver genes known as the JAK-STAT 28-gene signature (Fig. 2). We demonstrated prognosis of the signature in five cancer cohorts (n = 2976), where its performance was independent of TNM stage (Figs. 2 and 3). It is worth noting that JAK-STAT scores were calculated based on mean normalized log2 expression of signature genes as an indication of pathway activity following previous reports, which may be sensitive to potential biases (Ge et al., 2018; Zhang et al., 2017). Patients with aberrant JAK-STAT signaling exhibited interactions with other major oncogenic pathways, including MAPK, Wnt, TGF-β, PPAR and VEGF (Fig. 4). This suggests that co-regulation of intracellular signaling cascades could have direct functional effects and combinatorial therapies simultaneously targeting these pathways may improve treatment efficacy and overcome resistance.
We demonstrated that hyperactivation of JAK-STAT signaling promotes the loss of anti-tumor immunity in glioma and renal cancer patients (Fig. 6). In tumor cells, constitutive activation of STAT3 inhibits anti-tumor immune response by blocking the secretion of proinflammatory cytokines and suppressing dendritic cell function (Wang et al., 2004). Furthermore, hyperactivation of STAT3 is linked to abnormal differentiation of dendritic cells in colon cancer cells (Nefedova et al., 2004). STAT3 promotes interleukin-10-dependent Treg function (Hsu et al., 2015) while STAT5 promotes Treg differentiation (Wei et al., 2008). We demonstrated that in glioma and renal cancer, JAK-STAT scores were strongly correlated with Treg expression scores, suggesting that persistent activation of JAK-STAT could promote tumor immune evasion in these cancers (Fig. 6). Moreover, in tumors with high levels of JAK-STAT signaling and low IRF8 expression (a TF involved in regulating innate and adaptive immune responses), we observed a dramatic decrease in overall survival rates (Fig. 5). IRF8 is essential for dendritic cell development and proatherogenic immune responses (Clément et al., 2018). Moreover, IRF8 is crucial for NK-cell-mediated immunity against mouse cytomegalovirus infection (Adams et al., 2018). Promoter hypermethylation and gene silencing of IRF8 abrogates cellular response to interferon stimulation and overexpression of IRF8 in nasopharyngeal, esophageal and colon cancer cell lines could inhibit clonogenicity (Lee et al., 2008). IRF8 expression is also negatively correlated with metastatic potential by increasing tumor resistance to Fas-mediated apoptosis (Yang et al., 2007). Together, our results and those of others support the tumor suppressive roles of IRF8. Importantly, loss of IRF8 may further suppress tumor immunity in patients with hyperactive JAK-STAT signaling. A number of JAK inhibitors (tofacitinib, ruxolitinib and oclacitinib) have been FDA-approved and along with 2nd-generation JAKinibs and STAT inhibitors currently undergoing testing (O’Shea et al., 2015), our signature can be used for patient stratification before adjuvant treatment with these inhibitors to enable selective targeting of tumors that are likely to respond.
What started as an initiative to understand JAK-STAT pathway genes that are somatically altered in varying combinations and frequencies across diverse cancer types has now resulted in a framework that supports selective targeting of novel candidates in a broad spectrum of cancers. Our study also reveals important crosstalk between JAK-STAT and other oncogenic pathways and when targeted together, this could radically improve clinical outcomes. Other researchers can harness this novel set of data and the JAK-STAT gene signature in the future for prospective validations in clinical trials and functional studies involving JAK-STAT inhibitors and immune checkpoint blockade.
Availability of data and materials
The datasets supporting the conclusions of this article are included within the article and its additional files.
Area under the curve
Copy number variation
Differentially expressed genes
Janus kinase-signal transducer and activator of transcription
Kyoto Encyclopedia of Genes and Genomes
Permutational multivariate analysis of variance
Receiver operating characteristic
Tumor, node and metastasis
Regulatory T cell
Adams NM, Lau CM, Fan X, Rapp M, Geary CD, Weizman O-E, et al. Transcription factor IRF8 orchestrates the adaptive natural killer cell response. Immunity Elsevier. 2018;48:1172–82.
Arbuthnot P, Capovilla A, Kew M. Putative role of hepatitis B virus X protein in hepatocarcinogenesis: effects on apoptosis, DNA repair, mitogen-activated protein kinase and JAK/STAT pathways. J Gastroenterol Hepatol Wiley Online Library. 2000;15:357–68.
Chang WH, Forde D, Lai AG. Dual prognostic role for 2-oxoglutarate oxygenases in ten diverse cancer types: implications for cell cycle regulation and cell adhesion maintenance. Cancer Commun Cold Spring Harbor Laboratory. 2019b;39.
Chang WH, Lai AG. Aberrations in Notch-Hedgehog signalling reveal cancer stem cells harbouring conserved oncogenic properties associated with hypoxia and immunoevasion. Br J Cancer. 2019a. https://doi.org/10.1038/s41416-019-0572-9.
Chang WH, Lai AG. An integrative pan-cancer investigation reveals common genetic and transcriptional alterations of AMPK pathway genes as important predictors of clinical outcomes across major cancer types. bioRxiv [Internet]. Cold Spring Harbor Laboratory; 2019b; Available from: https://www.biorxiv.org/content/early/2019/08/15/735647
Chang WH, Lai AG. Timing gone awry: distinct tumour suppressive and oncogenic roles of the circadian clock and crosstalk with hypoxia signalling in diverse malignancies. J Transl Med [Internet]. Cold Spring Harbor Laboratory; 2019c;17:132. Available from: https://www.biorxiv.org/content/early/2019/02/21/556878
Chang WH, Lai AG. Transcriptional landscape of DNA repair genes underpins a pan-cancer prognostic signature associated with cell cycle dysregulation and tumor hypoxia. DNA repair (Amst) [internet]. Elsevier; 2019d;78:142–153. Available from: http://www.ncbi.nlm.nih.gov/pubmed/31054516.
Chang WH, Lai AG. Pan-cancer genomic amplifications underlie a Wnt hyperactivation phenotype associated with stem cell-like features leading to poor prognosis. Cold Spring Harbor Laboratory: Transl Res; 2019e.
Chang WH, Lai AG. The pan-cancer mutational landscape of the PPAR pathway reveals universal patterns of dysregulated metabolism and interactions with tumor immunity and hypoxia. Ann N Y Acad Sci [internet]. Cold Spring Harbor Laboratory; 2019f;1448:65–82. Available from: https://www.biorxiv.org/content/early/2019/02/28/563676
Chen EY, Tan CM, Kou Y, Duan Q, Wang Z, Meirelles GV, et al. Enrichr: interactive and collaborative HTML5 gene list enrichment analysis tool. BMC Bioinformatics. BioMed Central. 2013;14:128.
Clément M, Haddad Y, Raffort J, Lareyre F, Newland SA, Master L, et al. Deletion of IRF8 (interferon regulatory factor 8)-dependent dendritic cells abrogates proatherogenic adaptive immunity. Circ Res Am Heart Assoc. 2018;122:813–20.
De Simone M, Arrigoni A, Rossetti G, Gruarin P, Ranzani V, Politano C, et al. Transcriptional landscape of human tissue lymphocytes unveils uniqueness of tumor-infiltrating T regulatory cells. Immunity [Internet]. 2016;45:1135–47. Available from:. https://doi.org/10.1016/j.immuni.2016.10.021.
Ge Z, Leighton JS, Wang Y, Peng X, Chen Z, Chen H, et al. Integrated Genomic Analysis of the Ubiquitin Pathway across Cancer Types. Cell Rep. 2018;23:213–226.e3.
Grivennikov S, Karin E, Terzic J, Mucida D, Yu G-Y, Vallabhapurapu S, et al. IL-6 and Stat3 are required for survival of intestinal epithelial cells and development of colitis-associated cancer. Cancer Cell Elsevier. 2009;15:103–13.
Horiguchi A, Oya M, Shimada T, Uchida A, Marumo K, Murai M. Activation of signal transducer and activator of transcription 3 in renal cell carcinoma: a study of incidence and its association with pathological features and clinical outcome. J Urol Elsevier. 2002;168:762–5.
Hsu P, Santner-Nanan B, Hu M, Skarratt K, Lee CH, Stormon M, et al. IL-10 potentiates differentiation of human induced regulatory T cells via STAT3 and Foxo1. J Immunol Am Assoc Immnol. 2015;195:3665–74.
Kuleshov M V, Jones MR, Rouillard AD, Fernandez NF, Duan Q, Wang Z, et al. Enrichr: a comprehensive gene set enrichment analysis web server 2016 update. Nucleic acids res. Oxford University Press; 2016;44:W90--W97.
Lee KY, Geng H, Ng KM, Yu J, Van Hasselt A, Cao Y, et al. Epigenetic disruption of interferon-γ response through silencing the tumor suppressor interferon regulatory factor 8 in nasopharyngeal, esophageal and multiple other carcinomas. Oncogene. Nat Publ Group. 2008;27:5267.
Mattei F, Schiavoni G, Sestili P, Spadaro F, Fragale A, Sistigu A, et al. IRF-8 controls melanoma progression by regulating the cross talk between cancer and immune cells within the tumor microenvironment. Neoplasia Elsevier. 2012;14:1223–IN43.
Mermel CH, Schumacher SE, Hill B, Meyerson ML, Beroukhim R, Getz G. GISTIC2. 0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers. Genome Biol. BioMed Central; 2011;12:R41.
Monnien F, Zaki H, Borg C, Mougin C, Bosset J-F, Mercier M, et al. Prognostic value of phosphorylated STAT3 in advanced rectal cancer: a study from 104 French patients included in the EORTC 22921 trial. J Clin Pathol BMJ Publishing Group. 2010;63:873–8.
Nefedova Y, Huang M, Kusmartsev S, Bhattacharya R, Cheng P, Salup R, et al. Hyperactivation of STAT3 is involved in abnormal differentiation of dendritic cells in cancer. J Immunol. Am Assoc Immnol. 2004;172:464–74.
Nevalainen MT, Xie J, Torhorst J, Bubendorf L, Haas P, Kononen J, et al. Signal transducer and activator of transcription-5 activation and breast cancer prognosis. J Clin Oncol American Society of Clinical Oncology. 2004;22:2053–60.
Peck AR, Witkiewicz AK, Liu C, Stringer GA, Klimowicz AC, Pequignot E, et al. Loss of nuclear localized and tyrosine phosphorylated Stat5 in breast cancer predicts poor clinical outcome and increased risk of antiestrogen therapy failure. J Clin Oncol. Proc Am Soc Clin Oncol. 2011;29:2448.
Poschke I, Mougiakakos D, Hansson J, Masucci GV, Kiessling R. Immature immunosuppressive CD14+ HLA-DR−/low cells in melanoma patients are Stat3hi and overexpress CD80, CD83, and DC-sign. Cancer res. AACR. 2010;70:4335–45.
Sultan AS, Brim H, Sherif ZA. Co-overexpression of Janus kinase 2 and signal transducer and activator of transcription 5a promotes differentiation of mammary cancer cells through reversal of epithelial-mesenchymal transition. Cancer Sci Wiley Online Library. 2008;99:272–9.
Sultan AS, Xie J, LeBaron MJ, Ealley EL, Nevalainen MT, Rui H. Stat5 promotes homotypic adhesion and inhibits invasive characteristics of human breast cancer cells. Oncogene. Nat Publ Group. 2005;24:746.
Tabas-Madrid D, Nogales-Cadenas R, Pascual-Montano A. GeneCodis3: a non-redundant and modular enrichment analysis tool for functional genomics. Nucleic Acids Res Oxford University Press. 2012;40:W478–83.
Tartour E, Pere H, Maillere B, Terme M, Merillon N, Taieb J, et al. Angiogenesis and immunity: a bidirectional link potentially relevant for the monitoring of antiangiogenic therapy and the development of novel therapeutic combination with immunotherapy. Cancer Metastasis Rev Springer. 2011;30:83–95.
Tirosh I, Izar B, Prakadan SM, Wadsworth MH, Treacy D, Trombetta JJ, et al. Dissecting the multicellular ecosystem of metastatic melanoma by single-cell RNA-seq. Science (80- ). American Association for the Advancement of Science. 2016;352:189–96.
Yang D, Thangaraju M, Greeneltch K, Browning DD, Schoenlein PV, Tamura T, et al. Repression of IFN regulatory factor 8 by DNA methylation is a molecular determinant of apoptotic resistance and metastatic phenotype in metastatic tumor cells. Cancer res. AACR. 2007;67:3301–9.
Yoshida T, Hanada T, Tokuhisa T, Kosai K, Sata M, Kohara M, et al. Activation of STAT3 by the hepatitis C virus core protein leads to cellular transformation. J Exp Med Rockefeller University Press. 2002;196:641–53.
Zheng C, Zheng L, Yoo JK, Guo H, Zhang Y, Guo X, et al. Landscape of Infiltrating T Cells in Liver Cancer Revealed by Single-Cell Sequencing. Cell [Internet]. 2017;169:1342–1356.e16. Available from: https://doi.org/10.1016/j.cell.2017.05.035
AGL designed the study and supervised the research. WHC and AGL analyzed the data, interpreted the results and wrote the initial manuscript draft. AGL revised the manuscript draft. Both authors approved the final submitted manuscript.
Differentially expressed genes between patients separated by the 28-gene signature into 4th and 1st quartiles.
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.