The Prognostic Significance of Low-Frequency Somatic Mutations in Metastatic Cutaneous Melanoma

Background: Little is known about the prognostic significance of somatically mutated genes in metastatic melanoma (MM). We have employed a combined clinical and bioinformatics approach on tumor samples from cutaneous melanoma (SKCM) as part of The Cancer Genome Atlas project (TCGA) to identify mutated genes with potential clinical relevance.

Methods: After limiting our DNA sequencing analysis to MM samples (n = 356) and to the CANCER CENSUS gene list, we filtered out mutations with low functional significance (snpEFF). We performed Cox analysis on 53 genes that were mutated in ≥3% of samples, and had ≥50% difference in incidence of mutations in deceased subjects versus alive subjects.

Results: Four genes were potentially prognostic [RAC1, FGFR1, CARD11, CIITA; false discovery rate (FDR) < 0.2]. We identified 18 additional genes (e.g., SPEN, PDGFRB, GNAS, MAP2K1, EGFR, TSC2) that were less likely to have prognostic value (FDR < 0.4). Most somatic mutations in these 22 genes were infrequent (< 10%), associated with high somatic mutation burden, and were evenly distributed across all exons, except for RAC1 and MAP2K1. Mutations in only 9 of these 22 genes were also identified by RNA sequencing in >75% of the samples that exhibited corresponding DNA mutations. The low frequency, UV signature type and RNA expression of the 22 genes in MM samples were confirmed in a separate multi-institution validation cohort (n = 413). An underpowered analysis within a subset of this validation cohort with available patient follow-up (n = 224) showed that somatic mutations in SPEN and RAC1 reached borderline prognostic significance [log-rank favorable (p = 0.09) and adverse (p = 0.07), respectively]. Somatic mutations in SPEN, and to a lesser extent RAC1, were not associated with definite gene copy number or RNA expression alterations. High (>2+) nuclear plus cytoplasmic expression intensity for SPEN was associated with longer melanoma-specific overall survival (OS) compared to lower (≤ 2+) nuclear intensity (p = 0.048). We conclude that expressed somatic mutations in infrequently mutated genes beyond the well-characterized ones (e.g., BRAF, RAS, CDKN2A, PTEN, TP53), such as RAC1 and SPEN, may have prognostic significance in MM.

Introduction

The role of genetic aberrations in SKCM development and progression is unquestionable (1). Several studies in primary melanomas have shown the prognostic significance of frequent somatic mutations, such as NRAS and BRAF (2–4). Other studies have been more heterogeneous with respect to inclusion of both primary and metastatic melanoma samples (5, 6). Other genes have more complex genetic and epigenetic aberrations; therefore, non genetic-based assays have been used to assess prognostic significance (7). Even if studies are focused in metastatic melanoma (MM) specimens and for the most abundant BRAF and NRAS mutations, the prognostic significance in MM is less understood and, in some cases, controversial (8–10). Advances in next generation sequencing methodology have inspired the development of targeted sequencing mutation panels to assist clinicians toward personalized treatment decisions (11, 12). It is now commonplace that medical oncologists order targeted sequencing panels to identify genetic aberrations that may be predictive of response to anticancer therapies. With the exception of hematologic malignancies (13), however, sequencing panels have been infrequently used in solid tumors to assess prognosis. Is it possible that other, less well-characterized mutated genes are important for the understanding of genetics in MM as well as prognosis?

In this study, we explored the prognostic significance of somatic mutations in the TCGA SKCM cohort. Given the high somatic mutation burden in SKCM (8), we focused on mutations with intermediate/high functional impact. We have restricted our analysis to samples procured from MM for three reasons: (1) metastatic tumors are more likely to have mutations in driver genes (14); (2) passenger mutations that are associated with ultraviolet signature may be significantly less in MM (15); and (3) although not all patients with primary melanoma will succumb to their disease, patients with MM have a worse prognosis. Finally, we explored the clinical significance of our findings on genetic aberrations of MM in an independent cohort. We have identified that genetic aberrations in less characterized genes, such as RAC1 and SPEN, may have prognostic significance.

Patients and methods

The TCGA Cohort

Patients, Clinical Data, and Molecular Classification

Clinical, pathologic, and follow-up data were retrieved from the TCGA (Supplemental Material). Follow-up time from both the original diagnosis of melanoma (OS-original diagnosis) and from time-to-specimen collection to latest follow-up (OS-specimen collection) were recorded. Each specimen was subsequently classified according to the previously reported molecular and gene expression profiling classifications (Supplemental Material) (8).

Next Generation Sequencing

