Towards Personalized Medicine in Melanoma: Implementation of a Clinical Next-Generation Sequencing Panel

Molecular diagnostics are increasingly performed routinely in the diagnosis and management of patients with melanoma due to the development of novel therapies that target specific genetic mutations. The development of next-generation sequencing (NGS) technologies has enabled to sequence multiple cancer-driving genes in a single assay, with improved sensitivity in mutation detection. The main objective of this study was the design and implementation of a melanoma-specific sequencing panel, and the identification of the spectrum of somatic mutations in a series of primary melanoma samples. A custom panel was designed to cover the coding regions of 35 melanoma-related genes. Panel average coverage was 2,575.5 reads per amplicon, with 92,8% of targeted bases covered ≥500×. Deep coverage enabled sensitive discovery of mutations in as low as 0.5% mutant allele frequency. Eighty-five percent (85/100) of the melanomas had at least one somatic mutation. The most prevalent mutated genes were BRAF (50%;50/199), NRAS (15%;15/100), PREX2 (14%;14/100), GRIN2A (13%;13/100), and ERBB4 (12%;12/100). Turn-around-time and costs for NGS-based analysis was reduced in comparison to conventional molecular approaches. The results of this study demonstrate the cost-effectiveness and feasibility of a custom-designed targeted NGS panel, and suggest the implementation of targeted NGS into daily routine practice.

disease free survival and overall survival 4 . The selective BRAF kinase inhibitors (Vemurafenib and dabrafenib) are effective in BRAF mutant melanoma; MEK inhibitors (trametinib and cobimetinib) show efficacy against both BRAF-and KRAS/NRAS-driven tumors; KIT inhibitors (imatinib, dasatinib, sunitinib and nilotinib) have demonstrated clinical responses in melanoma arising from acral, mucosal, and chronic sun-damaged cutaneous sites; and additionally, there are novel therapeutic monoclonal antibodies targeted against immunosuppressive molecules such as CTLA4, PD-1 and PD-L1.
Therefore, molecular diagnostics are increasingly performed routinely in the diagnosis and management of patients with melanoma. Conventional molecular analyses for detecting cancer somatic alterations have relied on methods such as Sanger sequencing and real-time quantitative polymerase chain reaction (RQ-PCR). These approaches have the limitation that are performed separately for each gene and therefore consume a high turn-around time. Moreover, Sanger sequencing has a relatively low sensitivity, and sometimes it can be challenging to detect somatic mutations, especially when tumor material is mixed with normal tissue. In this context, the development of next-generation sequencing (NGS) technologies has enabled to massively analyze millions of DNA segments in parallel, thus allowing to sequence multiple cancer-driving genes in a single assay, with improved sensitivity in mutation detection. One of the developed NGS methodologies is the new Ion Torrent sequencing platform, based on the detection of hydrogen ions released on each cycle of DNA polymerization. It has been described to be cost and time effective 5 , and its applicability in formalin-fixed and paraffin-embedded (FFPE) specimens with small amounts of DNA has been proved in several reports 6,7 . The use of a multi-gene screening panel may potentially allow a more personalized approach to cancer therapy by identifying less common but potentially actionable mutations.
In the present study, we have used Ion Torrent sequencing technology with the Personal Genome Machine (PGM) and a custom AmpliSeq Panel including 35 genes. The main objective was the design and implementation of a melanoma-specific sequencing panel, based on the inclusion of relevant melanoma-genes with clinical diagnostic, prognostic or treatment value. Moreover, we aimed to describe the mutation profile in a series of primary melanoma samples in order to provide new insights into the molecular subclassification of melanoma. The translation of these study results may provide further understanding of the molecular alterations that lead to the development of melanoma, and therefore may contribute to the improvement of a personalized medicine.

Results
Next-generation sequencing quality metrics and sensitivity assay. Quality control analysis revealed excellent performance of the panel (Fig. 1). Of note, the average total reads per sample was 1,305,617 with an average coverage of 2,575.5 reads per amplicon ( Figure S1). Moreover, 97.8% of targeted bases covered ≥100x and 92.8% ≥500x. The coverage uniformity of amplicon sequencing was very high with an average of 91.9%. In addition, the sensitivity assay performed with serial dilutions in two independent experiments confirmed a high sensibility of the IonTorrent PGM (Table 1). BRAF mutations could be detected at an allele frequency as low as 0.05% with a rate of false reads of 0.1%. Therefore adjusting the sensitivity to our error-rate, the limit of detection was 0.5% ( Figure S2).

