Genetic Profiling of Advanced Melanoma: Candidate Mutations for Predicting Sensitivity and Resistance to Targeted Therapy.

BACKGROUND: Molecularly targeted therapy has revolutionized the treatment of advanced melanoma. However, despite its high efficiency, a majority of patients experience relapse within 1 year of treatment because of acquired resistance, and approximately 10-25% patients gain no benefit from these agents owing to intrinsic resistance. This is mainly caused by the genetic heterogeneity of melanoma cells.
OBJECTIVE: We aimed to validate the predictive significance of selected genes in advanced melanoma patients before treatment with BRAF/MEK inhibitors. PATIENTS AND METHODS: Archival DNA derived from 37 formalin-fixed paraffin-embedded pre-treatment advanced melanoma samples of patients treated with targeted therapy was used for next-generation sequencing analysis using the Ion Torrent platform. The AmpliSeq Custom Panel comprised coding sequences or hot spots of 23 melanoma genes: ATM, BRAF, CDK4, CDKN2A, CTNNB1, EGFR, HOXD8, HRAS, IDH1, KIT, KRAS, MAP3K8, MAP2K1, MAP2K2, MITF, MYC, NF1, NRAS, PAX5, PIK3R1, PTEN, RAC1, and RB1. The sequences were evaluated for genomic alterations and further validated using Sanger sequencing.
RESULTS: Our analysis revealed non-BRAF genetic alterations in 28 out of 37 samples (75.7%). Genetic changes were identified in PTEN, CDK4, CDKN2A, CTNNB1, EGFR, HOXD8, HRAS, KIT, MAP2K1, MAP2K2, MITF, MYC, NF1, PAX5, RAC1, and RB1. Fifteen known pathogenic mutations (single nucleotide variants or indels) and 11 variants of unknown significance were detected. Statistical analysis revealed an association between the presence of pathogenic mutations and time to progression during treatment with combination therapy.
CONCLUSIONS: Pathogenic mutations identified by gene panel sequencing have potential predictive value for targeted therapy of melanoma and are worth further validation in a larger series of cases. The role of some known mutations (e.g. CDK4R24, PTEN c.801 + 1G > A, CTNNB1S45F) as well as variants of unknown significance identified in this study (e.g. MITFR316K, KITG498S) in the generation of resistance to BRAF/MEK inhibitors should be further investigated.

Key points

Introduction

Melanoma treatment has dynamically developed in recent years. Currently, two main strategies are available for the treatment of patients with advanced melanoma: immunotherapy and targeted therapy. The latter comprises MAPK signaling pathway inhibitors that target either BRAF or MEK kinases. BRAF/MEK inhibitor therapy is administered to patients harboring BRAF mutations and results in response rates of more than 50% [1, 2]. The main limitations of this treatment are intrinsic (primary) and acquired drug resistance. Acquired resistance almost inevitably develops during therapy. Primary resistance affects approximately 20% of the patients, and these do not benefit from this treatment. The early identification of primary resistance to targeted therapy is particularly important since immunotherapy is also indicated in patients with BRAF mutated melanoma. Therefore, patients who are unlikely to benefit from BRAF/MEK-targeted therapies could be spared the side effects of these agents and treated with checkpoint inhibitors instead.

The main mechanism responsible for both primary and acquired resistance is reactivation of the MAPK signaling pathway and/or activation of the PI3K/AKT pathway caused by genetic alterations present before therapy or selected for during treatment. MAPK re-activation is triggered by activating mutations in genes involved in this signaling pathway such as NRAS, KRAS, MEK1, and MEK2, but also by BRAF amplification or alternative splicing of this oncogene [3]. Mutations in the RAS inhibitor, NF1 [4] and the hotspot mutation RAC1P29S [5] can also induce resistance to BRAFi/MEKi by activating this pathway. Loss of NF1 additionally causes activation of the PI3K/AKT signaling pathway and together with loss of PTEN and mutations in PI3K or AKT contribute to targeted therapy resistance by activating this alternative pro-survival pathway [6]. Deregulation of the aforementioned pathways can also be caused by overexpression or hyperactivation of RTKs (receptor tyrosine kinases) such as PDGFR, IGF1R, or EGFR, mainly due to epigenetic changes [6]. Signaling pathways regulating the melanocyte differentiation pathway, MITF [7], RB [8], or KIT pathway [9] are also potentially involved in the generation of resistance to targeted therapy in melanoma.

