ctDNA detected by ddPCR reveals changes in tumour load in metastatic malignant melanoma treated with bevacizumab

Bevacizumab is included in an increasing number of clinical trials. To find biomarkers to predict and monitor treatment response, cancer and angiogenesis relevant mutations in tumour and circulating tumour DNA (ctDNA) were investigated in 26 metastatic melanoma patients treated with bevacizumab. Patients with >1% BRAF/NRAS ctDNA at treatment start had significantly decreased progression free survival (PFS) and overall survival (OS) (PFS: p = 0.019, median 54 vs 774 days, OS: p = 0.026, median 209 vs 1064 days). Patients with >1% BRAF/NRAS ctDNA during treatment showed similar results (PFS: p = 0.002, OS: p = 0.003). ≤1% BRAF/NRAS ctDNA and normal lactate dehydrogenase (LDH) levels both significantly predicted increased response to treatment, but BRAF/NRAS ctDNA was better at predicting response compared to LDH at treatment start (OR 16.94, p = 0.032 vs OR 4.57, p = 0.190), and at predicting PFS (HR 6.76, p = 0.002) and OS (HR 6.78, p = 0.002) during therapy. ctDNA BRAF p.V600D/E/K and NRAS p.G12V/p.Q61K/L/R were better biomarkers for response prediction than TERT promoter mutations (OR 1.50, p = 0.657). Next generation sequencing showed that all patients with ≥2 mutations in angiogenesis-relevant genes had progressive disease, but did not reveal other biomarkers identifying responders. To conclude, ctDNA and LDH are useful biomarkers for both monitoring and predicting response to bevacizumab.

Tumour burden in malignant melanoma patients is in general followed using computer tomography (CT) to detect relapse or progression of disease. This method gives reliant results, but is costly and time consuming and can become a bottleneck in the hospital system. Liquid biopsies have therefore been proposed to be applicable in malignant melanoma to detect both response to therapy and early relapse 1 . These non-invasive biopsies use body fluids such as blood, urine and saliva to find circulating tumour cells (CTCs), circulating tumour DNA (ctDNA), tumour associated endothelial cells or tumour derived microvesicles. Sequencing of the tumour DNA permits detection of mutations classifying malignant cells in each patient. Tracking of these mutations in plasma can be performed using either PCR based methods or next generation sequencing (NGS) specifically targeting one mutation, or by defined gene panels that can be used on a larger patient subgroup.
For cutaneous malignant melanoma, mutations in BRAF and NRAS are found in approximately 50% and 20% of all patients 2 , respectively, and are regarded as early events in tumourigenesis. These mutations therefore represent attractive targets for monitoring tumour burden by ctDNA in blood samples. Another target is mutations in the TERT promoter. These are found in 60-70% of malignant melanomas 3,4 and are correlated with adverse outcome 5 , particularly when combined with BRAF or NRAS mutations 6,7 .
In this study 26 malignant melanoma patients with metastatic, non-resectable tumours were treated with bevacizumab, a monoclonal antibody specifically targeting VEGF-A 8 . The drug is currently investigated in several clinical trials, including melanoma, colorectal, ovarian and non-small cell lung cancer (ClinicalTrials.gov Identifiers NCT00790010, NCT03743428, NCT02884648, NCT03836066, respectively). In melanoma the drug gives significantly increased disease-free interval as monotherapy in an adjuvant setting 9 . Studies show that high serum concentration of Activin A 10 is associated with objective response to bevacizumab. Thus far, ctDNA has been detected in patients treated with bevacizumab 11 , but quantitative measurements have not been done. To identify biomarkers that are easy to measure we therefore performed mutational analysis using NGS on tumour biopsies and plasma samples, and digital droplet PCR (ddPCR) on plasma samples. We aimed to determine whether patients' mutational profile, ctDNA levels or lactate dehydrogenase (LDH) levels could serve as predictive markers for response to bevacizumab, prognostic markers for progression free survival (PFS) and overall survival (OS) or if changes in ctDNA or LDH during treatment could serve as pharmacodynamic markers. As care must be taken when using single gene analysis as measure for tumour burden, the investigated mutation should be a primary hit present in all tumour cells. We therefore compared the ctDNA fraction of TERT promoter mutations to BRAF/NRAS.