Concurrent molecular alterations in melanomas.
Several genes harbored concurrent mutations among the studied genes (Fig. 3).We observed mutual exclusivity between BRAF hotspot mutations (p.V600E, p.V600K and p.K601E) and NRAS/KRAS mutations in almost all samples; however, one sample concurrently harbored BRAF p.V600E and NRAS p.Q61R mutations. Moreover, we found concurrent mutations in BRAF non-hotspot mutations (p.L584F and p.L597R) with NRAS (p.Q61L and p.G12S). Strikingly, one sample showed two concurrent NRAS mutations (p.Q61R and p.E62Q), another sample showed concurrent mutations on KRAS (p.Q61R and p.A146V), and another one two concurrent mutations in KIT (p.Y553S and p.Y578C). The clinicopathological features of patients with concurrent mutations are listed in Table 3.
To verify the variant calling accuracy, BRAF, NRAS, KRAS and KIT mutated samples were subsequently analyzed by HRM, Sanger and RQ-PCR. BRAF mutations were detected by HRM and Sanger in 54% (27/50) of the melanomas; the average allele frequency of those wt samples was 4.4% with an average variant coverage of 229.5 reads. NRAS and KRAS mutations were confirmed by RQ-PCR in all the samples, and KIT mutations were confirmed in 80% (4/5) of the samples. Above mentioned concurrent mutations of NRAS, KRAS and KIT in the same sample were also confirmed by direct sequencing.
TERT promoter mutations. The prevalence of TERT promoter mutations was 33% (33/100), The most frequent somatic changes were −146C > T and −124C > T that were detected in 45% (15/33) and 30% (10/33) of the samples, respectively. Other recurrent changes identified were −124/−125CC > TT (6%; 2/33), −138/−139 CC > TT (15%; 5/33) and −57 A > C (3%; 1/33). We found mutual exclusivity between all the variants. Turn-around time and cost comparison. In order to compare NGS with conventional molecular testing under routine laboratory conditions, we calculated the turn-around time and costs for three commonly melanoma-related genes (BRAF, NRAS, KIT) in eight samples. Starting from FFPE samples, we were able to isolate DNA, prepare libraries and sequence eight samples within two sequencing runs on 318v2 chips in approximately three working days. The costs of consumables and laboratory personal for NGS was 415.5€ per sample. Grouping turn-around time analysis for BRAF (exon 15), NRAS (exons 2, 3 and 4) and KIT (exons 9, 11, 13 and 17) genes, the conventional methods resulted in approximately five working days. As expected, the mean turn-around time for NGS-based analysis was lower in comparison to routine methods. The total cost for conventional methods was more expensive in comparison to NGS (Table 4).