So far, the genes involved in resistance generation can be divided into two groups: those which seem to have a relatively well-established role in this process and those which were reported in individual studies but were not validated biologically. The first group comprises genes such as PTEN [10-12], NF1 [4], MEK [13], and RAC1 [14]. The second class consists of genes the role of which in intrinsic or acquired resistance is not fully elucidated; however, there is evidence for their potential involvement. For example, some of these genes have been identified in resistant samples by next generation sequencing (NGS). NGS is a technique that enables simultaneous sequencing of many genes in many samples and is widely used in high-throughput genomic analysis of cancerous cells/tissues. There are several NGS analyses of melanoma samples aimed at the identification of novel melanoma genes and potential therapeutic targets [15-19], investigation of epigenetic effects [20], and identification of biomarkers for immunotherapy [21, 22]. Deep sequencing of samples derived from patients who responded in various manners to targeted therapy enabled the identification of other mutated genes potentially involved in the generation of resistance to targeted therapy in melanoma. These are CDKN2A [11, 23], RB1 [10, 11] ], PIK3CA [12], AKT3, HOXD8 [12], PAX5 [23], MAP3K8 [24], and MITF [12]. All these genes are either involved in the regulation of the aforementioned MAPK and PI3K/AKT signaling pathways or they are tumor suppressors affecting drug resistance in other cancers (HOXD8, [25]; PAX5, [26]).The panel of selected genes and variants is of course not exhaustive and there are other putative genetic biomarkers of melanoma resistance to be discovered. However, they are so far the best candidates for developing a gene panel with potential predictive power.

Here, we have used NGS to analyze the coding sequence of 17 genes and hotspot sites of 6 genes of potential predictive significance in pre-treatment melanoma samples derived from 37 patients with different length of time to progression (TTP) on targeted therapy.

Patients, Materials, and Methods

Patients and Control DNA

Pre-treatment samples from 45 patients with BRAFV600-positive metastatic melanoma who were treated with targeted therapy at the Maria Sklodowska-Curie Institute, Oncology Center Gliwice Branch between 2014 and 2017, were included in this retrospective study. All patients were in good clinical condition (ECOG: 0–1) and had no serious comorbidities. The BRAF/MEK inhibitors (vemurafenib, dabrafenib/cobimetinib, trametinib) were administered daily. In Poland, combination therapy has only been reimbursed since March 2017, which is why the majority of patients (60%) in our project received monotherapy with maximal recommended doses of vemurafenib (960 mg, orally, twice daily) or dabrafenib (150 mg, orally, twice daily). Those treated with combined treatment additionally received trametinib (2 mg, orally, once a day) or cobimetinib (60 mg, orally, once a day). Patients with brain metastasis detected by computed tomography (CT) (29%; 8 out of 28 examined) received radiotherapy before the BRAF/MEK inhibitor treatment. Patient examination with blood tests was performed approximately every 28 days and CT was performed at least every 14 weeks of treatment.

Clinical outcome was assessed using TTP calculated from the date of initiation of the BRAF/MEK inhibitor therapy until progression documented by CT according to Response Evaluation Criteria in Solid Tumor (RECIST, ver 1.0 or 1.1) or clinical examination. Control samples were DNA isolated from a melanoma cell line (1205Lu) purchased from ATCC (American Type Culture Collection), and a whole blood sample of a healthy donor from whom informed consent was obtained. The study was approved by the Ethics Committee of the Maria Sklodowska-Curie Institute, Oncology Center Gliwice Branch (approval no. KB/430-29/17).

DNA Samples

We used archival DNA from formalin-fixed, paraffin-embedded (FFPE) samples isolated using cobas® DNA Sample Preparation Kit (Roche Molecular Systems, Pleasanton, CA, USA) in which BRAF mutation status was evaluated by the cobas® 4800 BRAF V600 Mutation Test (Roche Molecular Systems) according to the manufacturer’s instructions. All melanoma samples contained more than 50% of tumor cells and the median tumor cell percentage was 75%. The quantity of the DNA was measured using Qubit™ dsDNA HS Assay Kit and Qubit 3.0 Fluorimeter (ThermoFisher Scientific, Waltham, Massachusetts, USA) and the quality of selected samples was evaluated using Agilent D1000 ScreenTape Assay and 4200 TapeStation System (Agilent Technologies, Santa Clara, California, USA). When required, the DNA was purified and concentrated using the AgencourtAMPure XP solid-phase reversible immobilization kit (Beckman Coulter, Brea, CA, USA) according to the manufacturer’s protocol. Control DNA was isolated using Syngen DNA Mini Kit (Syngen Biotech, Wroclaw, Poland) and used as a quality control for NGS. The DNA from a healthy donor was additionally used as reference sample in CNV analysis.

Ion Torrent NGS