Results
Survival and response. Twenty-six patients with inoperable metastatic melanoma were treated with bevacizumab as monotherapy (Table 1, Methods section). Nine patients responded to treatment (34.6%); complete response was achieved in two patients (CR; 7.7%), partial response in five patients (PR; 19.2%) and stable disease longer than 6 months in two patients (SD; 7.7%). By June 2017, both complete responders were alive and disease free more than seven (P70) and 10 years (P31) after treatment was stopped (Fig. 1a). Seventeen patients progressed within 6 months of treatment (best overall response is shown in Fig. 1b). Median PFS was 64 days and median OS was 265 days. Details on the clinical trial are previously published 12 . Mutational landscape in tumour biopsies. To assess the mutation status across genes related to cancer and angiogenesis in tumour biopsies, we applied targeted sequencing using a panel of 419 genes in fresh frozen biopsies from metastatic lesions for 22 of the 26 patients ( Supplementary Fig. 1). The patients had a median of 10 mutations (range 0-21, Supplementary Table 1), and number of mutations did not correlate with treatment response. Somatic mutations were detected for all but one patient. Apart from BRAF and NRAS mutations, the only other recurrent mutations were KMT2C p.C1114R (P46, P57 and P61), DDX11 p.K272M and DDX11 p.K273R (both in P26 and P38). The two DDX11 variants were in close genomic proximity and had similar variant allele frequencies (VAFs) (P26: p.V272M; 10.8%, p.K273R; 12.8%, P38: p.V272M; 30.4%, p.K273R; 30.4%), but were found on separate reads. Based on the similarity of VAFs one can speculate that the DDX11 variants were biallelic in the tumor clone. These variants cannot be excluded as germline, however they do appear in the minor fraction of the cells. The three patients with lowest PFS (patient P14, P51 and P64) all harboured known cancer-related somatic hotspot mutations. Patient P14 (PFS 28 days, OS 66 days) had a loss-of-function mutation in CDKN2A (p.D108Y) 13 , a gain-of-function mutation in the catalytic domain of IDH1 (p.R132C) 14 and a stabilizing mutation in CTNNB1 (p.S45F) 15 . Patient P51 (PFS 9 days, OS 82 days) had a kinase activating mutation  16 , a loss-of-function mutation in the DNA-binding domain of TP53 (p.S241F) [17][18][19][20] , and an activating mutation of RAC1 (p.P29S) 21 known to increase PD-L1 expression significantly and giving resistance to RAF inhibitors in vitro 22,23 . Furthermore, all patients with ≥ 2 mutations in angiogenesis-associated genes had progressive disease (median PFS: 54 days vs 68 days, OS: 188 days vs 465 days (not statistically significant, Supplementary Table 2)).
BRAf and nRAS status in tumour biopsies. BRAF exon 15 and NRAS exon 1 and 2 were investigated by Sanger sequencing of formalin fixed paraffin embedded (FFPE) primary tumors from time of diagnosis, and metastatic tumors collected before inclusion in the clinical trial (Table 2). For BRAF, p.V600E was most frequent (n = 11), followed by p.V600K (n = 2) and p.V600D (n = 1). One patient was positive for BRAF p.V600E in the primary tumour, but positive for p.V600D in the metastasis. One patient had both BRAF codon 600 and NRAS codon 12 mutations. For NRAS, mutations were as follows; p.Q61R (n = 4), p.Q61K (n = 4), p.Q61L (n = 1) and p.G12C (n = 1) ( Table 2). NGS was performed on fresh frozen tumour metastasis biopsies for 24 of the patients. BRAF or NRAS hotspot mutations were found in 20 of 24 patients (BRAF p.V600E (n = 8) and p.V600D (n = 1); NRAS p.Q61K (n = 6), p.Q61R (n = 3), p.Q61L (n = 1) and p.G12V (n = 1)). The two methods identified in total 15 patients positive for BRAF hotspot mutations and 12 patients positive for NRAS hotspot mutations (Table 2). Discrepancies were . Two patients were complete responders (light green) and are still in remission seven (P70) and ten years (P31) after treatment start. Partial responders (dark green, n = 5), patients with stable disease (orange, n = 2) and progressive disease (red, n = 17) received 1-74 cycles of treatment before progression (red triangles). LDH and BRAF/NRAS ctDNA fractional abundances before treatment and during treatment (first sample point after treatment initiation) are indicated to the left of the Y-axis. (b) Best overall response to bevacizumab for 25 patients who had undergone at least one tumour assessment measured as the change from baseline in the sum of the largest diameters of each target lesion. P43 progressed clinically before first tumour assessment and is not shown. According to RECIST, progressive disease was defined by occurrence of new lesions in some patients in spite of stable target lesions. Stapled lines indicate cut-off for RECIST scores. PFS and BOR by August 2011 were previously published 12 , this figure depicts data updated per June 2017.
found for five patients when comparing Sanger and NGS (Table 2). This was assigned to four potential reasons; 1) the suspected primary tumor was collected 2-4 years before the metastatic lesion, and we cannot exclude that the patient had more than one primary tumour, 2) different tumour content between biopsies from primary and metastatic lesions, 3) different detection limits of NGS and Sanger sequencing, and 4) the fact that metastatic lesions in melanoma can be heterogeneous with respect to driver mutations 24 . BRAf and nRAS status in ctDnA. To evaluate whether ctDNA was detectable in plasma, cell free DNA (cfDNA) from Sanger or NGS-identified BRAF/NRAS mutation-positive patients was obtained at treatment start (n = 26 patients, range 0 to 310 days before treatment, median 3 days) and at follow up (n = 16 patients) with a total of 50 samples; Fig. 1a. Using mutation specific ddPCR-assays, ctDNA was detected in 23 of 26 patients (13 BRAF and 10 NRAS (Fig. 2, Table 2, Supplementary Table 3). This gives an 88.5% concordance with the sequencing results of the biopsied tumours. Diagnostic tumour biopsy FFPE DNA was available for one of the ctDNA negative patients (P67), and ddPCR analysis did not confirm the BRAF p.V600E detected by Sanger (802 wildtype copies detected, 0 mutant copies detected). P60 progressed 43 months after treatment initiation and the samples investigated were obtained 10 months prior to treatment and at treatment start. Low levels of wild-type copies were detected indicating low cfDNA concentration before and at inclusion (pre-inclusion; 450 copies, treatment start; 838 copies), however samples with similar cfDNA concentration and number of wild-type copies could be found to be ctDNA positive (e.g. P14, P31, P41, P53).

BRAF and NRAS ctDNA quantification.
To evaluate the response prediction potential of ctDNA, we calculated fractional abundance (percentage of mutation positive DNA fragments relative to total numbers of DNA fragments) for each sample. Patients with BRAF/NRAS ctDNA fractional abundance higher than 1% at treatment start had significantly lower PFS compared to patients with ctDNA fractional abundance ≤1% (median 54 vs 774 days, HR 2.52, p = 0.019) and OS (median 209 vs 1064 days, HR 2.44, p = 0.026) (Fig. 3a,b, Table 3,  Supplementary Table 3). Furthermore, patients with ctDNA fractional abundance ≤1% within 9 cycles of treatment (day 126) had significantly longer PFS compared to patients with continued elevated ctDNA (median 820 vs 56 days, HR 6.76, p = 0.002) (Fig. 3c). Corresponding numbers for OS were median 1064 vs 256 days (HR 6.78, p = 0.002) (Fig. 3d). Additionally, radiological responders (CR, PR and SD) had significantly lower BRAF/ NRAS ctDNA fractional abundance compared to non-responders (progressive disease; PD) both before treatment Patient Tumour -Sanger primary lesion Tumour -Sanger metastatic lesion Tumour -NGS § metastatic lesion cfDNA -ddPCR www.nature.com/scientificreports www.nature.com/scientificreports/ start (median 0.2% vs. 17.3%, PFS and OS p = 0.001) and during follow up (median 0% vs 12.5%, PFS and OS p = 0.049). Patients with BRAF/NRAS ctDNA fractional abundance ≤1% at treatment start also had higher chance of responding to therapy (achieving complete/partial response and stable disease) compared to patients with higher ctDNA fractional abundance (OR 16.94, p = 0.032) teRt promoter mutations. TERT promoter mutations are frequent in malignant melanoma 5 . In our study DNA from fresh frozen metastatic tumour biopsies were available for 19 of 26 patients. Sanger sequencing of the TERT promoter detected somatic mutations creating new Ets/TCF binding motifs in 15 patients (79%); 3 patients at −124 bp (C > T, chr5:1,295,228; rs1242535815) and 12 patients at −146 bp (C > T, chr5:1,295,250) ( Table 4) Table 4). TERT ctDNA fractional abundances were lower than BRAF/NRAS ctDNA fractional abundances (Fig. 4). Comparison of Sanger and ddPCR in biopsy samples confirmed that this was less likely due to differences in assay efficacy, however, the genomic region is repetitive and GC rich. TERT ctDNA-negative samples had little or no BRAF/NRAS ctDNA (Supplementary Table 3), or had little cfDNA available. Additionally, eight patients (42%) carried the likely benign c.−245T > C germline variant (rs2853669, chr5:1,295,349 hg19 coordinate) that has been shown to reduce promoter activity in presence of other promoter mutations. In our study, TERT promoter variants (c.−124C > T, c.−146C > T) did not predict treatment response (OR 1.50, p = 0.657), and no correlation was found between TERT promoter mutations and PFS or OS.