Discussion
In this study we analyzed the spectrum of mutations in 35 melanoma-related genes using a targeted NGS approach on 100 melanoma samples. This in-depth analysis provides important insights into the molecular alterations that contribute to the development and progression of melanoma. In this regard, NGS analysis using commercial pan-cancer panel has been performed in a large cohort of melanoma samples 9 ; however as it has been pointed out the majority of genes analyzed in this commercialized cancer panels are not relevant in melanoma pathogenesis, and other recently described genes are not included 10 . Furthermore, in our analysis the deep coverage enabled sensitive discovery of mutations in as low as 0.5% mutant allele frequency, which may be important in FFPE specimens where tumor content may be low and DNA may be degraded. The applicability of FFPE in NGS was previously demonstrated in a study conducted by Chen et al. 11 , that sequenced a pair of matched fresh-frozen and FFPE tumor samples and found high concordances between both types of samples. This is relevant as in routine practice FFPE specimens are usually available and fresh-frozen samples are difficult to achieve, especially in melanoma. In addition, as described above and consistent with other reports 12 , sample processing in our hands  took approximately three working days, and therefore in comparison to conventional molecular approaches, with NGS the overall time for sample mutation detection is reduced.
Recent WES studies have demonstrated that melanoma has one of the highest rates of somatic mutations among all cancers 13 . In our cohort, we found that the great majority of the analyzed samples (85%) showed at least one somatic mutation among the 35 genes of our customized panel. Moreover, the prevalence of mutations detected in the study among the different melanoma-related genes is in agreement with previous reports. In our primary melanoma series, BRAF mutations were the most prevalent (50%) followed by TERT promoter mutations (33%). The prevalence of mutations in RAS family was also similar to other studies, as NRAS mutations were present in 15% of the samples and KRAS in 3%. As previously described, mutations were detected at different   frequencies across all melanoma subtypes; BRAF and NRAS were more frequent in SSM and NM, KIT mutations were the most prevalent in ALM, and TERT promoter mutations were significantly associated with NM 2, 14 . Consistent with previous reports, we found mutual exclusivity between somatic mutations in BRAF and NRAS/KRAS in 97% (97/100) of the samples 15 . Interestingly, we also identified concurrent point mutations on the KRAS (Q61R and A146V) and NRAS (E62Q and Q61R) in two of the samples. To our knowledge these concurrent mutations have not been described before, as previously reported cases show coexisting mutations in codon 12 16 . It has been suggested that these co-mutations may cooperate to activate the MAPK pathway 9 , and therefore the identification could have prognostic and therapeutic implications.
Moreover, we aim to analyze the prevalence of mutations in new recurrently mutated genes in melanoma. Such is the case of PREX2 (a guanine nucleotide exchange factor and a PTEN regulating protein) and GRIN2A (a ionotropic glutamate receptor), both recently described in WES studies with mutations distributed along the length of the entire gene, and in frequencies around 14% and 33%, respectively 17,18 . In our cohort, PREX2 frequency is consistent with previous reports; however GRIN2A prevalence was lower than previously described. In this respect, it is worth mentioning that all variants with a benign in silico prediction were excluded. Additionally, we have analyzed the ERBB4 gene, that was initially found to be somatically mutated in 19% of patients with cutaneous melanoma 19 , although subsequent studies have reported lower prevalences, in coexistence with BRAF or NRAS mutations 20 . The prevalence of pathogenic variants in our series was 12%, of whom 50% (6/12) were concurrent with BRAF/NRAS. Combination strategies for targeted therapies with BRAF inhibitors and ERBB family kinases inhibitors have been pointed out as a promising therapeutic option in the future. In addition, among the newly identified cancer genes are RAC1 and PPP6C. The RAC1 p.P29S mutation was present in 7% of the samples and as previously described, it was significantly associated with sun-exposed melanomas 21 . This activating mutation may have clinical implications as it has been associated to aggressive melanoma features 22 , and it has been described to regulate PD-L1 expression 23 and to confer resistance to BRAF inhibitors 24 . Hotspot mutations in PPP6C have been found with an overall prevalence of 8% 22,25 , significantly associated with NRAS-mutated melanomas 26 . We have identified an overall prevalence of 6%, with the p.R264C PPP6C hotspot in 4% of our samples, mostly in concurrence with BRAF mutations rather than NRAS mutations. Although several studies have highlighted the relevance of these recently described genes, clinical and therapeutic impact of such genes will need to be determined in future studies.
As recently described by TCGA 3 , our panel provides a genomic classification of melanoma with four different subtypes. It should be noted that we have a higher prevalence in the triple-wt subtype. This could be explained in part due to the lower prevalence in NF1 gene, as we did not performed WES and therefore our panel did not  Table 4. Turn-around time and cost comparison between Next-generation sequencing and conventional molecular analysis. NGS: Next-generation sequencing; PGM: Personal Genome Machine PCR: polymerase chain reaction; RQ-PCR: real-time quantitative polymerase chain reaction; SS: Sanger sequencing; HRM: Highresolution melting; NA: Non applicable. *The duration of the second emulsion PCR required is not included as it is done at the same time as the first sequencing run. **Costs: Cost of consumables and laboratory personal (as cost calculated from the time that the technician/physician is required for the analysis).
covered the entire CDS. An alternative explanation of the different prevalence is that our cohort only included primary melanomas.
Considering the recent development of targeted therapies, it is of great importance to identify molecular alterations that contribute to the appearance of resistances in a high percentage of the patients 27 . In this regard, our panel design covers genes that have been related to such mechanisms of resistance, as MAP2K1, MAP2K2, MITF, PTEN, PIK3CA, NF1, RAC1 and HOXD8. Therefore, application of this melanoma-specific panel could be of interest in the management of metastatic patients, not only for the detection of actionable mutations, but also to identify those patients who may more likely benefit from those treatments.
In summary, the cost-effectiveness (reduced DNA input amount, increase in sensitivity, low turn-around time, and simultaneous analysis of multiple cancer-driving genes) of our new NGS approach based on Ampliseq libraries and Ion PGM sequencing, strongly suggest its implementation in routine diagnostics. In this study we have successfully performed NGS in FFPE melanoma and we have found a wide variety of somatic mutations that may contribute to the pathogenesis of melanoma. Since there are multiple mutations in individual tumors, and each tumor has a specific genetic profile, characterization of the molecular alterations of individual samples seems to be necessary in order to develop a personalized medicine. In addition, application of this panel may also provide further information about the genetic mechanisms of resistance to available therapies, so it would be especially valuable for clinicians in the management of patients with metastatic disease.

Materials and Methods
Patients and tumor samples. One  DNA preparation. Genomic DNA was isolated from two 10-μm thick FFPE sections using Deparaffinization Solution and the GeneRead DNA FFPE Kit (Qiagen, Hilden, Germany) according to manufacturer's protocol. This isolation kit contains Uracil-DNA-Glycosylase (UNG) treatment that prevents formalin-fixation induced artifacts, which may lead to false-positive mutation callings 28 . DNA concentration was quantified by Qubit dsDNA HS Assay Kit (ThermoFisher Scientific). Genomic DNA with at least 10 ng/μL was subjected to library preparation.  Table S2.