An Ion AmpliSeq Custom panel targeting coding regions of ATM, BRAF, CDK4, CDKN2A, CTNNB1, EGFR, HOXD8, HRAS, IDH1, KIT, KRAS, MAP3K8, MAP2K1, MAP2K2, MITF, MYC, NF1, NRAS, PAX5, PIK3R1, PTEN, RAC1, and RB1 genes was designed using the Ion AmpliSeq Designer tool (www.ampliseq.com). The panel size was 39.56 kb, contained 349 amplicons, with an in silico coverage of 78.79%. For each sample, 10 ng DNA (for one sample 8 ng) was used for automated library preparation with Ion AmpliSeq™ Kit for Chef DL8 and IonChef instrument with 21 cycles and by following the manufacturer’s protocol. The libraries were quantified with Ion Library TaqMan Quantitation Kit (Thermo Fisher Scientific). They were diluted to 50 pM and pooled for automated templating with Ion 510™ & Ion 520™ & Ion 530™ kit for IonChef Instrument. Sequencing was performed with the GeneStudio S5 system and Ion 530™ chips (24 samples/chip). The average total reads per sample was 935,377 with an average coverage of 2487 reads per amplicon. Moreover, 98.47% of targeted bases had coverage of ≥ 100× and 90.12% bases had ≥ 500×. The coverage uniformity of amplicon sequencing was high with an average of 93.4%. The average median depth after exclusion of low-quality samples (mean depth < 500) was 2400×.

Sequence Data Analysis and Selection of Samples

Sequence data were processed using the Torrent Suite 5.12.0 pipeline software optimized for Ion Torrent platform to perform raw data analysis and base calling, trim adapter sequences, remove low quality reads, and make alignments to human genome (hg19). Variant calling was performed with Ion Reporter Server 5.6 software and custom analysis pipeline (MelanSeq). The settings were as follows: Variant Types: SNV, INDEL, MNV, CNV; p ≤ 1.0E−5; allele frequency ≥ 5%; UCSC Common SNPs filtered out; alternate allele count ≥ 10. The final filtered variants were annotated using the Ion Reporter Server 5.6 software. Mutations were visually examined using .bam files generated with Torrent Suite—FileExplorer v.5.12.0.0 plugin together with Alamut Visual v.2.10 software (Sophia Genetics, Boston, Massachusetts, USA). The technical criteria of samples selection for further analysis were as follows: mean depth > 500, total read count > 200,000, and variant number < 150 (to avoid a risk for artefacts generated by improper formalin treatment). Based on these criteria, eight samples were excluded from the experiment after sequencing, and the final analysis was performed on 37 samples. To further reduce the risk of identification of unauthentic variants for samples with number of variants in the range of 10–150, the alternative allele frequency threshold was set to 10% for samples with 10–50 variants and 20% for those with > 50 variants. Copy number variants were detected with OncoCNV version 6.9, using the whole blood sample of a Healthy Donor as a control sample [27]. The analysis was performed only for genes having more than 3 amplicons sequenced (all but BRAF, CDKN2A, IDH1, NRAS, KRAS, and HRAS). ATM was covered by 17 amplicons and KIT by 9 amplicons. Coding sequences of the rest of the genes were covered in the range of 64–83% (mean 74.5%).

Interpretation of the Data

Variant annotation and significance prediction were performed with the aid of Alamut Visual v.2.10 software and web application VarSome [28]. The variants were grouped into three classes according to their biological significance: Group 1 containing mutations of known pathogenic significance and nonsense or frameshift mutations; Group 2 containing variants of unknown significance (VUS) but likely to impact protein function as assessed by the six prediction tools Provean, Mutation Taster, Sorting Intolerant From Tolerant (SIFT), Polymorphism Phenotyping v2 (PolyPhen 2), DANN [29], and Functional Analysis through Hidden Markov Models (FATHMM-MKL); and Group 3 containing VUS with inconsistent in silico results of significance prediction.

Sanger Sequencing

DNA fragments containing mutations of interest were amplified using Gold Taq DNA Polymerase (Syngen) and primers listed in Table S1 (Electronic Supplementary Material, ESM). The Sanger sequencing was performed by Genomed S.A. (Warszawa, Poland).

Statistics

For each of the categorical/qualitative variables, we carried out a univariate survival analysis. Survival curves were estimated using the Kaplan–Meier method and compared using the log-rank test [30]. For multivariate analysis (performed for two variables: therapy type and presence of Group 1 mutations) we used Cox regression on right-censored data with likelihood statistical test. To assess the robustness of the results, bootstrapping analysis was conducted using the censboot function [31] (from the bootstrap package) for bootstrapping survival models with censored observations. In the bootstrap procedure, we constructed 1000 bootstrap samples in which the underlying survival models were recalculated for each sample. A p value of less than 0.05 was considered to be statistically significant. Analyses were performed using R versions 3.5.1 and 3.6.0 and survival analysis package (R Core Team [32]).

Results

Characteristics of Detected non-BRAF Mutations (SNV, indels)

Out of 45 FFPE pre-treatment melanoma samples derived from primary or metastatic tumors, 37 passed our technical criteria (see Sect. 2.4) and were taken for analysis of coding sequences of 16 genes and hotspots of 7 melanoma genes (349 amplicons) using the Ion Torrent NGS platform. Baseline patient characteristics, therapies, and best responses are listed in Table 1. The patients treated with either monotherapy or combination therapy (see Sect. 2.1) were characterized with various responses to the treatment. The TTP ranged from 1 to more than 27 months and the median was 5.3 months, which was lower than median reported in clinical trials (approximately 7 months for monotherapy, [33] and 12.7 for combined inhibitors, [34]). We divided the patients into three groups according to their TTP: low responders representing primary/early resistance (TTP < 3 months as previously classified by Van Allen et al. [12]), high responders with TTP higher than the median for combination therapy in clinical trials (> 12.7 months), and medium responders with TTP between those two cut-off points (TTP = 3–12.7 months).