LDH.
To compare our BRAF/NRAS ctDNA findings with a more established marker for disease progression in melanoma we evaluated the levels of LDH in serum (n = 2-21 time-points per patient, Supplementary Table 5). Patients with LDH lower than upper normal limits (ULN) at therapy start showed significantly longer PFS www.nature.com/scientificreports www.nature.com/scientificreports/  Table 3). Samples taken during treatment showed similar results (PFS: median 774 vs 54 days, HR 5.24, p = 0.002; OS: median 1064 vs 221 days, HR 3.76, p = 0.003) (Fig. 3g,h). Furthermore, clinical responders had significantly lower LDH compared to non-responders both before treatment (PFS: p = 0.022, OS: p = 0.019), and during treatment (PFS: p = 0.008, OS: p = 0.017). Patients with normal LDH at treatment start had higher chance of responding to therapy (achieving complete/partial response and stable disease) compared to patients with LDH above ULN, but this was not statistically significant (OR 4.57, p = 0.190). ctDNA and LDH levels were positively correlated in both pre-treatment samples (r(24) = 0.805, p = 6.97 × 10 -7 ) and samples taken during treatment (r(14) = 0.517, p = 0.040). Furthermore, having both ≤1% ctDNA and normal LDH at treatment start (OR 12.01, p = 0.002) was better at predicting response compared to normal levels of LDH alone, but not ≤1% ctDNA alone. However, not all patients could be discriminated by pre-treatment LDH measurements. Patient P31 (Fig. 5) is a complete responder presenting with high pre-treatment ctDNA fractional abundance (23.3%) and high tumour load according to RECIST     www.nature.com/scientificreports www.nature.com/scientificreports/ indicating benefit from treatment, further supported by a decrease in size of target lesions in the right ovary and both lungs measured by CT (SLD 62 mm). Whereas LDH levels fluctuate in correspondence with tumour size, it is noticeable that the levels are consistently within normal range (blue, stapled line). The following increase in ctDNA (9.03%) two years after treatment initiation was mirrored by tumour growth in the ovary (Fig. 5). The patient was operated to remove this lesion, and at follow-up 10 years after inclusion the patient is still tumour free as measured by both CT and ctDNA (0%). ctDnA mutation load. To evaluate the potential of an extended gene panel for determining mutational load in cfDNA samples, a 900 cancer-gene panel designed by the Norwegian Cancer Genomics Consortium 25 was used for five selected patients (Fig. 6). As genetic variants in corresponding normal DNA were not available, cancer specific mutations were defined as present in less than 1% allele frequency in population databases. Samples from two patients (P19 and P43) were evaluated to investigate patients having short PFS and high ctDNA fractional abundance, one patient (P61) was evaluated due to lack of measurable cfDNA NRAS mutation, and two patients (P31 and P68) were evaluated due to large decrease and increase, respectively, in CT measured tumour burden. Patients with higher mutational burden had lower PFS compared to patients with fewer mutations (Supplementary Table 6). Number of mutations detected and ctDNA BRAF/NRAS fractional abundances were correlated (r(3) = 0.988, p = 0.002). Plasma mutation frequencies were comparable by both NGS and ddPCR (Supplementary Table 7), but ddPCR showed higher sensitivity.
Comparison of mutations in plasma and mutations detected by the 419-gene panel in tumour biopsies was possible for three of the above patients. For P19 and P31, plasma samples showed higher frequency of BRAF p.V600E compared to the biopsy (56.6% vs 28.8%, and 25.6% vs 16.9%, respectively), suggesting invasion of normal cells in the biopsy samples. All mutations found in the biopsy of P19 were verified in plasma, including a missense mutation of the putative tumour suppressor PPP6C (p.R264C) 26 that is sensitive to Aurora Kinase inhibition 27 (Supplementary Table 6, Fig. 7). For patient P31, two low frequent variants were found exclusively in the tumour biopsy (CEBPA: p.188_189del, VAF 10.3%; MYO5B: p.V1703A, VAF 10.1%).
For patient P61, no variants present in the biopsy were found in plasma (Fig. 6). Of the three variants found in plasma (Fig. 7, Supplementary Table 6), all were confirmed somatic in COSMIC (EP400, MCC, PTPRN2), and the EP400 variant was also seen in P68. A selection of 10 infrequent single nucleotide polymorphisms (SNPs) were compared between the biopsy and plasma sample for patient P61, confirming that the samples were from the same patient. The sample was below detection levels for NRAS p.Q61R ctDNA (0.05%), explaining the lack of tumour specific mutations found in plasma by the 900-gene panel. This patient did not progress until 820 days after the time-point for plasma NGS analysis, and all continuous blood sample time-points (n = 4, Fig. 1) were negative for NRAS p.Q61R ctDNA by ddPCR (Supplementary Table 3).
Patient P43 with the highest ctDNA fractional abundance in this study (Fig. 7, Supplementary Table 6), presented with the highest mutational load (n = 33), including missense mutations in GRM3 (p.D257N) and PIK3R5 (p.R741Q). Contribution of germline variants cannot be excluded, however 14 variants were verified as somatic as they are not present in any population databases. The two evaluated clinical responders, P31 (CR) and P61 (PR), both had mutations in TEK exon 10. Interestingly, P31 that had a significant reduction in ctDNA (from 23.4% to 1.7%) and tumour load (tumour volume reduction of 62%), also showed a reduction from 12 mutations www.nature.com/scientificreports www.nature.com/scientificreports/ to no mutations above VAF cut off-limits at corresponding time-points (Fig. 7, first and second time-point). An opposite pattern was observed for patient P68 that had an increase in ctDNA (from 1.8% to 39.8%) and tumour load (sum of largest diameters (SLD) from 135 to 158 mm with two new lesions) from first to second sample time-point. In plasma all four variants found in the first sample had an increase in VAF at progression of disease, and three new variants were identified at this point.

