Comprehensive Genetic Landscape of Uveal Melanoma by Whole-Genome Sequencing (original) (raw)

Abstract

Uveal melanoma (UM) is a rare intraocular tumor that, similar to cutaneous melanoma, originates from melanocytes. To gain insights into its genetics, we performed whole-genome sequencing at very deep coverage of tumor-control pairs in 33 samples (24 primary and 9 metastases). Genome-wide, the number of coding mutations was rather low (only 17 variants per tumor on average; range 7–28), thus radically different from cutaneous melanoma, where hundreds of exonic DNA insults are usually detected. Furthermore, no UV light-induced mutational signature was identified. Recurrent coding mutations were found in the known UM drivers GNAQ, GNA11, BAP1, EIF1AX, and SF3B1. Other genes, i.e., TP53BP1, CSMD1, TTC28, DLK2, and KTN1, were also found to harbor somatic mutations in more than one individual, possibly indicating a previously undescribed association with UM pathogenesis. De novo assembly of unmatched reads from non-coding DNA revealed peculiar copy-number variations defining specific UM subtypes, which in turn could be associated with metastatic transformation. Mutational-driven comparison with other tumor types showed that UM is very similar to pediatric tumors, characterized by very few somatic insults and, possibly, important epigenetic changes. Through the analysis of whole-genome sequencing data, our findings shed new light on the molecular genetics of uveal melanoma, delineating it as an atypical tumor of the adult for which somatic events other than mutations in exonic DNA shape its genetic landscape and define its metastatic potential.

Main Text

Despite having the very rare incidence of 5–8 new cases per million per year,1, 2 uveal melanoma (UM [MIM: 155720]) is the most common primary intraocular tumor of the adult. It develops from melanocytes in the choroid, the ciliary body, or the iris (collectively called the “uvea,” one of the inner layers of the eye) and usually metastasizes through the blood stream to the liver.3, 4 Symptoms include variable and painless visual disturbances, often presenting when the tumor has already reached a considerable mass. Survival and potential therapeutic options depend, among other things, on the presence of specific genetic alterations.5 Although population studies suggest ethnic predisposition,6, 7 environmental factors that are directly involved in the transformation process have not been clearly delineated. For instance, a possible association with UV light exposure has been suggested,8, 9, 10, 11 but questioned recently by molecular data.12 Research on UM molecular genetics has been performed mostly by investigating coding mutations or copy-number variations detectable by direct sequencing of target genes, karyotype, array-CGH, MLPA, or SNP-array analyses.13, 14, 15, 16, 17 As a result, mutations at codon 209 of the paralogous oncogenes GNAQ (MIM: 600998) and GNA11 (MIM: 139313)18, 19 and in the tumor suppressor BAP1 (MIM: 603089)20 have been identified in the majority of UMs, whereas insults in EIF1AX (MIM: 300186) and SF3B1 (MIM: 605590) or other genes seem to be less frequent, accounting for at most 20% of cases.12, 21, 22, 23, 24, 25, 26 Moreover, copy gains and losses are common events in this tumor, typically involving chromosome 3 monosomy, 6p gain, and 8q gain.14, 17 After whole-genome sequencing of a series of tumor-control pairs, we present here an analytically unbiased and comprehensive assessment of the genetic landscape of UM.

We screened 33 UM samples (24 primary tumors and 9 unrelated metastases; Table S1) and corresponding normal tissue pairs by deep-coverage whole-genome sequencing (WGS), using the processing platform by Complete Genomics.27 Written informed consent was obtained from all individuals enrolled in this study, and approval for human subject research was obtained from the Institutional Review Boards of all participating Institutions. Sequence reads were mapped to the human reference genome, assembly GRch37, and somatic variants in tumors were called by comparison with the matched normal genomes, as previously described.28 All samples were surgically collected from eye enucleations or resected from liver metastases, allowing very clear post-surgery macroscopic isolation of tumor tissue from the surrounding environment. None of the 33 affected individuals received any treatment prior to primary tumor removal. Average coverage was 112× (range 102–118) for both tumor and control samples (>96% of the genome was covered 40× or more times), with minimal inter-individual variations (not shown). Mutation calls were performed genome-wide and included single-nucleotide variants (SNVs), copy-number variations (CNVs), as well as structural variations (SVs) such as chromosomal rearrangements. CNVs and SVs were assessed by comparing sequence coverage and especially de novo assembly of reads defining novel genetic junctions.28 Data were extracted from MasterVar files and other relevant matrices by ad hoc Perl, bash, and R scripts, available upon request. Overall, we detected 37,321 SNVs (average per sample 1,166; range 576–2,131), 1,584 SVs (average per sample 50; range 13–182), and a number of CNVs corresponding to an average of 13.6% of the genome (range 0.03%–33.9%) (Table S2).

