Introduction
Melanomas arising within the iris, ciliary body, and choroid of the eye make up the rare subset of uveal melanomas with approximately 2500 cases annually in the United States [1]. Despite aggressive management of primary uveal melanoma (pUM) with radiation or enucleation, approximately 50% of patients will develop metastasis and this primary disease phenotype is associated with specific gene expression patterns [2]. Metastasis of uveal melanoma preferentially develops in the liver and is thought to be driven by high levels of growth factors such as epidermal-, insulin-like-, and hepatic growth factors [1]. Survival of patients following metastasis is poor at 6–12 months and the treatment of uveal melanoma has not benefited by the recent improvements for melanoma more broadly [1,3]. Currently, no systemic therapies are associated with an improvement in survival. Therapies such as immune-checkpoint inhibitors [4,5] and targeted therapies [6] have demonstrated little benefit such that hepatic tumor embolization remains a default approach [1].
In contrast to cutaneous melanoma where mutations in BRAF, NRAS, NF1, and KIT genes are hallmarks of disease [3], pUM rarely if ever harbors these and instead commonly carries mutations in the G-protein α proteins q (GNAQ) or 11 (GNA11) [7,8]. Activated signaling via GNAQ/GNA11 is mediated through phospholipase C, eventually leading to outputs of the Yes-associated protein (YAP) [9,10] and mitogen-activated protein kinase (MAPK) pathway [7]. In addition to Gα genes, CYSTLR2 [11] and PLCB4 [12] act in the same pathway but are mutually exclusive with Gα. Additionally, several other genes have also been identified as recurrently mutated, dysregulated, or over-expressed such as BAP1, SF3B1, and EIF1AX [13].
The molecular biology of uveal melanoma, and associations with survival or treatment outcomes, have been detailed almost exclusively from pUM samples and cell lines derived from them. The most comprehensive description of uveal melanoma biology to date has been provided from The Cancer Genome Atlas (TCGA) in which analysis of 80 pUM suggested four molecular groups defined the presence or absence of chromosomal monosomy 3 (M3) vs. disomy 3 (D3) [13]. More specifically poor-prognosis M3-UM is associated with BAP1 loss, aberrant DNA methylation, somatic copy number aberrations, and RNA alterations relative to good-prognosis D3-UM, where alternative genomic and transcriptional changes are observed. Protein-based analyses have also suggested expression patterns associated with survival in pUM. For example, in a retrospective review of pUM samples, c-MET expression was associated with a higher risk of death from metastatic disease [14]. This led to the rationale for targeting MET kinase therapeutically, which despite preclinical support [15] was unsuccessful in a randomized study of cabozantinib as compared with chemotherapy (A091201) [16]. This experience highlights the need for molecular data associated with outcomes from tumors in the metastatic setting.
An inadequate literature describes the biology of human uveal melanoma metastases, although some studies have begun to report somatic mutational landscapes similar to that described in primary lesions. From a study of 52 metastatic samples obtained in a prospective clinical trial of metastatic disease [17], massively parallel sequencing of clinically relevant cancer genes identified mutations known from pUM within TCGA [13] and other databases [12] including GNAQ, GNA11, BAP1, PLCB4, and amplification of chromosome arm 8q. An analysis of 33 uveal melanomas (including nine metastases) described similar findings with low mutation frequency among similar genes, although also describing novel mutated genes such as CSMD1, TTC28, DLK2, and KTN1 at low frequency [18]. More recently, targeted sequencing of 500 genes in pUM and matched metastases from 35 patients revealed new driver mutations in genes such as CDKN2A, PBRM1, EZH2, PIK3R2, PIK3CA, PTEN, and MED12 with clonal and subclonal events [19]. Minimal data exists surrounding gene expression patterns in metastatic UM (mUM); however, analysis by TCGA primary tumors suggested four-clusters in pUM that separated M3- and D3-UM [13]. Pathway profiling of these data suggested the importance of signatures including DNA damage response, hypoxia, MYC signaling, and MAPK/AKT programs differentiated subgroups within both the M3-UM and D3-UM subgroups.
Although current systemic treatments are limited, there is some suggestion that emerging immune insights [20] and immuno-therapeutics may have utility [21–23]. Recently, single-cell sequencing of mUM samples suggested a complex landscape of tumor and immune cell interactions [24,25]. In this context, a better understanding of mUM biology may nominate molecular targets worth further exploration. Here, we investigated gene expression patterns in metastatic tumor specimens from patients with uveal melanoma (Supplementary Table 1, Supplemental digital content 1, https://links.lww.com/MR/A245) correlating these with survival and performing differential gene expression based on patient survival. These data nominate an important role for epithelial-mesenchymal transition (EMT) impacting outcomes for metastatic disease and nominate neuropilin-1 (NRP1) as a novel therapeutic target for future translational and potentially clinical investigation.
Methods
Tumor sample collection
Nineteen formalin-fixed paraffin-embedded (FFPE) metastatic tumor samples were collected and analyzed (Supplementary Table 1, Supplemental digital content 1, https://links.lww.com/MR/A245). Tumor percentage per sample was graded by a pathologist (T.C.) for 15 out of 19 samples, where sufficient tumor tissue remained after nucleic acid extraction. On average, 66% tumor content was observed per sample.
RNA sequencing of tumor samples
Total RNA was isolated from tumor FFPE tissue sections using QIAGEN Allprep DNA/RNA FFPE kit (Qiagen, Inc.) according to manufacturer’s instruction at the Human Immunologic Monitoring Facility at the University of Chicago. The quality and quantity of RNA was measured on an Agilent 2100 Bio-analyzer using Agilent reagents and protocols (Agilent Technologies, Santa Clara, USA). RNAseq libraries were generated using Illumina TruSEQ Total RNA stranded library making kits using Illumina protocols (Illumina, San Diego, USA). The quality and quantity of the library was determined using an Agilent 2100 Bio-analyzer using Agilent reagents and protocols. RNAseq libraries were sequenced on an Illumina HiSeq 4000 instrument using Illumina reagents and protocols at the University of Chicago Genomics Core Facility.
Gene expression quantification and differential gene expression detection
The quantification of gene-level expression and detection of differentially expressed genes (DEGs) were performed using a protocol similar as previously published [26]. In brief, the quality of raw reads was assessed by FastQC (v0.11.5) [27]. Transcript-level read counts were quantified by Kallisto (v0.44.0) [28] in a strand-specific mode with GENCODE [29] annotation of human reference transcriptome (v28, GRCh38). Kallisto implements a kmer-based pseudoalignment algorithm to accurately quantify transcripts from RNASeq data while robustly detecting errors in the reads. Transcript-level abundance was summarized into gene level by tximport (v1.4.0) [30], normalized by trimmed mean of M values (TMM) method, and log2-transformed. Out of 19 883 protein-coding genes in total, 14 307 genes that are expressed [defined as, counts per million of reads > 3] in at least six samples were kept for further analysis. Genes differentially expressed between groups were identified using limma voom algorithm with precision weights (v3.38.3) [31], followed by Benjamini-Hochberg (BH)-FDR correction for multiple testing [32].
Gene set enrichment analysis
Gene set enrichment analysis (GSEA) was performed by fGSEA (v1.8.0) [33]. The human Hallmark (H) gene sets from the Molecular Signatures Database (MSigDB) (v6.2) [34] was used for analyzing the full list of 14 307 genes from the differential gene expression analysis. Genes were preordered by the log2-transformed expression fold change metrics (log2FC). Enrichment nominal P-values were calculated by permutation test (10 000 permutations), with BH-FDR correction for multiple testing [32]. Enrichment score was computed same as in Broad GSEA implementation [34], and normalized enrichment score was computed by normalizing enrichment score to mean enrichment of random samples of the same size [33].
The Cancer Genome Atlas primary cancer cohorts
Raw sequencing FastQ files of tumor RNAseq data were downloaded from TCGA pUM cohort (primary tumors, n = 80) using the Genomic Data Commons data portal (https://portal.gdc.cancer.gov/) (accessed 9/17/2018). Demographic and clinical data of pUM were obtained from supplementary tables of the published TCGA study [13]. All RNAseq data were harmonized using the same analysis pipelines described above.
Survival analysis
The significance of association between molecular markers and patients’ overall survival (OS) was tested by Cox proportional hazard univariable or multivariable model using coxph function from R library survival (v2.42-3), which also estimates hazard ratio between patient groups. Patients alive upon last follow-up were censored.
Cell culture
Mel290 and OCM3 were received from Robert Folberg in 2009 (University of Illinois, Chicago, IL). Mel285 and OMM1.3 were kindly provided by Boris Bastian in 2010 (Memorial Sloan-Kettering Cancer Center, New York, NY). Mel290, Mel285, and OCM3 were established from primary tumors by Bruce Ksander (Schepens Eye Research Institute, Boston, MA) [35]. OMM1.3 was established from liver metastases also by Bruce Ksander [36]. OCM3 uveal melanoma cell lines have been sequenced for the presence of activating mutations in codons 209 (exon 5) and 183 (exon 4) of GNAQ and GNA11. OMM1.3 had GNAQ mutation while Mel290, Mel285, and OCM3 are wild-type for the gene. A karyotype test was also performed for each cell line in 2012. Cells were cultured in RPMI medium supplemented with 10% FBS, 100 units/ml penicillin, and 100 mg/ml streptomycin and maintained at 37°C in 5% CO2.
Gene silencing
The small interfering RNA (siRNA) sequences for scramble control (sc-37007) and NRP1 (sc-36038) were purchased from Santa Cruz Biotechnology (Dallas, TX). NRP1 siRNA was a pool of three siRNA duplexes:
5′-GAAGGCAGACAGAGAUGAA-3′
5′-CGAGACAGUCCAUUCAUCU-3′
5′-GAAGGAAAGCACUAAGAAA-3′
Mel290, Mel285, OMM1.3, and OCM3 cells were plated on 60-mm plates, and transfected with scramble control or NRP1 siRNA using Lipofectamine RNAiMAX (Invitrogen, Carlsbad, CA) according to the manufacturer’s protocol. The transfections were performed twice, each time in overnight incubations with a recovery phase of 6 h in between transfections.
Quantitative real-time PCR
Total RNA was harvested from cells transfected with either scramble control or NRP1 siRNA using PureLink RNA Mini Kit (Invitrogen, Carlsbad, CA) according to manufacturer’s instructions. cDNA was then synthesized from 1 µg total RNA using SuperScript IV First-Strand Synthesis System (Invitrogen, Carlsbad, CA). TaqMan mutation detection assays specific for NRP1 primer-probe set (Applied Biosystems, Foster City, CA) were used to quantitate mRNA expression of target genes normalized to GAPDH. Quantitative real-time PCR (qRT-PCR) assays were done using the 7500 Real-Time PCR System (Applied Biosystems, Foster City, CA). The relative quantity of genes was determined by the ΔΔCT method.
Immunoblotting
Cells transfected with either scramble control or NRP1 siRNA were lysed with radioimmunoprecipitation assay (RIPA) buffer supplemented with protease inhibitor cocktail tablets (Roche Diagnostics) and 1 mmol/l Na3VO4. Equal amounts of protein were loaded on 4–12% PAGE gels (Invitrogen, Carlsbad, CA). Polyvinylidene difluoride (PVDF) membranes were blocked with 5% nonfat dried milk and probed with NRP1, GAPDH (Cell Signaling Technology, Danvers, MA), and p27Kip1 (Santa Cruz Biotechnology, Dallas, TX) antibodies.
Flow cytometry
Cells transfected with either scramble control or NRP1 siRNA were washed and fixed in ice-cold 70% ethanol before staining with propidium iodide (50 μg/ml) containing RNase (5 μg/ml) to measure DNA content. Samples were sorted using LSR II Cell Analyzer (BD Biosciences, Franklin Lakes, NJ) for cell-cycle distribution and analyzed using the FCS Express 6 software. A total of 10 000 events were examined per sample.
Cell viability assays
Cells transfected with either scramble control or NRP1 siRNA were plated in 96-well plates in triplicates. Viability was assessed after 72 h of treatment using the Cell Counting Kit 8 (CCK8) from Dojindo Molecular Technologies (Rockville, MD) according to the manufacturer’s instructions. Survival is expressed as a percentage of untreated cells.
Cell invasion assays
Cells transfected with either scramble control or NRP1 siRNA were seeded in triplicates for 24 h in media with 0.1% serum on BioCoat Matrigel Invasion Chambers (BD Biosciences, San Jose, CA) according to manufacturer’s instructions. RPMI medium with 10% serum was used as chemoattractant. Noninvading cells were then removed from the matrigel and cells on the other side of the matrix were fixed with 100% methanol and stained with 1% Toluidine Blue. Images of stained cells were obtained from three random sections of each matrigel to account for cell distribution. Invading cells were then quantified by adding cells from the three sections and calculating the mean of each triplicate. Cell invasion is expressed as the number of cells migrated.
Statistics
For the analysis of gene signature expression between groups, two-sided Student’s t-test was used. Differential gene expression comparison between groups was performed using linear regression models implemented in limma voom with precision weights [31]. For multiple comparisons, P-value was adjusted using BH-FDR correction for multiple testing [32]. Statistical analysis was performed using R (v3.5.1) and Bioconductor (release 3.8). qPCR experiment was performed twice, and all in vitro experiments were carried out three times; P-values were calculated using two-sided Student’s t-test, standard error was calculated as the SD divided by the square root of the number of samples.
Study approval
Tumor samples were collected within the context of Alliance for Clinical Trials in Oncology A091201 [16], which was reviewed and approved by the NCI Central Institutional Review Board (CIRB) or the IRB of each participating site (ClinicalTrials.gov Identifier: NCT01835145). Each participant signed an IRB-approved, protocol-specific informed consent document in accordance with federal and institutional guidelines.
Results
Association of transcriptional programs with survival in uveal melanoma
Uveal melanoma has previously been described as a relatively simple genomic disease characterized by four overarching molecular subtypes [13]. Given the lack of therapeutic options for treatment of advanced disease, we were interested in investigating pathways and potential therapeutic targets that might be relevant to survival outcomes. To pursue this, we analyzed the association between cancer-associated signaling pathways and survival in the pUM cohort from TCGA. To reduce false discovery rate, we focused on 50 curated gene sets representative of cancer-associated pathways defined by human Hallmark (H) category from the Molecular Signatures Database (MSigDB; v6.2) (Supplementary Table 2, Supplemental digital content 1, https://links.lww.com/MR/A245). One score was assigned to each gene set by taking the average of expression across all genes involved in that gene set, and then testing for significant associations with OS. Using Cox proportional hazard univariable models, we observed higher expression score of 18 gene sets significantly associating with longer OS (FDR-adjusted P < 0.05; Supplementary Table 3, Supplemental digital content 1, https://links.lww.com/MR/A245). Using multivariable models, after adjusting for covariates including demographic factors (age and gender) and previously reported significant clinical variables from TCGA (histology and pigmentation) [13], we observed top gene sets (ranked by P-value from low to high) of MTORC1 signaling, glycolysis, PI3K/AKT/MTOR signaling, and IL6/JAK/STAT3 signaling (Supplementary Table 3, Supplemental digital content 1, https://links.lww.com/MR/A245). As the first three gene sets share some degree of functional redundancy, Fig. 1a conveys the impact of MTORC1 signaling as an example as well as the other distinct pathway of IL6/JAK/STAT3 signaling. The adjusted P-values and HRs for gene sets of interest and covariates are shown in Fig. 1b. Results of glycolysis and PI3K/AKT/MTOR signaling are provided in Supplementary Figure 1, Supplemental digital content 1, https://links.lww.com/MR/A246. As previous studies have suggested that pUM can be categorized into low-risk (Disomy 3, D3) and high-risk (Monosomy, M3) groups based on genomic signatures [13], we repeated the same analysis within the high-risk M3 group (42/80 patients), and observed no gene sets demonstrating significant association with OS by cox proportional hazard univariable or multivariable models at FDR-adjusted P < 0.05 (Supplementary Table 3, Supplemental digital content 1, https://links.lww.com/MR/A245).
;)
Fig. 1: Transcriptional programs associated with overall survival in primary uveal melanoma (pUM) from The Cancer Genome Atlas (TCGA). (a) Kaplan–Meier survival curves of MTORC1 signaling and IL6/JAK/STAT3 signaling gene expression in pUM, split by median expression of each signature (high vs. low). Survival risk table is shown below the Kaplan–Meier plot in each panel. (b) Forest plots showing the hazard ratio and P-values in Cox proportional hazards (PHs) multivariable model of the signaling pathways with demographic and clinical covariates. n = 80 patients in the TCGA primary UM cohort were shown for (a) and (b). Log-rank test was used in (a), and Cox PH multivariable model was used in (b).
Epithelial-mesenchymal transition signature separates 1-year overall survival in metastatic uveal melanoma
We recently completed a prospective phase II study in metastatic uveal melanoma demonstrating no difference in outcome between cabozantinib and chemotherapy (Alliance for Clinical Trials in Oncology A091201) [16]. From this study, we obtained adequate pretreatment tissue from a group of patients (n = 19) and pursued RNAseq analysis. To identify the molecular factors associated with survival in mUM, we investigated the 50 transcriptional programs outlined above comparing OS, as the time from the start of treatment in the trial until death from any cause, in A091201 [16]. Using Cox proportional hazard univariable model, EMT was the top gene program associated with OS ranked by P-values (P = 0.02), with higher EMT expression score predicting worse outcome (Supplementary Table 4, Supplemental digital content 1, https://links.lww.com/MR/A245). This association lost significance after adjusting for multiple testing correction (FDR-adjusted P = 0.61), potentially due to small size. Considering these findings from a clinical context, we then split patients into two groups centered on 1 year of survival from starting treatment on the clinical trial (Supplementary Figure 2, Supplemental digital content 1, https://links.lww.com/MR/A246), observing that patients with OS ≤ 1 year showed higher expression of the EMT gene expression signature (P = 0.15, Student’s t-test, two-sided) (Fig. 2a). Considering the OS is a time-to-event endpoint, we next tested the association between EMT gene signature expression as a continuous variable in Cox proportional hazard univariable model of OS and observed a significant association (P = 0.046, HR = 2.45) (Fig. 2b). Taking together, these data suggested that EMT may be a prognosis predictor in mUM.
Fig. 2: Epithelial-mesenchymal transition (EMT) signature is significantly associated with overall survival (OS) in patients with metastatic uveal melanoma. (a) Expression of the EMT gene signature in patient survival groups split by 1-year OS; n = 14 in patients with OS ≤ 1 year, n = 5 in patients with OS > 1 year. (b) Kaplan–Meier survival analysis of OS in patients with tumors EMThigh and EMTlow, split by median expression of the EMT gene signature. Nineteen patients were shown in (a) and (b), split by two different metrics [1-year OS in (a), and median EMT expression in (b)]. Survival risk table is shown below the plot. Two-sided Student’s t-test was used in (a). Cox proportional hazards univariable model was used in (b).
Patients with overall survival of less than 1 year demonstrated differential gene expression including expression of NRP1
Cognizant of the lack of effective systemic therapeutic options for patients facing mUM, we were interested to investigate possible targets associated with differential survival outcomes observed in A091201. To pursue this, we initially compared the whole transcriptome profiles of patients with OS ≤ 1 year to those with OS > 1 year. We ranked individual genes by P-value smaller to larger and selected top DEGs at P < 0.005 (unadjusted) and fold change ≥2.0 or ≤ –2.0 (Fig. 3a). A total of 76 genes passed the threshold, with 22 upregulated in patients who lived less than 1 year (Fig. 3b; Supplementary Table 5, Supplemental digital content 1, https://links.lww.com/MR/A245). Among those, NRP1, encoding a transmembrane glycoprotein neuropilin-1 that binds various vascular endothelial growth factor isoforms and TGFβ, has been shown to be associated with EMT in multiple tumor types [37–40]. Additionally, previous [41] and on-going clinical trials (clinicaltrials.gov: NCT03565445) are investigating monoclonal antibodies targeting neuropilin-1, raising therapeutic interest. Considering the overall lower tumor purity in the metastatic samples from A091201 compared to primary tumors from TCGA, we sought to confirm that differences observed in NRP1 gene expression were not due to tumor purity differences between groups and observed no significant tumor purity differences between mUM patient OS groups (P = 0.182, Student’s t-test, two-sided). Collectively, these analyses suggest expression of NRP1 in patients with OS > 1 year may represent disease-specific biological context.
Fig. 3: Differentially expressed genes of metastatic uveal melanoma in patient survival groups split by 1-year overall survival (OS). (a) Expression heatmap of 76 differentially expressed gens (DEGs) comparing tumors from patients who lived less than 1 year to those who lived longer. Genes were filtered by P < 0.005 (unadjusted) and fold change ≥2.0 or ≤–2.0. Samples were clustered on the column with dendrogram shown above the heatmap. Annotation bar labels patient groups with OS ≤ or >1 year. Genes are shown in the boxes to the right of the heatmap, following the same order as the gene dendrogram to the left side of the heatmap. (b) Expression of 22 DEGs upregulated in tumors from patients with OS ≤ 1 year relative to those with OS > 1 year. FC = expression fold of change calculated by comparing patients with OS >1 year to patients with OS ≤1 year. 1 yr = one year. The limma voom regression model with precision weights was used in (a) and (b).
To confirm our analysis by an orthogonal approach, we additionally investigated the overall pattern of signaling pathway changes between the two OS groups utilizing the full list of genes without P-value or fold change cutoffs. We calculated enrichment scores and significance level of 50 HALLMARK gene sets from MSigDB database. EMT genes and IL6/JAK/STAT3 signaling genes were identified to be significantly enriched and activated in patients with OS ≤ 1 year relative to >1 year (FDR-adjusted P < 0.05) (Supplementary Figure 3, Supplemental digital content 1, https://links.lww.com/MR/A246). This result supports our findings above investigating association between gene sets and survival using Cox proportional hazard analysis. Despite no genes passing FDR-adjusted P < 0.05 in multivariable models previously, we subsequently demonstrate that using the GSEA approach, the same transcriptional programs are significantly enriched in tumors from patients who lived less than 1 year.
Knockdown of NRP1 gene expression in vitro induces G1 arrest and inhibits cell proliferation and invasion
Observing association between NRP1 expression and OS < 1 year, we were interested to assess the functional impact of blocking NRP1 in uveal melanoma cells. To pursue this, we initially evaluated the gene expression of NRP1 in several uveal melanoma cell lines and identified at least two NRP1 expressing cell lines (Mel290, Mel285), and two with relatively low NRP1 mRNA expression (OMM1.3, OCM3) (Fig. 4a). Western blot analysis of uveal melanoma also showed greater protein expression of NRP1 in Mel290 and Mel285 relative to OMM1.3 and OCM3 (Fig. 4b). NRP1 has been shown to stimulate GIPC1 and Syx complex formation and activating RhoA thus inducing p27Kip1 degradation and that suppression of NRP1 expression by siRNA results in increased p27Kip1 expression [42]. In our study, siRNA knockdown of NRP1 markedly induced p27Kip1 expression compared to the scramble transfected control. Since p27Kip1 is induced by NRP1 suppression, we conducted a cell cycle analysis of uveal melanoma cells transfected with scramble control or NRP1 siRNA (Fig. 4c). We found that NRP1 suppression for 72 h significantly induced G1 arrest in Mel290 and Mel285 cells (P < 0.0001 and P < 0.05, respectively) while cells that do not express NRP1, OMM1.3 and OCM3 were not affected, when compared to scramble controls. Similarly, in viability assays, NRP1 siRNA for 72 h significantly inhibited cell proliferation in NRP1 expressing cell lines Mel290 (P = 0.009) and Mel285 (P = 0.003) (Fig. 4d). In contrast, NRP1 downregulation induced lesser effects in cells with low expression of NRP1, as cell viability was slightly reduced in OCM3 cells (P < 0.001) or increased in OMM1.3 cells (P = 0.18) (Fig. 4d). In cell invasion assays through matrigel, NRP1 siRNA knockdown significantly suppressed cell migration of Mel290 (P < 0.05) and Mel285 (P < 0.01) but not OMM1.3 and OCM3 after 24 h (Fig. 4e and f), conditions under which cell viability is not affected (data not shown).
Fig. 4: Knockdown of neuropillin-1 (NRP1) gene expression in vitro induces G1 arrest and inhibits cell proliferation and invasion. (a) NRP1 mRNA expression in uveal melanoma cell lines Mel290, Mel285, OMM1.3, and OCM3. n = 2 in each of the four cell lines, quantitative PCR was performed in triplicates, experiment was repeated two times. (b) Western blots of Mel290, Mel285, OMM1.3, and OCM3 cells transfected with NRP1 small interfering RNA (siRNA) shows inhibition of NRP1 expression and induction of p27Kip1 expression. (c) Cell cycle analysis of uveal melanoma cells transfected with NRP1 siRNA for 72 h shows a significant increase in G1 population in Mel290 (P < 0.0001) and Mel285 (P = 0.035) cells but not OMM1.3 and OCM3 cells; n = 3 in each of the four cell lines, with two groups each, experiment was repeated three times. (d and e) The siRNA knockdown of NRP1 expression significantly inhibits cell viability (d) after 72 h (Mel290 P = 0.009, Mel285 P = 0.003, OMM1.3 P = 0.18, OCM3 P < 0.001) and cell invasion (e) after 24 h of uveal melanoma cells that highly express NRP1 (Mel290 and Mel285). In (d), n = 3 in each of the four cell lines, with two groups each, experiment was repeated three times. (f) Quantitation of migrated cells shows the selective inhibition of invasion by NRP1 siRNA knockdown in Mel290 (P = 0.046) and Mel285 (P = 0.007) cells but not OMM1.3 (P = 0.40) and OCM3 (P = 0.42) cells. n = 3 in each of the four cell lines, with two groups each, experiment was repeated three times. In (c), (d), and (f), each bar is shown as mean ± S.E.M, with standard error calculated as the SD divided by the square root of the number of samples. Two-sided Student’s t-test was used in (c), (d), and (f); ****P < 0.0001, ***P < 0.001, **P < 0.01, *P < 0.05.
Discussion
Treatment options for patients with mUM are limited and outcomes are poor emphasizing a high, unmet need for patients facing this rare disease [1]. While the biology of pUM is increasingly understood through efforts such as TCGA, comparatively less has been known about metastatic disease. In our transcriptional analysis of patient samples derived from a prospective clinical trial, we identify transcriptional signatures associated with patient survival in metastatic disease and a number of DEGs associated with poor outcomes. Of these, NRP1 may be especially interesting as knockdown experiments in uveal melanoma cells suggest a functional role in proliferation, migration, and other processes relevant to metastasis and cancer progression in mUM.
We initiated our analysis focusing on pUM via analysis of TCGA observing growth signals through PI3K/mTOR and IL6 signaling being associated with progression and death from primary lesions. Multiple studies have described pUM as paradoxically demonstrating a tumor microenvironment densely packed by immune cells [1]; however, large-scale sequencing efforts also suggest the majority of pUM lack an active tumor microenvironment or T cell-inflamed phenotype [43]. Relevant to this, signaling through the PI3K cascade has been associated with immuno-suppression [44] and despite the lack of genomic alterations observed in pUM, transcriptional activation through PI3K may be driving this phenotype as has been described for other oncogenes such as β-catenin [45]. Additionally, the tumor microenvironment of pUM has been described as macrophage rich [1] and our finding of IL6-associated pathway activation is consistent with this. Considering these findings, it may be interesting to speculate on the utility of PI3K or mTOR inhibitor having a potential adjuvant role for patients with pUM after definitive therapy to limit metastasis and eventual death from disease.
Our analysis of transcriptional programs from patients with mUM suggests a potentially important role for EMT programs in dictating survival outcomes for this population of patients. An impact of EMT on survival in patients with metastatic cancer has been described across multiple tumor settings [46], although it previously had not been described in tumors from patients with mUM. Identification of the EMT process in uveal melanoma is potentially of relevance toward drug development and clinical trials, as approaches targeting EMT-related molecules such as β-catenin, TGF-β, and MYC have not been previously prioritized. It is also potentially noteworthy that EMT or stem-like states have been linked to the non-T-cell-inflamed tumor microenvironment and resistance to cancer immunotherapy [20,45,47,48]. In this setting, inhibitors of EMT may be useful as combination partners for immunotherapy, especially in mUM.
Our analysis of DEGs, in the context of patient survival, nominated a group of potentially intriguing genes that might be explored as therapeutic targets. While it should be noted that the sample size of the analysis was relatively small, mUM is an orphan disease and this is one of the first reports of transcriptional data from human tumor samples collected in a clinical trial. In our analysis, we observed that expression of NRP1 was associated with survival of patients of less than 1 year. NRP1 is particularly interesting as it has well-described functions relating to angiogenesis [49], EMT [37] as well as immune-modulation of Treg cells [50] and M2 macrophages [51]. Here, we particularly observed a cell-intrinsic role for NRP1 in regulating uveal melanoma proliferation, viability, and invasion, which could be reversed with NRP1 siRNA exposure in cell line experiments. A previous generation of anti-NRP1 monoclonal antibodies was brought unsuccessfully into cancer clinical trials in attempts to target angiogenesis in advanced solid tumors [41]. More recently, a new generation of NRP1 antibodies, including ASP1948 and others, has come forward with intent to be used as checkpoint blocking antibodies either alone or with anti-PD1 agents. These data suggest that consideration of a clinical trial in mUM might be prioritized.
We acknowledge that there are limitations to our report. As a rare cancer, investigation of patient-derived biospecimens of mUM is inherently limited by sample size. The tumor samples analyzed here were obtained in the context of a National Clinical Trials Network clinical trial (A091201) in which diagnostic tumor samples were obtained predominantly from the community practice setting. While extremely valuable, this approach has the potential to limit sample quality and limits the discovery power for differential gene expression analysis (e.g. multiple testing correction). This approach also eliminates the availability of further tissue to do confirmatory protein-based studies (e.g. immunohistochemistry). Obtaining biospecimens in this manner does have the utility of facilitating prospective clinical annotation and limiting the biases surrounding analysis of samples collected and analyzed in retrospective cohorts.
In summary, we have performed a transcriptional analysis of human mUM and observed novel associations with survival in the metastatic setting. We identified from the prospective clinical trial A091201 an impact of EMT on survival and both nominated and performed preliminary functional validation of NRP1 as a potential therapeutic target. These results immediately suggest avenues for novel investigation and clinical trials for patients with metastatic uveal melanoma. ClinicalTrials.gov Identifier: NCT01835145.
Acknowledgements
We thank M. Jarsulic for technical assistance on the high-performance computing clusters.
Research reported in this publication was supported by the National Cancer Institute of the National Institutes of Health under Award Numbers U10CA180821, U10CA180882, U24CA196171 (to the Alliance for Clinical Trials in Oncology), U10CA180836, UG1CA189960, and U10CA180863 (CCTG); Also supported in part by CCSG award P30CA047904, and funds from Exelixis. J.J.L. acknowledges Department of Defense Career Development Award (W81XWH-17-1-0265), National Cancer Institute Cancer Clinical Investigator Team Leadership Award (NIH P30CA014599), the Arthur J Schreiner Family Melanoma Research Fund, the J. Edward Mahoney Foundation Research Fund, Brush Family Immunotherapy Research Fund, Buffet Fund for Cancer Immunotherapy, the Sy Holzer Endowed Immunotherapy Research Award. J.J.L. and R.B. acknowledge funding from the Hillman Fellows for Innovative Cancer Research Program. D.J.O. acknowledges the Clinical Therapeutics Training Grant (NIH/NIGMS T32GM007019). The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health. https://acknowledgments.alliancefound.org.
Author contributions: J.L. conducted the study. J.L., B.B., M.B., D.H., P.M., and G.S. contributed patients to the clinical trial where tumor samples were collected. R.B. and J.L. developed the methodology. R.B. performed the bioinformatics analysis of sequencing data. O.S. performed the cell line experiments and analyzed the biological data. Y.Z. acquired the samples and processed samples for sequencing. T.C. performed the pathology review of the H&E slides. R.B., J.L., D.O., B.L., O.S., and G.S. interpreted the results and wrote the article. All authors reviewed the article.
Conflicts of interest
J.J.L. declares Scientific Advisory Board: (no stock) 7 Hills, Spring bank (stock) Actym, Alphamab Oncology, Arch Oncology, Kanaph, Mavu, Onc.AI, Pyxis, Tempest. Consultancy with compensation: Abbvie, Aligos, Array, Bayer, Bristol-Myers Squibb, Checkmate, Cstone, Eisai, EMD Serono, KSQ, Janssen, Merck, Mersana, Novartis, Partner, Pfizer, RefleXion, Regeneron, Ribon, Rubius, Silicon, Tesaro, Werewolf, Xilio, Xencor. Research Support: (all to institution for clinical trials unless noted) AbbVie, Agios (IIT), Array (IIT), Astellas, Bristol-Myers Squibb (IIT & industry), Corvus, EMD Serono, Immatics, Incyte, Kadmon, Macrogenics, Merck, Spring bank, Tizona, Xencor. Travel: Bristol-Myers Squibb, Janssen, Mersana, Pyxis. Patents: (both provisional) Serial #15/612,657 (Cancer Immunotherapy), PCT/US18/36052 (Microbiome Biomarkers for Anti-PD-1/PD-L1 Responsiveness: Diagnostic, Prognostic and Therapeutic Uses Thereof). RB: None declared, Patents: (all provisional) PCT/US15/612657 (Cancer Immunotherapy), PCT/US18/36052 (Microbiome Biomarkers for Anti-PD-1/PD-L1 Responsiveness: Diagnostic, Prognostic and Therapeutic Uses Thereof), PCT/US63/055227 (Methods and Compositions for Treating Autoimmune and Allergic Disorders). For the remaining authors, there are no conflicts of interest.
References
1. Luke JJ, Triozzi PL, McKenna KC, Van Meir EG, Gershenwald JE, Bastian BC, et al. Biology of advanced
uveal melanoma and next steps for clinical therapeutics. Pigment Cell Melanoma Res. 2015; 28:135–147
2. Onken MD, Worley LA, Ehlers JP, Harbour JW. Gene expression profiling in
uveal melanoma reveals two molecular classes and predicts metastatic death. Cancer Res. 2004; 64:7205–7209
3. Luke JJ, Flaherty KT, Ribas A, Long GV. Targeted agents and immunotherapies: optimizing outcomes in melanoma. Nat Rev Clin Oncol. 2017; 14:463–482
4. Luke JJ, Callahan MK, Postow MA, Romano E, Ramaiya N, Bluth M, et al. Clinical activity of ipilimumab for metastatic
uveal melanoma: a retrospective review of the Dana-Farber Cancer Institute, Massachusetts General Hospital, Memorial Sloan-Kettering Cancer Center, and University Hospital of Lausanne experience. Cancer. 2013; 119:3687–3695
5. Algazi AP, Tsai KK, Shoushtari AN, Munhoz RR, Eroglu Z, Piulats JM, et al. Clinical outcomes in metastatic
uveal melanoma treated with PD-1 and PD-L1 antibodies. Cancer. 2016; 122:3344–3353
6. Carvajal RD, Piperno-Neumann S, Kapiteijn E, Chapman PB, Frank S, Joshua AM, et al. Selumetinib in combination with dacarbazine in patients with metastatic
uveal melanoma: a phase III, multicenter, randomized trial (SUMIT). J Clin Oncol. 2018; 36:1232–1239
7. Van Raamsdonk CD, Griewank KG, Crosby MB, Garrido MC, Vemula S, Wiesner T, et al. Mutations in GNA11 in
uveal melanoma. N Engl J Med. 2010; 363:2191–2199
8. Van Raamsdonk CD, Bezrookove V, Green G, Bauer J, Gaugler L, O’Brien JM, et al. Frequent somatic mutations of GNAQ in
uveal melanoma and blue naevi. Nature. 2009; 457:599–602
9. Yu FX, Luo J, Mo JS, Liu G, Kim YC, Meng Z, et al. Mutant Gq/11 promote
uveal melanoma tumorigenesis by activating YAP. Cancer Cell. 2014; 25:822–830
10. Feng X, Degese MS, Iglesias-Bartolome R, Vaque JP, Molinolo AA, Rodrigues M, et al. Hippo-independent activation of YAP by the GNAQ
uveal melanoma oncogene through a trio-regulated rho GTPase signaling circuitry. Cancer Cell. 2014; 25:831–845
11. Moore AR, Ceraudo E, Sher JJ, Guan Y, Shoushtari AN, Chang MT, et al. Recurrent activating mutations of G-protein-coupled receptor CYSLTR2 in
uveal melanoma. Nat Genet. 2016; 48:675–680
12. Johansson P, Aoude LG, Wadt K, Glasson WJ, Warrier SK, Hewitt AW, et al. Deep sequencing of
uveal melanoma identifies a recurrent mutation in PLCB4. Oncotarget. 2016; 7:4624–4631
13. Robertson AG, Shih J, Yau C, Gibb EA, Oba J, Mungall KL, et al.; TCGA Research Network. Integrative analysis identifies four molecular and clinical subsets in
uveal melanoma. Cancer Cell. 2017; 32:204–220.e15
14. Mallikarjuna K, Pushparaj V, Biswas J, Krishnakumar S. Expression of epidermal growth factor receptor, ezrin, hepatocyte growth factor, and c-Met in
uveal melanoma: an immunohistochemical study. Curr Eye Res. 2007; 32:281–290
15. Wu X, Zhou J, Rogers AM, Jänne PA, Benedettini E, Loda M, Hodi FS. c-Met, epidermal growth factor receptor, and insulin-like growth factor-1 receptor are important for growth in
uveal melanoma and independently contribute to migration and metastatic potential. Melanoma Res. 2012; 22:123–132
16. Luke JJ, Olson DJ, Allred JB, Strand CA, Bao R, Zha Y, et al. Randomized phase II trial and tumor mutational spectrum analysis from cabozantinib versus chemotherapy in metastatic
uveal melanoma (alliance A091201). Clin Cancer Res. 2020; 26:804–811
17. Piperno-Neumann S, Kapiteijn E, Larkin JMG, Carvajal RD, Luke JJ, Seifert H, et al. Landscape of genetic alterations in patients with metastatic
uveal melanoma. J Clin Oncol. 2014; 3215_Suppl9043–9043
18. Royer-Bertrand B, Torsello M, Rimoldi D, El Zaoui I, Cisarova K, Pescini-Gobert R, et al. Comprehensive genetic landscape of
uveal melanoma by whole-genome sequencing. Am J Hum Genet. 2016; 99:1190–1198
19. Shain AH, Bagger MM, Yu R, Chang D, Liu S, Vemula S, et al. The genetic evolution of metastatic
uveal melanoma. Nat Genet. 2019; 51:1123–1130
20. Johnson DB, Bao R, Ancell KK, Daniels AB, Wallace D, Sosman JA, Luke JJ. Response to anti-PD-1 in
uveal melanoma without high-volume liver metastasis. J Natl Compr Canc Netw. 2019; 17:114–117
21. Sato T, Nathan PD, Hernandez-Aya LF, Sacco JJ, Orloff MM, Truscello J, et al. Intra-patient escalation dosing strategy with IMCgp100 results in mitigation of T-cell based toxicity and preliminary efficacy in advanced
uveal melanoma. J Clin Oncol. 2017; 3515_Suppl9531–9531
22. Chandran SS, Somerville RPT, Yang JC, Sherry RM, Klebanoff CA, Goff SL, et al. Treatment of metastatic
uveal melanoma with adoptive transfer of tumour-infiltrating lymphocytes: a single-centre, two-stage, single-arm, phase 2 study. Lancet Oncol. 2017; 18:792–802
23. Heppt MV, Amaral T, Kähler KC, Heinzerling L, Hassel JC, Meissner M, et al. Combined immune checkpoint blockade for metastatic
uveal melanoma: a retrospective, multi-center study. J Immunother Cancer. 2019; 7:299
24. Durante MA, Rodriguez DA, Kurtenbach S, Kuznetsov JN, Sanchez MI, Decatur CL, et al. Single-cell analysis reveals new evolutionary complexity in
uveal melanoma. Nat Commun. 2020; 11:496
25. Karlsson J, Nilsson LM, Mitra S, Alsén S, Shelke GV, Sah VR, et al. Molecular profiling of driver events in metastatic
uveal melanoma. Nat Commun. 2020; 11:1894
26. Godfrey J, Tumuluru S, Bao R, Leukam M, Venkataraman G, Phillip J, et al. PD-L1 gene alterations identify a subset of diffuse large B-cell lymphoma harboring a T-cell-inflamed phenotype. Blood. 2019; 133:2279–2290
27. Andrews S. FastQC: A quality control application for high throughput sequence data. Babraham Institute Project page.
http://wwwbioinformaticsbabrahamacuk/projects/fastqc 2016. [Accessed 13 September 2018]
28. Bray NL, Pimentel H, Melsted P, Pachter L. Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol. 2016; 34:525–527
29. Harrow J, Frankish A, Gonzalez JM, Tapanari E, Diekhans M, Kokocinski F, et al. GENCODE: the reference human genome annotation for The ENCODE Project. Genome Res. 2012; 22:1760–1774
30. Soneson C, Love MI, Robinson MD. Differential analyses for RNA-seq: transcript-level estimates improve gene-level inferences. F1000Res. 2015; 4:1521
31. Law CW, Chen Y, Shi W, Smyth GK. voom: precision weights unlock linear model analysis tools for RNA-seq read counts. Genome Biol. 2014; 15:R29
32. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Statist Soc. 1995; B:289–300
33. Sergushichev A. An algorithm for fast preranked gene set enrichment analysis using cumulative statistic calculation. bioRxiv. 2016060012
34. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005; 102:15545–15550
35. Zuidervaart W, van Nieuwpoort F, Stark M, Dijkman R, Packer L, Borgstein AM, et al. Activation of the MAPK pathway is a common event in uveal melanomas although it rarely occurs through mutation of BRAF or RAS. Br J Cancer. 2005; 92:2032–2038
36. Luyten GP, Mooy CM, De Jong PT, Hoogeveen AT, Luider TM. A chicken embryo model to study the growth of human
uveal melanoma. Biochem Biophys Res Commun. 1993; 192:22–29
37. Chu W, Song X, Yang X, Ma L, Zhu J, He M, et al. Neuropilin-1 promotes epithelial-to-mesenchymal transition by stimulating nuclear factor-kappa B and is associated with poor prognosis in human oral squamous cell carcinoma. PLoS One. 2014; 9:e101931
38. Luo M, Hou L, Li J, Shao S, Huang S, Meng D, et al. VEGF/NRP-1axis promotes progression of breast cancer via enhancement of epithelial-mesenchymal transition and activation of NF-κB and β-catenin. Cancer Lett. 2016; 373:1–11
39. Matkar PN, Singh KK, Rudenko D, Kim YJ, Kuliszewski MA, Prud’homme GJ, et al. Novel regulatory role of neuropilin-1 in endothelial-to-mesenchymal transition and fibrosis in pancreatic ductal adenocarcinoma. Oncotarget. 2016; 7:69489–69506
40. Ma L, Zhai B, Zhu H, Li W, Jiang W, Lei L, et al. The miR-141/neuropilin-1 axis is associated with the clinicopathology and contributes to the growth and metastasis of pancreatic cancer. Cancer Cell Int. 2019; 19:248
41. Weekes CD, Beeram M, Tolcher AW, Papadopoulos KP, Gore L, Hegde P, et al. A phase I study of the human monoclonal anti-NRP1 antibody MNRP1685A in patients with advanced solid tumors. Invest New Drugs. 2014; 32:653–660
42. Yoshida A, Shimizu A, Asano H, Kadonosono T, Kondoh SK, Geretti E, et al. VEGF-A/NRP1 stimulates GIPC1 and Syx complex formation to promote RhoA activation and proliferation in skin cancer cells. Biol Open. 2015; 4:1063–1076
43. Spranger S, Luke JJ, Bao R, Zha Y, Hernandez KM, Li Y, et al. Density of immunogenic antigens does not explain the presence or absence of the T-cell-inflamed tumor microenvironment in melanoma. Proc Natl Acad Sci U S A. 2016; 113:E7759–E7768
44. Peng W, Chen JQ, Liu C, Malu S, Creasy C, Tetzlaff MT, et al. Loss of PTEN promotes resistance to T cell-mediated immunotherapy. Cancer Discov. 2016; 6:202–216
45. Luke JJ, Bao R, Sweis RF, Spranger S, Gajewski TF. WNT/β-catenin pathway activation correlates with immune exclusion across human cancers. Clin Cancer Res. 2019; 25:3074–3083
46. Dongre A, Weinberg RA. New insights into the mechanisms of epithelial-mesenchymal transition and implications for cancer. Nat Rev Mol Cell Biol. 2019; 20:69–84
47. Spranger S, Bao R, Gajewski TF. Melanoma-intrinsic β-catenin signalling prevents anti-tumour immunity. Nature. 2015; 523:231–235
48. Hugo W, Zaretsky JM, Sun L, Song C, Moreno BH, Hu-Lieskovan S, et al. Genomic and transcriptomic features of response to anti-PD-1 therapy in metastatic melanoma. Cell. 2016; 165:35–44
49. Soker S, Takashima S, Miao HQ, Neufeld G, Klagsbrun M. Neuropilin-1 is expressed by endothelial and tumor cells as an isoform-specific receptor for vascular endothelial growth factor. Cell. 1998; 92:735–745
50. Delgoffe GM, Woo SR, Turnis ME, Gravano DM, Guy C, Overacre AE, et al. Stability and function of regulatory T cells is maintained by a neuropilin-1-semaphorin-4a axis. Nature. 2013; 501:252–256
51. Roy S, Bag AK, Singh RK, Talmadge JE, Batra SK, Datta K. Multifaceted role of neuropilins in the immune system: potential targets for immunotherapy. Front Immunol. 2017; 8:1228