Discussion
ctDNA has the potential to play an important part in cancer precision medicine 28 . Particularly for cancers where new biopsies are difficult to obtain, the analysis of ctDNA can provide a tool for identifying treatment targets and monitoring treatment response, minimal residual disease and relapse. With the increasing numbers of clinical trials investigating the role of bevacizumab in cancer (ClinicalTrials.gov Identifiers NCT00790010, NCT03743428, NCT02884648, NCT03836066), commonly in combination with other drugs [29][30][31][32][33] , the discovery of an early biomarker for treatment response, using a method that allows repetitive, minimal invasive testing at low cost, would be of high value.
CtDNA investigations can range in scale from analysis of single mutations to whole-genome analyses. Monitoring treatment response by multigene analysis using massive parallel sequencing of gene-panels or exomes has the advantage of being less dependent on identifying early mutational events. It also permits detection of subclonal tumour response and progression, and can identify treatment targets and mutational load. It is however costly, time consuming and has lower sensitivity compared to ddPCR. As approximately 70% of malignant melanoma patients have a hotspot mutation in BRAF or NRAS 34 , and the mutation being an early event in most cases, single gene analysis is suitable for monitoring these patients during the treatment course. ddPCR provides the necessary sensitivity 35 , is fast, robust and cheap, and therefore well suited for future diagnostics.
A slight discrepancy was found between mutation profiles obtained by Sanger sequencing and NGS ( Table 2). As NGS is more sensitive than Sanger, BRAF p.V600E for P26 and NRAS p.Q61K for P53 were only detected by NGS. Additionally, high degree of normal cell-infiltration of some biopsies cannot be excluded as Sanger results for P38 and P42 could not be validated by NGS despite sufficient read depth (312 and 184 reads at BRAF codon 600, respectively). Furthermore, BRAF p.V600E for P41 could not be verified by NGS nor ddPCR (read depth 259). www.nature.com/scientificreports www.nature.com/scientificreports/ In our study we used ddPCR to investigate the clinical validity of ctDNA measurements in melanoma patients treated with bevacizumab monotherapy by quantifying BRAF, NRAS and TERT promoter mutations in plasma. We showed that both normal LDH and ≤1% BRAF/NRAS positive ctDNA before and during treatment reflected response to therapy by increased PFS and OS, and that the levels of ctDNA and LDH correlated, the latter finding being supported by research by Wong et al. 36 and Gray et al. 37 . Pre-treatment ctDNA has previously been found to not predict survival in malignant melanoma 38 , but we found that ≤1% ctDNA at inclusion was a positive predictor for response to bevacizumab and that ctDNA was a better predictor than LDH. Diem and colleagues investigated LDH as prognostic marker for metastatic melanoma patients treated with anti-PD-1 therapy, and concluded in concordance with our study that elevated levels at treatment start resulted in significantly shorter OS 39 . An association between LDH and PFS is also found in patients treated with mitogen-activated protein kinase inhibitors, suggesting this to be a universal mechanism independent by drug class 36 . Herbreteau et al. showed that detection of ctDNA after two weeks of treatment was 100% correlated with progressive disease in melanoma patients treated with anti-PD1 immunotherapy 40 . Felix et al. showed that low LDH at 3 and 6 weeks after ipilimumab treatment is significantly associated with increased OS 41 . However, normal LDH levels resulted in false prediction of treatment response in patients with progressive disease in our study, and did not reflect the change in tumour burden as accurately as ctDNA. This finding is in concordance with previous studies investigating BRAF inhibition and immune checkpoint blockade therapy 42,43 .
Recently, McEvoy and colleagues showed that detection of TERT mutation positive ctDNA was predictive for PFS (n = 15) 44 . A trend for increased melanoma-specific survival in TERT wild-type patients has been reported in larger studies 4 , and others have found that patients with mutations have shorter disease-free survival 7 . In our study ctDNA TERT promoter mutations were not correlated to PFS or OS. McEvoy also showed that for a substantial fraction of the patients with mutation positive tumours, TERT mutations could not be detected in plasma samples 44 , reflecting our findings. In addition, we observed that the fractional abundance of TERT was consistently lower that BRAF/NRAS, a finding supported by Calapre et al. 45 and Wong et al. 36 . This could indicate that only a subclone of the tumour contains TERT promoter mutations, or that false low VAFs due to technical issues caused by GC rich gene regions makes detection of this variant difficult in plasma.
In ctDNA detected by the expanded 900 cancer-gene panel 25 , we observed that patients with higher mutational burden had lower PFS compared to patients with fewer mutations. This was however only based on five patients, and the observation therefore needs to be supported in a larger cohort to conclude whether our results reflects mutational burden or is a consequence of higher levels of cfDNA and ctDNA in patients with low PFS.
Whereas high mutation load is associated with treatment response in melanoma patients treated with adoptive T-cell therapy 46 , correlating with improved OS after anti-PD-1 therapy (pembrolizumab or nivolumab) 47 , no such association has been explored for anti-angiogenesis therapies. Relevant to the angiogenic effects of bevacizumab, www.nature.com/scientificreports www.nature.com/scientificreports/ all patients harbouring two or more mutations in angiogenesis-related genes had progressive disease (not statistically significant). Additionally, the patients with shortest PFS had several known pathogenic somatic mutations in metastatic lesions. Particularly IDH1 p.R132C and CTNNB1 p.S45F, as well as TP53 p.S241F and RAC1 p.P29S are likely pathogenic variants. RAC1 p.P29S is an oncogenic hotspot in melanoma 48 , resulting in a spontaneous activation of this GTPase 21 giving resistance to vemurafenib and dabrafenib 23 . Patients with this mutation are however found to have upregulated PD-L1, suggesting they might benefit from immunotherapy 22 . IDH1 p.R132C is a hotspot mutation commonly found in acute myeloid leukaemia, but has more recently been identified in melanoma 49 , and frequently co-occurs with NRAS mutations 50 . The mutation alters cancer cell differentiation and results in metabolic reprogramming (reviewed in 51 ). Several IDH1 inhibitors have been developed, and are undergoing phase I clinical trials reducing 2HG levels produced by the mutant protein to normal physiological levels (reviewed in 52 ).
In treatment monitoring, single gene analysis require identification of an early event in the tumour and for melanoma BRAF/NRAS hotspot mutations are good candidates. This is supported by a high concordance between BRAF/NRAS mutations in plasma and biopsy (88.5%, Table 2), between NGS VAFs and ddPCR fractional abundances (Supplementary Table 7), and is in agreement with previous studies 53 . TERT promoter mutations are on the other hand present at lower fractional abundance (Table 4, Supplementary Table 4) and show lower concordance to tumour biopsy (73.3%), both in our and other studies 45 . Due to high GC content of the TERT promoter, and this being a chronologically later mutational hit, we believe BRAF p.V600 and NRAS p.Q61/p.G12 variants to be the best pharmacodynamic biomarkers to monitor treatment response in mutation positive malignant melanoma patients. Furthermore, the level of ctDNA detected by ddPCR correlates to tumour burden 54 , and studies show that undetectable ctDNA at baseline or within 8 weeks of treatment is an independent marker of response and prolonged survival in melanoma patients 38 . Collectively, this could indicate that ctDNA measurements may not only be a predictive marker, but also a prognostic marker indicating which patients would benefit from closer follow-up by CT, and potentially contribute to reducing number of CT analyses needed for patients with undetectable ctDNA.
To conclude, ddPCR of BRAF and NRAS hotspot mutations in plasma showed high concordance with NGS and Sanger sequencing of biopsy samples, and ctDNA and LDH levels reflected disease progression in the majority of patients. Moreover, low level of detected ctDNA at inclusion (≤1%) was a positive predictor for response to bevacizumab (odds ratio 16.94) and was a better predictor for response compared to LDH. NGS analyses www.nature.com/scientificreports www.nature.com/scientificreports/ showed that patients having mutations in two or more angiogenesis-relevant genes progressed earlier compared to patients with fewer mutations, however this was not statistically significant. No other novel genetic biomarkers for response to this drug were found. Our data support the clinical validity for ctDNA in monitoring treatment response and indicate that number of CT scans might be reduced for sub-group of patients. To confirm these findings, a randomized clinical trial with paired CT and ctDNA plasma sampling is ongoing (NCT02872259). For other clinical trials investigating the utility of bevacizumab in cancer treatment, we suggest ddPCR analyses monitoring ctDNA to be of higher clinical benefit compared to mapping the tumour mutational profile.