Unsupervised hierarchical clustering of our samples on CNVs revealed four major subgroups associated with mutational and metastatic status, branched two by two (Figure 1). Classes A and B (first branch) involved samples carrying chromosome 3 monosomy (by Fisher test, p value = 4.4 × 10−6), chr 8q gain (p value = 2.8 × 10−9), and in some instances chr 8p loss (p value = 3.0 × 10−2). In addition, class B tumors also had loss of chr 6q (p value = 1.0 × 10−3). Conversely, classes C and D represented more distinct subtypes with relatively few chromosomal rearrangements; class C tumor had no major aneuploidies, whereas class D reported gains of the distal part of chr 8q (p value = 2.0 × 10−3). Seven samples presenting chr 1q gain were scattered across all classes, whereas loss of chr 1p was typical of class A tumors (p value = 5.0 × 10−3).

Figure 1.

Figure 1

Unsupervised Hierarchical Clustering and Global Genetic Landscape of All Tumors Analyzed in This Study

Sample IDs are indicated on the right. CNV events are depicted in blue (copy losses) or in shades of red (copy gains) and ploidy is indicated in the legend provided at the top. SNVs in genes found to carry mutations in six or more individuals are shown on the left, with mutation classes provided within gray or blue boxes. Clustering identifies four classes—A, B, C, and D (see text)—indicated within the dendrogram.

Samples with monosomy of chr 3 were also associated (77% of cases) with somatic mutations in the tumor suppressor BAP1, lying on chr 3p21.1, in accordance with Knudson’s two-hit model of tumorigenesis.29 Indeed, BAP1 SNVs included all kinds of somatic events, but mostly mutations leading to premature stop codons and therefore to functional protein knockout (Figure S1, Table S3). Hallmark driver mutations in the GNAQ and GNA11 paralogs, encoding the components of the alpha subunit of the Gq protein heterotrimer, were present in 100% of the samples examined. They occurred in a perfectly mutually exclusive pattern and involved only four specific missenses—c.626A>T (p.Gln209Leu) (GenBank: NM_002072.3; 11 samples) and c.626A>C (p.Gln209Pro) (GenBank: NM_002072.3; 8 samples) in GNAQ, and c.626A>T (p.Gln209Leu) (GenBank: NM_002067.2; 13 samples) and c.626A>C (p.Gln209Pro) (GenBank: NM_002067.2; 1 sample) in _GNA11_—affecting the same functional amino acid residue and conferring oncogenic potential to this G protein.18, 19 Six tumors (18%) had missense mutations in SF3B1 (Splicing Factor 3B, subunit 1), affecting codon 625 (5 cases) and codon 626 (1 case) (Figure S1, Table S3), a previously described hotspot region.21 Finally, 7 other tumors had mutations impacting the first codons of EIF1AX (Eukaryotic Translation Initiation Factor 1A, X-Linked) (Figure S1, Table S3).22 Mutations in SF3B1 and EIF1AX seemed to occur in a mutually exclusive pattern and to be enriched with classes C and D (p value = 1.6 × 10−4), with SF3B1 preferentially mutated in class D. Except for one sample, BAP1 mutations were never observed in case subjects carrying mutations in SF3B1 or EIF1AX (p value = 1.4 × 10−5), in agreement with findings from previous literature.23, 30 Also, consistent with the fact that all of the tumors analyzed harbored variants affecting either GNAQ or GNA11 Gln209, no somatic SNVs were observed in PLCB4 or CYSLTR2, two genes that have been found to be mutated in a mutually exclusive pattern with respect to GNAQ or GNA11 variants.25, 26