Table 1

Baseline patient demographic and clinical characteristics

Variable	No. of patients (n = 37)	%
Age, years (Median)	58
Sex
Female	16	43
Male	21	57
Location of primary tumor
Head and neck	2	5
Trunk	18	49
Extremities	12	32
Multiple melanoma	1	3
Lack of data	4	11
Histologic type
Superficial spreading type	7	19
Nodular type	10	27
Lack of data	20	54
Primary tumor thickness
T1	2	5
T2	2	5
T3	7	19
T4	9	25
Lack of data	17	46
TNM Stage IV
M1a	10	27
M1b	12	33
M1c	6	16
M1d	9	24
Therapy
BRAFi monotherapy	22	60
BRAFi/MEKi	15	40
Best response
CR	1	3
PR	17	46
SD	11	30
PD	8	21

CR complete response, PR partial response, SD stable disease, PD progressive disease

The gene panel consisted of genes that have been shown to influence the response of melanoma cells to targeted therapy in functional studies (e.g. MAP2K1, [13]; NF1, [4]) or have been identified in samples resistant to BRAF/MEK inhibitors (e.g. HOXD8, [12]; PAX5, [23]).

Using the Ion Reporter software and custom analysis pipeline, we filtered variants that met the required criteria (see Sect. 2.4). All samples had the BRAFV600 mutation identified previously with COBAS mutation testing. Three samples had BRAF non-V600E mutations: two with V600K and one with V600E and S602T. All samples differed in the frequency of mutated alleles within the cells, which ranged from 8 to 95% (mean 39%) and from 10 to 152% when values were normalized to the amount of tumor cells in the samples. The mean normalized allele frequency was 55%. The frequency of the mutant BRAF allele above 50% suggests amplification of the BRAF oncogene, or may by caused by miscalculation of tumor cell percentage in tumor samples. For technical reasons described in Sect. 2.4 we were not able to evaluate the copy number of BRAF.

As for the non-BRAF mutation analysis, 26 genetic variants (SNV, indels) in total were identified in 21 samples. Eighteen of these samples had one non-BRAF mutation, two samples presented with two additional changes each, and one sample had four alterations. The vast majority of these mutations (19) were of missense type. There were also two cases each of nonsense, frameshift, and 3 splicing mutations. We decided to include the pathogenic missense mutation in CTNNB1 in the list of genetic variants. This gene was identified in our pilot study using the Ion AmpliSeq™ Cancer Hotspot Panel v2 (unpublished data). This variant was not identified in this study since the target sequence was not covered by the primers designed for CTNNB1 sequencing by the Ion AmpliSeq Designer program. In total fifteen mutations were classified into Group 1, seven into Group 2, and four into Group 3 (see Sect. 2.5 and Table 2 for Group’s criteria). All identified variants are described in Table 2.

Table 2

Genetic variants identified in cutaneous melanoma samples

Group 1 (known pathogenic or likely pathogenic mutations)	Group 2 (VUS, likely pathogenic by in silico analysis)	Group 3 (VUS, inconsistent in silico results)
Oncogenes
MAP2K1 (MEK1)	MITF  chr3:70005615; NM_198159.2:c.947G > A (missense); p.Arg316Lys MITF  chr3:70014341; NM_198159.2:c.1505C > T (missense); p.Ser502Phe KIT  chr4:55592168; NM_000222.2:c.1492G > A (missense); p.Gly498Ser EGFR  chr7:55227884; NM_005228.3:c.1351C > T (missense); p.Arg451Cys	MYC
chr15:66729163; NM_002755.3:c.371C > T (missense); p.Pro124Leu		chr8:128751239; NM_002467.4:c.776C > T (missense); p.Thr259Ile
MAP2K2 (MEK2)		MITF
achr19:4099365; NM_030662.3:c.753G > A; (nonsense); p.Trp251Ter		chr3:69788764; NM_198159.2:c.16G > A (missense); p.Gly6Arg
CDK4		MAP2K2 (MEK2)
chr12:58145431; NM_000075.3:c.70C > T (missense); p.Arg24Cys		chr19:4101262; NM_030662.3:c.545C > T (missense); p.Ala182Val
bchr12:58145430; NM_000075.3:c71G > T (missense); p.Arg24Leu
chr12:58145430; NM_000075.3:c71G > T (missense); p.Arg24Leu
RAC1
chr7:6426892; NM_018890.3:c.85C > T (missense); p.Pro29Ser HRAS
chr11:533881; NM_001130442.2:c.175G > A (missense); p.Ala59Thr
CTNNB1
chr3:41266137; NM_001904.3:c.134C > T (missense); p.Ser45Phe
Tumor suppressors
PTEN	PTEN  bchr10:89717777; NM_000314.6:c.801 + 1G > T, (splicing site) CDKN2A  chr9:21971200; NM_001195132.1:c.158T > C (missense); p.Met53Thr PAX5  chr9:36966553; NM_016734.2:c.773C > T (missense); p.Pro258Leu	NF1
chr10:89717777; NM_000314.6:c.801 + 1G > A, (splicing site)		chr17:29661894; NM_001042492.2:c.5851C > T (missense); p.Pro1951Ser
chr10:89717777; NM_000314.6:c.801 + 1G > A, (splicing site)
achr10:89717703; NM_000314.6:c.731_737del (frameshift); p.Pro244fs
PTEN
chr10:89653782; NM_000314.6:c.80A > G (missense); p.Tyr27Cys
RB1
achr13:49027167; NM_000321.2:c.1734_1735delinsTT (nonsense); p.Arg579Ter
CDKN2A
achr9:21971148 NM_001195132.1:c.210del (frameshift); p.Asn71fs
CDKN2A
chr9:21971108; NM_001195132.1:c.250G > A (missense); p.Asp84Asn

