Abstract
Background/Aim: The FOXC2 transcription factor promotes the progression of several cancer types, but has not been investigated in the context of melanoma cells. To study FOXC2's influence on melanoma progression, we generated a FOXC2-deficient murine melanoma cell line and evaluated The Cancer Genome Atlas (TCGA) patient datasets. Materials and Methods: We compared tumor growth kinetics and RNA-seq/qRT-PCR gene expression profiles from wild-type versus FOXC2-deficient murine melanomas. We also performed Kaplan–Meier survival analysis of TCGA data to assess the influence of FOXC2 gene expression on melanoma patients' response to chemotherapy and immunotherapy. Results: FOXC2 promotes melanoma progression and regulates the expression of genes associated with multiple oncogenic pathways, including the oxidative stress response, xenobiotic metabolism, and interferon responsiveness. FOXC2 expression in melanoma correlates negatively with patient response to chemotherapy and immunotherapy. Conclusion: FOXC2 drives a tumor-promoting gene expression program in melanoma and is a prognostic indicator of patient response to multiple cancer therapies.
Melanoma, a highly aggressive form of cancer arising from pigment-producing melanocytes, is responsible for the majority of skin cancer-related mortality, accounting for ~60,000 annual deaths worldwide (1). Importantly, the incidence of melanoma has risen substantially over the last 40 years (2), a trend that is expected to continue to at least 2031 (3). Although surgical removal of primary lesions is typically successful in the treatment of early-stage disease, many melanoma patients are not diagnosed until later stages of metastatic disease in which surgery is either not possible or largely ineffective. Unfortunately, malignant melanoma is highly resistant to radiation and chemotherapy, and the only FDA-approved chemotherapeutic for the treatment of melanoma, dacarbazine (DTIC), has a minimal impact on patient survival (4). While advances in targeted therapy and immunotherapy have improved the prognosis for melanoma in recent years, there are still patients who do not respond to these regimens, and relapse of therapy-resistant tumors remains an ongoing challenge in many patients who do achieve clinical benefit (5, 6). Therefore, in order to improve the clinical outcome of melanoma patients going forward, it is necessary to gain additional insight into factors that promote melanoma progression and resistance to these therapeutic modalities.
FOXC2 is a member of the forkhead box family of transcription factors that control a variety of cellular processes in embryonic and adult tissues. In addition to its normal regulation of development, growth, and metabolism in various tissue types, FOXC2 has recently emerged as a driver of several hallmarks of cancer progression as well. Within vascular endothelial cells, FOXC2 promotes expression of multiple genes that enhance angiogenesis (7-9). FOXC2 can also become overexpressed or dysregulated in tumor cells themselves, where it is associated with aggressive characteristics and/or poor survival in patients with colorectal, ovarian, cervical, prostate, esophageal, oral, and basal-like breast cancers (10-16). In these and other studies with murine and human cancer cell lines, FOXC2 expression in tumor cells has been shown to promote proliferation, epithelial-mesenchymal transition (EMT), metastatic behavior, and drug resistance (10, 11, 15-19). Despite the significance of these oncogenic functions of FOXC2, however, much remains to be learned about the mechanisms underlying their induction. Moreover, FOXC2 has not been previously investigated in the context of melanoma cells, and little is known about whether this transcription factor plays a role in the progression of cancers arising from non-epithelial origin.
In this study, we describe a novel FOXC2-deficient murine melanoma cell line that we engineered through CRISPR-Cas9 gene editing. Using this model, we show for the first time that FOXC2 promotes the progression of melanoma, and we highlight data from our FOXC2-associated gene expression studies suggesting a role for this transcription factor in the regulation of melanoma cell responsiveness to oxidative stress, xenobiotics, and interferons (IFN). Finally, we provide clinical evidence for the prognostic value of FOXC2 gene expression in predicting melanoma patient response to chemotherapy and immunotherapy. Together, our work provides a strong rationale for evaluating FOXC2's influence on melanoma progression in larger patient cohorts, and this study highlights the utility of our murine model for investigating the biology of FOXC2's oncogenic functions in melanoma.
Materials and Methods
Mice. Female C57Bl/6 mice were purchased from Taconic Biosciences (Germantown, NY, USA) and used for experiments between the ages of 8-12 weeks. All experiments were approved by the Hampden-Sydney College Animal Care and Use Committee (Approval #113) and were performed in accordance with Guiding Principles in the Care and Use of Animals approved by the Council of the American Physiological Society.
Cell lines. B16-F1 murine melanoma cells purchased from the American Type Culture Collection (Manassas, VA, USA) were grown in RPMI-1640 medium supplemented with 2 mM L-glutamine, 2 g/l glucose, and 2 g/l sodium bicarbonate (Thermo Fisher Scientific, Waltham, MA, USA), as well as 10% fetal bovine serum (Premium Select, Atlanta Biologicals, Norcross, GA, USA). B16-F1ΔFOXC2 cells were generated as described below and maintained in the same growth media as the parental cell line. All cultures were grown at 37°C in a 5% CO2 incubator and passaged at 80-90% confluence.
Generation of B16-F1ΔFOXC2 melanoma cell line by CRISPR-Cas9 gene editing. Genome editing of the Foxc2 gene in B16-F1 melanoma cells was performed by transfecting cells with a pSpCas9 BB-2A-Puro (PX459) V2.0 all-in-one vector (GenScript, Piscataway, NJ, USA) encoding Cas9 and a gRNA sequence (5’ CTTCTTCTCTGGCGCGTTC 3’) targeting a region within the murine Foxc2 gene that encodes a portion of the N-terminal forkhead-binding domain of FOXC2. Transfection was performed with Attractene Transfection Reagent (Qiagen, Germantown, MD, USA) according to the manufacturer's protocol. Transfected cells were cultured for 24 h in the absence of antibiotics and then grown under puromycin selection (2 μg/ml) for 48 h, at which time puromycin-induced death of all cells in an untransfected control group was observed. The surviving polyclonal cell transfectants were then removed from puromycin selection and cloned by limiting dilution. Following expansion of clones, DNA was isolated with a DNeasy Blood and Tissue Kit (Qiagen) to screen for successful editing of the Foxc2 gene. Q5 Hot Start High Fidelity 2X Master Mix (New England BioLabs, Ipswich, MA, USA) and PCR primers (forward primer=5’ CATGCAGGCGCGTTACTC 3’; reverse primer=5’ ATAGCCCGCATACTGCACTGGTAG 3’) specific to regions of the Foxc2 gene flanking the gRNA target site were used to amplify a PCR product that was gel extracted using a QIAquick Gel Extraction Kit (Qiagen). Sanger sequencing was then performed on the purified PCR products using GenScript's DNA Sequencing Service in order to assess Foxc2 gene editing.
Tumor challenge. For tumor challenge experiments, 4×105 B16-F1 or B16-F1ΔFOXC2 melanoma cells were injected subcutaneously at the nape of the neck in a volume of 0.2 ml sterile, endotoxin-free 1X PBS (Teknova, Hollister, CA, USA). Mice were monitored daily for tumor formation, at which time tumor area was determined every 1-2 days using digital calipers to take perpendicular diameter measurements of the tumor. Mice were euthanized once tumor area reached >250 mm2 or tumors became necrotic.
Western blotting. To validate functional disruption of the Foxc2 gene in successfully edited clones, tumor cell pellets were flash frozen in a 95% ethanol/dry ice bath and shipped to Zyagen (San Diego, CA, USA) for their protein extraction and western blot analysis services. In brief, 10 μg of extracted protein was fractionated through an SDS-PAGE gel and transferred to a PVDF membrane. After blocking, membranes were incubated overnight at 4°C with a primary anti-FOXC2 antibody (sc-515472; Santa Cruz Biotechnology, Dallas, TX, USA) and then with an HRP-conjugated mouse IgG kappa binding protein (sc-516102, Santa Cruz Biotechnology) for 30 min, followed by chemiluminescent detection. Blots were then stripped and reprobed with an antibody optimized by Zyagen for detection of murine GAPDH.
RNA-sequencing. B16-F1 or B16-F1ΔFOXC2 melanoma cells were plated and cultured for 24 h to ~90% confluence before isolating RNA with a RNeasy Mini Kit (Qiagen) according to the manufacturer's recommendations. During extraction, on-column DNase-digestion was performed using Qiagen's RNase-free DNase Set. RNA was quantified with an Epoch Spectrophotometer (BioTek, Winooski, VT, USA) and shipped overnight on dry ice to Arraystar, Inc. (Rockville, MD, USA) for analysis using the company's Illumina Hi-seq 6G RNA-sequencing service. A260/280 and A260/230 ratios were both ≥2.0 for all samples. RNA integrity and genomic DNA contamination were examined by standard denaturing agarose gel electrophoresis, and all samples (5 independent replicates per group) passed quality control assessment.
To prepare sequencing libraries, mRNA was isolated from total RNA with oligo (dT) magnetic beads using the NEBNext® Poly(A) mRNA Magnetic Isolation Module (New England BioLabs). RNA fragmentation, random hexamer primed 1st strand cDNA synthesis, 2nd strand synthesis to incorporate dUTP into strand-specific libraries, end-repairing, A-tailing, adaptor ligation, and library PCR amplification were all performed using the KAPA Stranded RNA-Seq Library Prep Kit (Illumina). Completed libraries were qualified with an Agilent 2100 Bioanalyzer and quantified by absolute quantification qPCR. Barcoded libraries were mixed in equal amounts, denatured to single stranded DNA with 0.1M NaOH, loaded onto channels of the flow cell at 8 pM concentration, and amplified in situ using a TruSeq SR Cluster Kit v3-cBot-HS (Illumina). Sequencing was performed by running 150 cycles for both ends on an Illumina HiSeq 4000 instrument.
Generation of a FOXC2-deficient melanoma cell line. (A) CRISPR-Cas9 gene editing strategy for disruption of the Foxc2 gene in B16-F1 melanoma. Protein-coding domains of the murine Foxc2 gene were adapted from (7). (B) Sanger sequencing chromatograms of purified PCR products from an amplified region of the Foxc2 gene in wild-type B16-F1 melanoma and a clone of gene-edited cells derived from B16-F1. CRISP-ID deconvolution of overlapping chromatograms reveals biallelic out-of-frame indels in a clone of gene-edited cells. (C) Western blot analysis of FOXC2 and GAPDH proteins in wild-type B16-F1 and the gene-edited clone shown in (B) that verifies functional knockout of FOXC2 protein in the clone we designate B16-F1ΔFOXC2.
RNA-seq data processing and analysis. Image analysis and base calling were performed using Solexa pipeline v1.8 (Off-Line Base Caller software, v1.8). Sequence quality was examined using FastQC software, and raw sequencing data that passed Illumina chastity filtering were used for analysis. Fragments were 5’, 3’-adaptor trimmed and filtered ≤20 bp reads with cutadapt software. The trimmed reads were aligned to reference genome GRCm38 using Hisat 2 software. Transcript abundances for each sample were estimated with StringTie, and the normalized expression level (FPKM value) of known genes was calculated with the R package ballgown. The number of identified genes per group was calculated based on an FPKM mean ≥0.5 in that group. Differential gene expression analysis was performed with ballgown using the following cutoffs to filter differentially expressed genes: fold change ≥1.5, p-value ≤0.05, and FPKM ≥0.5 mean in at least one group. Gene ontology (GO) enrichment analysis of differentially expressed genes was performed using standard GO Terms from the Gene Ontology Resource (http://www.geneontology.org) and a Fisher's exact test to estimate statistical significance of the enrichment of terms between the B16-F1 and B16-F1ΔFOXC2 cell lines.
FOXC2 promotes the outgrowth of B16-F1 melanoma. C57Bl/6 mice were challenged subcutaneously with 4×105 B16-F1 or B16-F1ΔFOXC2 melanoma cells and monitored for tumor growth and progression. Images of representative mice at Day 12 post-challenge are shown in (A), and pooled data from 3 separate experiments, each with 2-4 mice per group, are shown in (B) and (C). The rate of tumor outgrowth was calculated as the size of tumor at the time of death divided by the number of days from the appearance of a measurable tumor to death. Data are graphed as the average of all replicates with error bars that represent standard deviation of the mean. ***p<0.001.
Data deposition. RNA-seq data discussed in this publication have been deposited in NCBI's Gene Expression Omnibus (20) and are accessible through GEO Series accession number GSE134296 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE134296). In addition to raw data in .fastq format, a matrix table of processed data (.xlsx format) with the normalized FPKM expression values for known genes from each sample is freely available via this accession number.
Quantitative RT-PCR. RNA was isolated as described above and reverse transcribed with the Qiagen RT2 First Strand Kit according to the manufacturer's protocol. The resulting cDNA was used for SYBR Green-based quantitative PCR (qPCR) with the following RT2 Profiler PCR Arrays from Qiagen: Mouse Oxidative Stress, Mouse Drug Metabolism, and Mouse Type I IFN Response (the latter of which also includes many genes associated with the IFNγ pathway). All reactions were performed on an Applied Biosystems StepOnePlus instrument (Thermo Fisher Scientific) under the following conditions: a 95°C pre-incubation stage for 10 min, followed by 40 cycles of a 95°C denaturation for 15 sec and a 60°C annealing/extension stage for 60 sec. SYBR Green incorporation into amplified DNA was recorded at the end of each 60°C annealing/extension step of the reaction, and a fluorescence threshold of 0.1 was set for determination of Ct values for each gene. Relative gene expression was calculated with Qiagen's GeneGlobe Data Analysis Center resource using the ΔΔCt method and normalization against the Gapdh housekeeping gene. Differentially expressed genes were considered significant if they passed the thresholds of fold change ≥1.5 and p≤0.01.
Kaplan–Meier survival analysis of The Cancer Genome Atlas (TCGA) datasets. The KM Plotter data analysis tool (21) was used to assess the prognostic value of FOXC2 gene expression in melanoma patients treated with either dacarbazine or ipilimumab. Beeswarm plots of RNA-seq expression data from TCGA were used to stratify melanoma patients from each treatment group into FOXC2low and FOXC2high cohorts for a univariate Cox regression analysis. Kaplan–Meier survival plots and hazard ratios with 95% confidence intervals as well as log-rank p-values were generated using the KM Plotter online platform (http://www.kmplot.com).
Results
Generation of a novel FOXC2-deficient murine melanoma cell line as a tool to study the role of FOXC2 in melanoma progression. The FOXC2 transcription factor has been implicated in the progression of multiple cancer types, but has not been previously investigated in the context of melanoma. In order to study whether FOXC2 plays a role in melanoma progression, we used CRISPR-Cas9 gene editing to knock out FOXC2 in the B16-F1 murine melanoma cell line. Cells were transfected with a plasmid encoding Cas9 and a gRNA designed to target a sequence within the Foxc2 gene that encodes a portion of the N-terminal forkhead-binding domain of FOXC2 (Figure 1A). Following transfection, clones were generated and screened for indels through Sanger sequencing, and we used CRISP-ID to deconvolute overlapping spectra from clones with biallelic gene edits resulting in indels too similar in size to be resolved by gel electrophoresis. As shown in Figure 1B, we identified a clone (Clone 43) with two out-of-frame indels in the targeted region of the Foxc2 gene, with one allele carrying a 2 bp insertion and the other carrying a 5 bp deletion. We confirmed functional disruption to the coding sequence of the Foxc2 gene in this clone by western blot analysis (Figure 1C), and we have named this FOXC2 protein-deficient cell line B16-F1ΔFOXC2.
To assess whether FOXC2 plays a role in the progression of melanoma, we challenged mice with wild-type B16-F1 or B16-F1ΔFOXC2 melanoma cells and monitored tumor outgrowth over time. B16-F1ΔFOXC2 tumor cells exhibited a slight delay in the onset of outgrowth, and tumors from mice challenged with these cells remained smaller than those in mice challenged with their wild-type counterparts at all time points analyzed (Figure 2A and B). Though both tumors ultimately progressed, the rate of B16-F1ΔFOXC2 progression was significantly slower than that for wild-type B16-F1 tumors. These data therefore demonstrate a critical role for the FOXC2 transcription factor in the progression of B16-F1 melanoma and highlight the utility of our novel B16-F1ΔFOXC2 model as a tool for investigating FOXC2's contribution to melanoma progression.
RNA-seq analysis of differential gene expression in B16-F1 versus B16-F1ΔFOXC2 melanomas. To gain insight into FOXC2's regulation of melanoma progression, we isolated RNA from B16-F1 versus B16-F1ΔFOXC2 melanoma cells and performed RNA-seq differential gene expression analysis. In order for FOXC2's contribution to the induction or repression of gene expression in B16-F1 melanoma to be most easily assessed, B16-F1ΔFOXC2 melanoma cells lacking this transcription factor served as the reference sample for this analysis. Volcano plot filtering of RNA-seq data identified 598 genes that are differentially expressed (fold-change ≥1.5, p≤0.05) by these cell lines: of these, 254 genes were up-regulated in B16-F1, implicating a role for FOXC2 in the expression of these genes, and 344 genes were down-regulated in B16-F1, reflecting a role for FOXC2 in the suppression of these genes (Figure 3A). While the complete details of this analysis will be published elsewhere, here we bring focus to select Gene Ontology (GO) Terms that were found to be significantly enriched with differentially expressed genes from these tumor cell lines. Interestingly, several Biologic Process-related GO Terms related to tumor progression, including those associated with the oxidative stress response, metabolism of xenobiotics, and IFN responsiveness, were enriched with genes from our differentially expressed cohort. In order to highlight the most significant genes associated with these particular GO Terms, we have presented only those relevant GO Terms that passed a statistical value threshold of p≤0.01 (Figure 3B), and we report the differentially expressed genes within these GO Terms that meet the statistical criteria of p≤0.01 and q≤0.05 (Table I). Of note, with respect to genes associated with oxidative stress- and xenobiotic metabolism-related GO Terms, several with direct functions in reactive oxygen species (ROS) and drug detoxification were up-regulated in B16-F1 melanoma. Of the differentially expressed genes associated with IFN response-related GO Terms, many IFN-stimulated genes as well as genes encoding interfering signaling pathway components were significantly down-regulated in B16-F1. Together, these data suggest that FOXC2 may positively and negatively regulate diverse pathways to ultimately achieve its tumor-promoting functions in melanoma.
Validation of RNA-seq data by qRT-PCR. To validate our RNA-seq data with a more sensitive qRT-PCR-based method for assessing gene expression, we used cDNA prepared from RNA isolates of the B16-F1 and B16-F1ΔFOXC2 melanoma cell lines for analysis with Pathway-Focused RT2 Profiler PCR Arrays relevant to oxidative stress, drug metabolism, and IFN responsiveness. While these arrays did not provide coverage of all the differentially expressed genes from our RNA-seq data presented in Table I, we did validate differential expression of many of these genes and also discovered additional differentially expressed genes from these pathways that were not identified by our RNA-seq approach. These data confirm significant up-regulation of several genes with antioxidant functions in B16-F1 (Table II), including the oxidoreductase Nqo1 (28.05-fold higher expression than in B16-F1ΔFOXC2) and various contributors to the glutathione detoxification system (Gstm4, Gstm5, Gstk1, Gstt1, and Gsr). These data also validate the down-regulation of many genes in B16-F1 involved in the type I IFN and IFNγ pathways (note that genes from both pathways are included in the Type I IFN Response RT2 Profiler PCR Array). The most highly down-regulated genes in B16-F1 identified by this array were IFN-stimulated genes (Oas1a=14.22-fold down-regulation, Oas1b=7.54-fold down-regulation, Isg15=5.96-fold down-regulation), though we also observed down-regulation of several genes involved in IFN signaling, including the Ddx58 gene that encodes RIG-I as well as the Stat1/Stat2/Stat3 and Irf7/Irf9 transcription factor genes. Collectively, these data provide additional support that the FOXC2 transcription factor is a key regulator of genes that control melanoma cell responses to oxidative stress, xenobiotics, and IFN.
FOXC2 expression predicts melanoma patient response to chemotherapy and immunotherapy. Based on FOXC2's regulation of genes associated with the response of melanoma cells to xenobiotics and IFN in our murine model, we wished to assess whether FOXC2 expression levels might have any influence on melanoma patient response to either chemotherapy or immunotherapy. Using TCGA datasets, we evaluated how FOXC2 gene expression in melanoma biopsies influenced progression-free survival (PFS) and overall survival (OS) in a small cohort of patients treated with dacarbazine (n=28 for PFS, n=29 for OS) or ipilimumab (n=22 for PFS and OS). As shown in the beeswarm plots in Figure 4, while little to no FOXC2 expression was observed in many of the patients from these cohorts, FOXC2 was expressed at higher levels in some patients, and we therefore manually set the expression cutoffs in each cohort to stratify patients into FOXC2low (black dots) and FOXC2high (red dots) populations. While there was no significant difference in PFS of these populations under dacarbazine treatment (Figure 4B), there was a statistically significant difference in OS, with a median OS of 67.67 months in the FOXC2low group and 47.93 months in the FOXC2high group (Figure 4C). In the ipilimumab-treated cohort, FOXC2 expression levels did correlate negatively with PFS (Figure 4D), with a median PFS of 24.17 months versus 8.27 months in the FOXC2low and FOXC2high groups, respectively (p=0.00016). Though the difference in OS between these groups did not reach statistical significance (p=0.18) in this small cohort (Figure 4F), a similar trend in clinical benefit was observed in the FOXC2low patient population, with a median OS of 62 months in this group compared to 34.4 months in the FOXC2high group. Collectively, these analyses provide rationale for assessing FOXC2 expression and clinical outcome in larger cohorts of melanoma patients undergoing these and other therapeutic regimens, as they highlight the potential prognostic value of FOXC2 expression in predicting patient response to widely used treatments for this cancer.
Discussion
Using a novel murine melanoma cell line engineered to lacking expression of the FOXC2 transcription factor, we showed for the first time that FOXC2 plays a critical role in the progression of melanoma. We also provide key insights into FOXC2's regulation of gene expression in melanoma cells, as our RNA-seq and qRT-PCR studies highlight FOXC2's influence over several genes associated with important tumor-promoting pathways, including the oxidative stress response, xenobiotic metabolism, and IFN responsiveness. The clinical significance of these findings is underscored by our analysis of TCGA RNA-seq data from melanoma patients, which demonstrate a negative correlation between FOXC2 expression in tumor biopsies and patient response to both chemotherapeutic and immunotherapeutic agents.
FOXC2 has previously been identified as a tumor-promoting transcription factor in several cancers of epithelial origin (22). Its oncogenic potential was first described in breast cancer, where up-regulation of FOXC2 was shown to correlate with high tumor grade and the aggressive, basal-like subtype of invasive ductal breast carcinoma (15). This seminal study and others to follow identified FOXC2 as a key regulator of EMT, a pathway key to tumor cell metastasis that is driven by FOXC2's direct repression of p120ctn and direct induction of Zeb1 (23, 24), both of which lead to down-regulation of E-cadherin, a major hallmark of EMT. In another study, FOXC2 has been shown to promote non-small cell lung cancer cell migration and invasion through up-regulation of Itgb1 gene expression (25). It has also been shown to drive osteosarcoma metastasis to lung tissue via a CXCR4-dependent mechanism (26). A large body of evidence has also recently demonstrated that FOXC2 promotes drug resistance in various cancer types, including osteosarcoma, nasopharyngeal carcinoma, and non-small cell lung, prostate, ovarian, colorectal, and basal-like breast cancers (16, 18, 19, 27-30).
To date, the only study investigating FOXC2 function in the context of melanoma comes from Sano et al., who evaluated the role of FOXC2 in tumor-associated vascular endothelial cells rather than in melanoma cells themselves (31). In their study, B16-F10 tumor growth was significantly reduced in Foxc2+/− haplodeficient mice as compared to wild-type hosts, and this was attributed to the diminished neovascularization of tumors observed in haplodeficient animals, suggesting a role for FOXC2 in tumor angiogenesis. Our study is the first to evaluate FOXC2 activity within melanoma cells, where we show that it also displays functions that contribute to tumor progression. Importantly, through comparative gene expression studies between the wild-type B16-F1 melanoma cell line and our novel B16-F1ΔFOXC2 model, we have identified a FOXC2-associated gene signature in B16-F1 that provides potentially key insights into this transcription factor's promotion of melanoma progression. Though we cannot distinguish direct versus indirect regulation of the differentially expressed genes that are linked to FOXC2 expression in our current study, we note that our B16-F1ΔFOXC2 cell line will be a powerful tool for future chromatin immunoprecipitation (ChIP) studies that aim to make this distinction, as it will serve as an ideal negative control for defining the background signal against which we can determine true enrichment of FOXC2-bound DNA in chromatin preparations from B16-F1 melanoma. Such studies will ultimately provide mechanistic insight into FOXC2's regulation of gene expression in melanoma. Regardless of which genes are found to be direct versus indirect targets of this transcription factor, though, the panel of differentially expressed genes identified in our current work reveals important new insights into FOXC2's tumor-promoting functions.
Between our RNA-seq and qRT-PCR analyses, we found that FOXC2 expression is associated with up-regulation of several genes involved in the oxidative stress response and xenobiotic metabolism. There was considerable overlap in the differentially expressed genes in B16-F1 versus B16-F1ΔFOXC2 identified from these pathways by our two approaches. Included among these genes were several members of the glutathione detoxification system and Nqo1 (28.05-fold higher expression in B16-F1 by qRT-PCR analysis), the latter of which encodes a two-electron oxidoreductase that is frequently up-regulated in aggressive cancers and confers protection against xenobiotics and toxic free radical/non-radical derivatives of ROS (32-34). Based on our findings, it will be interesting to investigate in our model how FOXC2 influences melanoma cell response to ROS-inducing stressors that are naturally encountered by progressing tumors, such as hypoxia and nutrient deprivation, as the ability to adapt to these harsh conditions in the tumor microenvironment might explain the enhanced rate of tumor progression observed in mice challenged with B16-F1 versus B16-F1ΔFOXC2. Additionally, melanoma is well-known for its resistance to many chemotherapeutic agents, some of which function in part by inducing ROS accumulation. With FOXC2's known role in chemoresistance in other cancer types and our TCGA analysis demonstrating that FOXC2 expression correlates negatively with melanoma patient response to dacarbazine, we are eager to explore in our model how FOXC2 contributes to chemotherapy resistance in melanoma. Our findings suggest that interfering with FOXC2 in melanoma might increase tumor cell susceptibility to chemotherapeutics, which has important implications not only for direct treatment of patients but also for developing strategies that could utilize such agents to promote immunogenic cell death as part of cancer vaccine regimens (35), a possibility we are currently exploring in our laboratory.
RNA-seq analysis of differential gene expression in B16-F1 versus B16-F1ΔFOXC2. (A) RNA was isolated from B16-F1 and B16-F1ΔFOXC2 melanoma cells for RNA-seq analysis. Volcano plot filtering was performed to identify differentially expressed genes (≥1.5-fold up- or down-regulation) that passed a statistical significance of p≤0.05. B16-F1ΔFOXC2 served as the reference sample so that genes up-regulated in B16-F1 could be interpreted as being positively regulated by FOXC2 and genes down-regulated in B16-F1 could be interpreted as being negatively regulated by FOXC2. Data are derived from sequencing of 5 separate RNA preparations for each melanoma cell line. (B) Select GO Terms related to the oxidative stress response, xenobiotic metabolism, and IFN responsiveness are presented from a GO enrichment analysis of differentially expressed genes in the B16-F1 versus B16-F1ΔFOXC2 cell lines.
RNA-seq differentially expressed genes from select GO terms in B16-F1 vs B16-F1ΔFOXC2 melanoma.
We are also particularly interested in the FOXC2-associated down-regulation of several IFN-related genes in our model. Consistent with our RNA-seq data, the most heavily down-regulated gene in B16-F1 that was identified by our IFN pathway-focused qRT-PCR array was Oas1a, whose human homolog is also down-regulated in aggressive breast and prostate cancers (36). Among other IFN-stimulated genes differentially expressed between B16-F1 and B16-F1ΔFOXC2, we also observed down-regulation in B16-F1 of several guanylate-binding protein-encoding genes, many of which correlate positively with survival when expressed at high levels in melanoma patients (37). Therefore, FOXC2-associated down-regulation of these and other IFN pathway genes in B16-F1 is likely to contribute to the rapid progression of this tumor as compared to its FOXC2-deficient counterpart. Interestingly, though IFN adjuvant therapy is an FDA-approved treatment for melanoma, its efficacy in this setting is limited, and it has minimal impact on OS of melanoma patients (38), possibly due to acquired resistance by tumor cells. The significance of such resistance is underscored by recent work showing that type I IFN signaling in tumor cells does control melanoma progression and enhances the therapeutic efficacy of both BRAF and PD-1 inhibitors (39, 40). Likewise, IFNγ responsiveness in melanoma is also key to the success of anti-PD-1 immunotherapy (41). IFNγ signaling in melanoma cells contributes to the efficacy of anti-CTLA-4 immunotherapy as well, as patients who fail to respond to checkpoint blockade with ipilimumab are defined by tumors with an accumulation of genomic defects in IFNγ signaling pathway genes (42), including many of those down-regulated by FOXC2 in our murine model. These and other studies highlight various mechanisms by which IFN signaling may be disrupted in tumor cells, including copy number alterations and single nucleotide variant/missense mutations or frameshift mutations of IFN pathway genes (42, 43), metabolic stressors such as hypoxia and nutrient deprivation that interfere with posttranslational activation of IFN pathway signaling components (44), elevated expression of IFN signaling inhibitors (45, 46), and loss-of-function mutations in proteins that directly interact with kinases essential for IFN signaling (47). In addition to these mechanisms, our current work suggests that epigenetic regulation of IFN pathway genes by transcriptional regulators might also contribute to dysregulation of IFN signaling in tumor cells, and we highlight the FOXC2 transcription factor as a novel suppressor of several IFN pathway genes in B16-F1 melanoma. The importance of these findings is highlighted by our demonstration that FOXC2 expression correlates negatively with PFS of melanoma patients treated with ipilimumab, and we are therefore eager to employ our murine model to gain insights into the FOXC2-mediated resistance to this and possibly other immunotherapies.
qRT-PCR analysis of differential gene expression from select pathways in B16-F1 vs. B16-F1ΔFOXC2.
FOXC2 expression correlates negatively with melanoma patient response to chemotherapy and immunotherapy. RNA-seq data from TCGA melanoma cases were used for Kaplan–Meier survival analyses of patients treated with dacarbazine or ipilimumab. Patients from these treatment groups were divided into FOXC2low and FOXC2high populations as shown in the beeswarm plots (A, D) for univariate Cox regression analysis, and plots of PFS (B, E) and OS (C, F) for each treatment cohort are shown.
Collectively, the data reported herein offer previously unappreciated insight into the role of the FOXC2 transcription factor in melanoma progression and describe a novel murine model for investigating the oncogenic functions of FOXC2 in melanoma. Importantly, though our clinical analyses were restricted to relatively small patient populations for which FOXC2 gene expression data was available, they provide strong rationale for examining FOXC2's influence on clinical outcome in larger patient cohorts in the future. The prognostic value of such data may not only prove useful in stratifying patient populations for appropriate treatment, but it may also lead to the design of combinatorial regimens that aim to improve the efficacy of existing therapies for melanoma by concurrent interference with FOXC2-mediated resistance mechanisms. Indeed, combinatorial strategies that interfere with oncogenic pathways in melanoma are already in use or currently being explored (48), and the introduction of new inhibitors that target FOXC2 or other components of FOXC2-associated pathways may expand the repertoire of therapeutic options for melanoma patients and enable further personalization of therapy for this cancer in the future.
Acknowledgements
Clinical analyses shown herein are based upon data generated by the TCGA Research Network: https://www.cancer.gov/tcga. This research was supported by funding from Virginia's Commonwealth Health Research Board, a Jeffress Trust Awards Program in Interdisciplinary Research Grant from the Thomas F. and Kate Miller Jeffress Memorial Trust (Bank of America, N.A., Trustee), a Virginia Foundation for Independent Colleges (VFIC) Maurice L. Mednick Memorial Fellowship, and a Hampden-Sydney College Research Grant from the Arthur Vining Davis endowment (to K. Hargadon). This work was also supported by VFIC Undergraduate Science Research Fellowships to J. Thompson and C. Williams and a Virginia Academy of Science Undergraduate Research Grant to D. Bushhouse. The Authors also thank Mr. Michael Hargadon and Mrs. Patricia Hargadon for generous donations to support the involvement of Hampden-Sydney College undergraduate students in this research.
Footnotes
Authors' Contributions
K.M. Hargadon designed the study, made contributions to all experimental work and interpretation of data presented herein, and wrote the manuscript. B. Győrffy developed the KM Plotter data analysis tool and assisted with Kaplan-Meier survival analyses of TCGA datasets. J. Thompson and D. Bushhouse contributed to the Foxc2 genome editing experiments. C. Williams and C. Johnson assisted with PCR amplification and DNA sequencing studies to verify editing of the Foxc2 gene in B16-F1 melanoma. E. Strong and B. Tarnai contributed to the qRT-PCR studies to validate RNA-seq data. All Authors agreed to the submission of this manuscript.
This article is freely accessible online.
Conflicts of Interest
The Authors declare no conflicts of interest regarding this study.
- Received August 5, 2019.
- Revision received August 18, 2019.
- Accepted August 19, 2019.
- Copyright© 2019, International Institute of Anticancer Research (Dr. George J. Delinasios), All rights reserved