Five genes (TP53BP1 [MIM: 605230], CSMD1 [MIM: 608397], TTC28 [MIM: 615098], DLK2, and KTN1 [MIM: 600381]) harbored somatic missense or truncating mutations in at least two samples, across all tumor classes (Table S3). TP53BP1 is a partner of the tumor suppressor protein p53, known to play a crucial role in maintaining genomic integrity as a mediator and effector of homologous recombination in response to double-strand breaks. This protein acts as a molecular scaffold that recruits responsive proteins, in order to repair damaged chromatin,31 and its depletion has been associated with increased cell proliferation.32 CSMD1 (Cub and Sushi Multiple Domains-1) is a candidate tumor suppressor gene, the hyper-expression of which increases survival in mice with xenografted tumors.33 Loss of CSMD1 was detected in a large set of cancers, including head and neck, lung, breast, and skin primary tumors,34 and associated with high tumor grade in invasive ductal breast carcinoma.35 TTC28 (Tetratricopeptide Repeat Domain 28) is a ubiquitous protein, associated with diverse biological functions. Of note, TCC28 plays a critical role in the progress of mitosis and cytokinesis during mammalian cell cycle and its dysfunction was described as a potential component of tumorigenesis and tumor progression.36, 37 DLK2 (Delta-Like 2 homolog) is a transmembrane epidermal growth factor-like protein. It is highly homologous to DLK1, a protein that was found to be present at high levels in gliomas and involved in cell proliferation.38 Similar to DLK1, DLK2 can bind to NOTCH1,39 modulating the oncogenic potential of cultured melanoma cells.40 Finally, KTN1 (Kinectin 1) is a protein of the endoplasmic reticulum membrane that interacts with kinesins.41, 42 Its role in cancer may be linked to dysregulation of cytoskeletal activity and mitosis. Two of these five genes were previously found to be mutated in UM: a missense in TTC28 was detected in 1 out 35 samples from a WES screen,22 and the cBioPortal repository reports a missense mutation in KTN1 in one out of 80 tumors profiled by TCGA.43, 44

Taken together, our clustering analysis indicates that initial events involve GNAQ or GNA11 mutations, followed by a major branching determined by the functional loss of BAP1 and copy gains of chr 8q versus cases with a relatively normal chromosomal ploidy. These latter samples have conversely mutations in EIF1AX or in SF3B1 (classes A and B versus C and D, respectively). In tumors with BAP1 mutations, the long arm of chr 6 could eventually be lost, differentiating class B from class A. In tumors that are negative for BAP1 mutations, part of chr 8q could undergo amplification, differentiating class D from class C (Figure 2).

Figure 2.

Figure 2

Inferred Somatic Events Defining Tumor Classes, as Identified by Clustering

Colors are the same as those shown in Figure 1. All steps determining branching are statistically significant.

Aggregate analyses on genetic data showed no significant differences between primary UM (PUM) and metastatic UM (MUM) samples, in terms of number of somatic coding SNVs, non-coding SNVs, CNVs, and SVs, indicating that the extent of genomic instability was here not associated with metastatic potential (Figures 3B–3E). Although singularly none of the main somatic drivers (chr 3, 6q, and 8q aneuploidies, as well as SNVs in BAP1, GNAQ, GNA11, SF3B1, EIF1AX) were computed as being statistically different, enrichment in PUMs versus MUMs showed very clear association trends (Figure 3A).

Figure 3.

Figure 3

Genetic Features in Primary and Metastatic Tumors

(A) Overview of all major somatic events with respect to PUMs and MUMs. Each circle indicates a specific genetic event; its center corresponds to the percentage of samples carrying this feature in PUMs versus MUMs, whereas its diameter indicates the total number of such samples. Asterisks indicate statistical significance. The gray area depicts the surface of the plot for which there is an enrichment in MUMs.

(B–E) Boxplots of different types of genetic alterations, at the genome-wide scale.

