University of Groningen Genomic heterogeneity of clear cell renal cell carcinoma Ferronika, Paranita

(1)

Genomic heterogeneity of clear cell renal cell carcinoma

Ferronika, Paranita

DOI:

10.33612/diss.101437783

IMPORTANT NOTE: You are advised to consult the publisher's version (publisher's PDF) if you wish to cite from it. Please check the document version below.

Document Version

Publisher's PDF, also known as Version of record

Publication date: 2019

Link to publication in University of Groningen/UMCG research database

Citation for published version (APA):

Ferronika, P. (2019). Genomic heterogeneity of clear cell renal cell carcinoma. Rijksuniversiteit Groningen. https://doi.org/10.33612/diss.101437783

Copyright

Other than for strictly personal use, it is not permitted to download or to forward/distribute the text or part of it without the consent of the author(s) and/or copyright holder(s), unless the work is under an open content license (like Creative Commons).

Take-down policy

If you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediately and investigate your claim.

Downloaded from the University of Groningen/UMCG research database (Pure): http://www.rug.nl/research/portal. For technical reasons the number of authors shown on this cover page is limited to 10 maximum.

(2)

537042-L-sub01-bw-Ferronika 537042-L-sub01-bw-Ferronika 537042-L-sub01-bw-Ferronika 537042-L-sub01-bw-Ferronika Processed on: 24-10-2019 Processed on: 24-10-2019 Processed on: 24-10-2019

Processed on: 24-10-2019 PDF page: 51PDF page: 51PDF page: 51PDF page: 51

51

Mutational heterogeneity between

different regional tumour grades

of clear cell renal cell carcinoma

Paranita Ferronika, Gursah Kats-Ugurlu, Helga Westers, Sofia M. Haryana, Totok Utoro, Hanggoro Tri Rinonce, Raden Danarto, Kim de Lange, Martijn M. Terpstra, Rolf H. Sijmons, Klaas Kok

Manuscript in preparation

3

(3)

52

Abstract

Background

The histomorphological features of clear cell renal cell carcinoma (ccRCC) are clinically classified using the WHO/ISUP tumour grading system. Different grades can be observed within and between ccRCCs and have been associated with different clinical outcomes. Inter- and intra-tumour heterogeneity has also been observed at the genomic level. One of the questions in trying to understand the development of tumour features is if, and which, mutated genes drive those features. Only a limited number of studies have explored the possible associations between tumour grade and mutated genes in ccRCC, and we set out to investigate this further using a multiple sampling and next generation sequencing (NGS) approach in a series of ccRCCs.

Methods

Multiple regions were sampled from formalin-fixed paraffin-embedded ccRCC tumour blocks from seven patients. In 27 samples from six patients, we performed targeted NGS using a custom 42-gene panel based on the most frequently mutated 42-genes in ccRCC reported in public databases. In four samples from a seventh patient, we performed whole exome sequencing (WES) and array comparative genomic hybridisation for detection of copy number variants (CNVs). Mutated genes and the tumour grades of the samples in which they had been identified were compared both within and between all individual tumours. CNVs were compared across all samples from patient 7.

Results

We identified clear genetic heterogeneity within and across tumours. Within individual tumours, we did not observe any mutated genes that uniformly marked tumour grades, i.e. that were present in all samples of a particular grade from a single tumour but absent in samples with a different tumour grade. Looking across all samples, we did identify, in total, eleven genes that were only mutated in samples with a particular tumour grade. However, these were never observed in all samples with that particular grade. Increasing chromosomal instability corresponded with increasing tumour grade, but we observed minimal association between tumour grade and total mutational load in the WES data.

Conclusion

Our study confirms the genetic heterogeneity and tumour grade heterogeneity of ccRCC. Although a relatively small number of samples was analysed, genes were identified that could potentially be specific, though insensitive, markers of higher ccRCC tumour grades.

Keywords

Tumour heterogeneity, tumour grade, clear cell renal cell carcinoma, next generation sequencing, whole exome sequencing, target sequencing, gene variants, copy number variation.

(4)

53

1. Introduction

Kidney cancer is the 14th most common cancer worldwide, with clear cell renal cell carcinoma (ccRCC) accounting for 70% of total cases in adult patients [1, 2]. The global incidence of ccRCC was over 403,262 in 2018 with a mortality rate of 1.8 per population of 100,000 [1]. Conventionally, tumour grading is used to guide clinical decisions on ccRCC patient management and is based on the highest dedifferentiated morphological features in the tumour [3, 4]. Due to the rapid developments in molecular diagnostics, mutational profiling is being introduced clinically for a range of tumour types. Although tumour morphology has been shown to correlate with both tumour progression and prognosis [3], the prognostic value of mutational profiling in ccRCC is less clear.

In the next generation sequencing era, the most frequently mutated genes in ccRCC have been described in several large studies such as The Cancer Genome Atlas project [5]. Subsequently, the mutation rank system was introduced, in which mutations were scored and ranked by their functional impact based on in silico analysis using gene variant predictor tools [6]. Using this ranking, some genes were found to harbour more frequent mutations with higher predicted functional impact compared to other genes, and these were referred to as cancer driver mutations. For ccRCC, several cancer driver genes were identified by these studies, including VHL, PBMR1, SETD2, BAP1, and MTOR. However, knowledge on the associations between these driver mutations and the prognosis for ccRCC patients is still very incomplete. So far, in a few studies, inactivating mutations in BAP1, SETD2, and TP53 have been shown to be correlated with the prognosis of ccRCC patients [7, 8].

One of the challenges in studying the relationship between tumour characteristics, including the somatic driver gene variants, and clinical outcome, is intratumour heterogeneity. This heterogeneity exists both at the level of histological tumour grading and at the level of the mutational profile [9]. Although still limited in number, some studies have reported that variations in morphological features between different regions of an individual ccRCC, as indicated by regional tumour grades, are correlated with the mutation profile [9, 10]. In the last decade, several studies have shown intratumour heterogeneity in ccRCC at the genomic level. However, these studies did not examine this mutational heterogeneity in the context of regional tumour grading. As it is of interest to know if regional tumour grading, which plays an important role in current tumour classification and ccRCC prognosis, matches certain regional mutational profiles, we set out to study this in a series of six ccRCC patients using targeted sequencing and in one patient using whole exome sequencing (WES). In particular, we wished to know whether mutational profiles would be distinct between regions with distinct grades within individual tumours and if, in the total series of tumours, particular mutated genes would uniquely be associated with higher tumour grades.

(5)

54

2. Materials and Methods

2.1. Patient Materials and Tumour Grading

Formalin-fixed paraffin-embedded (FFPE) tissues were collected from nephrectomy specimens of seven patients with ccRCC. Hematoxylin and eosin (H&E) stained tissue sections were examined by two pathologists (P.F. and G.K-U.). Each section within those tumours with a homogenous histomorphology was assigned a regional tumour grade using the 4-grade system of the World Health Organization/International Society of Urological Pathology (WHO/ISUP grading) [2, 3]. For each grade in a patient, homogeneous tissue sections were collected using macrodissection from all available blocks without degraded tissue. In total, DNA was isolated from 31 tumour regions and seven matched normal kidney FFPE samples from seven patients (Table 1, Figure 1). The use of human material and clinical data from patients in this study has been approved by the Board Medical and Health Research Ethics Committee at the Faculty of Medicine, Universitas Gadjah Mada (UGM), Indonesia, as recognized by the Forum for Ethical Review Committee in Asia and the Western Pacific (FERCAP). Approval was given by the committee on September 4, 2013 (Reference: KE/FK/795/EC). The study was also performed in accordance with the University Medical Centre Groningen Medical Ethical Review board guidelines (project number 20190251, filed January 4, 2016) and Dutch ethical guidelines and laws, and complied with the regulations stated in the Declaration of Helsinki.