In bold are variants taken for Sanger validation

aVariants not described in ClinVar, but with deleterious effect on the protein

bVariants not validated by Sanger sequencing

The most frequently mutated genes were PTEN, CDK4, and CDKN2A, which were mutated in five, three, and three samples, respectively. Two out of four different PTEN mutations were identified in the 3′ splicing site of exon 7. Mutation NM_000314.6(PTEN):c.801 + 1G > A affects splicing and was originally discovered as germline mutation associated with Cowden Syndrome, PTEN Hamartoma Tumor Syndrome, and other cancer syndromes. This mutation was present in two samples. The other variant at the same site (NM_000314.6(PTEN):c.801 + 1G > T) has not been described so far; however, in silico analysis suggests that this alteration affects the splicing mechanism in the same way as that of the aforementioned mutation. The third variant involved the deletion of seven nucleotides located in exon 7, which causes a frameshift (PTENP244fs). This mutation has not been described in ClinVar so far. The fourth PTEN mutation in our cohort was located at codon 27 (PTENY27C). This variant is a known pathogenic mutation associated with PTEN Hamartoma Tumor Syndrome and is also present in somatic samples.

We identified three pathogenic or likely pathogenic, recurrent mutations at codon 24 of the CDK4 gene: CDK4R24L (2 samples) and CDK4R24C. These mutations have been identified previously in both germline and somatic samples.

The third most commonly mutated gene was CDKN2A. CDKN2AD84N mutation is a likely pathogenic variant detected mainly as germline mutation (according to ClinVar), but also present in various cancers as a somatic mutation (according to Cosmic). The second mutation (CDKN2AM53T) located in the same exon is a variant of unknown significance; however, it is predicted to exert a deleterious effect on protein function. The third variant (CDKN2AN71fs) is a frameshift mutation caused by the deletion of a C nucleotide in the melanoma-associated commonly mutated region.

Other mutations with deleterious effects on the coded proteins were detected in the following genes: CTNNB1 (CTNNB1S45F), MAP2K1 (MAP2K1P123L), MAP2K2 (MAP2K2W251Ter), RB1 (RB1R579Ter), RAC1 (RAC1P29S), and HRAS (HRASA59T). Four of these mutations are known pathogenic or likely pathogenic mutations already identified in melanoma (CTNNB1, RAC1, MAP2K1, HRAS). The other two were considered of deleterious effect because of their nonsense (MAP2K2) or stop gain (RB1) mutations. Variants of unknown significance were identified in MITF (three variants), NF1, MAP2K2, KIT, CDKN2A, PTEN, EGFR, MYC, and PAX5. Most of them (7 out of 11) were considered deleterious in silico analyses (see Table 2). Among the 26 genetic variants identified, 14 were taken for verification and 11 were positively validated with Sanger sequencing (see Table 2). For one mutation, the results were ambiguous (CTNNB1), and two mutations were not validated. All three had alternative allele frequency on the border of Sanger sensitivity (around 10%).

Copy Number Analysis

We used OncoCNV package for copy number analysis. For technical reasons mentioned in Sect. 2.4, only the following genes were taken for analysis: ATM, CDK4, CTNNB1, EGFR, HOXD8, KIT, MAP3K8, MAP2K1, MAP2K2, MITF, MYC, NF1, PAX5, PIK3R1, PTEN, RAC1, and RB1. Due to the small number of reference samples (one) and possible amplification bias caused by technical reasons (AmpliSeq technology), only changes in the copy number of ≥ 1 and ≥ 3 for loss and gain, respectively, and those previously identified in melanoma were considered significant. The number of changes per sample ranged from 0 to 5 (median 1) and were mainly losses. Twelve samples showed copy number changes of potential biological significance. Eight samples had deletions in PTEN and four in PAX5. KIT and MITF were amplified in two samples each. One sample had 17.5 copies of MITF gene. Remaining samples showed either HOXD8 or ATM deletions or MYC or EGFR copy gain. Neither the number of CNV nor presence of any aberrations associated with TTP (data not shown).