Remarkably, when considering specific levels of 8q amplification detectable by algorithms querying non-coding WGS data for CNVs and aneuploidies, we found a very clear association between metastatic potential and 8q ploidy of five copies or more (p value = 8.6 × 10−4, Figure 3A). In addition, single-copy amplification of 8q (ploidy = 3) was indeed associated with primary tumors (p value = 4.7 × 10−3, Figure 3A). Similarly, when mutational sets defining tumor sub-classes were considered, a significant association between sub-class B and metastases was identified (Figure 1, p value = 2.0 × 10−2).

With respect to metastatic samples, the aneuploidies identified correlated well with those of 66 liver metastases from UM investigated previously, detecting chr 3 monosomy (73%), 8q gain (89%), 6q loss (64%), 1p loss (47%), 8p loss (45%), 1q gain (35%), and 16q loss (32%).14 Similarly, the identified SNVs matched those on another study on five liver metastases.45 Finally, the identification of a SF3B1 mutation in one metastatic sample from our series is also in line with late-onset metastasis occurring in individuals with SF3B1 mutations.23

Mutations targeting BAP1 are one of the genetic landmarks of UM20 and were found here to be associated with classes A and B (p value = 1.2 × 10−4), classes that are in fact defined mostly by the presence of chr 3 monosomy. To test for functional inactivation of the BAP1 protein, we assessed its nuclear staining in histological sections of all tumors (Figure 4). Of the 33 samples (60%), 20 displayed loss of nuclear localization (Figure 4). Of these, 17 (85%) presented chr 3 monosomy and a coding SNV (a truncating SNV in 14 cases and a missense or an in-frame deletion in 3 cases), accounting for loss of heterozygosity and protein delocalization.

Figure 4.

Figure 4

Landscape of Genetic Alterations Involving BAP1 and Immunohistochemistry of BAP1 Protein

(A) Samples are ordered with respect to absence/presence of BAP1 nuclear localization, indicating loss or preservation of protein function, respectively. In general, loss of BAP1 function correlates with presence of a somatic SNV and loss of heterozygosity.

(B and C) Representative micrographs of paraffin-embedded UM samples showing absence (B) and presence (C) of BAP1 protein. Magnification: 252×.

The number of somatic SNVs involving coding and non-coding regions was strikingly low (Figure S2). Globally, the average load of coding mutations was 0.24/Mb (range 0.08–0.42/Mb, Table S4), one of the lowest detected so far in tumors. Comparison with other cancer types revealed that UM mutational load for coding regions was closer to that of pediatric tumors such as rhabdoid tumors, medulloblastoma (MIM: 155255), neuroblastoma, etc.,46 rather than that of adult cancers (Figure S3). Pediatric tumors typically develop over a shorter period than most adult malignancies, frequently harbor few driver mutations, and may thereby have fewer sources of heterogeneity, facilitating the assessment of both the genetic and epigenetic determinants underlying their pathogenesis. Our data seem to suggest that, similarly to pediatric tumors, UM development may rely more on epigenetic drivers of transformation and tumor progression, rather than the classical accumulation of genetic events observed in the vast majority of adult malignancies.

The number of non-coding SNVs was also relatively low (736 per tumor on average, range 371–1,347) and mostly proportional to the number of coding SNVs (17 per tumor on average, range 7–28) (Figure S2), confirming that, compared to both cutaneous and conjunctival melanomas, which also originate from melanocytes, UM follows a different oncogenic pathway, characterized by significantly fewer mutations.47, 48 In addition, we failed to identify any statistically relevant non-coding SNVs for tumor-specific sites that were present in four samples or more, suggesting the absence of common regulatory variants in the landscape of these tumors, at least in our cohort.