Integrative molecular analysis by
Ion Torrent Library preparation. Library preparation was performed using the Ion Ampliseq Library kits 2.0 (Thermo Fisher Scientific). Multiplex PCR amplification of 10 ng of DNA was performed using the custom Ion AmpliSeq Primer Pool and the Ion AmpliSeq Hifi Mix. (Thermo Fisher Scientific) according to manufacturer's protocol. Primer sequences were partially digested using FuPa Reagent, and Ion Torrent adapters and Ion Xpress Barcodes were ligated with DNA ligase. Following adapter's ligation, amplicons were purified with Agencourt ® AMPure ® XP (Beckam Coulter), and subsequently quantification of the amplified library was performed using Ion Library Equalizer Kit (ThermoFisher Scientific) according to manufacturer's protocol.
Emulsion PCR and DNA sequencing. The library pool was clonally amplified in an emulsion PCR reaction using Ion Sphere Particles (ISPs) on the One Touch 2 Instrument, and subsequently template-positive ISPs were enriched on the Ion One Touch ES (Thermo Fisher Scientific) as described by the manufacturer. Enriched template-positive ISPs were subjected to sequencing on the Ion Torrent PGM on a 318v2 Chip (four samples per chip) (Thermo Fisher Scientific).
Sensitivity assay. Sensitivity was assessed by sequencing serially diluted DNA isolated from two cell lines: DNA isolated from HT-29 (ACC-299), a human colorectal adenocarcinoma cell line with a known heterozygous mutation in BRAF (p.V600E) and a polymorphism in MC1R (p.R160W), was diluted into DNA obtained from another colorectal adenocarcinoma cell line Caco-2 (ACC-169), known to be wild-type for the mentioned variations, in ratios of 1:1, 1:3, 1:24, 1:49, 1:99, 1:999 resulting in 25%, 12.5%, 2%, 1%, 0.5% and 0.05% dilutions of the mutated allele, respectively. To asses reproducibility sequencing runs of the diluted DNA were performed in two independent experiments. Sun related site (%) Non exposed 21 Occasionally exposed 63  Variant Calling and experimental validation. Data from sequencing runs were transferred to the Torrent Server, and Ion Torrent Suite Software was used to generate initial variant calling. Filtered variants were annotated using Ion Reporter software. Exclusion of non-exonic variants and synonymous mutations was carried out. Subsequently, several steps were used in order to filter out variants with low read numbers: a minimum depth of total coverage ≥500 reads, an each variant coverage of ≥20 reads, and P-value < 0.01. Moreover, mutations were visually examined using Integrative Genomics Viewer (IGV) software. In addition, detected missense mutations in BRAF, NRAS, KRAS and KIT were subsequently confirmed by different methods. BRAF mutations were validated by High Resolution Melting (HRM) as previously described 29 , and Sanger's sequencing; NRAS and KRAS mutations were confirmed by RQ-PCR using AmoyDx ® KRAS/NRAS Mutations Detection Kit, according to manufacturer's instructions; and KIT mutations were validated by direct sequencing.
Prediction tools analysis. Databases such as COSMIC, TCGA and dbSNP were used to assess recurrent known mutations and to exclude reported germline polymorphisms. Furthermore, in all the variants of unknown significance we aimed to identify those ones likely to impact protein function using five prediction tools such as Provean, SIFT (Sorting Intolerant From Tolerant), PolyPhen-2 (Polymorphism Phenotyping v2), SNPS&GO, and Condel, which use algorithms that predict the effect of amino acid substitution on the protein structure and function. We excluded the variants predicted as "benign" by at least three of the five applied prediction tools.
Detection of TERT promoter mutations by direct sequencing. Mutational status of the TERT promoter region (from position −27 to −286 from ATG start site) was determined by PCR and Sanger sequencing. A 260 bp region was amplified using a pair of primers previously described 30 . The PCR was performed in a conventional thermal cycler using the following cycling conditions: initial heating at 95 °C for 5 min followed by 40 cycles of 45 s denaturation at 95 °C, 45 s annealing at 59 °C, 36 s extension at 72 °C, and, finally, 72 °C for 10 min. PCR product was subsequently sequenced on the ABIprism 3130 (Applied Biosystems).
Statistical analysis. Quantitative variables were summarized by their mean and standard deviation, and categorical variables by relative and absolute frequencies. The relationship between mutations and clinicopathological features was evaluated using logistic regression analysis with estimation of OR and 95% CI. P values < 0.05 were considered statistically significant. Computations were performed using the SPSSv21.statistical package (Chicago, IL).