The Variant Call Format (VCF) files were retrieved from the TCGA access-controlled data portal (https://cbiit.cancer.gov) using the reference GRCh37/hg19 genome assembly. Somatic mutation calls were annotated by the UNCseq™ pipeline (version 2016.07): snpEFF (v4.3) (16) was used to annotate the variant calls to ascertain gene information [gene id, HGVS_C for Variant in Human Genome Variation Society (HGVS) DNA notation; HGVS_P for Variant in HGVS protein notation] and functional impact [IMPACT for HIGH, MODERATE, LOW, MODIFIER, EFFECT for Effect in Sequence Ontology terms, http://snpeff.sourceforge.net/SnpSift.html]. The variants were also annotated with the Catalog of Somatic Mutations in Cancer (COSMIC, v77) and ExAC (v0.3) to obtain the allele frequency of the variants.

Identification of Mutations With Potentially High Biological Impact

Multi-step filters for somatic mutation calling were applied to identify non-synonymous mutations in potentially important genes (Supplemental Material). To ensure a high quality mutation calling, only mutations with mutant allele frequency (MAF) ≥5% and mutated allele count (MAC) ≥5 were kept (Figure S1). To further refine the gene list, the following clinically driven decisions were made: first, multiple mutations per gene per given specimen were counted once; second, only genes mutated in ≥3% of the total number of patient samples were considered [n = 128 (OS-original diagnosis analysis), n = 133 (OS-specimen collection analysis)]; third, only genes that were more or less frequently mutated in deceased vs. living subjects at the time of last follow-up were considered. By convention, these genes were defined by the fold of the mean frequency of samples mutated in deceased vs. living patients ($\text{x}=\frac{\frac{\text{No pts with a given mutated gene in deceased pts}}{\text{ Deceased pts}}}{\frac{\text{No pts with a given mutated gene in living pts}}{\text{Living pts}}}$) which should be either ≥1.5 or ≤ 0.67. Time to last follow-up was not taken into consideration at this stage of gene list refinement.

Univariate Cox regression analysis was subsequently performed on the final list of mutated genes to investigate the impact of the mutation frequency for each gene on OS. The null hypothesis was that the mutation frequency for each gene does not correlate with OS. Mutated genes with per gene-adjusted FDR ≤ 0.4 were considered promising.

Both DNA and RNA sequencing reads [per sample binary version of SAM files (BAM) against the hg19 assembly] were obtained to perform UNCeqR (v0.05) analysis (17), and confirm the somatic mutations of prognostically important genes.

The Multi-Institution Validation Cohort for the TCGA Data

This multi-institution cohort consists of published datasets whose MM tumors were subjected to whole exome (18–20), genome (21), or targeted panel sequencing (22). Updated follow-up data were obtained from each dataset, if possible (18, 19, 21). We also added the MM cohort from the University of North Carolina at Chapel Hill (UNC-CH), whose tumors underwent targeted panel sequencing (IRB protocol 11-1115) (23), and clinical data were collected (IRB protocol 16-2959; the UNC-CH UNCseq™ MM cohort). Only non-synonymous somatic mutations identified in tumor samples procured from MM of cutaneous or unknown primary (i.e., no cell lines established from tumors) were included.

The UNC-CH Metastatic Melanoma Cohort to Investigate SPEN Protein Expression in Metastatic Melanoma

To investigate protein expression of SPEN in relation to several clinical and histopathologic features in MM, we performed single-color immunohistochemical analysis of the SPEN protein expression in MM tissue samples from patients who were treated at the UNC-CH Melanoma Program. Collection of MM tumor tissues, clinicopathologic data, and clinical follow-up were allowed under the UNC-CH IRB approved protocol 09-1737. These MM tissue samples were different from those included in the UNC-CH UNCseq™ MM cohort, which was described above. Representative 0.7-mm diameter tissue cores from each MM specimen were spotted in duplicates or triplicates in the UNC-CH MM tissue microarray (TMA).

A commercially available anti-SPEN antibody (polyclonal rabbit, product #HPA-15825, Lot# A105440, Sigma-Aldrich, St. Louis, MO) (24) was optimized in normal human tissues and a nevus-primary melanoma TMA prior to staining the UNC-CH MM TMA. Briefly, 5 μm-thick tissue sections from the UNC-CH MM TMA were baked at 60°C for 90 min followed by heat-induced epitope retrieval using HIER Buffer L (Thermo Scientific, TA-135-HBL, Waltham, MA). Tissues were blocked using 10% normal goat serum for 1 h at room temperature, then incubated with the anti-SPEN antibody overnight at 4°C. Following incubation with a biotinylated goat anti-rabbit IgG (1:500; Jackson ImmunoResearch Laboratories Inc., 111-065-144, West Grove, PA) for 60 min at room temperature, tissues were treated with ABC-AP (Vector Laboratories, AK-5000, Burlingame, CA) and Impact Vector Red (Vector Laboratories, SK-5105). Finally, the tissues were counterstained with hematoxylin (Thermo Scientific, 6765003), dehydrated, cleared, and coverslipped using DPX (Electron Microscopy Sciences, 13512, Hatfield, PA).

Histopathologic analysis was performed using the 0, 1+ (< 25% of melanoma cells with nuclear stain), 2+ (25–80% of melanoma cells with nuclear stain), and 3+ (>80% of melanoma cells with strong stain) semiquantitative scale (PBG). Each tissue core was also evaluated for the presence or absence of any small-diameter mononuclear cells with stippled chromatin and indistinctive cytoplasm suggestive of lymphocytes. Contingency tables and the chi-square test were used to assess the association between expression intensity of SPEN protein by melanoma cells and the presence vs. absence of tumor-infiltrating lymphocytes (R project, www.r-project.org). Kaplan-Meier curves were constructed using GraphPad Prism (v 8.0, GraphPad, La Jolla, CA) to estimate the melanoma-specific overall survival (OS) in patients with high (>2+) nuclear plus cytoplasmic nuclear-only signal vs. high (>2+) nuclear-only signal vs. low (≤ 2+) expression of SPEN protein by melanoma cells. Given the multiple replicates per patient, the average of the expression intensity was considered. In addition, a specimen was classified as having both nuclear and cytoplasmic signal if at least one of the replicate cores exhibited nuclear plus cytoplasmic localization for the SPEN protein.

Results

Bioinformatics Analysis of Somatic Mutations in TCGA MM Samples

Figures 1, 2 show details of our filtering strategy. 25,102,889 mutations, both somatic and presumed germline, were identified in 474 tumors. For the OS-original diagnosis survival analysis, we excluded samples from stage 0 melanoma, whereas for the OS-specimen collection analysis, we excluded patients with a discrepant biopsy date. For the OS-original diagnosis analysis, 5,351 mutations (0.22% total; 5,285 single nucleotide variations, 37 dinucleotide variants, and 29 insertions/deletions) in 537 genes were found in 356 tumors [median, 10 mutations, ~95% confidence interval (~95CI) 9–11 mutations; range, 0–132 mutations]. At a median follow-up of 53.0 months (~95CI, 47.2–59.1 months), 185 (51.8%) patients were deceased. We observed similar findings in the OS-specimen collection analysis [5,481 mutations in 541 genes were found in 363 tumor specimens; median follow-up 53.0 months (~95CI, 47.8–59.2 months); 190 (52.3%) patients were deceased]. Table 1 shows demographics, molecular, and gene expression classification of the 356 MM patients.

FIGURE 1

Figure 1. Strategy to identify important somatic mutations in the TCGA SKCM cases. Flow diagram shows important clinical steps, such as removal of cases procured from patients with primary melanoma and lack of long-term follow-up. For the OS-original diagnosis analysis, we excluded patients who were originally diagnosed with AJCC stage 0. For the OS-specimen collection analysis, time was calculated from time of specimen collection.

FIGURE 2

Figure 2. Flow diagram shows number of mutations that are filtered out during each step for the OS-original diagnosis analysis. Similar results were seen for the OS-specimen collection analysis. Abbreviations: pts, patients; mts, mutations; mdn, median; AJCC, American Joint Committee on Cancer; dx, diagnosis; AF, allele frequency; OxoG, oxidative artifact from oxidation of guanine to 8-oxoguanine during DNA library preparation; COSMIC, catalog of somatic mutations in cancer.

TABLE 1

Table 1. Patient characteristics of the TCGA SKCM MM tissue cohort (n = 356).

To assess whether our filtering strategy retained known somatic mutations, we tested the impact of this strategy on the hotspot mutations in five known cancer-associated genes: BRAF, RAS family (HRAS, NRAS, and KRAS), and stop-gain NF1 gene mutations. We found only 4/303 total hotspot mutations fell outside the quality filtering criteria (Supplemental Material, Figure S2). We also compared the incidence of ultraviolet signature mutations in the unfiltered primary melanomas (n = 113), the unfiltered MM (n = 357), and the filtered MM (n = 356). Incidence of ultraviolet signature mutations (C>T substitutions) is significantly higher in unfiltered primary melanomas compared to MM samples (Figure 3, left panel, p < 2.2 × 10⁻¹⁶, Wilcoxon rank-sum test). Overall, C>T transitions, C>T transitions at a dipyrimidine site, and CC → TT accounted for 79.4%, 77.8%, and 0.11% of the total 5,322 mutations, respectively. Our filtering strategy found a significantly increased incidence of ultraviolet signature in retained vs. filtered mutations in the MM subgroup (p < 2.2 × 10⁻¹⁶, Wilcoxon). Nevertheless, mutations in other previously reported genes (e.g. PTEN, RAF1, CTNNB1, PBRM1, and KIT) were also retained (Figure 3, right panel). Thus, our filtering strategy retained high-quality mutation calls within known oncogenic driver genes.

FIGURE 3

Figure 3. Frequency of ultraviolet-signature mutations across different sample subgroups before and after our filtering strategy. Results are shown as boxplots with median and 25/75 percentiles. Asterisks show significant differences at a p-value ≤ 0.05.

Table S1 shows that 128 (OS-original diagnosis analysis) and 133 genes (OS-specimen collection analysis) were mutated in >11 patients (3%). To further limit the number of multiple comparisons testing for OS analysis, we conventionally used the survival status at the time of last follow-up. Although survival status (dead vs. alive) is inherently dependent on the duration of follow-up, the majority of surviving subjects (n = 171) had sufficient follow-up; only 7 (4%) vs. 33 (19%) of the censored subjects had follow-up for < 1-year vs. < 2-years, respectively. Of these, 53 (OS-original diagnosis analysis) and 57 (OS-specimen collection analysis) genes were mutated ≥1.5 times or ≤ 0.67 times in deceased subjects compared to living subjects. Gene lists from each of the two OS analyses were similar. Although BRAF, NRAS, NF1, TP53, PTEN, PPP6C, CTNNB1, PDGFRA, IDH1, and KIT were significantly mutated, they did not exhibit >50% change in deceased subjects compared to living subjects. Additionally, no change was seen when analysis was restricted to hotspot mutations previously reported for BRAF, NRAS, NF1, TP53, PTEN, ATM, CTNNB1, and PPP6C genes (25).

Cox Analysis of Somatically Mutated Genes in MM TCGA SKCM Samples

We performed the age-adjusted univariate Cox regression analysis for each of the 53 genes (Supplemental Material, Table S1, panel A; OS-original diagnosis analysis). Using a FDR cutoff of 0.2, we found that 4/53 genes had prognostic potential: fibroblast growth factor receptor 1 (FGFR1), RAS-related C3 botulinum toxin substrate 1 (RAC1), caspase recruitment domain-containing protein 11 (CARD11), and class II major histocompatibility complex transactivator (CIITA). An additional 10 genes had less significance (FDR 0.2–0.4). When OS time from specimen collection was used for Cox analysis, 26 genes looked promising (FDR < 0.4). With the exception of NTRK1, CUX1, TSC2, and AFF3, all 10 genes that were identified in the OS-original diagnosis analysis were also found to have considerable prognostic significance (FDR < 0.4) in the OS-specimen collection analysis. Eight additional genes (CNTRL, GNAS, AKAP9, KMT2C, PRDM16, KDM5A, PDGFRB, and STAG2) also had prognostic value (FDR < 0.4, Cox p-value ≤ 0.10). Figures 4, 5 show hazard ratio point estimates with 95CI for these genes that were identified out of the OS from original diagnosis and OS from time to specimen collection analyses, respectively. The previously reported CDKN2A, ARID2, and FBXW7 genes (26) exhibit at least 50% (fold change ≥1.5/ ≤ 0.67) difference in deceased subjects vs. living subjects. However, none of these genes met the FDR cutoff, even when analysis was focused on hotspot mutations (25).

FIGURE 4

Figure 4. Identifying prognostically significant mutated genes using univariate Cox regression analysis (OS-original diagnosis) in the TCGA stage III/IV cohort. Following identification of clinically significant, high-quality mutations, genes bearing somatic mutations in more than 3% of cases were selected for age-adjusted univariate Cox analysis of OS. Hazard ratio [exp(coef)] point estimates with 95% confidence intervals for genes using the false discovery cutoff of 0.2 (red) or 0.2–0.4 (black) are shown.

FIGURE 5

Figure 5. Identifying prognostically significant mutated genes using univariate Cox regression analysis (OS-specimen collection diagnosis) in the TCGA stage III/IV cohort. Following identification of clinically significant, high-quality mutations, genes bearing somatic mutations in more than 3% of cases were selected for age-adjusted univariate Cox analysis of OS. Hazard ratio [exp(coef)] point estimates with 95% confidence intervals for genes using the false discovery cutoff of 0.2 (red) or 0.2–0.4 (black) are shown.

Properties of the Somatic Mutations in the 22 Genes in MM Samples

Twenty-two genes were mutated in 190 (53.4%) specimens. These genes include the 14 identified to be promising for prognostic significance from the OS analysis from original diagnosis of melanoma (Figure 4) plus the eight genes that were uniquely found to be promising for prognostic significance from the OS analysis from time of specimen collection, since CIITA, RAC1, CARD11, and FGFR1 were common in both datasets (Figure 5). As shown in Figure 3, right panel, the majority of the mutations bear the ultraviolet signature (median frequency, 83.3%, ~95CI 78.6–88.0%). There was significant difference in the frequency of ultraviolet signature mutations seen in the 22 genes compared to the remaining 111 genes that were mutated in ≥3% of patients in either OS analyses (p = 0.04, Wilcoxon). Compared with the MAF of BRAFV600 mutations, the MAFs for most genes from the 22-gene list were significantly lower (Supplemental Material, Figure S2, panel A). However, while the MAF of BRAFV600 mutations was higher in primary as opposed to MM samples, the corresponding MAF for most of the 22 genes was the opposite—namely, higher in MM compared to primary melanoma specimens (Supplemental Material, Figure S3, panel B). The higher frequency of BRAFV600 mutations in SKCM suggests their founder status in this disease, while supporting the notion that mutations in the other 22 genes are potentially clonal or subclonal events that occur more frequently in MM samples.

Figure 6 shows the landscape of the 22 genes that demonstrated prognostic significance in the 190 specimens. As expected, vital status was not associated with total somatic mutation burden, but was associated with “immune-high” RNA-seq signatures (p = 0.873 and p = 1.78 × 10⁻⁷, respectively; Chi-square test). Overall, the frequency of mutations in any of these 22 genes was < 10% in the MM samples. Mutated RAC1, FGFR1, AFF3, and NTRK1 genes were more frequently seen in MM samples from patients who died within 5 years from the original diagnosis. In contrast, mutated CIITA was more frequently seen in MM samples from patients who have lived >5 years after original diagnosis (Figure 6). All 22 genes, except MAP2K1, were significantly associated with high somatic mutation burden (FDR = 0.09 for MAP2K1, FDR < 0.001 for the rest). Somatic mutations in the haploinsufficient homeodomain transcription factor gene, cut-like homeobox 1 (CUX1) (27), were differentially distributed in the “immune-high” vs. non “immune-high” specimens (p = 0.027).

FIGURE 6

Figure 6. Landscape of somatically mutated genes that are promising for being prognostic (n = 22, FDR < 0.4) in the TCGA stage III/IV cohort following OS-original diagnosis and OS-time from specimen collection analyses. Only tumors with mutations in at least one of the 22 genes are shown in the center of the matrix (190 out of a total of 356). Mutated genes are shown as a binary outcome (presence/absence) for each sample and are color-coded according to the presence/absence of mutations as identified by DNA-seq (dark gray or light gray) and RNA-seq analysis [filled circles (∙) or open circles (°)]. Total number of somatic mutations (# Mutations, top), vital status, age at specimen collection (≤ 65, >65 years old), OS follow-up (in years), AJCC stage at specimen procurement, molecular classification [BRAF hotspot, RAS hotspot, NF1 mutations (□ indicates nonsense mutation), triple wild-type (TWT)], RNA-seq classification, mutation spectra (left panel), and mutated gene frequency (right panel) for each of the 14 mutated genes (left) are shown for each sample.

UNCeqR analysis revealed that 271 of these 423 mutation events (64.1%) were sufficiently covered with the exception of NTRK1, EGFR, and GNAS (6, 10, and 23%, respectively), implying low-to-absent gene expression (Figure 6). Of the remaining 19 genes, somatic mutations in only 9 genes were confirmed by the RNA-seq data in more than 75% of samples (RAC1, AFF4, TSC2, CNTRL, CUX1 MAP2K1, KDM5A, SPEN, AKAP9), although corresponding genes were expressed based on the RNA-Seq by Expectation-Maximization (RSEM) data. Figure S4 (Supplemental Material) further characterizes somatic mutations in these 9 genes that were confirmed by RNA-seq in stage III/IV SKCM.

Validation of the 22 Genes in an Independent Cohort

To validate our findings related to the 22 genes from the TCGA cohort in a separate MM cohort, we analyzed the presence of these mutations in a cohort of 33 stage III/IV patients with SKCM who were followed at the UNC-CH Melanoma Program and had consented to the UNCseq™ project (Supplemental Material, Table S2), combined with 6 other previously published melanoma datasets (18–20, 22, 26, 28). Table 2 shows the non-synonymous somatic mutations from MM tumor samples of cutaneous or unknown primaries that were subjected to next generation sequencing analysis and grouped according to primary (n = 107) vs. metastatic (n = 417) status. As was the case with the TCGA SKCM cohort, the frequency for the 22 somatically mutated genes in MM was low, with few exceptions (KMT2C, PRDM16, CARD11, PTPRC). In contrast with the TCGA SKCM cohort, however, more mutated genes were seen in primaries compared to metastases (2/22 in the TCGA cohort vs. 10/22 in the validation cohort). Of the 9 mutated genes whose mutations were confirmed by RNA-seq, all but one (i.e., AFF4) were equivocally found to be expressed in the other two datasets that reported both RNA and DNA sequencing analysis (19, 26).

TABLE 2

Table 2. Frequency of nonsynonymous somatic mutations in the validation (non-TCGA) cohort for the 22 genes that were considered promising for prognostic significance in melanoma from the TCGA cohort analysis.

Of the seven studies that comprised our validation cohort, we were able to obtain reliable follow-up data from 224 patients across five studies (UNC-CH plus 18, 19, 21, 22. The median follow-up of this patient cohort was 20.0 months (~95CI, 16.9–27.0 months), during which time 30.4% of deaths occurred. Cox analysis for each of the 8 mutated genes whose RNA expression was confirmed across all studies showed trends of mutations in RAC1 in MM with worse prognosis (HR = 2.1, 95CI 0.66–6.63, log-rank p = 0.07), whereas mutations in SPEN showed trends of mutations with better prognosis (HR = 0.50, range 0.27–0.93, log-rank p = 0.09) (Figure 7).

FIGURE 7

Figure 7. Overall survival analysis of patients with MM in the validation cohort that had follow-up according to the mutation status of SPEN and RAC1 using the Kaplan-Meier method. Although 225 patients had available follow-up across 5 studies, only 3 studies had available mutation data for these two genes (18, 19, 21).

Beyond SPEN and RAC1 Somatic Mutations; Copy Number Alterations and Protein Expression

In general, mutations are frequently associated with copy number alterations and can influence gene expression (29–31). Furthermore, genetic aberrations may have an impact on the function and/or stability of RNA which may lead to corresponding changes in protein function and/or abundance. Figure 8 shows integrated analysis of somatic mutations, copy number alterations, and gene expression alterations for RAC1 and SPEN for the 357 TCGA MM samples. All, but two, somatic mutations were not associated with definite copy number alterations (GISTIC score 2+ or 2–) and all, but one, somatic mutations were not associated with changes in gene expression. However, 22/34 patients with any genetic aberrations in RAC1 (somatic mutations and/or copy number amplifications) were deceased at the time of analysis. The unfavorable outcome of patients with any genetic aberrations in RAC1 is in line with a report on the adverse prognostic significance of high RAC1 protein expression by immunohistochemistry in primary cutaneous melanoma samples (32).

FIGURE 8

Figure 8. Integration of somatic mutations (red), gene copy number (blue), and gene expression (green) for RAC1 and SPEN in the MM TCGA dataset. Results are shown as RNA expression (log₂-transformed RSEM data) over copy number alterations that have been generated using the GISTIC algorithm. For better visualization, samples with no somatic mutations are shown as small gray dots or squares.

Somatic mutations in SPEN were not associated with any definite copy number alterations (GISTIC score 2+ or 2–) or clear changes in gene expression. To investigate the protein expression pattern in SPEN, we performed immunohistochemical analysis in MM samples procured from 87 patients who were treated at UNC-CH Melanoma Program between 2000 and 2011. Clinicopathologic characteristics of this patient cohort are shown in Table S3. Protein expression of SPEN by melanoma cells was predominantly nuclear (76%); 23.5% of cores expressed both strong (2+ and 3+) nuclear and cytoplasmic stain (23.5%) and 0.5% of cores expressed predominantly cytoplasmic stain (Figures 9A,B). SPEN expression was overall higher by melanoma cells compared to stromal cells (Figure 9B). Using the semiquantitative 0, 1+, 2+, 3+ scale, 0%, 9.3%, 28.3%, 62.3% of tissue cores exhibited absent, low, intermediate, or strong expression of SPEN by melanoma cells. Strong expression staining intensity (2+ and 3+) for SPEN protein, if present in both nuclear and cytoplasmic compartments in melanoma cells, trended to correlate with the presence of tumor-infiltrating lymphocytes (2-way contingency table, χ² test p = 0.08). At a median follow-up of 15 months (range 1–68 months) and for the patients who had melanoma-specific overall survival data (n = 82) patients with high (>2+) nuclear plus cytoplasmic expression of SPEN by melanoma cells had longer melanoma-specific OS compared with patients with low expression (≤ 2+) (HR and 95CI, 0.56, 0.29–1.09, log-rank p-value 0.048) (Figure 9C). No significant differences were observed between the other group comparisons. Given that SPEN somatic mutations do not associate with significant differences in RNA expression, we speculate that distinct yet-to-be identified SPEN mutations may regulate SPEN localization.

FIGURE 9

Figure 9. Protein expression of SPEN in MM. (A,B). Images of representative tissue sections from the MM UNC-CH tissue microarray (10X magnification) that were stained with an anti-SPEN antibody (see materials and methods for details). Differences in SPEN staining intensity (A, 3+; B, 2+), expression of SPEN in melanoma vs. stromal cells (B), and localization (A, nuclear plus cytoplasmic; B, nuclear only) are shown. (C). Melanoma-specific OS analysis from time of specimen collection to the last follow-up in the 87 MM UNC-CH dataset. Of these patients only 82 had available melanoma-specific OS data.

Discussion

We have analyzed the largest to-date sample cohort for somatic mutations in MM. With available follow-up information from the TCGA, we validated our findings in an independent dataset. In contrast to previous work that defined the landscape of somatic mutations in melanoma irrespective of stage and subtype (8, 19–21, 26, 28), our focus on MM was unique in that it attempted to understand the prognostic significance of mutations in established metastases specifically. We have selected for high-quality, functional, statistically significant mutations. The large number of tissue specimens, the focus on MM samples and somatic mutations only, and our filtering strategy are a few of the critical factors that differentiate our study results from previous efforts to identify prognostic factors in melanoma using a combination of clinicopathologic and multiple—omics platforms (33, 34). Our study has several important conclusions. Only a handful of genes that had been previously described as recurrently mutated may have prognostic significance (e.g., RAC1). Other well-characterized and more frequently mutated genes (e.g., BRAF, RAS, NF1, TP53, PTEN, CDKN2A) were not prognostic in MM, perhaps because such mutations are already found in nevi and primary melanomas (1). Less well-known genes in melanoma biology not only may be significantly mutated due to the random effect of ultraviolet radiation, but their mutated status may have potential prognostic significance [e.g., the Spen homolog transcriptional regulator (SPEN)]. Not all somatic mutations are expressed at the RNA level. Validation of expressed somatic mutations detected at the DNA level by RNA sequencing maybe the future in clinical oncology; current commercially available molecular tests that perform targeted panel sequencing identify mutations at the DNA level only (e.g., Foundation Medicine).

To date, there is limited knowledge about the significance of low-frequency somatic mutations in melanoma. Studies have established that ultraviolet radiation induces mutations, which in turn may potentially yield neoantigens. The neoantigens provoke an immune response, which may account for the immunogenicity of cutaneous melanoma and other cancers (35). The majority of the 22 genes that were identified from the TCGA and validation cohort were infrequently mutated, more frequently found in specimens with high somatic mutation burden, are evenly distributed across the entire gene. Also, they are of the missense type (except for RAC1 and perhaps MAP2K1). The lack of mutations with high functional impact does not mitigate the importance of such missense mutations in protein function, as in the case of KEAP1 in lung cancer (36). Instead, it may provide an explanation about the lack of association between high somatic mutation burden with host immune response and OS in melanoma, which is complex: somatic mutations in genes associated with immune surveillance (e.g., PTPRC/CD45, FCRL4, CARD11) may be associated with favorable prognosis because potentially damaging mutations are not ultimately expressed. In contrast, somatic mutations in tumor-promoting genes (FGFR1, RAC1, AFF3, AFF4, TSC2) may be intuitively associated with unfavorable prognosis.

Our validation cohort was comprised of patient samples with significant heterogeneity with respect to selection for BRAFV600 mutations and, systemic treatment type (i.e., FDA approved treatments before or after 2011), which may have influenced OS (37, 38), and the geopolitical origin of patients (39). Despite the challenges in data interpretation, we were able to confirm the low somatic mutation frequency and the UV signature type of most mutations. Use of different bioinformatics algorithms may have accounted for the discrepancy in validation of expressed mutations by RNA-seq between the TCGA dataset and others. Nevertheless, certain mutated genes were confirmed by RNA-seq across all datasets (e.g., RAC1, MAP2K1, SPEN, CUX1, TSC2, CNTRL, AKAP9, and STAG2) whereas others were found not to be expressed, irrespective of the RNA-seq validation algorithm used (e.g., FCRL4, CARD11, PDGFRB).

Our survival analysis of the validation cohort was significantly underpowered to confirm the prognostic significance of the eight mutated genes that were also found to be expressed. Furthermore, the validation cohort had shorter follow-up and less events. For example, assuming that a given gene whose incidence of mutation is 5.3% in the study population (e.g., NTRK1) and is associated with worse OS (HR = 1.93), the power to detect significant prognostic difference if n = 245 and 70% of patients had an event (e.g., death) is only 0.64. Nevertheless, we were able to show that RAC1 and SPEN may be promising prognostic factors for OS in MM. RAC1 been previously reported in SKCM and other datasets (8, 20) and was associated with melanoma progression, suppression of host immune response, and drug resistance (40). A recent study has shown that high expression of RAC1 protein in primary cutaneous melanoma samples was associated with thinner melanomas, BRAFV600 mutation and with RAC1 mutation (32). Little is known about the role of SPEN, a transcriptional regulator of NOTCH1 and hormone receptor signaling, in melanoma and other cancers (41). In our work we have shown that high (>2+) nuclear and cytoplasmic SPEN expression trends to associate with melanomas that have present tumor-infiltrating lymphocytes. This association may account for our finding regarding the favorable prognostic significance of patients with high melanoma nuclear and cytoplasmic SPEN expression compared to those with lower (≤ 2+) protein expression. The results are comparable with a recent report that investigated the prognostic significance of high vs. low SPEN RNA expression in luminal A breast cancer (24). Irrespective of the limitations of our validation cohort, our study suggests that the prognostic impact for most of these genes may be small. Rather, multiple genetic, and epigenetic aberrations within the cancer itself (42), the host (43), and the environment (44) may play an even larger role in established metastases as opposed to primary melanoma.

Ethics Statement

This study was carried out in accordance with the guidelines and recommendations from the Genomic Data Commons (GDC) Data Portal of the Cancer Genome Atlas Project. GDC Data Portal is an interactive data system for researchers to search, download, upload, and analyze harmonized cancer genomic data sets, including TCGA.

Author Contributions

GP, MH, AS, MP, DB, AE, JF, WH, JW, CQ, and SM involved in data collection. XZ, Al, PL, AH, JP, DM, MS, HJ, SJ, KH, DO, NT, PG, FC, CL, HE, NS, WH, JW, CQ, NW, PJ, JT, NH, GM, RL, DJ, RS, DH, and SM involved in data analysis. XZ, GP, KH, CQ, RS, DH, and SM involved in manuscript preparation.

Funding

Research reported in this publication was supported by the National Cancer Institute Cancer Clinical Investigator Team Leadership Award (5P30CA016086-38, NS and SM), the National Cancer Institute of the National Institutes of Health under award number U10CA181009 (DH), the University of North Carolina at Chapel Hill University Cancer Research Fund (SM), the National Health and Medical Research Council of Australia, Cancer Council and Cancer Institute New South Wales, Australia, Australian Cancer Research Foundation (CQ, NW, PJ, JT, NH, GM, RS).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank the TCGA Data Access Committee. The results shown here are in whole or part based upon data generated by the TCGA Research Network (http://cancergenome.nih.gov). We thank Mses. Colleen MacDonagh and Erin Laurie for assistance in editing the manuscript. We thank Drs. Ian R. Watson, Philippe Lefrancois, Mathieu Lajoie, and Rached Alkallas for assistance with providing updated TCGA clinical data file. We also thank Drs. Michael Krauthammer and Ruth Halaban for providing critical review of the manuscript as well as updated follow-up data from the Yale Melanoma Cohort (Yale SPORE in Skin Cancer, P50CA121974; Ruth Halaban, PI).

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2018.00584/full#supplementary-material

Abbreviations

SKCM, Skin Cutaneous Melanoma; MM, metastatic melanoma; TCGA, The Cancer Genome Atlas Project; VCF, variant call format; UNC-CH, the University of North Carolina at Chapel Hill; OS, overall survival; RNA-seq, RNA sequencing; 95CI, 95% confidence intervals; FDR, false discovery rate; MAF, mutant allele frequency; MAC, mutant allele count.

References

1. Shain AH, Yeh I, Kovalyshyn I, Sriharan A, Talevich E, Gagnon A, et al. The genetic evolution of melanoma from precursor lesions. N Engl J Med. (2015) 373:1926–36. doi: 10.1056/NEJMoa1502583

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Thomas NE, Edmiston SN, Alexander A, Groben PA, Parrish E, Kricker A, et al. Association between NRAS and BRAF mutational status and melanoma-specific survival among patients with higher-risk primary melanoma. JAMA Oncol. (2015) 1:359–68. doi: 10.1001/jamaoncol.2015.0493

CrossRef Full Text | Google Scholar

3. Hosler GA, Davoli T, Mender I, Litzner B, Choi J, Kapur P, et al. A primary melanoma and its asynchronous metastasis highlight the role of BRAF, CDKN2A, and TERT. J Cutan Pathol. (2015) 42:108–17. doi: 10.1111/cup.12444

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Nagore E, Heidenreich B, Rachakonda S, Garcia-Casado Z, Requena C, Soriano V, et al. TERT promoter mutations in melanoma survival. Int J Cancer (2016) 139:75–84. doi: 10.1002/ijc.30042

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Bai X, Kong Y, Chi Z, Sheng X, Cui C, Wang X, et al. MAPK pathway and TERT promoter gene mutation pattern and its prognostic value in melanoma patients: a retrospective study of 2,793 cases. Clin Cancer Res. (2017) 23:6120–27. doi: 10.1158/1078-0432.CCR-17-0980

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Griewank KG, Murali R, Puig-Butille JA, Schilling B, Livingstone E, Potrony M, et al. TERT promoter mutation status as an independent prognostic factor in cutaneous melanoma. J Natl Cancer Inst. (2014) 106:dju246. doi: 10.1093/jnci/dju246

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Roh MR, Gupta S, Park KH, Chung KY, Lauss M, Flaherty KT, et al. Promoter methylation of PTEN is a significant prognostic factor in melanoma survival. J Invest Dermatol. (2016) 136:1002–11. doi: 10.1016/j.jid.2016.01.024

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Cancer Genome Atlas N. Genomic Classification of Cutaneous Melanoma. Cell (2015) 161:1681–96. doi: 10.1016/j.cell.2015.05.044

CrossRef Full Text | Google Scholar

9. Jakob JA, Bassett RL, Ng CS, Curry JL, Joseph RW, Alvarado GC, et al. NRAS mutation status is an independent prognostic factor in metastatic melanoma. Cancer (2012) 118:4014–23. doi: 10.1002/cncr.26724

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Long GV, Menzies AM, Nagrial AM, Haydu LE, Hamilton AL, Mann GJ, et al. Prognostic and clinicopathologic associations of oncogenic BRAF in metastatic melanoma. J Clin Oncol. (2011) 29:1239–46. doi: 10.1200/JCO.2010.32.4327

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Frampton GM, Fichtenholtz A, Otto GA, Wang K, Downing SR, He J, et al. Development and validation of a clinical cancer genomic profiling test based on massively parallel DNA sequencing. Nat Biotechnol. (2013) 31:1023–31. doi: 10.1038/nbt.2696

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Conley BA, Doroshow JH. Molecular analysis for therapy choice: NCI MATCH. Semin Oncol. (2014) 41:297–9. doi: 10.1053/j.seminoncol.2014.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Schlenk RF, Döhner K, Krauter J, Fröhling S, Corbacioglu A, Bullinger L, et al. Mutations and treatment outcome in cytogenetically normal acute myeloid leukemia. N Engl J Med. (2008) 358:1909–18. doi: 10.1056/NEJMoa074306

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Sanborn JZ, Chung J, Purdom E, Wang NJ, Kakavand H, Wilmott JS, et al. Phylogenetic analyses of melanoma reveal complex patterns of metastatic dissemination. Proc Natl Acad Sci USA. (2015) 112:10995–10000. doi: 10.1073/pnas.1508074112

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Harbst K, Lauss M, Cirenajwis H, Isaksson K, Rosengren F, Törngren T, et al. Multiregion whole-exome sequencing uncovers the genetic evolution and mutational heterogeneity of early-stage metastatic melanoma. Cancer Res. (2016) 76:4765–74. doi: 10.1158/0008-5472.CAN-15-3476

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Cingolani P, Platts A, Wang le L, Coon M, Nguyen T, Wang L, et al. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome of drosophila melanogaster strain w1118; iso-2; iso-3. Fly (2012) 6:80–92. doi: 10.4161/fly.19695

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Wilkerson MD, Cabanski CR, Sun W, Hoadley KA, Walter V, Mose LE, et al. Integrated RNA and DNA sequencing improves mutation detection in low purity tumors. Nucleic Acids Res. (2014) 42:e107. doi: 10.1093/nar/gku489

PubMed Abstract | CrossRef Full Text | Google Scholar

18. 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. doi: 10.1016/j.cell.2016.02.065

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Krauthammer M, Kong Y, Bacchiocchi A, Evans P, Pornputtapong N, Wu C, et al. Exome sequencing identifies recurrent mutations in NF1 and RASopathy genes in sun-exposed melanomas. Nat Genet. (2015) 47:996–1002. doi: 10.1038/ng.3361

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Krauthammer M, Kong Y, Ha BH, Evans P, Bacchiocchi A, McCusker JP, et al. Exome sequencing identifies recurrent somatic RAC1 mutations in melanoma. Nat Genet. (2012) 44:1006–14. doi: 10.1038/ng.2359

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Hayward NK, Wilmott JS, Waddell N, Johansson PA, Field MA, Nones K, et al. Whole-genome landscapes of major melanoma subtypes. Nature (2017) 545:175–80. doi: 10.1038/nature22071

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Catalanotti F, Cheng DT, Shoushtari AN, Johnson DB, Panageas KS, Momtaz P, et al. PTEN loss-of-function alterations are associated with intrinsic resistance to BRAF inhibitors in metastatic melanoma. JCO Precis Oncol. (2017). doi: 10.1200/PO.16.00054

CrossRef Full Text | Google Scholar

23. Zhao X, Wang A, Walter V, Patel NM, Eberhard DA, Hayward MC, et al. Combined targeted DNA sequencing in non-small cell lung cancer (NSCLC) using UNCseq and NGScopy, and RNA sequencing using UNCqeR for the detection of genetic aberrations in NSCLC. PLoS ONE (2015) 10:e0129280. doi: 10.1371/journal.pone.0129280

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Légaré S, Cavallone L, Mamo A, Chabot C, Sirois I, Magliocco A, et al. The estrogen receptor cofactor SPEN functions as a tumor suppressor and candidate biomarker of drug responsiveness in hormone-dependent breast cancers. Cancer Res. (2015) 75:4351–63. doi: 10.1158/0008-5472.CAN-14-3475

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Chang MT, Asthana S, Gao SP, Lee BH, Chapman JS, Kandoth C, et al. Identifying recurrent mutations in cancer reveals widespread lineage diversity and mutational specificity. Nat Biotechnol. (2016) 34:155–63. doi: 10.1038/nbt.3391

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Hodis E, Watson IR, Kryukov GV, Arold ST, Imielinski M, Theurillat JP, et al. A landscape of driver mutations in melanoma. Cell (2012) 150:251–63. doi: 10.1016/j.cell.2012.06.024

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Wong CC, Martincorena I, Rust AG, Rashid M, Alifrangis C, Alexandrov LB, et al. Inactivating CUX1 mutations promote tumorigenesis. Nat Genet. (2014) 46:33–8. doi: 10.1038/ng.2846

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Berger MF, Hodis E, Heffernan TP, Deribe YL, Lawrence MS, Protopopov A, et al. Melanoma genome sequencing reveals frequent PREX2 mutations. Nature (2012) 485:502–6. doi: 10.1038/nature11071

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Ding J, McConechy MK, Horlings HM, Ha G, Chun Chan F, Funnell T, et al. Systematic analysis of somatic mutations impacting gene expression in 12 tumour types. Nat Commun. (2015) 6:8554. doi: 10.1038/ncomms9554

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Zack TI, Schumacher SE, Carter SL, Cherniack AD, Saksena G, Tabak B, et al. Pan-cancer patterns of somatic copy number alteration. Nat Genet. (2013) 45:1134–40. doi: 10.1038/ng.2760

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Roszik J, Wu CJ, Siroy AE, Lazar AJ, Davies MA, Woodman SE, et al. Somatic copy number alterations at oncogenic loci show diverse correlations with gene expression. Sci Rep. (2016) 6:19649. doi: 10.1038/srep19649

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Mar VJ, Wong SQ, Logan A, Nguyen T, Cebon J, Kelly JW, et al. Clinical and pathological associations of the activating RAC1 P29S mutation in primary cutaneous melanoma. Pigment Cell Melanoma Res. (2014) 27:1117–25. doi: 10.1111/pcmr.12295

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Jiang Y, Shi X, Zhao Q, Krauthammer M, Rothberg BE, Ma S. Integrated analysis of multidimensional omics data on cutaneous melanoma prognosis. Genomics (2016) 107:223–30. doi: 10.1016/j.ygeno.2016.04.005

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Jayawardana K, Schramm SJ, Haydu L, Thompson JF, Scolyer RA, Mann GJ, et al. Determination of prognosis in metastatic melanoma through integration of clinico-pathologic, mutation, mRNA, microRNA, and protein information. Int J Cancer (2015) 136:863–74. doi: 10.1002/ijc.29047

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Brown SD, Warren RL, Gibb EA, Martin SD, Spinelli JJ, Nelson BH, et al. Neo-antigens predicted by tumor genome meta-analysis correlate with increased patient survival. Genome Res. (2014) 24:743–50. doi: 10.1101/gr.165985.113

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Hast BE, Cloer EW, Goldfarb D, Li H, Siesser PF, Yan F, et al. Cancer-derived mutations in KEAP1 impair NRF2 degradation but not ubiquitination. Cancer Res. (2014) 74:808–17. doi: 10.1158/0008-5472.CAN-13-1655

CrossRef Full Text | Google Scholar

37. Wolchok JD, Chiarion-Sileni V, Gonzalez R, Rutkowski P, Grob JJ, Cowey CL, et al. Overall survival with combined nivolumab and ipilimumab in advanced melanoma. N Engl J Med. (2017) 377:1345–56. doi: 10.1056/NEJMoa1709684

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Long GV, Eroglu Z, Infante J, Patel S, Daud A, Johnson DB, et al. Long-term outcomes in patients with BRAF V600-mutant metastatic melanoma who received dabrafenib combined with trametinib. J Clin Oncol. (2018) 36:667–73. doi: 10.1200/JCO.2017.74.1025

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Karimkhani C, Green AC, Nijsten T, Weinstock MA, Dellavalle RP, Naghavi M, et al. The global burden of melanoma: results from the global burden of disease study 2015. Br J Dermatol. (2017) 177:134–40. doi: 10.1111/bjd.15510

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Watson IR, Li L, Cabeceiras PK, Mahdavi M, Gutschner T, Genovese G, et al. The RAC1 P29S hotspot mutation in melanoma confers resistance to pharmacological inhibition of RAF. Cancer Res. (2014) 74:4845–52. doi: 10.1158/0008-5472.CAN-14-1232-T

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Stephens PJ, Davies HR, Mitani Y, Van Loo P, Shlien A, Tarpey PS, et al. Whole exome sequencing of adenoid cystic carcinoma. J Clin Invest. (2013) 123:2965–8. doi: 10.1172/JCI67201

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Mo Q, Wang S, Seshan VE, Olshen AB, Schultz N, Sander C, et al. Pattern discovery and cancer gene identification in integrated cancer genomic data. Proc Natl Acad Sci USA. (2013) 110:4245–50. doi: 10.1073/pnas.1208949110

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Blank CU, Haanen JB, Ribas A, Schumacher TN. CANCER IMMUNOLOGY. The “cancer immunogram”. Science (2016) 352:658–60. doi: 10.1126/science.aaf2834

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Zitvogel L, Daillère R, Roberti MP, Routy B, Kroemer G. Anticancer effects of the microbiome and its products. Nat Rev Microbiol. (2017) 15:465–78. doi: 10.1038/nrmicro.2017.44

PubMed Abstract | CrossRef Full Text | Google Scholar