Another difference between UM and cutaneous and conjunctival melanomas involved its mutational spectrum (Figure 5). Analysis of all coding and non-coding somatic single-nucleotide substitutions (SNSs) from our series showed the clear absence of an UV-induced mutation signature. This particular spectrum results from sunlight-driven formation of pyrimidine dimers on the DNA49 and is found in both cutaneous and conjunctival melanomas (Figure 6A).47, 48 Direct analysis of genes known to be involved in cutaneous melanoma, such as BRAF (MIM: 164757), NRAS (MIM: 164790), and NF1 (MIM: 162200) revealed no somatic mutations in UM, supporting again the notion that uveal and cutaneous melanomas have a different molecular etiology.50, 51 Conversely, the UM mutational spectrum was remarkably similar across all PUMs and MUMs and resembled that of apparently unrelated tumors, such as clear cell renal carcinoma, thyroid tumor, and glioblastoma (Figure 6). Notably, despite a different cellular origin, UM shares with these tumors recurrent genetic modifications; BAP1 mutations and chr 3 monosomy are frequently seen in clear cell renal carcinoma,52 while hotspot mutations in the first codons of EIF1AX are recurrent in papillary thyroid carcinoma.53

Figure 5.

Figure 5

Mutational Signature of Our Samples, for SNSs, Genome-wide

(A) Main graph: comparison of mutational load of the UM samples studied with respect to two melanomas of the conjunctiva (CM), sequenced and processed according to the same methods.48 The number of mutations is radically different in UM versus CM. Inset: mutational spectrum of each UM sample, in percentage, showing a relatively homogeneous pattern.

(B) Results of the analysis of mutational events according to the methods and the classification proposed by Alexandrov et al.54, 55 Three main signatures are detected in our samples, evocative of Alexandrov’s signatures 12/16, 1B, and 6. The different peaks indicate specific genetic contexts of the altered nucleotide and are ordered according to the original article.54

Figure 6.

Figure 6

Analysis of the Mutational Spectrum in Our Samples versus Other Cancer Types, in Coding Regions Only

(A) The spectrum from UM is dissimilar from those from cutaneous and conjunctival melanomas, which are dominated by UV light-induced events (C to T transitions) and is conversely closer to that of thyroid and renal papillary cancer.

(B) Principal-component analysis (PCA) of the same data, showing the relatedness of UM with a few cancers but, again, not with other melanomas. Dimensions of the PCA are indicated by the arrowed axes. Primary data other than UM are from previously published sources.48, 62 Cancer types are indicated by the following abbreviations. aLung: lung, adenocarcinoma; scLung: lung, squamous cell; pRenal: renal, papillary; ccRenal: renal, clear cell; GBM: glioblastoma multiforme; AML: acute myeloid leukemia; Cut_mel: cutaneous melanoma; CM: conjunctival melanoma.

Analysis of all specific SNS types along with composition of the flanking bases allowed us to determine specific mutational signatures for UM, according to the classification of Alexandrov et al.54 Our samples appeared to be enriched for signatures 12 or 16 (55% of the score), signature 1B (25%), and signature 6 (20%) (Figure 5B). Signature 1B corresponds to a rather ubiquitous pattern in cancer, resulting from the spontaneous deamination of 5-methyl-cytosine, which in turn is thought to correlate with the process of aging.54, 55 Conversely, signatures 12/16 and 6 are associated with defects in nucleotide excision repair and DNA mismatch repair, respectively.

A more global approach, considering the intersection between the SNVs detected in our series and the most frequently mutated genes in cancer (TCGA PANCAN list),56 also revealed a minimal overlap, limited to BAP1 and SF3B1 (Figure S4).

A non-negligible number of structural variants (SVs) such as deletions, duplications, inversions, or inter- and intra-chromosomal rearrangements were also observed (Figure S5). Only a few of these events were recurrent, indicating the absence of major common drivers constituted by genetic events involving large parts of the genome. Among these, however, there were three inter-chromosomal events that were present in at least two individuals (Figure 7). Three samples (PUM20, PUM18, and PUM5) had a translocation involving chr 6 and chr 8 disrupting UBE2W (MIM: 614277) and MYO6 (MIM: 600970) for the two first samples, respectively. The third event occurred in intergenic regions. Translocations between chr 13 and chr 17 (no genes involved) were also present in MUM9 and PUM5 and between chr 3 and chr 12 in PUM17 (no genes involved) and in PUM5 (disrupting KDM2B [MIM: 609078]) (Table S5). Although these translocations impacted roughly the same genomic areas, highlighting possible hotspot regions in UM genome, they neither targeted the same genes nor defined a specific tumor sub-category. Of note, one individual (PUM5) appeared to harbor a higher number of interchromosomal events and SVs than the average value of the other cases (Figure S5). Notably, this individual was also an outlier of our clustering analysis (Figure 1). However, neither the medical history nor tumor pathology revealed any uncommon feature, compared to the rest of the cohort.