Methods
Study design. All participating patients provided signed informed consent before enrolment. Fifty-two patients were screened between April 2005 and August 2009 for this phase II, open-label, single-arm, single institution clinical trial for bevacizumab treatment (ClinicalTrials.gov Identifier: NCT00139360, 31/08/2005), performed at the Haukeland University Hospital, Bergen, Norway. Thirty-five patients met the inclusion criteria and were enrolled as described previously 12 . In short, each treatment cycle consisted of bevacizumab 10 mg/kg IV on day 1 in a 2-weekly schedule. Drug toxicity was assessed after each cycle; response rate was evaluated after every 4 cycles. Patients with disease progression or unmanageable toxicity were discontinued and offered further melanoma treatment at the clinician's discretion. For the current study, the 26 patients with tumour biopsies positive for BRAF p.V600 or NRAS p.Q61 mutations were included (Tables 1 and 2). PFS and response rates per August 2011 were previously used as endpoints in the study by Schuster et al. 12 , but updated PFS, OS and response ( Fig. 1) in addition to not previously published data on LDH is included in this study. Blood was collected in EDTA-tubes at treatment start, follow-up, and/or relapse for 26 patients, and immediately centrifuged at 4 °C, 1600 g. cfDNA was extracted from 2-3 mL plasma using the QIAamp DSP Circulating Nucleic Acid Kit (Qiagen) as recommended by the manufacturer. The study was conducted in accordance with the Declaration of Helsinki and the International Conference on Harmonization of Good Clinical Practice, and the protocol was approved by the Regional Ethics Committee of Western Norway and the Norwegian Medicines Agency (REK2012/910).
Response assessment. All patients were in clinical and/or radiological progression at the time of inclusion (0 to 19 days before treatment start, median 3.5 days). The primary endpoint was objective response (OR) defined as complete response (CR) or partial response (PR) according to RECIST 55 as well as disease control (DC) defined as CR + PR and including stable disease (SD) for more than 6 months. Briefly, response definitions were based on measurements of the longest diameter of target lesions; for CR, the target lesion was completely disappeared, for PR the sum of target lesion diameters decreased by ≥30% and for progressive disease (PD) the sum of target lesion diameters increased by ≥20%. OR and DC were calculated on the basis of investigator assessment. P43 had and P51 clinical disease progression before first radiological progression and was defined as having PD. Thus, best overall response (BOR) was not available for these patients. Image I and II of Fig. 5a have previously been published as part of a figure by Schuster et al. 12 , but is updated and revised in our study. Time to progression was defined as the time from enrolment to disease progression or death due to melanoma.
Tissue sampling and sanger sequencing. Tumour tissue was acquired as core biopsies or surgical biopsies, and either snap frozen in liquid nitrogen or fixed in formalin and cast in paraffin immediately upon dissection. DNA was extracted from tumour tissue manually dissected from three paraffin sections (10 μm) using the E.Z.N.A Tissue DNA Kit (Omega Bio-Tek, Inc.,). Due to limited metastatic tissue, Sanger was not performed on the same biopsies as NGS (exceptions seen in Table 2). The DNA was screened for mutations in BRAF (NM_004333) exon 15, as well as NRAS (NM_002524) exon 1 and 2, by direct Sanger sequencing after PCR amplification (primers are described previously [56][57][58] ). Mutations in the TERT promoter (NM_198253) were investigated in DNA from fresh frozen tumour metastasis biopsies by direct PCR amplification and Sanger sequencing. PCR conditions and primers are described in Supplementary Materials. nGS. Genomic DNA was isolated from 22 fresh frozen tumour metastasis biopsies and matched normal controls from whole blood using QIAmp DNA Mini kit (Qiagen) according to the manufacturer's instructions. Tissue from the metastasis occurring closest to time of inclusion was used for NGS. Targeted massive parallel sequencing of 419 genes was performed for tumour and matched normal DNA for each patient. These 419 genes represented a panel of 360 cancer related genes previously described 59 , extended with 59 angiogenesis related genes (Supplementary Table 8). Gene panel specific sequencing libraries were generated by hybridisation to custom RNA baits (Agilent Technologies) according to the Agilent SureSelect protocol, and sequenced paired-end (2 × 75 bp) on an Illumina MiSeq instrument. To explore the mutational profile present in plasma derived cell free DNA (cfDNA) from 5 patients, an extended, custom SureSelect (Agilent Technologies) in-solution capture panel (NCGC 900) developed by the Norwegian Cancer Genomics Consortium was used to enrich for exons of 900 cancer-related genes, selected promoters and introns frequently involved in fusions as described previously 25 . Libraries were constructed using the ThruPLEX Plasma-Seq kit from Takara, captured using the Agilent SureSelect protocol and sequenced paired-end (2 × 100 bp) on an Illumina HiSeq. 2500 (Illumina Inc.) using TruSeq SBS v3 chemistry.
Mutation calling. 419-cancer and angiogenesis-gene panel. Sequence reads were aligned to the human reference genome (NCBI build 37) using BWA 60 . Subsequent mutation calling was performed by Illumina's MiSeq reporter software v2.2 using default settings. To exclude germline SNPs, all called tumour variants with a variant allele frequency in matched normal DNA > 0.05 were removed. Further, a set of additional post processing filters were applied in order to remove false positive calls. These included; sequencing depth at < 1000 × , total read depth > 10x and absolute number of reads carrying the variant ≥ 3.
Droplet digital pcR. ctDNA was assessed using ddPCR ™ Mutation Detection Assay: NRAS p.Q61R, NRAS p.Q61L, NRAS p.G12V, BRAF p.V600K, BRAF p.V600D; PrimePCR ™ Mutation Assay: NRAS p.Q61K, NRAS WT for p.Q61K, BRAF p.V600E, BRAF WT for p.V600E; ddPCR ™ Expert Design Assay: TERT C288T_88, TERT C250T_88 (all from Bio-Rad). All samples were run in duplicate. Samples generating a total of ≤ 2 droplets positive for the mutation assay were defined as having ctDNA not detected. Droplets were generated using the QX200 Droplet Generator (Bio-Rad), and the PCR reaction is described in Supplementary Material. Droplets were read using the QX200 Droplet Reader (Bio-Rad) and data was analysed using QuantaSoft version 1.7.4. Results are presented as % ctDNA (number of droplets positive for mutant-assay/total number of mutant-and wt-assay positive droplets) defined as fractional abundance.
Statistical methods. PFS and OS related to ctDNA and LDH were assessed using the Log-rank (Mantel-Cox) test (GraphPad Prism Version 6.0d). Hazard ratios were calculated by Cox Multivariate regression and odds ratios by Logistic Regression (IBM SPSS Statistics Version 24.0.0.1); Contingency between responding groups was analysed by the Pearson's Chi-squared test describing two-tailed p-values (GraphPad Prism) and correlation between ctDNA and LDH was performed by Pearson's bivariate analysis describing two-tailed p-values (IBM SPSS Statistics).