Genetic Variants in Samples from Patients with Various Responses to Targeted Therapy

Twelve patients in our cohort presented with primary resistance (low responders). In seven of them, we identified additional pathogenic mutations in the following genes MAP2K1, MAP2K2, CDK4, and PTEN as well as variants of unknown significance in MYC, MITF and NF1. Four of thse patients did not respond to targeted therapy at all (best response as progressive disease). The patient with the shortest TTP (< 1 month) had a mutation in the MAP2K1 gene (MAP2K1P124L), previously shown to be linked with intrinsic resistance [13]. The other two patients with 1.1 and 1.4 months of TTP had two VUSs (MITFG6R and MYCT259I) or likely pathogenic VUS in the MITF gene (MITFR316K), respectively. The fourth patient showed nonsense mutation in the MAP2K2 (MAP2K2W251Ter) gene.

In the group of medium responders the vast majority of identified variants were present in samples derived from patients treated with combination therapy (11 out of 15) and all those patients but one harbored pathogenic or likely pathogenic mutation in the following genes: CTNNB1, CDK4, PTEN (4 cases), RB1, RAC1, CDKN2A (3 cases). Patient no. 19 presented with four variants in: HRAS (HRASA59T), CDKN2A (CDKN2AD84N), PTENSplice, and KITG498S. Genetic alterations in HRAS, PTEN, and CDKN2A are known pathogenic mutations, while the one in KIT is a VUS, but strong in silico evidence for its pathogenicity (Group 2) is available. Likely pathogenic VUS were also identified in PTEN and CDKN2A in this group.

The group of high responders consisted of 8 patients with an exceptionally good response to targeted therapy. All of them had TTP > 18 months. Three variants of unknown significance (Group 2) in PAX5, EGFR and MITF were identified with relatively low frequency (6.6–13%). No pathogenic mutation (Group 1) was detected in these samples. The summary of all the results along with selected clinical characteristics of the patients is presented in Fig. 1 and Table S2 (ESM).

Fig. 1

Summary of sequencing results and basic clinical characteristic of the patients. The identified SNV and indels as well as copy number variations are presented. No alteration were identified in the following genes of the panel: IDH1, NRAS, KRAS, PIK3R1, MAP3K8. The results are sorted by time to progression (months). C cobimetinib, D dabrafenib, T trametinib, V vemurafenib; Astrerisk signifies genes which were not taken for CNV analysis or presence of two mutations: V600E and S602T in sample no. 7; up arrow signifies patients still on therapy when the data was censored

Factors Associated with Worse Response to Targeted Therapy Evaluated Using Statistical Analysis

To evaluate the association of clinical and genetic factors with time to progression during targeted therapy, we performed Kaplan–Meier analysis and log-rank test. The analysis revealed two factors that significantly affected the TTP of the patients: treatment scheme, and presence of Group 1 variants. The combination therapy was more effective than monotherapy (p = 0.04 in log-rank test, Fig. 2) with vemurafenib + cobimetinib being the best treatment option in our group of patients (Fig. S1, ESM). Further, those patients who had additional known pathogenic or likely pathogenic mutations (Group 1) showed worse response to targeted therapy (p = 0.032 in log-rank test; Fig. 3). Since our cohort was heterogeneous in therapy schemas, we performed multivariate analysis which showed that both therapy type and pathogenic mutation (Group 1) were associated with TTP. The multivariate model proved strongly significant with both covariates significantly predictive of the outcome (Cox regression based likelihood test p = 0.003; Figure S2, ESM). Association of the presence of Group 1 mutations with TTP was not validated via Kaplan–Meier analysis with bootstrap percentile estimation. We additionally analyzed the influence of this genetic factor on TTP within more homogenous groups (monotherapy and combination therapy). The results revealed that the presence of Group 1 mutations significantly affected TTP in patients treated with combination therapy (Fig. S3, ESM). Those patients with no Group 1 pathogenic mutations had significantly longer TTP (p = 0.009 in bootstrap-based log-rank test). The results of selected comparisons are presented in Table 3.

Fig. 2

Kaplan-Meier curves for time to progression in patients treated with monotherapy or combination therapy. The scheme shows the results of association analysis of TTP with treatment schemas. The statistical significance was evaluated by log-rank test and the p value < 0.05 was regarded significant

Fig. 3

Kaplan–Meier curves for time to progression in patients harboring or not harboring mutation from Group 1. The scheme shows the results of association analysis of TTP with the presence of pathogenic mutations (Group 1). The statistical significance was evaluated by log-rank test and the p value < 0.05 was regarded significant

Table 3

Time to progression analysis