Figure 7.

Figure 7

Circos Plot of All Somatic Interchromosomal Events, in All UM Samples

Red lines indicate events involving the same chromosomal regions in more than two individuals.

Amplification of chr 8q is a well-known and important feature of UM.17, 57, 58 Levels of chr 8q amplification seem to define prognostic status and metastatic potential in UM and differentiate class C from D (Figures 1 and 2). However, the molecular bases for this phenomenon are not known. One possible explanation is that the amplification is driven by the MYC oncogene (MIM: 190080), which lies in this region.59 By comparing the pattern of chr 8q amplification in our samples, we determined the minimal region of overlap, involving a 2.3 Mb fragment toward its telomeric site (chr 8: 126,404,000–128,682,000).

Surprisingly, this region was very close to MYC (chr 8: 128,748,314–128,753,680) but did not include it. Conversely, it harbored six other genes (POU5F1B [MIM:615739], FAM84B [MIM: 609483], TRIB1 [MIM: 609461], LOC100130231/LINC00861, LOC100507056/CCAT1, and LOC727677/CASC8). The most interesting of them was POU5F1B (POU Class 5 Homeobox 1B), a pseudogene of the POU5F1/OCT4 family, recently shown to be involved in prostatic and gastric cancer.60, 61 Real-time quantitative PCR experiments showed that only POU5F1B, TRIB1, LINC00861, and CCAT1 were expressed in our UM samples, but no statistically significant correlation between their expression levels and 8q amplification or tumor class could be detected. The same held true for the MYC transcript, suggesting that none of these genes play a key role in UM pathogenesis (Figures S6 and S7).

By using a WGS-based, untargeted approach to investigate the genetic components of UM, we had the unique opportunity of assessing its genomic landscape as a whole, from single-nucleotide variants to interchromosomal rearrangements, providing the basis for future functional studies that go beyond the scope of our analysis. The global picture emerging from our work indicates that, genetically, UM is a relatively atypical tumor, mostly in virtue of the paucity of somatic events that characterize it. Driver mutations are very few and are confined to a relatively low number of genes, such as BAP1, GNAQ, and GNA11. Other genes, including those that were identified in this study, may have a role in tumorigenesis, but they are nonetheless present in a small fraction of the samples studied. Conversely, larger events such as extended copy-number and structural variations seem to shape UM’s genome in a much more relevant way, possibly determining tumor progression and fate. Taken together, our results point to a critical role for non-canonical mechanisms of cellular transformation in UM development, where chromosomal rearrangements and non-coding SNVs potentially affecting distal regulatory elements may collaborate in the establishment of a permissive oncogenic landscape.

Acknowledgments

This work was supported by the “Fond'action contre le cancer” foundation, Switzerland, the Swiss National Science Foundation (Div III grant# 156260 to C.R. and SNF Professorship grant# PP00P3-157468/1 to N.R.), and the MEDIC Foundation. We wish to thank Complete Genomics (Mountain View, CA) for providing the whole-genome sequencing and Patricia Martin for her invaluable help. We thank as well the Lausanne Institutional Biobank and Danielle Minaidis for help in collecting biological material from affected individuals.

Published: October 13, 2016

Footnotes

Web Resources

Supplemental Data

Document S1. Figures S1–S7 and Tables S1–S3 and S5

Table S4. Somatic Coding SNVs Detected

Document S2. Article plus Supplemental Data

References

Associated Data

This section collects any data citations, data availability statements, or supplementary materials included in this article.

Supplementary Materials

Document S1. Figures S1–S7 and Tables S1–S3 and S5

Table S4. Somatic Coding SNVs Detected

Document S2. Article plus Supplemental Data