2.2. Targeted Sequencing and Whole Exome Sequencing

The targeted sequencing protocol is based on the Single Primer Enrichment Technique as implemented in the OvationTM Target Enrichment System (NuGEN, San Carlos, CA, USA). Our in-house designed set of landing probe areas covers the entire consensus coding region of the 42 most frequently mutated genes in ccRCC [5, 11] (https://cancer.sanger.ac.uk/cosmic, database accessed 23 June 2014) and genes associated with the VHL/HIF pathway and the (PI3K)-AKT-MTOR pathway in ccRCC with a frequency ≥1% [12] (see Supplementary materials and methods, Supplementary Table S1). This set includes 31 acknowledged cancer driver genes, 12 of them for ccRCC (https://www.intogen.org/search, release 2014.12) [6]. Library preparation was done according to the manufacturer’s protocol, starting from 500 ng of DNA. Enriched libraries were sequenced with the Illumina HISEQ 2500TM_{(Illumina, San Diego, CA, USA) using single-end} sequencing with 100bp reads.

WES was performed using the ThruPLEXTM-Fd Prep Kit (Rubicon, Ann Arbor, MI, USA) followed by exome capturing using SureSelectAll exon V2TM_{(Agilent, Santa Clara, CA, USA), starting from} 100 ng DNA sheared into 200-bp fragments. DNA integrity was measured using CaliperLabChipTM GX (Perkin Elmer, Waltham, MA, USA). Paired-end sequencing with 100-bp reads was done on the Illumina HISEQ 2500TM_{(Illumina, San Diego, CA, USA).}

(6)

55

Age (y) Sex TNM Tumour size

(max. diameter)

Tumor

region WHO/ ISUP Grade

Patient 1 53 F pT3aN0M0 10 cm P1T1 P1T2 P1T3 P1T4 P1T5 G2 G2 G2 G3 G2 Patient 2 62 F pT3bN1M1 11.5 cm P3T1 P3T2 P3T3 P3T4 P3T5 P3T6 P3T7 G4 (rhabdoid) G3 G2 G4 (sarcomatoid) G4 (sarcomatoid) G2 G4 (rhabdoid) Patient 3 57 M pT3bN0M1 8 cm P3T1 P3T2 P3T3 G3 G3 G3 Patent 4 46 F pT1bN0M1 7 cm P4T1 P4T2 P4T3 P4T4 P4T5 G4 G3 G4 G3 G4 Patient 5 59 M pT1N0M1 3.8 cm P5T1 P5T2 P5T3 P5T4 G2 G2 G3 G3 Patient 6 69 F pT3N0M1 4.2 cm P6T1 P6T2 P6T3 G2 G2 G2 Patient 7 55 M T2bN0M1b 12 cm P7T1 P7T2 P7T3 P7T4 G2 G4 G1 G4

Abbreviations: M, male; F, female; G, grade; TNM, tumour-node-metastasis classification of malignant tumors; FFPE, formalin fixed paraffin-embedded; WHO/ ISUP, World Health Organization/ International Society of Urological Pathology.

Table 1. Patient and sample characteristics

Abbreviations

M, male; F, female; G, grade; TNM, tumour-node-metastasis classification of malignant tumors; FFPE, formalin-fixed paraffin-embedded; WHO/ ISUP, World Health Organization/ International Society of Urological Pathology.

(7)

56

Figure 1. Sample selection and DNA isolation. Multiple FFPE blocks containing tumour and normal tissue were collected from each patient (A). Each FFPE block was serially cut in 10 μm sections (B). The first and the last section (3 μm in thickness) were stained with H&E and used as reference in identifying tumour nests with different WHO/ISUP tumour grade (C). DNA was isolated from each tumour region with different WHO/ISUP grade (G1-4) in multiple block sections.

(8)

57

2.3. Sequence Data Analysis and Somatic Mutation Identification

Data analysis was done using a pipeline based on the Genome Analysis Toolkit (GATK) best practice recommendation [13]. For both targeted sequencing and WES, two different variant calling algorithms have been used: HaplotypeCaller from GATK and FreeBayes [13, 14]. Called variants were annotated and filtered to identify true somatic mutations (see Supplementary materials and methods).

For each patient, the mutations present in each tumour region were classified into major clonal, minor clonal, absent, or inconclusive using the following approach. For each sample, we first determined the somatic variant with the highest mutant read frequency (MRF). Variants with a total number of mutant reads ≥5 and an MRF ≥50% of the highest MRF seen for that sample, were considered to be major clonal variants likely to be present in the majority of the tumour cells. Mutations with a total number of mutant reads equal to 3 or 4, or with a total number of mutant reads ≥5 and an MRF <50% of the highest MRF, were defined as minor clonal variants likely to be present in a minority of the tumour cells. If none of the above criteria were met and the total read count was ≥10, the mutation was defined as absent. If none of the above criteria were met and the total read count was <10, the mutation was considered inconclusive. Only variants with a major clone in at least one sample from a patient were included in the final analysis. Matched normal kidney samples were included to remove personal variants in the targeted sequencing and WES data. The Integrative Genomic Viewer was used to confirm the authenticity of the identified somatic mutations [15]. Sequencing data will be available in the European Nucleotide Archive repository.

2.4. Copy Number Variation Analysis

ArrayCGH was carried out using 500 ng genomic DNA from four tumour samples from Patient 7 using the Complete Genomic SureTag DNA Enzymatic Labelling Kit protocol and an OligoaCGH/ ChIP-on-Chip Hybridization kit (Agilent, Santa Clara, CA, USA). Details of the protocol have been reported previously [16].

2.5. Data interpretation

First, we analysed the distribution of mutations between different tumour regions per patient. We then grouped all 31 samples based on their regional tumour grade and analysed the distribution of mutations between the four different tumour grade groups. A Venn diagram was made to visualize the association of the mutated genes with specific tumour grades.

3. Results

3.1. Mutations Identified by Targeted Sequencing

The genes were selected based on their high mutational frequency in ccRCC, as described in several studies and in The Cancer Genome Atlas Network database [5, 11, 12, 17-19], and from genes associated with the VHL/HIF pathway and the (PI3K)-AKT-MTOR pathway in ccRCC [12].

(9)

58

The probe designs for the 42 genes are given in Supplementary Table S1. The average distance between different landing probes was optimized for FFPE based on preliminary data obtained with a catalogue assay from NuGEN (data not shown).

With the targeted sequencing assay, we analysed 27 tumour regions and six matched normal kidney tissues from six patients. The average tumour cell percentage within the samples is 78%, as established by the two pathologists. Sequencing resulted in a mean target coverage of 45x (Supplementary Table S2). An overview of the mutations detected is given in Figure 2 and Supplementary Table S3. For each patient, we focused on those variants that were regarded as a major mutation in at least one of the sections from that patient.

(10)

59

Figure 2. The mutated genes identified by targeted sequencing in Patient 1 (A), Patient 2 (B), Patient 3 (C), Patient 4 (D), Patient 5 (E), and Patient 6 (F). The samples are arranged in order from left to right based on tumour grade (G1- G4). Well-known cancer driver genes in ccRCC and other cancer(s) are highlighted in colors. The classification of the mutations is indicated by colors as shown in the legend at the bottom of the Figure. Mutations are included only if they are a major mutation in at least one of the subregions of a patient’s tumour. Reads counts and minor allele frequencies for the mutations are listed in Supplementary Table S3.

(11)

60

We analysed the distribution of mutations among different tumour regions in each patient. For Patient 1, four grade 2 regions and one grade 3 region were analysed from the primary tumour (Table 1). Although the samples were from five different locations within the primary tumour, we observed little variation in their mutational spectrum for our selection of 42 genes (Figure 2A, Supplementary Table S3). All five regions shared identical mutations in two cancer driver genes: VHL and ERBB3. The VHL mutation was present as a major clone in all regions, whereas the ERBB3 mutation was present as a major clone in one grade 2 region and as a minor clone in all other regions.

In Patient 2, we assessed two grade 2 regions, one grade 3 region, and four grade 4 regions from the available FFPE blocks. We identified rhabdoid and sarcomatoid dedifferentiation in the grade 4 regions. Mutations in VHL, BAP1, and ROS1 were shared by all regions, and were major clonal in the majority of them (Figure 2B, Supplementary Table S3). The VHL and ROS1 mutations were each identified as a minor clone in a single region. Mutations in TSC1 and KDM5C were present in some of the regions, but were never confined to one grade.

For Patient 3, all three FFPE blocks were of tumour grade 3. Three cancer driver genes, ARID1A, VHL, and PBRM1, were mutated in all tumour regions, with ARID1A and VHL being major clonal. The PBRM1 mutation was major clonal in two samples and appeared to be a minor clone in the third sample. Six additional genes, including the cancer driver genes MTOR, KMT2C, and ZFHX3, were mutated in only a single region (Figure 2C, Supplementary Table S3).

From Patient 4, two grade 3 regions and three grade 4 regions were analysed. Surprisingly no single mutation was shared between all regions. Interestingly, a VHL mutation was shared between only one of the two grade 3 regions and all grade 4 regions (Figure 2D, Supplementary Table S3). A PBRM1 mutation was shared by one grade 3 and two grade 4 regions.

For Patient 5, we analysed two grade 2 regions and two grade 3 regions. Mutations were detected for two ccRCC driver genes from our panel: VHL and KDM5C. The mutations were shared by all samples as a major or minor clone (Figure 2E, Supplementary Table S3).

For Patient 6, all available FFPE blocks were of grade 2, and we analysed three tumour regions. A VHL mutation was detected in all samples. A PBRM1 mutation was detected in two samples, but was inconclusive in the third due to low coverage. Three additional genes were mutated in only a subset of the samples (Figure 2F, Supplementary Table S3).

The VHL gene was mutated in all the patients’ tumours as a major clone, and was thus clearly a trunk mutation in all cases. For patients 3 and 6, PBRM1 appeared to be a second trunk mutation. BAP1, ARID1A, and KDM5C were second trunk mutations in one case each.

(12)

61

3.2. Mutational Profiles Identified by Whole Exome Sequencing

To look beyond the 42 panel genes, we carried out WES on tumour material of one additional patient. For this patient, we analysed four tumour regions and matched normal kidney cortex as control. We analysed one grade 1 region, one grade 2 region, and two grade 4 regions. Sequencing resulted in a mean target coverage of 41x in which 93% of the target region had >10x coverage (Supplementary Table S2).

We identified 39 nonsynonymous somatic mutations located in 38 genes, which were variously distributed among tumour regions and responsible for the intratumour heterogeneity (Figure 3 and Supplementary Table S4). For the four regions, we observed 14, 11, 15 and 19 major clonal mutations, respectively. Thus the mutational load of the grade 4 samples is slightly higher compared to the low-grade samples. This also holds true when the minor clonal mutations are taken into consideration. Out of the 39 somatic mutations found in this patient, 12 were shared by all regions, including two as a major mutation. Fifteen mutations were shared by multiple regions, either as major or minor clone, but were not present in all of them. Eight of the 38 mutated genes are known cancer driver genes, of which five are included in our custom panel for targeted sequencing. Three of these cancer driver genes (VHL, PBRM1, and ARID1A) are related to ccRCC development. Whereas all regions harboured a VHL mutation, the two grade 4 regions shared mutations in ARID1A, TP53, and PIK3CA. A PBRM1 mutation was only observed in one grade 2 region. The two grade 4 tumour regions shared the highest number of mutations in comparison to the grade 1 and grade 2 regions, including major clonal mutations in MYO7B, ARID1A, and FMN1.

(13)

62

Figure 3. Mutated genes identified in Patient 7 detected by whole exome sequencing. The mutations are classified as major or minor clonal. The classification of identified mutations is indicated by colors as shown in the legend at the bottom of the figure. The samples are arranged from left to right based on increasing tumour grade (WHO/ISUP G1-4). The identified somatic mutations involved 38 genes, including 8 cancer driver genes, as highlighted by colors. Reads counts and minor allele frequencies of the mutations are listed in Supplementary Table 4.

(14)

63 Intratumour heterogeneity was also observed for the pattern of copy number variations (CNVs) of patient 7 (Figure 4). All samples shared loss of chromosomal arm 3p, which is characteristic for ccRCC. The grade 1 and grade 2 tumour samples, P7T3 and P7T1, are near diploid and differ from each other at five genomic segments. P7T3 has a 9q gain, and appears to be mosaic for del18. For P7T1, we observed loss of 9p, 10q, and 16q. A small fraction of the cells of P7T1 appear to have a gain of chromosome 9q and a loss of chr20. The two grade 4 tumour samples, P7T2 and P7T4, had virtually identical and complex copy number profiles that showed much more CNVs than the lower grade tumour samples. The copy number profile in these two samples suggests an aneuploidic nature and includes a characteristic chromotripsis-like structure for chromosome 5.

3.3. Tumour-grade specificity of mutated genes

To obtain more insight into possible grade-related mutation profiles, we pooled the observed mutated genes of all six patients per tumour grade (Figure 5A). We included the mutations observed for patient 7, but limited ourselves to the genes present in our targeted sequencing panel. We also restricted ourselves to the major clones of each tumour sample. To represent the overlap between grades and genes uniquely mutated per grade, we constructed a Venn diagram (Figure 5B). VHL mutations were present in all tumour grades. The mutated genes that were present in multiple tumour grades, except tumour grade 1, were BAP1, ROS1, KDM5C, and PBRM1 (Figure 5B). PBRM1 was present in multiple tumour grade regions from four patients, whereas BAP1 and ROS1 were present in multiple tumour regions of only one patient (Figure 2). A number of genes were mutated in regions with only one specific tumour grade (Figure 5B), namely ERBB3, LRRK and DROSHA in grade 2; six genes including MTOR in grade 3; and TP53 and PIK3CA in grade 4. It is important to recognise that mutations in these genes were not identified in all samples with a particular tumour grade, so, for example, not all samples with grade 2 showed mutations in ERBB3.

(15)

64

Figure 4. ArrayCGH analysis of four tumour samples of Patient 7. The samples are arranged in order of WHO/ ISUP tumour grade (G1-4); P7T3 with tumour grade 1 (A), P7T1 with tumour grade 2 (B), P7T2 with tumour grade 4 (C), P7T4 with tumour grade 4 (D). The X-axis of the plots represents the genomic position starting from 1pter until Xqter, while the Y-axis indicates log2 intensity ratio between tumour and reference. Samples P7T2 and P7T4 which were classified as tumour grade 4, show a virtual identical, and highly aberrant CNV pattern, including a chromotripsic-like pattern for chromosome 5.

(16)

65

Figure 5. Distribution of major clonal mutated genes in different tumour grades among all seven patients, (A) shown by a bar chart for each tumour grade group and (B) by a Venn diagram indicating the overlap between different grades.

(17)

66

4. Discussion

For the six patients subjected to targeted sequencing, the majority of the mutations were indeed observed in the most frequently mutated genes in ccRCC [5, 20]. Most of these belong to the genes that encode proteins involved in transcription regulation [21]: PBRM1, ARID1A, KDM5C, BAP1, ZFHX3, and HUWE1. The first four genes have been recognized as chromatin modifiers, regulators of genomic architecture and DNA accessibility which are crucial for the cellular processes including gene expression programmes and DNA damage repair [22].

If tumour phenotypes in ccRCC were consistently driven by mutations in specific genes, these genes would have very likely ended up high in the ranks of the mutation databases and thus be captured by our panel, which is based on published gene mutation frequencies and cancer driver-gene rankings. Indeed, we observed genomic heterogeneity testing those genes. Eighteen genes of our 42 panel genes, including five ccRCC driver genes, contributed to the intratumour heterogeneity in our group of seven patients. Perhaps surprisingly, given its known role as an early driver of ccRCC development, this included VHL. VHL was a qualified as a minor clone in two samples, and found not to be mutated in one grade 3 sample from patient 4. In the last case, functional loss of VHL cannot be ruled out because it could have been inactivated by promoter methylation. This mechanism is known to occur frequently for VHL in ccRCC [12], and would not have been detected by our panel sequencing.

Although genetic heterogeneity was present in our series, at the level of individual patients we could not identify mutated genes that clearly marked one grade as opposed to others. As shown in Figures 2 and 3, for individual patients there were no major clone mutated genes that were present in all samples with a particular grade but not present in those with another grade. We acknowledge that this study is limited by the number of tumours and the overall homogenous grade of the primary tumours of two patients.

Using a targeted panel of the 42 most frequently mutated genes may not be sufficient to determine the full spectrum of subclonal mutations. WES did reveal a number of unique mutations for each sample from patient 7. Moreover, we observed a small increase in the number of major clone mutations for the samples with a high-grade. A previous study also reported an increase of the mutational load with the tumour grade of ccRCC tumours [23], although they did not make an intratumour comparison. An increase in the mutational load became more apparent in our WES patient when looking at large-scale genomic changes, i.e. CNVs (Figure 4). We saw a marked increase in chromosomal instability from the low-grade to high-grade samples in this patient, including a characteristic chromotripsis-like event for the grade 4 samples that carry a TP53 mutation. Chromotripsis events have been shown to have a strong association with TP53 mutations in medulloblastoma and acute myeloid leukemia [24]. Chromotripsis for chromosome 5 has previously been found in renal cancer cell line TK10 [25]. These observations are also in line with the intratumour heterogeneity that we observed in another patient [16]. If confirmed in a larger series, the total mutational load at the CNV level might be a phenomenon that could be used to distinguish individual tumour grades.

(18)

67 Our study identified 11 mutated genes as potential markers for tumour grade, but the value of these findings is limited by the small size of our series. Moreover, these 11 genes were found to be wild type in several samples with the same tumour grade, so there is no way that they can fully reflect the genomic changes that drive ccRCC evolution and heterogeneity at tumour-grade level. Apparently, variations in tumour grade are not accomplished by single mutations in one of the 42 most-frequently- mutated genes in ccRCC, or even in groups of genes representing specific pathways. Instead, the underlying mechanism may be an expressional shift in pathways caused by an interplay between rare mutations and epigenetic changes. In addition to DNA sequencing, gene expression (and gene methylation) data might help in understanding the biology behind the tumour grades. The ultimate goal would be help improve ccRCC patient prognostics, and this requires larger series of patients and, ideally, RNA as well as DNA analysis in their tumours.. Our study confirms the phenomenon of intratumour heterogeneity of ccRCC. Through multi-region sampling we were able to detect a higher number of mutations, as a subset of the mutations were present in local tumour regions only. The identification of these local mutations, including CNVs, may be crucial in the future when these alterations become predictive for prognosis as markers for metastatic potential and indicators for targeted therapy. Using data from multiregional samples would in fact help to distinguish trunk (driver) mutations from local mutations, which could help in selecting the optimal targeted therapy. Thus, either way, multiregion sampling may guide selection of the treatment that can most effectively target all, or the most threatening, tumour populations.

Authors’ contributions

PF and GK-U carried out the molecular genetic studies and drafted the manuscript. SMH, TU, HTR, and RD collected patient material and participated in the design of the study and drafting of the manuscript. JB and MMT participated in the study design and data analysis. KK, RHS, and HW participated in study design, analysis and drafting the manuscript. All authors read and approved the final manuscript.

Competing interests

None of the authors have conflicts of interest to declare. Acknowledgements

This work was supported by a Netherlands Fellowship Program Grant, Netherlands Organization for International Cooperation in Higher Education, NUFFIC [NFP-PhD.13/ 119 to PF] and International Research Collaboration Grant, Kementerian Pendidikan dan Kebudayaan, Direktorat Jenderal Pendidikan Tinggi, Indonesia [0094/E5.1/PE/2015 to SMH] . We thank Kate McIntyre for editorial advice.

(19)

68

References

Bray F, Ferlay J, Soerjomataram I et al. Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin 2018; 68: 394-424.

Moch H, Cubilla AL, Humphrey PA et al. The 2016 WHO Classification of Tumours of the Urinary System and Male Genital Organs-Part A: Renal, Penile, and Testicular Tumours. Eur Urol 2016; 70: 93- 105. Dagher J, Delahunt B, Rioux-Leclercq N et al. Clear cell renal cell carcinoma: validation of World Health Organization/International Society of Urological Pathology grading. Histopathology 2017; 71: 918- 925. Delahunt B, Srigley JR, Egevad L et al. International Society of Urological Pathology grading and other prognostic factors for renal neoplasia. Eur Urol 2014; 66: 795-798.

Network CGAR. Comprehensive molecular characterization of clear cell renal cell carcinoma. Nature 2013; 499: 43-49.

Gonzalez-Perez A, Perez-Llamas C, Deu-Pons J et al. IntOGen-mutations identifies cancer drivers across tumor types. Nat Methods 2013; 10: 1081-1082.

Hakimi AA, Chen YB, Wren J et al. Clinical and pathologic impact of select chromatin-modulating tumor suppressors in clear cell renal cell carcinoma. Eur Urol 2013; 63: 848-854.

Manley BJ, Zabor EC, Casuscelli J et al. Integration of Recurrent Somatic Mutations with Clinical Outcomes: A Pooled Analysis of 1049 Patients with Clear Cell Renal Cell Carcinoma. Eur Urol Focus 2017; 3: 421-427. Singh RR, Murugan P, Patel LR et al. Intratumoral morphologic and molecular heterogeneity of rhabdoid renal cell carcinoma: challenges for personalized therapy. Mod Pathol 2015; 28: 1225-1235.

Bi M, Zhao S, Said JW et al. Genomic characterization of sarcomatoid transformation in clear cell renal cell carcinoma. Proc Natl Acad Sci U S A 2016; 113: 2170-2175.

Forbes SA, Bindal N, Bamford S et al. COSMIC: mining complete cancer genomes in the Catalogue of Somatic Mutations in Cancer. Nucleic Acids Res 2011; 39: D945-950.

Sato Y, Yoshizato T, Shiraishi Y et al. Integrated molecular analysis of clear-cell renal cell carcinoma. Nat Genet 2013; 45: 860-867.

McKenna A, Hanna M, Banks E et al. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res 2010; 20: 1297-1303.

Garrison E, Gabor M. Haplotype-based variant detection from short-read sequencing. arXiv preprint arXiv 2012; 1207.3907

Robinson JT, Thorvaldsdottir H, Winckler W et al. Integrative genomics viewer. Nat Biotechnol 2011; 29: 24-26.

Ferronika P, Hof J, Kats-Ugurlu G et al. Comprehensive Profiling of Primary and Metastatic ccRCC Reveals a High Homology of the Metastases to a Subregion of the Primary Tumour. Cancers (Basel) 2019; 11. Dalgliesh GL, Furge K, Greenman C et al. Systematic sequencing of renal carcinoma reveals inactivation of histone modifying genes. Nature 2010; 463: 360-363.

Duns G, Hofstra RM, Sietzema JG et al. Targeted exome sequencing in clear cell renal cell carcinoma tumors suggests aberrant chromatin regulation as a crucial step in ccRCC development. Hum Mutat 2012; 33: 1059-1062.

Pena-Llopis S, Vega-Rubin-de-Celis S, Liao A et al. BAP1 loss defines a new class of renal cell carcinoma. Nat Genet 2012; 44: 751-759.

Forbes SA, Beare D, Boutselakis H et al. COSMIC: somatic cancer genetics at high-resolution. Nucleic Acids Res 2017; 45: D777-D783.

Huang da W, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc 2009; 4: 44-57.

de Cubas AA, Rathmell WK. Epigenetic modifiers: activities in renal cell carcinoma. Nat Rev Urol 2018; 15: 599-614. 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19. 20. 21 22

(20)

69 23.

24. 25.

Turajlic S, Xu H, Litchfield K et al. Deterministic Evolutionary Trajectories Influence Primary Tumor Growth: TRACERx Renal. Cell 2018; 173: 595-610 e511.

Rausch T, Jones DT, Zapatka M et al. Genome sequencing of pediatric medulloblastoma links catastrophic DNA rearrangements with TP53 mutations. Cell 2012; 148: 59-71.

Stephens PJ, Greenman CD, Fu B et al. Massive genomic rearrangement acquired in a single catastrophic event during cancer development. Cell 2011; 144: 27-40.

(21)

70

Supplementary materials

Supplementary Table S1. Custom gene panel for targeted sequencing

References

1. Network CGAR. Comprehensive molecular characterization of clear cell renal cell carcinoma. Nature 2013; 499: 43-49.

2. Gonzalez-Perez A, Perez-Llamas C, Deu-Pons J et al. IntOGen-mutations identifies cancer drivers across tumor types. Nat Methods 2013; 10: 1081-1082

No. Genes Mutation frequency in 512

ccRCCs of TCGA projects [1] Cancer driver status based on https://www.intogen.org/search [2] Number of target Number of target covered >90% 1 VHL 48.0% cancer driver gene in ccRCC 3 3 2 PBRM1 46.0% cancer driver gene in ccRCC 30 30 3 SETD2 16.0% cancer driver gene in ccRCC 21 21 4 BAP1 13.0% cancer driver gene in ccRCC 17 17 5 MTOR 9.0% cancer driver gene in ccRCC 57 57 6 NSD1 9.0% cancer driver gene in other cancer(s) 22 22 7 KDM5C 6.0% cancer driver gene in ccRCC 27 27 8 KMT2C 5.0% cancer driver gene in ccRCC 59 59

9 CSMD3 5.0% - 72 72

10 ARID1A 4.0% cancer driver gene in ccRCC 20 20 11 PTEN 4.0% cancer driver gene in ccRCC 9 9 12 ATM 4.0% cancer driver gene in ccRCC 62 62

13 LRP1B 4.0% - 91 91

14 PIK3CA 4.0% cancer driver gene in other cancer(s) 20 20 15 AKAP9 4.0% cancer driver gene in other cancer(s) 50 50 16 TP53 3.0% cancer driver gene in other cancer(s) 12 12 17 SMARCA4 3.0% cancer driver gene in other cancer(s) 34 34 18 SRGAP3 3.0% cancer driver gene in ccRCC 23 23

19 ROS1 3.0% - 43 43

20 KMT2D 2.8% cancer driver gene in other cancer(s) 54 54 21 CDKN2A 2.5% cancer driver gene in other cancer(s) 4 4

22 LRRK2 2.3% - 51 51

23 MLLT4 2.3% cancer driver gene in other cancer(s) 33 33

24 HUWE1 2.3% - 81 81

25 CUL3 2.0% cancer driver gene in other cancer(s) 17 17 26 NF2 2.0% cancer driver gene in other cancer(s) 17 17 27 ZFHX3 1.7% cancer driver gene in other cancer(s) 9 9 28 NFE2L2 1.7% cancer driver gene in other cancer(s) 5 5

29 PIK3CG 1.7% - 10 10

30 TCEB1 1.7% - 3 3

31 UBR5 1.7% - 59 59

32 MLL 1.4% cancer driver gene in other cancer(s) 36 36

33 RNF213 1.4% - 67 67

34 SMARCA2 1.4% - 33 33

35 TSC1 1.1% cancer driver gene in other cancer(s) 21 21 36 ARID2 0.8% cancer driver gene in other cancer(s) 21 21

37 IGF1R 0.8% - 21 21

38 KEAP1 0.8% cancer driver gene in other cancer(s) 5 5 39 ERBB3 0.6% cancer driver gene in other cancer(s) 28 28 40 SMARCB1 0.6% cancer driver gene in other cancer(s) 9 9 41 EGFR 0.6% cancer driver gene in ccRCC 30 30 42 CHEK2 0.3% cancer driver gene in other cancer(s) 14 14 References

1. Network CGAR. Comprehensive molecular characterization of clear cell renal cell carcinoma. Nature 2013; 499: 43-49.

2. Gonzalez-Perez A, Perez-Llamas C, Deu-Pons J et al. IntOGen-mutations identifies cancer drivers across tumor types. Nat Methods 2013; 10: 1081-1082.

(22)

71 Samp le Tum ou r pe rc ent ag e Se que nc ing appr oa ch Re ad s r aw Fr acti on rea ds alig ne d M ea n rea d len gt h Fr acti on du pl ica te Uni que rea ds alig ne d Uni que ba ses alig ne d On ta rg et ba ses Fr acti on us abl e ba ses M ea n ta rg et co ve ra ge Fr acti on ba ses ta rg et 10 X P1 T1 70% Ta rg et 1084222 0. 97 62. 40 0. 64 357213 2. 00E +07 1. 50E +07 0. 23 61. 00 0. 76 P1 T2 70% Ta rg et 1052453 1. 00 88. 00 0. 71 291497 2. 30E +07 1. 60E +07 0. 18 65. 00 0. 83 P1 T3 70% Ta rg et 4350925 0. 96 62. 90 0. 92 300547 1. 70E +07 1. 20E +07 0. 05 48. 00 0. 77 P1 T4 70% Ta rg et 518718 0. 99 82. 70 0. 76 117689 8. 00E +06 5. 00E +06 0. 14 23. 00 0. 59 P1 T5 70% Ta rg et 1896395 0. 96 60. 80 0. 94 104650 5. 00E +06 3. 00E +06 0. 03 15. 00 0. 43 P2 T1 80% Ta rg et 1307754 0. 96 60. 60 0. 77 264621 1. 30E +07 9. 00E +06 0. 13 39. 00 0. 57 P2 T2 90% Ta rg et 1407958 1. 00 83. 60 0. 75 331220 2. 40E +07 1. 70E +07 0. 15 71. 00 0. 76 P2 T3 85% Ta rg et 740511 0. 96 61. 30 0. 64 246110 1. 30E +07 1. 00E +07 0. 23 42. 00 0. 64 P2 T4 85% Ta rg et 1359297 0. 99 83. 10 0. 76 310272 2. 20E +07 1. 60E +07 0. 15 66. 00 0. 74 P2 T5 80% Ta rg et 1675835 0. 97 62. 10 0. 79 316249 1. 60E +07 1. 10E +07 0. 12 47. 00 0. 63 P2 T6 85% Ta rg et 1122030 0. 96 60. 70 0. 77 229173 1. 10E +07 8. 00E +06 0. 13 33. 00 0. 56 P2 T7 90% Ta rg et 2716742 0. 97 65. 70 0. 85 379723 2. 10E +07 1. 50E +07 0. 09 62. 00 0. 77 P3 T1 85% Ta rg et 25813870 0. 96 95. 40 0. 99 222683 1. 78E +07 1. 02E +07 0. 00 21. 93 0. 65 P3 T2 80% Ta rg et 9826301 0. 97 90. 87 0. 97 239000 1. 74E +07 1. 16E +07 0. 01 25. 21 0. 66 P3 T3 85% Ta rg et 16855881 0. 97 93. 75 0. 96 457899 3. 62E +07 2. 22E +07 0. 02 47. 20 0. 85 P4 T1 75% Ta rg et 16778946 0. 98 98. 00 0. 96 524908 4. 27E +07 2. 60E +07 0. 02 55. 26 0. 83 P4 T2 80% Ta rg et 19455326 0. 98 99. 65 0. 96 527354 4. 41E +07 2. 60E +07 0. 01 55. 30 0. 85 P4 T3 80% Ta rg et 28634520 0. 97 99. 92 0. 97 700515 5. 96E +07 3. 59E +07 0. 01 75. 93 0. 88 P4 T4 75% Ta rg et 24196114 0. 98 97. 70 0. 97 475184 3. 87E +07 2. 34E +07 0. 01 49. 72 0. 82 P4 T5 80% Ta rg et 15084913 0. 96 97. 81 0. 96 457959 3. 69E +07 2. 20E +07 0. 02 46. 79 0. 79 P5 T1 85% Ta rg et 7365836 0. 98 92. 71 0. 95 276953 1. 98E +07 1. 20E +07 0. 02 25. 70 0. 71 P5 T2 90% Ta rg et 7203330 0. 97 91. 22 0. 95 302660 2. 14E +07 1. 31E +07 0. 02 27. 91 0. 72 P5 T3 70% Ta rg et 6908896 0. 97 89. 42 0. 94 303636 2. 23E +07 1. 36E +07 0. 02 29. 05 0. 76 P5 T4 85% Ta rg et 10531575 0. 98 95. 17 0. 95 368850 2. 81E +07 1. 74E +07 0. 02 37. 24 0. 78 P6 T1 70% Ta rg et 13877303 0. 97 93. 36 0. 95 554940 4. 41E +07 2. 82E +07 0. 02 60. 03 0. 83 P6 T2 70% Ta rg et 11487580 0. 96 93. 11 0. 92 697097 5. 62E +07 3. 50E +07 0. 04 74. 17 0. 87 P6 T3 70% Ta rg et 18599565 0. 97 96. 72 0. 96 627004 5. 34E +07 3. 34E +07 0. 02 70. 86 0. 85 P7 T1 80% WE S 94079930 0. 99 101. 00 0. 39 53138216 4. 92E +09 2. 35E +09 0. 25 47. 00 0. 94 P7 T2 70% WE S 88393336 0. 99 101. 00 0. 51 39544100 3. 73E +09 1. 65E +09 0. 19 33. 00 0. 94 P7 T3 80% WE S 86079748 0. 99 101. 00 0. 55 35293871 3. 29E +09 1. 44E +09 0. 17 29. 00 0. 90 P7 T4 70% WE S 80609054 0. 99 101. 00 0. 57 31491644 2. 98E +09 1. 27E +09 0. 16 25. 00 0. 89 Supplementar y T able S2. Q ualit y r epor

t of sequencing data gener

at ed b y tar get ed sequencing (255855-bp tar get t errit or y) and whole ex ome sequencing (50390601-bp tar get t errit or y) for the se ven patients

(23)

72 Supplementar y T able S3. Somatic m utations identified b y tar get ed sequencing in six c cR CC patients Pat ient 1 Gene Pos iti on Ref Al t T1 T2 T3 T5 T4 Grad e 2 Grad e 2 Grad e 2 Grad e 2 Grad e 3 VH L Ch r3 :101838 00 A T 19/ 38/ 0. 50* 22/ 70/ 0. 31 27/ 56/ 0. 48 5/ 18/ 0. 28 5/ 13/ 0. 38 ERBB3 Ch r1 2:56479 077 T A 16/ 42/ 0. 38 11/ 75/ 0. 15 8/ 47/ 0. 17 3/ 10/ 0. 30 2/ 9/ 0. 22 Pat ient 2 Gene P os iti on Ref Al t T3 T6 T2 T1 T7 T4 T5 Grad e 2 Grad e 2 Grad e 3 Grad e 4 Grad e 4 Grad e 4 Grad e 4 rha bd oi d rha bd oi d sarc om at oi d sarc om at oi d BA P1 Ch r3 :524368 87 C A 19/ 39/ 0. 49 11/ 26/ 0. 42 24/ 89/ 0. 27 13/ 40/ 0. 33 25/ 43/ 0. 58 33/ 63/ 0. 52 36/ 54/ 0. 67 VH L Ch r3 :101837 62 C A 13/ 20/ 0. 65 6/ 19/ 0. 32 18/ 56/ 0. 65 4/ 19/ 0. 21 16/ 48/ 0. 33 27/ 55/ 0. 49 19/ 26/ 0. 73 ROS 1 Ch r6 :117683 890 T C 9/ 22/ 0. 41 3/ 24/ 0. 13 9/ 36/ 0. 25 6/ 20/ 0. 30 16/ 37/ 0. 43 14/ 26/ 0. 54 9/ 20/ 0. 45 KD M5C Ch rX :53222480 C A 2/ 39/ 0. 05 4/ 21/ 0. 19 6/ 55/ 0. 11 5/ 24/ 0. 21 5/ 56/ 0. 09 4/ 56/ 0. 07 12/ 28/ 0. 43 TS C1 Ch r9 :135781 248 G A 0/ 99/ 0. 00 0/ 84/ 0. 00 45/ 193/ 0. 23 24/ 84/ 0. 29 60/ 160/ 0. 38 0/ 130/ 0. 00 0/ 115/ 0. 00 Pat ient 3 Gene Pos iti on Ref Al t T2 T3 T1 Grad e 3 Grad e 3 Grad e 3 AR ID 1A Ch r1 :270927 46 A G 8/ 21/ 0. 38 9/ 22/ 0. 41 5/ 10/ 0. 50 VH L Ch r3 :101838 72 G A 17/ 38/ 0. 45 37/ 73/ 0. 51 10/ 20/ 0. 50 HU W E1 Ch rX :53574757 C T 6/ 12/ 0. 50 0/ 13/ 0. 00 0/ 8/ 0. 00 PBRM 1 Ch r3 :526686 92 C G 4/ 14/ 0. 29 6/ 14/ 0. 43 7/ 13/ 0. 54 M TO R Ch r1 :113077 09 G A 0/ 18/ 0. 00 0/ 27/ 0 6/ 10/ 0. 60 KM T2 C Ch r7 :151851 403 G A 0/ 10/ 0. 00 0/ 28/ 0. 00 6/ 10/ 0. 60 ZFH X3 Ch r1 6:72845 872 C T 0/ 26/ 0. 00 0/ 33/ 0. 00 5/ 14/ 0. 36 CS MD 3 Ch r8 :113933 976 G A 0/ 3/ 0. 00 0/ 28/ 0. 00 7/ 18/ 0. 39

(24)

73 Pat ient 4 Gene P os iti on Ref Al t T4 T2 T1 T3 T5 Grad e 3 Grad e 3 Grad e 4 Grad e 4 Grad e 4 VH L Ch r3 :101915 06 C G 2/ 15/ 0. 13 7/ 15/ 0. 47 8/ 13/ 0. 62 10/ 24/ 0. 42 9/ 14/ 0. 64 PBRM 1 Ch r3 :526206 65 G T 2/ 25/ 0. 08 4/ 25/ 0. 16 11/ 23/ 0. 48 1/ 35/ 0. 03 4/ 10/ 0. 40 N F2 Ch r2 2:30051 598 C T 5/ 11/ 0. 45 0/ 14/ 0. 00 3/ 17/ 0. 18 2/ 14/ 0. 14 2/ 13/ 0. 15 LR P1B Ch r2 :141359 029 G A 5/ 17/ 0. 29 0/ 15/ 0. 00 0/ 20/ 0. 00 0/ 28/ 0. 00 0/ 20/ 0. 00 CS MD 3 Ch r8 :113353 856 C T 6/ 23/ 0. 26 1/ 24/ 0. 04 1/ 16/ 0. 06 0/ 40/ 0. 00 0/ 17/ 0. 00 Pat ient 5 Gene Pos iti on Ref Al t T1 T2 T3 T4 Grad e 2 Grad e 2 Grad e 3 Grad e 3 VH L Ch r3 :101837 65 T G 7/ 14/ 0. 50 13/ 22/ 0. 59 6/ 24/ 0. 25 13/ 22/ 0. 59 KD M5C Ch rX :53224239 de lC 3/ 6/ 0. 50 7/ 16/ 0. 43 9/ 15/ 0. 60 11/ 20/ 0. 55 Pat ient 6 Gene Pos iti on Ref Al t T1 T2 T3 Grad e 2 Grad e 2 Grad e 2 VH L Ch r3 :101883 16 de l GC C AG insC CCA 8/ 22/ 0. 36 18/ 40/ 0. 45 21/ 48/ 0. 43 PBRM 1 Ch r3 :526760 53 de lA A 0/ 4/ 0. 00 8/ 26/ 0. 31 2/ 5/ 0. 40 LRRK2 Ch r1 2:40699 647 G T 14/ 48/ 0. 29 19/ 84/ 0. 23 0/ 62/ 0. 00 DROS HA Ch r5 :314015 68 C A 6/ 24/ 0. 25 0/ 26/ 0. 00 0/ 26/ 0. 00 HU W E1 Ch rX :53616722 G A 6/ 33/ 0. 18 0/ 37/ 0. 00 0/ 36/ 0. 00 *nu m bers in di cat e AL T read s, read d ept h an d M AF of in di cat ed m ut at io n fo r eac h an al ysed reg io n, respec tiv el y. Layou t o f sam pl es is id ent ical to Fi gu re 2. *n umbers indic at e AL T r eads , r

ead depth and MAF of indic

at

ed m

utation for each analysed r

egion, r

espectiv

ely

. L

ay

out of samples is identic

al t

o F

igur

(25)

74

Supplementary Table S4. Somatic mutations identified by whole exome sequencing in Patient 7

*numbers indicate alternate (ALT) reads, read depth and minor allele frequency (MAF) of indicated mutation for each analysed region, respectively (**) gene included in the targeted sequencing panel of 42 genes used in patients 1-6.

Gene Position Reference Alternate T3 T1 T2 T4

grade 1 grade 2 grade 4 grade 4

SPHKAP chr2:228881750 C G 29/80/0.36* 13/35/0.37 11/31/0.35 12/31/0.39 SRPX2 chrX:99922376 A T 9/16/0.56 5/15/0.33 9/11/0.82 5/11/0.45 VHL (**) chr3:10183694 GAGGCC G 9/28/0.32 4/19/0.21 9/22/0.41 6/13/0.46 CRTC1 chr19:18888079 G A 54/169/0.32 12/86/0.14 36/85/0.42 8/48/0.17 LIMS1 chr2:109300348 T C 10/26/0.38 2/21/0.10 5/14/0.36 0/12/0.00 CAMK1D chr10:12870789 G T 17/77/0.22 14/47/0.30 18/59/0.31 7/27/0.26 RGSL1 chr1:182458305 A T 8/26/0.31 10/23/0.43 9/13/0.69 3/10/0.30 PABPC3 chr13:25672098 C G 3/34/0.09 5/15/0.33 5/19/0.26 6/33/0.18 MMRN1 chr4:90857007 G T 4/15/0.27 4/22/0.18 8/17/0.47 3/8/0.38 PDZD2 chr5:32074097 T A 24/104/0.23 5/71/0.07 20/68/0.29 5/85/0.06 INTS1 chr7:1519120 C A 14/89/0.16 1/59/0.02 18/73/0.25 0/41/0.00 GDPD2 chrX:69650918 C G 5/21/0.24 1/12/0.08 5/17/0.29 0/10/0.00 RFX4 chr12:107002655 G A 0/49/0.00 0/50/0.00 8/31/0.26 0/19/0.00 RAMP2 chr17:40914468 G A 0/37/0.00 0/10/0.00 7/22/0.32 0/19/0.00 TRPA1 chr8:72951105 T C 5/9/0.56 11/12/0.92 4/7/0.57 7/9/0.78 SHANK3 chr22:51117347 C G 9/53/0.17 9/49/0.18 11/53/0.21 6/20/0.30 SLC4A8 chr12:51796889 ATTAGGTAATAGGTATAA ACTGCTCAGAACAGTGCT CAGAACATAATAAATGC A 8/29/0.28 4/13/0.31 4/9/0.44 6/14/0.43 CCR10 chr17:40831655 T G 4/19/0.21 1/14/0.07 4/19/0.21 5/16/0.31 BAI1 chr8:143614742 C T 0/87/0.00 22/73/0.30 4/50/0.08 16/43/0.37 CXCR2 chr2:219000402 A G 36/110/0.33 0/67/0.00 0/60/0.00 0/43/0.00 PBRM1 (**) chr3:52643768 G A 9/17/0.53 0/16/0.00 0/13/0.00 0/11/0.00 TDG chr12:104376913 G T 6/19/0.32 0/13/0.00 0/22/0.00 0/13/0.00 CAPN3 chr15:42652022 G A 11/29/0.38 0/24/0.00 0/38/0.00 0/17/0.00 TANC2 chr17:61489025 C T 5/16/0.31 0/37/0.00 0/22/0.00 0/16/0.00 GNA12 chr7:2771104 A G 4/100/0.04 12/56/0.21 0/60/0.00 17/69/0.25 RQCD1 chr2:219452081 A G 3/30/0.10 0/15/0.00 0/19/0.00 6/22/0.27 MYO7B chr2:128367139 C T 0/62/0.00 12/33/0.36 0/32/0.00 5/19/0.26 ARID1A(**) chr1:27101195 C T 1/97/0.01 14/55/0.25 1/52/0.02 18/41/0.44 FMN1 chr15:33192175 A T 0/14/0.00 7/21/0.33 0/20/0.00 5/16/0.31 LARP7 chr4:113568502 C T 0/57/0.00 9/46/0.20 2/49/0.04 10/38/0.26 PIK3CA (**) chr3:178936094 C A 0/43/0.00 8/35/0.23 2/20/0.10 3/20/0.15 NRG2 chr5:139267002 C A 1/103/0.01 12/55/0.22 0/32/0.00 3/37/0.08 OR2B2 chr6:27879979 A G 0/32/0.00 5/20/0.25 0/21/0.00 3/14/0.21 TP53 (**) chr17:7578403 C A 0/72/0.00 15/68/0.22 0/58/0.00 7/40/0.18 KRT14 chr17:39740149 C T 0/64/0.00 0/39/0.00 0/33/0.00 10/29/0.34 AHNAK (**) chr11:62294618 CC CAC 0/21/0.00 0/17/0.00 0/21/0.00 5/19/0.26 AHNAK (**) chr11:62294626 A T 0/25/0.00 0/17/0.00 0/23/0.00 5/18/0.28 RBM5 (**) chr3:50155882 T C 0/31/0.00 6/23/0.26 0/8/0.00 0/19/0.00 DST chr6:56473952 C T 0/24/0.00 0/21/0.00 0/8/0.00 5/15/0.33

*numbers indicate alternate (ALT) reads, read depth and minor allele frequency (MAF) of indicated mutation for each analysed region, respectively (**) gene included in the targeted sequencing panel of 42 genes used in patients 1-6.

(26)

75

Supplementary Materials and Methods

Material Transfer Agreement

A material transfer agreement (MTA) with respect to sample mobility from Indonesia to The Netherlands has been made on 13 May 2013 and has been signed by both parties; Universitas Gadjah Mada (UGM) and University Medical Center Groningen (UMCG).

DNA Isolation

The FFPE blocks were serially cut in 10-μm slides, in which the first and the last slides (3 μm in thickness) were stained with H&E and used as reference in identifying the tumor regions and their WHO/ISUP grade (Figure 1). DNA isolation was done following the protocol Adaptive Focused AcousticsTM_{-based DNA extraction of FFPE of the truXTRAC}TM _{FFPE DNA kit (Covaris, Woburn, MA,} USA). The concentration of isolated DNA was measured using the QubitTM_{dsDNA Broad Range} Assay Kit (Life Technologies, Carlsbad, CA, USA).

Targeted Sequencing (TS) using Single Primer Enrichment Technology (SPET)

For TS, 42 genes were selected based on publications focusing on the most frequent mutated genes and therapy-related genes in ccRCC [1-6]. For all genes the consensus coding regions were enriched using Single Primer Enrichment Technology (SPET), in which the amplification step is based on a single target-specific primer. Primers designed up- and downstream of a target region work independently on both strands of this region to allow variant identification from both strands with only single-end sequencing. A custom panel kit (OvationTM_{Custom Target Enrichment} System, NuGEN, San Carlos, CA, USA) was designed to cover the entire coding region of the selected genes (in total 1308 target regions), resulting in 4003 landing probes (Supplementary Table 1). The distance between landing probes was optimized for sequencing the fragmented DNA of FFPE samples (average target length 195 nt, designed peak distance = 80 nt). The library preparation steps were done according to the protocol from the manufacturer. Briefly, aliquots of 500 ng isolated DNA were sheared into 500-bp fragments using truXTRAC DNA FFPE Kits and Adaptive Focused AcousticsTM_{technology (Covaris, Woburn, MA, USA). The DNA fragments were} end-repaired, and ligated with indexed adapters, which contain an 8-base sample index and a 6-base random sequence (N6 unique molecular identifier (UMI)) to create a unique identification of each original DNA fragment. This unique identification by N6 molecular barcodes was applied to remove duplicate reads resulting from PCR amplification during data analysis. After indexed adapter ligation, the process was continued with annealing of the landing probes to the indexed DNA fragments, an extension step, and a subsequent amplification step, in which the PCR cycle was determined by the optional qPCR step. The optimal PCR cycles were between 18-21 cycles, depended on the quality of each sample pool.

Sequence Data Analysis Workflow and Variant Filtering

The WES data analysis was done using a pipeline which is implemented as a molgenis compute workflow at https://github.com/mmterpstra/molgenis-c5-TumorNormal/ tree/459417cc9553fae8c3040953970938860dafdfea [7].

(27)

76

Briefly, BCL files generated by Illumina sequencing systems, were de-multiplexed and converted into FASTQ files by bcl2fastq conversion software. The reads were adapter-trimmed, and then aligned to the human genome reference (GRCh37) using Burrows-Wheeler Aligner version 0.7.12 [8]. The quality of aligned reads was improved by INDELs realignment and base quality score filtering using GATK. Picard tool was used to convert BAM to VCF format and mark duplicates. For TS, the following changes were made to the analysis. The N6 UMI’s were added to the read descriptors. The reads were trimmed for linker sequence ‘GAGAGCGATCCTTGC’, and poly G artifacts using Trim Galore. The reads were then aligned to the human genome reference (GRCh37) using Burrows-Wheeler Aligner version 0.7.12 [8]. Using the alignment the reads were trimmed for the landing probe sequence when the end intersected with one of the landing probes (https://github.com/mmterpstra/DigitalBarcodeReadgroups). Next, Picard version 2.10.0 UmiAwareMarkDuplicatesWithMateCigar tool was used to mark PCR duplicates. This was done using the N6 barcodes, in which the hamming distance of 1 was allowed.

For both WES and TS, two different variant calling algorithms have been used; the HaplotypeCaller from GATK and FreeBayes [9, 10]. Called variants were annotated using snpeff/snpsift 3.5, with the Ensembl release 75 gene annotations, 1000 genome phase 1, dbNSFP2.7, and ExAC 0.3 databases [11-14]. The annotated variants were filtered using the following exclusion criteria: mutant allele frequency >2% in 1000 genome project phase 1 or > 0.01% in ExAC database, the possibility of error >1/100 in calling (QUAL<20), low quality by depth (QD<2 and QD/AF < 8.0), present with strand bias (FS >60 for SNVs and >200 for Indels), present in tandem repeat units (RPA>8), no alterations found in tumor compared to normal samples, putative non harmful variant e.g. synonymous variants, and variants located in non-coding regions.

(28)

77

References

Dalgliesh GL, Furge K, Greenman C et al. Systematic sequencing of renal carcinoma reveals inactivation of histone modifying genes. Nature 2010; 463: 360-363.

Duns G, Hofstra RM, Sietzema JG et al. Targeted exome sequencing in clear cell renal cell carcinoma tumors suggests aberrant chromatin regulation as a crucial step in ccRCC development. Hum Mutat 2012; 33: 1059-1062.

Sato Y, Yoshizato T, Shiraishi Y et al. Integrated molecular analysis of clear-cell renal cell carcinoma. Nat Genet 2013; 45: 860-867.

Gerlinger M, Rowan AJ, Horswell S et al. Intratumor heterogeneity and branched evolution revealed by multiregion sequencing. N Engl J Med 2012; 366: 883-892.

Pena-Llopis S, Vega-Rubin-de-Celis S, Liao A et al. BAP1 loss defines a new class of renal cell carcinoma. Nat Genet 2012; 44: 751-759.

Network CGAR. Comprehensive molecular characterization of clear cell renal cell carcinoma. Nature 2013; 499: 43-49.

Byelas H, Dijkstra M, Neerincx PB et al. Scaling Bio-Analyses from Computational Clusters to Grids. CEUR-WS.org 2013; CEUR Workshop Proceedings 993. Link:http://ceur-ws.org/Vol- 993/paper992.pdf.

Li H, Durbin R. Fast and accurate long-read alignment with Burrows-Wheeler transform. Bioinformatics 2010; 26: 589-595.

McKenna A, Hanna M, Banks E et al. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res 2010; 20: 1297-1303.

Garrison E, Gabor M. Haplotype-based variant detection from short-read sequencing. arXiv preprint arXiv 2012; 1207.3907

Cingolani P, Platts A, Wang le L et al. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome of Drosophila melanogaster strain w1118; iso-2; iso-3. Fly (Austin) 2012; 6: 80-92.

Consortium GP, Abecasis GR, Auton A et al. An integrated map of genetic variation from 1,092 human genomes. Nature 2012; 491: 56-65.

Liu X, Jian X, Boerwinkle E. dbNSFP v2.0: a database of human non-synonymous SNVs and their functional predictions and annotations. Hum Mutat 2013; 34: E2393-2402.

Lek M, Karczewski KJ, Minikel EV et al. Analysis of protein-coding genetic variation in 60,706 humans. Nature 2016; 536: 285-291. 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14.

(29)