Variable	No. of patients	Median time to progression (months)	p (log-rank test)	p (bootstrap-based log-rank test)
Therapy
Monotherapy	22	4.55	0.04	0.13
Combination	15	9.7
Brain involvement before therapy
No	28	5.5	0.16	0.13
Yes	9	2.3
TNM Stage IV
M1a	10	5.6	0.42	0.70
M1b	12	7.6
M1c	6	4.4
M1d	9	2.3
PAX5 deletion
No	33	5.4	0.28	0.07
Yes	4	1.9
PTEN deletion
No	30	5.3	0.74	0.84
Yes	7	6.5
PTEN mutation
No	32	5.3	0.53	1.0
Yes	5	5.3
Group 1 mutations
No	24	5.6	0.032	0.58
Yes	13	4.2
Group 2 mutations
No	30	4.5	0.14	0.38
Yes	7	10.7
No. of variants
< 4	16	4.7	0.79	0.48
≥ 4	21	5.3
BRAF allele frequency (normalized)
< 50	15	5.6	0.35	0.52
≥ 50	22	4.7
Group 1 mutations (Monotherapy only)
No	17	4.7	0.31	0.13
Yes	5	2.3
Group 1 mutations (Combination therapy only)
No	7	19.6	0.002	0.009
Yes	8	4.7

C cobimetinib, D dabrafenib, T trametinib, V vemurafenib

Significant p values are in bold

Discussion

To improve our understanding of the genetic background of primary resistance to targeted therapy, we analyzed coding sequences or hotspots of 23 genes using 37 archival, pre-treatment advanced melanoma samples.

Our analysis revealed non-BRAF genetic alterations in 28 out of 37 samples (75.7%) and in 9 out of 12 patients (75%) who showed primary/intrinsic resistance (TTP of less than 3 months). We divided the selected variants into three groups according to their status of pathogenicity and performed association analysis of clinical and genetic factors with TTP. Statistically significant results were obtained for therapy schemas, and the presence of Group 1 mutations (pathogenic or likely pathogenic) but only in standard log-rank test. Lack of significantly better treatment effects of combination therapy over monotherapy (in bootstrap-based log-rank test) in this group of patients is probably caused by low efficiency of combined therapy of trametinib and dabrafenib observed in our cohort. On the contrary to clinical data [35, 36], in our analysis the combination of cobimetinib and vemurafenib appears to be a better treatment option than dabrafenib and trametinib. As the analyzed group was small, heterogeneous, and not representative, the results should be interpreted with caution. Of the genetic factors being analyzed, the presence of known pathogenic or likely pathogenic mutations correlated with shorter TTP. The P value was on the border of significance and the results were not confirmed by Kaplan–Meier analysis with bootstrap percentile estimation. The association was statistically significant (in both tests) when we analyzed only samples treated with combination therapy. Using multivariate analysis we discovered that both therapy type and the presence of pathogenic mutation were associated with TTP. High level of hazard ratio in both variables shows worse response to monotherapy in case of the presence of Group 1 mutations and better response to combination therapy in case of the absence of these mutations, respectively. The obtained results suggest that the analyzed gene panel is worth further validation in a larger, more homogenous cohort of patients. The genes of interest that comprised this group are: PTEN, CDK4, CDKN2A, MAP2K1, MAP2K2, RB1, HRAS, CTNNB1, and RAC1.

Mutations in PTEN have previously been detected in the samples of patients who showed primary resistance [11, 12]. Further, patients with inactivated PTEN showed a trend toward shorter TTP [37]. As PTEN is a negative regulator of PI3K/AKT signaling pathway, any alteration in PTEN that results in its suppression before treatment and is present at high frequency may lead to primary resistance to MAPK inhibitors. Nathanson et al. [37] observed a trend for shorter median progression-free survival in patients with baseline PTEN loss/mutation. Similarly, Trunzer et al. [38] observed lower PTEN expression in non-responders than in responders. In our study, all PTEN pathogenic variants were identified in patients with significantly shorter TTP than the median for the corresponding treatment group in clinical trials.

Another recurrent variant in our cohort was mutation in the CDK4 gene at codon 24. This is a known pathogenic mutation identified in both somatic and germline samples [39] as well as cell lines [8]. The identified missense variants CDK4R24C and CDK4R24L have not been associated with resistance to targeted therapy so far. Some in vitro studies suggest that mutations in CDK4 do not confer resistance to targeted therapy on their own, but they may potentiate the inhibitory effect of overexpression of cyclin D1 [40]. In our study, however, all identified mutations were present in pre-treatment samples of patients with TTP less than median for our cohort (5.3 months). Also, the frequency of CDK4 mutations (8%) in our cohort was higher than in a representative population (1.9%; [41]). This may stem from the fact that for the aim of the project, the cohort was deliberately enriched with patients with short TTP.

Another variant worth discussing is MAP2K1P124L. This mutation has been described as a resistance-generating mutation [13]. Mutations in this codon are recurrent variants detected in NGS studies [42], including those aimed at the identification of resistance genes [11, 12]. In our study, this mutation was identified in a patient with the lowest TTP (less than a month) and its allelic frequency was 34%. This supports its role in resistance to MAPK inhibitors. On the contrary, Shi et al. [43] showed that the MAP2K1P124S mutation does not confer resistance to BRAF inhibitors in vitro. Similar observations were presented by Trunzer et al. [38] in their genetic analysis of 132 samples. Indeed, not all MAP2K1 mutations, for example, C121S, generate resistance to MAPK inhibitors [11]. The discrepancies between the aforementioned two in vitro studies and genetic analysis of clinical samples are worth further investigation. The discrepancies may be because of differences in the molecular context in which the mutation exists; for example, coexistence with V600K mutation rather than V600E. In the studies in which BRAF mutation was already defined [11, 43], most mutations in the MAP2K1 gene coexisted with BRAF V600K, which is also true in our study. Moreover, both abovementioned in vitro analyses of MAP2K1P124L which showed an influence on resistance to targeted therapy [13, 43] were done on BRAFV600E/MAP2K1 double mutant melanoma cell lines, and in fact, the influence of this mutation on V600K melanoma cells was not investigated.

Another variant that has a potential influence on resistance is a recurrent hotspot mutation in codon P29 of RAC1 [19, 44]. This gain-of-function oncogenic event was detected in three pre-treatment tumor biopsies of low responders by Van Allen et al. [12] and medium responders by Long et al. [45]. Additionally, it was proved to confer resistance to RAF inhibitors in in vitro studies [5]. In this study, we detected this mutation in one sample derived from a patient treated with combination therapy. The alternative allele frequency was 10.5%, which is possibly too low to exert primary resistance. Still the patient had TTP shorter than the median for combination therapy in clinical trials (9.8 vs 12.7 months; [34]).

Novel genes, which have potential influence on resistance are CTNNB1 and HRAS. The variants were both identified in patients treated with combination therapy and short TTP (4 and 5.3 months respectively). In melanoma, mutations in CTNNB1 are relatively infrequent (4–7%; [15, 19]) and so far have not been associated with resistance to targeted therapy. This variant (CTNNB1S45F) is a well-described somatic mutation identified in many cancers including melanoma but with very low frequency [46]. Previous data suggest that Wnt/β-catenin signaling regulated by CTNNB1 may influence response to immunotherapy [47]. In vitro studies are required to elucidate whether this mutation may also generate functional resistance to targeted therapy. Similarly to CTNNB1, HRAS is also rarely mutated in melanoma and not linked to resistance. Its role in intrinsic or acquired resistance requires further validation.

A very interesting and original observation was the presence of PAX5 deletions in three low responders. This gene, which encodes paired box transcription factor, regulates B cell differentiation. Alterations in its expression contribute to neoplastic transformation [48]. To our knowledge, only two reports have associated it with melanoma biology [23, 26]. Thus, relatively high frequency of its alteration in our cohort should be verified. Our CNV analysis lacks robustness and the results were not validated with another technique. For these reasons we interpret the CNV results with caution. The only exception is the identification of 17.5 copies of MITF in sample no 24 derived from a patient with medium response to monotherapy (6.5 months). MITF amplification in pre-treatment samples of patients with less than median TTP was detected by Wheler et al. [23], and Van Allen et al. [12]. MITF is a melanocyte-specific transcription factor whose role in melanoma resistance to targeted therapy is not fully elucidated. In one study, its overexpression increased melanoma cells’ resistance to PLX4720 [12], whereas in another study, the opposite effect was observed [49]. Consistent with the observations of the latter study, Konieczkowski et al. [50] proved that low MITF expression and high NFĸB level may confer melanoma intrinsic resistance to MAPK inhibitors. In our study, apart from two cases of amplification in medium responders, we also detected two MITF VUS in patients who progressed very early (1.1 and 1.4 months). It would be interesting to evaluate the effects of these variants on the sensitivity of melanoma cells to BRAF/MEK inhibitors.

Conclusions

Despite some limitations, our study not only validated previous observations regarding potential predictive genes or mutations (MAP2K1P124L, RAC1P29S, PTEN) but also pointed out novel variants that are worth further investigation in the context of intrinsic or acquired resistance to targeted therapy (CTNNB1S45F, CDK4R24, HRASA59T, MITFR316K). The set of pathogenic mutations detected by our AmpliSeq Custom Panel has potential predictive value for treatment of melanoma patients with targeted therapy. Thus, the proposed genetic panel is worth further validation in a larger series of cases.

Below is the link to the electronic supplementary material.

Supplementary material 1 (PDF 432 kb)

Our sequencing gene panel has potential predictive power and following its further optimization may constitute part of a predictive test implemented in the clinics for the selection of melanoma patients for targeted therapy.

We have identified new variants potentially linked to melanoma resistance. Evaluation of these genetic changes could increase our knowledge of the genetic background of melanoma resistance to targeted therapy.