Determinants of anti-PD1 response and resistance in clear cell renal cell carcinoma

The TRACERx Renal Consortium

The members of TRACERx Renal Consortium are: Lewis Au, Ben Challacombe, Ashish Chandra, Simon Chowdhury, William Drake, Archana Fernando, Nicos Fotiadis, Andrew Furness, Emine Hatipoglu, Karen Harrison-Phipps, Steve Hazell, Peter Hill, Catherine Horsfield, James Larkin, Jose I. Lopez, Teresa Marafioti, David Nicol, Tim O’Brien, Jonathon Olsburgh, Lisa Pickering, Alexander Polson, Sergio Quezada, Sarah Rudman, Scott Shepherd, Charles Swanton, Samra Turajlic, Mary Varia, Hema Verma.

The PEACE Consortium

The members of PEACE Consortium are: Chris Abbosh, Kai-Keen Shiu, John Bridgewater, Daniel Hochhauser, Martin Forster, Siow-Ming Lee, Tanya Ahmad, Dionysis Papadatos-Pastos, Sam Janes, Peter Van Loo, Katey Enfield, Nicholas McGranahan, Ariana Huebner, Sergio Quezada, Stephan Beck, Peter Parker, Henning Walczak, Tariq Enver, Rob Hynds, Mary Falzon, Ian Proctor, Ron Sinclair, Chi-wah Lok, Zoe Rhodes, David Moore, Teresa Marafioti, Elaine Borg, Miriam Mitchison, Reena Khiroya, Giorgia Trevisan, Peter Ellery, Mark Linch, Sebastian Brandner, Crispin Hiley, Selvaraju Veeriah, Maryam Razaq, Heather Shaw, Gert Attard, Mita Afroza Akther, Cristina Naceur-Lombardelli, Lizi Manzano, Maise Al-Bakir, Simranpreet Summan, Nnenna Kanu, Sophie Ward, Uzma Asghar, Emilia Lim, Faye Gishen, Adrian Tookman, Paddy Stone, Caroline Stirling, Lewis Au, Andrew Furness, Kim Edmonds, Nikki Hunter, Sarah Sarker, Sarah Vaughan, Mary Mangwende, Karla Lingard, Lavinia Spain, Scott Shepherd, Haixi Yan, Ben Shum, Eleanor Carlyle, Steve Hazell, Annika Fendler, Fiona Byrne, Nadia Yousaf, Sanjay Popat, Olivia Curtis, Gordon Stamp, Antonia Toncheva, Emma Nye, Aida Murra, Justine Korteweg, Nahid Sheikh, Debra Josephs, Ashish Chandra, James Spicer, Ula Mahadeva, Anna Green, Ruby Stewart, Lara-Rose Iredale, Tina Mackay, Ben Deakin, Debra Enting, Sarah Rudman, Sharmistha Ghosh, Lena Karapagniotou, Elias Pintus, Andrew Tutt, Sarah Howlett, Vasiliki Michalarea, James Brenton, Carlos Caldas, Rebecca Fitzgerald, Merche Jimenez-Linan, Elena Provenzano, Alison Cluroe, Grant Stewart, Colin Watts, Richard Gilbertson, Ultan McDermott, Simon Tavare, Emma Beddowes, Patricia Roxburgh, Andrew Biankin, Anthony Chalmers, Sioban Fraser, Karin Oien, Andrew Kidd, Kevin Blyth, Matt Krebs, Fiona Blackhall, Yvonne Summers, Caroline Dive, Richard Marais, Fabio Gomes, Mat Carter, Jo Dransfield, John Le Quesne, Dean Fennell, Jacqui Shaw, Babu Naidu, Shobhit Baijal, Bruce Tanchel, Gerald Langman, Andrew Robinson, Martin Collard, Peter Cockcroft, Charlotte Ferris, Hollie Bancroft, Amy Kerr, Gary Middleton, Joanne Webb, Salma Kadiri, Peter Colloby, Bernard Olisemeke, Rodelaine Wilson, Ian Tomlinson, Sanjay Jogai, Christian Ottensmeier, David Harrison, Massimo Loda, Adrienne Flanagan, Mairead McKenzie, Allan Hackshaw, Jonathan Ledermann, Kitty Chan, Abby Sharp, Laura Farrelly, and Hayley Bridger.

Lead Contact

Further information and requests for resources and reagents should be directed to and will be fulfilled by the Lead Contact, Samra Turajlic (samra.turajlic{at}crick.ac.uk).

Materials Availability

This⍰ study did not generate new unique reagents.

Data and code availability

Clinical studies

ADAPTeR (NCT02446860) is a single-arm, open-label, phase 2 study of nivolumab therapy as pre-operative therapy in metastatic ccRCC. Planned interim analysis took place after six months after the last patient enrolled had their first Response Evaluation Criteria in Solid Tumours (RECIST version 1.1) defined objective response assessment. ADAPTeR was initially approved by NRES Committee London Fulham on 01/12/2014. ADAPTeR is performed in accordance with the ethical principles in the Declaration of Helsinki, Good Clinical Practice and applicable regulatory requirements.

Nivolumab was administered at a dose of 3mg per kilogram of body weight as a 60-minute intravenous infusion every 2-weeks. Eligible patients were 18 years of age or older, had histologic confirmation of advanced or metastatic renal-cell carcinoma (RCC) with predominantly clear cell component with at least one site of disease outside the kidney measurable according to the RECIST version 1.1, with no prior systemic therapy for ccRCC. All patients had an Eastern Cooperative Oncology Group (ECOG) performance status of 0 or 1. Key exclusion criteria were need for immediate nephrectomy, any active, known or suspected autoimmune disease or another condition requiring systemic treatment with either corticosteroids (>10mg daily prednisolone equivalent) or other immunosuppressive medications within 14-days of study drug administration (excluding vitiligo, Type 1 diabetes mellitus, residual hypothyroidism due to autoimmune condition only requiring hormone replacement, psoriasis not requiring systemic treatment or conditions not expected to recur in the absence of an external trigger). During the course of the study, inclusion expanded to those who have had a prior nephrectomy but are suitable for on treatment biopsies. The prognostic factors assessed for the risk categorisation are as per the published IMDC criteria³⁷: time to systemic therapy (<1 year), performance status, anaemia, hypercalcaemia, neutrophilia and thrombocytosis. Presence of zero (favourable-risk), one (intermediate-risk), and two or three (poor-risk) factors provides the categorisation.

The primary endpoint was the safety profile of nivolumab given pre- and post-operatively to patients with metastatic ccRCC undergoing nephrectomy. Secondary endpoints were overall response rate (ORR), progression free survival (PFS), and overall survival (OS). Exploratory endpoints pertain to biomarker analyses. Patients deemed clinically suitable for nephrectomy at baseline were scheduled for surgery after the fourth cycle of treatment. Patients not deemed clinically suitable for nephrectomy at baseline would undergo surgery if an excellent clinical response is observed and if surgery was clinically appropriate. Nivolumab treatment was recommenced post-operatively upon sufficient recovery, and until disease progression. Patients who remained clinically unsuitable for nephrectomy continued nivolumab treatment until disease progression.

For translational study sample collection, baseline tumour biopsy via appropriate guidance (ultrasound or computer tomography [CT]) at least 3 days and up to 14 days prior to starting nivolumab was obtained. Tumour multiple regions of nephrectomy specimen were sampled, as well as image guided biopsy of regressing lesions or at disease progression either at site of progression or, if not possible, percutaneous primary renal tumour biopsy, prior to commencement of any subsequent treatment. Blood samples were collected at each tumour sampling timepoint.

Autopsy samples from ADR001, ADR005, and ADR015 were obtained through the PEACE Study (NIHR 18422; NCT03004755), where samples are harvested within ∼48 hours from death. All patients were co-recruited to the TRACERx Renal study (NCT03226886; see secondary author list for the full list of TRACERx Renal consortium investigators). Patient and sample metadata (i.e. age a diagnosis, sex, clinical response, biopsy site) are provided as Table S1 and Supplemental Data Table 1. All the patients provided written informed consent. The protocols, amendments and informed consent forms were approved by the institutional review board or independent ethics committee at each trial site for each trial.

Sample collection

Tumour and normal tissue were collected via image-guided percutaneous biopsies, ex vivo sampling at nephrectomy, and at autopsy. Multiregion samples were obtained with all modalities. For samples obtained at nephrectomy, resected specimens were reviewed macroscopically by a pathologist to guide multiregion sampling for this study and to avoid compromising diagnostic requirements. Spatially separated regions sampled from the ‘‘tumour slice’’ using a 6mm punch biopsy needle. The punch was changed between samples to avoid contamination. The total number of samples obtained reflects the tumour size with a minimum of 3 biopsies that are non-overlapping and equally spaced. Areas which are obviously fibrotic or haemorrhagic are avoided during sampling and every attempt is made to reflect macroscopically heterogeneous tumour areas. Primary tumour regions are labelled as R1, R2, R3.Rn and locations are recorded. Normal kidney tissue was sampled from areas distant to the primary tumour and labelled N1. For all samples collected, each were split into two for snap freezing and formalin fixing respectively, such that the fresh frozen sample has its mirror image in the formalin-fixed sample which is subsequently paraffin embedded. Fresh samples were placed in a 1.8 ml cryotube and immediately snap frozen in liquid nitrogen for >30 seconds and transferred to −80 C for storage. Peripheral blood was collected at the time of surgery and processed to separate buffy coat and peripheral blood mononuclear cells (PBMCs).

Nucleic acid extraction, DNA and RNA library preparation and sequencing

DNA and RNA were co-extracted from fresh-frozen tumour tissue using AllPrep DNA/RNA mini kit (Qiagen). RNA from peripheral blood mononuclear cells (PBMC) were extracted from blood stored in Tempus tubes using the Tempus™ Spin RNA Isolation Kit (Invitrogen). Germline DNA was isolated from whole blood using the DNeasy Blood and Tissue kit (Qiagen). DNA yield and quality were assessed on TapeStation4200 (Agilent) and Qubit Fluorometric quantification (ThermoFisher Scientific). Samples were normalised to either 3 ug or 200ng and sheared to 150-200bp using a Covaris-E220 or LE220-plus. Agilent SureSelectXT enriched libraries were constructed following the manufacturer’s manual or automated (using the Agilent Bravo liquid handling platform) SureSelectXT Target Enrichment System for Illumina Paired-end Multiplexed Sequencing Library protocol. Hybridisation and capture were performed using the Agilent SureSelectXT Human All Exon v5 capture library. Final libraries were sequenced to a target coverage of 250x with 101bp paired-end reads multiplexed on the Illumina HiSeq4000 sequencing platform. The extracted RNA was normalised to 100ng for library construction using RNA-Ribozero (ribodeplete) Library Preparation Kits. The prepared libraries were multiplexed and QC’ed before paired-end sequencing with target coverage of 50 million reads per sample on HiSeq4000 sequencing platforms (Illumina). RNA was extracted from blood for TCR sequencing from the following cases and timepoints: all cases (n=15) pre- and post-treatment.

SNV, and INDEL calling from multiregion WE sequencing

Paired-end reads (2×100bp) in FastQ format sequenced by Hiseq were aligned to the reference human genome (build hg19), using the Burrows-Wheeler Aligner (BWA) v0.7.15. with seed recurrences (-c flag) set to 10000 ⁹⁰. Intermediate processing of Sam/Bam files was performed using Samtools v1.3.1 and deduplication was performed using Picard 1.81 (http://broadinstitute.github.io/picard/).Single Nucleotide Variant (SNV) calling was performed using Mutect v1.1.7 and small scale insertion-and-deletions (INDELs) were called running VarScan v2.4.1 in somatic mode with a minimum variant frequency (--min-var-freq) of 0.005, a tumour purity estimate (--tumour-purity) of 0.75 and then validated using Scalpel v0.5.3 (scalpel-discovery in --somatic mode) (intersection between two callers taken) ^91-93. SNVs called by Mutect were further filtered using the following criteria: i) ≤5 alternative reads supporting the variant and variant allele frequency (VAF) of ≤1% in the corresponding germline sample, ii) variants falling into mitochondrial chromosome, haplotype chromosome, HLA genes or any intergenic region were not considered, iii) presence of both forward and reverse strand reads supporting the variant, iv) >5 reads supporting the variant in at least one sample, v) variants were required to have a VAF of 0.01 in at least one sample, vi) sequencing depth need to be ≥20 and ≤3000 across all samples. Dinucleotide substitutions (DNV) were identified when two adjacent SNVs were called and their VAFs were consistently balanced (based on proportion test, P≥0.05). In such cases the start and stop positions were corrected to represent a DNV and frequency related values were recalculated to represent the mean of the SNVs. Variants were annotated using Annovar⁹⁴. Individual tumour biopsy regions were judged to have failed quality control and excluded from analysis based on the following criteria: i) sequencing coverage depth below 100X, ii) low tumour purity such that copy number calling failed. Driver variants are manually reviewed and predicted for variant effect, and variant annotations on the heatmap are only for confident driver events.

Methylation specific PCR

Methylation of the VHL promoter was detected after bisulphite treatment of 500ng of patient DNA using the EZ DNA Methylation-Direct kit (Zymo Research). Bisulphite treated DNA was amplified in the PCR using methylation specific oligonucleotides followed by Big Dye terminator Sanger sequencing. Methylation was confirmed by comparing and contrasting patient tumour and normal renal tissue for methylation protected CpG sequences.

Neoantigen calling

Neoantigen predictions were derived by first determining the 4-digit HLA type for each patient, along with mutations in class I HLA genes, using POLYSOLVER⁹⁵. Next, all possible 9, 10 and 11-mer mutant peptides were computed, based on the detected somatic non-synonymous SNV and INDEL mutations in each sample. Binding affinities of mutant and corresponding wildtype peptides, relevant to the corresponding POLYSOLVER-inferred HLA alleles, were predicted using NetMHCpan (v3.0)⁹⁶ and NetMHC (v4.0)⁹⁷. Neoantigen binders were defined as strong binders if their %rank was below <0.5 for the mutant and >0.5 for the wildtype protein.

TMB, fsINDEL burden, neonatigen burden, wGII

To calculate sample-level TMB, fsINDEL burden and neoantigen burden, variants were filtered with the same criteria as described above for multi-region WE sequencing but on a per sample-basis. Variants were restricted to positions falling inside the targeted capture range (±50bp padding). Tumour mutational burden (TMB) was calculated as the number of exonic non-synonymous SNVs per mega base. The frameshift INDEL (fsINDEL) burden was calculated as the total number of exonic frameshift INDELs per sample. The neoantigen burden was calculated as the total number of predicted strong binders per sample. The average proportion of the genome with aberrant copy number, weighted on each of the 22 autosomal chromosomes, was estimated as the weighted genome instability index (wGII).

SNP calling

Single nucleotide polymorphisms (SNPs) were called in the germline sample using Platypus v0.8.1 with default parameters apart from --genIndels=0 and --minMapQual=40. Tumour regions were genotyped based on the variants identified in the germline (parameters set to --minPosterior=0 -- getVariantsFromBAMs=0). SNPs with a minimum coverage of 50x in the germline and the tumour sample were used for allele-specific copy number segmentation.

Copy number analysis

CNVkit v0.7.3 was used with default parameters on paired tumour-normal sequencing data⁹⁸. Outliers of the derived log2-ratio (logR) calls from CNVkit were detected and modified using Median Absolute Deviation Winsorization before case-specific joint segmentation of fresh-frozen samples to identify genomic segments of constant logR⁹⁹. Formalin-fixed and paraffin-embedded (FFPE) samples were segmented separately while leveraging the segment information from the fresh-frozen samples. Copy number alterations were called as losses or gains relative to overall sample wide estimated ploidy. Driver copy number was identified by overlapping the called somatic copy number segments with putative driver copy number regions identified by Beroukhim et al.¹⁰⁰. Allele-specific segmentation was performed using the paired PSCBS method after removal of single-locus outliers (R package PSCBS v0.61.0 ¹⁰¹).

Purity and ploidy estimate

Tumour sample purity, average ploidy and absolute allelic copy number per segment were estimated using ABSOLUTE v1.2 in allelic mode¹⁰². In line with recommended best practice all ABSOLUTE solutions were reviewed by 3 researchers, with solutions selected based on majority vote. Purity assigned 0.1 for samples below ABSOLUTE estimate thresholds for comparison analysis of samples between responders and non-responders.

Subclonal deconstruction

To estimate the CCF of a mutation, we used the following formula: Where VAF is the variant allele frequency of the mutation, p the estimated tumour purity, CN_mut the number of copies carrying the mutation and CN_t the local copy number in the tumour cells. CN_n is the local copy number in the non-tumour proportion of the sample which was assumed to be 2. The CN_mut and CCF were estimated through iteration of all possible combinations of CCF (range 0.01 to 1, by 0.01) and CN_mut (range 1 to CN_t, by 1) using the formula above to identify the best fit CCF.

Genomic contraction, persistence and expansion analysis

For each patient with matched pre- and post-treatment WES data (N=8 patients), the pre-treatment CCF was compared to the CCF in the post-treatment samples. In patients with multiple pre-treatment samples, the median pre-treatment CCF was used. In each post-treatment sample, all nonsynonymous SNVs and fsINDELs were assigned to one of the following categories: “Genomic contraction”, “Genomic persistence”, “Genomic expansion”³⁹. A mutation was defined to have undergone genomic contraction if the CCF decreased by ≥10% from pre- to post-treatment or if the mutation was present in the pre-treatment but not the post-treatment sample; genomic expansion if the CCF increased by ≥10% from pre- to post-treatment or if the mutation was present in the post-treatment but not the pre-treatment sample; genomic persistence if the CCF in the post-treatment sample was within the range of ±10% of the pre-treatment CCF. The proportion of mutations falling into each category was calculated for nonsynonymous SNVs and fsINDELs and repeated with only neoantigen encoding mutations. An enrichment test (Fisher’s exact test) was performed to determine whether mutations which encode neoantigens were more likely to undergo genomic contraction than the remaining nonsynonymous SNVs and fsINDELs.

Mutational signature analysis

Mutational signatures were estimated using the deconstructSigs package in R¹⁰³. Sample specific mutational signature analysis was restricted to samples with at least 50 mutations.

Analysis for mismatch repair deficiency

Analysis for mutations in the following nominated genes was performed: POLD3, MLH3, MSH6, RPA4, LIG1, MLH1, MSH2, MSH3, PCNA, PMS2, POLD1, POLD2, POLD4, RFC1, RFC2, RFC3, RFC4, RFC5, RPA1, RPA2, RPA3, SSBP1, EXO1.

Analysis for mutations associated with defective antigen presentation

Analysis for mutations in the following nominated genes was performed: B2M, CIITA, IRF1, PSME1, PSME2, PSME3, ERAP1, ERAP2, HSPA, PSMA7, HSPC, HSPBP1, TAP1, TAP2, TAPBP, CALR, CNX, CANX, PDIA3.

Detection of B2M mutations by Sanger sequencing

Validation of the B2M mutation was performed using PCR followed by Big Dye Terminator Sanger sequencing on the ABI 3700. 20ng of patient DNA was amplified for exon 1 of B2M, to enable detection of B2M:c.42_45delTCTT:p.S14fs. PCR conditions involved 35 cycles of denaturation at 950C, followed by oligonucleotide primer annealing at 55oC and sequence extension at 720C using Qiagen Taq polymerase and reagents. Oligonucleotide sequences used are: Forward: aacgggaaagtccctctctc; Reverse: agatccagccctggactagc.

Bulk RNAseq data processing

RNAseq data were mapped to the hg19 reference human genome using the STAR¹⁰⁴ algorithm, and transcript and gene abundance were estimated by RSEM¹⁰⁵ with default parameters. Samples were excluded if they had less than 15,000 genes detected.

Human endogenous retrovirus (HERV) analysis

Expression of previously annotated HERVs^{28, 50, 51} was analysed. HERV loci used in these three studies^{28, 50, 51} were taken from Mayer et al. ⁵³ and Vargiu et al.⁵⁴ with 66 and 3173 loci respectively. BLASTn was used to match example sequences from HERVs in Mayer et al. to GRCh38, chromosome coordinates with the greatest homology over the greatest length were taken as the best match. The Lift Genome Annotations tool from UCSC (https://genome.ucsc.edu/cgi-bin/hgLiftOver) was used to convert annotated GRCh37 HERV loci coordinates from Vargui et al. to GRCh38 coordinates. Comparing the new coordinates, 47 of the 66 HERVs from Mayer et al. were present in the list of 3173. Coordinates of all the unique elements were then compared to a custom repeat region annotation previously built using the Dfam 2.0 library (v150923) for GRCh38⁵⁷. For this custom annotation, different regions of the same provirus (e.g. the LTR and internal genes) were annotated separately, these regions were merged to allow accurate quantitation of reads from the same provirus⁵⁷. LTR-containing repeat regions from the custom annotation had to begin, end, or be fully contained within previously annotated loci to be considered a match, a buffer of 5 bases either end of the locus was included. Previously annotated HERV loci from Mayer et al. and Vargiu et al. were found to overlap multiple repeat regions per locus in our custom annotations, or were found to overlap no repeat regions at all. Some loci also overlapped other endogenous retroelement types such as LINEs and SINEs, as well as overlapping canonical gene exons. For this analysis, only expression of matching LTR-containing elements was considered rather than expression of all repeats and genes overlapping previously annotated loci. Expression was measured using read counts calculated by the featureCounts function from the Subread package¹⁰⁶ (with parameters -p -C -B -f -T 2), multi-mapping reads were not counted. Analysis for purified immune cell subset expression were performed on publicly available datasets from Linsley et al. ¹⁰⁷ (accession no. GSE60424 (GEO)) and Kazachenka et al.¹⁰⁸ (E-MTAB-8208 (EMBL-EBI)). LTR-overlapping transcripts expressed highly specifically in ccRCC were previously described⁵⁷. These transcripts were identified through de novo transcriptome assembly and their expression quantified in by transcript per million calculations, as previously described⁵⁷.

Differential gene expression analysis, pathway analysis and gene set enrichment

DESeq2¹⁰⁹ was used for differential expression analysis, using the binomial Wald test after estimation of size factors and estimation of dispersion. To identify genes differentially expressed between responders and non-responders, we considered only transcripts with normalized count number >5 in at least 5 patients. Pathway analysis was performed using the R package XGR¹⁰⁸ using the gene ontology biological process (GOBP) databases. Induced and suppressed transcripts were analysed separately against the background of all tested transcripts. The “lea” ontology algorithm was used.

T-cell subset gene signature

Gene signature or single gene enrichment was evaluated using RSEM abundance, z score scaled across all samples for which RNA-Seq was available. Signature analysis was performed using 22 immune-related signatures listed below: i) the Danaher immune score is a 60-marker gene signature derived from pan-cancer RNAseq analysis for 14 immune cell populations, where marker genes have been benchmarked against histological tumour-infiltrating lymphocyte (TIL) estimates and flow cytometry data^{58, 111}; ii) IMmotion150 ²⁰; iii) Javelin101²¹.

1. Danaher Tcells : CD3D, CD3E, CD3G, CD6, SH2D1A, TRAT1

2. Danaher CD8 : CD8A, CD8B

3. Danaher Cytotoxic : CTSW, GNLY, GZMA, GZMB, GZMH, KLRB1, KLRD1, KLRK1, PRF1, NKG7

4. Danaher Bcells : BLK, CD19, MS4A1, TNFRSF17, FCRL2, KIAA0125, PNOC, SPIB, TCL1A

5. Danaher NKcells : NCR1, XCL2, XCL1

6. Danaher CD45 : PTPRC

7. Danaher DC : CCL13, CD209, HSD11B1

8. Danaher CD8Ex : CD244, EOMES, LAG3, PTGER4

9. Danaher Mac : CD163, CD68, CD84, MS4A4A

10. Danaher Mast : MS4A2,TPSAB1,CPA3,HDC,TPSB2

11. Danaher Neut : CSF3R, S100A12, CEACAM3, FCAR, FCGR3B, FPR1, SIGLEC5

12. Danaher NKCD56 : IL21R, KIR2DL3, KIR3DL1, KIR3DL2

13. Danaher Th1 : TBX21

14. Danaher Treg : FOXP3

15. IMmotion150 Angio : VEGFA, KDR, ESM1, PECAM1, ANGPTL4, CD34

16. IMmotion150 Teff : CD8A, IFNG, PRF1, EOMES, CD274

17. IMmotion150 Myeloid : CXCL1, CXCL2, CXCL3, CXCL8, IL6, PTGS2

18. Javelin101 TCR : CD3G, CD3E, CD8B, THEMIS, TRAT1, GRAP2, CD247

19. Javelin101 TCell : CD2, CD96, PRF1, CD6, IL7R, ITK, GPR18, EOMES, SIT1, NLRC3

20. Javelin101 NK : CD2, CD96, PRF1, CD244, KLRD1, SH2D1A

21. Javelin101 chemo : CCL5, XCL2

22. Javelin101 other : CST7, GFI1, KCNA3, PSTPIP1

The signature score was calculated as the arithmetic mean of z score scaled expression of all genes in that signature for each sample.

TCR sequencing

TCR β-chain sequencing was performed by utilizing whole RNA extracted from tissue samples or from cryopreserved PBMC samples, by using a quantitative experimental and computational TCR sequencing pipeline described previously^112-115. An important feature of this protocol is the incorporation of a UMI attached to each cDNA TCR molecule that enables correction for PCR and sequencing errors, which allows higher quantitative precision compared to alternate protocols in the TCR sequences retrieved^{113, 116}. The suite of tools used for TCR identification, error correction and CDR3 extraction is freely available at https://github.com/innate2adaptive/Decombinator.

For each TCR, we computed the abundance as the count of UMIs mapping to this TCR divided by the total number of UMIs in the sample. If several samples were available at a given patient-timepoint pair, the resulting abundance was calculated as the sum of counts for this TCR across the available samples divided by the sum of total counts across these samples.

Repertoire similarity measure

The similarity between two TCR repertoires was assessed with the normalised dot product (also known as the cosine similarity) between the vectors of TCR abundance. This measure is a well-established metric widely used in machine learning to compare numerical vectors and gives a value between 0 (no similarity, that is, orthogonal vectors) and 1 (complete similarity, from vectors with an identical magnitude and direction in the feature space). Each pair of repertoires is represented as two vectors of equal length, indexed by the union of TCRs found in both repertoires and containing the number of times each TCR is detected in each of the two repertoires (each position contains an integer ≥0). The similarity between the two vectors is given as where and are the abundance vectors, represents the vector product and paired vertical bars represent the Euclidean norm of the vector.

For longitudinal similarity (Figures 4, S14, and S15), the similarity measure was performed on the TCR abundance vectors derived from (patient, timepoint) pairs.

For spatial similarity (Figure S13), the similarity measure was performed on the TCR abundance vectors derived from each sample within a (patient, timepoint) pair. For this analysis, samples from different timepoints were not compared.

Repertoire clonality index

The clonality index was estimated for each sample by using the command entropy from the entropy R package, on the basis of the observed frequency of the TCRs in that sample where p_i is the frequency of the ith TCR in the repertoire and is the number of TCRs in that repertoire.

Classification of expanded, contracted and persistent TCRs

The difference in abundance between Pre-treatment and On-treatment was calculated with the poisson.test function in R, as the data were counts. TCRs with P values above 0.01 were labelled as persistent.

Classification of expanded TCRs

We counted the number of TCRs detected with frequencies above a range of frequency thresholds in the tumour repertoires. To measure how such defined expanded TCRs were representative of the shape of the TCR distribution captured by the clonality score, we computed the prevalence of the expanded population amongst the entire repertoire, for each threshold. To do so, we took the sum of counts for expanded TCRs and divided it by the sum of all counts in the sample. The proportion obtained was then correlated to the matched clonality score with the Spearman’s rank correlation.

To focus on the most expanded TCRs (Figures 4E-4F and S15-S16), we examined those present above a threshold frequency of 2/1,000 (corresponding to the top 1% of the empirical TCR frequency distribution). At this threshold, which we already described in previously published work ⁶¹, the correlation between clonality and proportion of repertoire occupied by expanded TCRs is very strong and the number of TCRs labelled as expanded is greater than for higher thresholds for which this correlation is also significant, which enables to keep the greatest amount of data whilst still applying a stringent filtering step.

CDR3 amino acid clustering

The pairwise similarity between pairs of TCRs was measured on the basis of amino acid triplet sharing. Sharing was quantified using the normalized string kernel function stringdot (with parameters stringdot (type = ‘spectrum’, length = 3, normalized = TRUE) from the Kernlab package. The kernel is calculated as the number of amino acid triplets (sets of three consecutive amino acids) shared by two CDR3s, normalized by the number of triplets in each CDR3 being compared. The TCR similarity matrix was converted into a network diagram by using the iGraph package in R. Two TCRs were considered connected if the similarity index was >0.82 (threshold previously optimised in a separate study).

Per (patient, timepoint) pair, we counted the number of clusters containing an expanded CDR3. To normalize the counts of clusters obtained (Nreal)for the input size, for each sample, we randomly selected, outside of the real clustering structure, a number of CDR3s equal to the number of expanded CDR3s in that sample and looked for clusters around those. This control step was repeated 10 times for each (patient, timepoint) pair and we computed the average number of clusters obtained for those control (Ncon) and used Nreal/ Ncon as the normalised cluster count value (Figure 4F, Figure S15F).

For Figure S15D, as depicted in Figure S15E, we used the clustering structure built as described above for Pre-treatment samples and retrospectively labelled expanded clones at that time-point as maintained if they were also expanded post-treatment or as replaced if they were not. By doing so, we could derive the number of Pre-treatment clusters containing maintained (resp. replaced) expanded clones which was then divided by the initial count of maintained (resp. replaced) expanded clones present in that sample to obtain the proportion displayed in FigureS15D.

Frequency ratio

We wanted to capture the rate of clonal replacement that occurs in the tumour repertoires. To do so, for each expanded TCR at baseline that could also be detected after treatment, we computed the ratio of the observed frequency at baseline divided by the observed frequency after-treatment. To derive a metric for each patient, we computed the average of ratio scores obtained for all expanded TCRs at baseline (those that could not be detected after treatment were excluded).

Multiplex Immunofluorescence Staining and Image Analysis

FFPE blocks were cut in 2 micron thick slides. The slides were baked for 60 minutes and stained using the antibodies listed below and opal fluorophores. Leica Bond III machine was used for the immunofluorescence staining. Images of the stained slides were acquired by using the Vectra 3 automated quantitative pathology imaging system (Akoya Biosciences). Matching H&E image of each slide was reviewed by a pathologist and areas to annotate on the immunofluorescent images for analysis were identified. Necrotic and stromal areas as well as non-tumour areas were excluded and tumour areas were scored. Slides for patient ADR009 were not evaluable due to necrosis. Total of 61 samples (41 pre-treatment and 20 on treatment samples) for the first mIF panel and 60 samples (40 pre-treatment and 20 on treatment samples) were for the second IF panel were used for analysis. The following antibodies were used for mIF staining: CD3 (Mouse monoclonal, LN10, 1:100 dilution on Opal 520 in 1:50 dilution), CD4 (Mouse monoclonal, 4B12, 1:50 dilution on Opal 540 in 1:100 dilution), CD8 (Mouse monoclonal, 4B11, 1:100 dilution on Opal 540 in 1:150 dilution and on Opal 620 in 1:150 dilution), FoxP3 (Mouse Monoclonal, 236A/E7, 1:80 dilution on Opal 570 in 1:150 dilution), CD163 (Mouse monoclonal, 10D6, 1:100 dilution on Opal 690 in 1:50 dilution), Granzyme B (Mouse monoclonal, 11F1, 1:80 dilution, on Opal 620 in 1:150 dilution)

Up to 25 multispectral images (MSI) were acquired per slide depending on the size of the tumour to include all representative areas of the tumour. Representative MSIs from different slides were used while training the algorithms for each marker. Scoring of each slide was performed using the inForm software on Vectra. The quality and accuracy of the scoring was checked by two clinicians one of whom was a histopathologist. MSIs with poor tissue quality were excluded from the analysis. Merged data obtained by using the inForm software was analysed using the phenoptrReports tool (Akoya Biosciences) on R. T-cells subsets (CD8⁺, CD4⁺ effectors, Tregs and CD8⁺CD4⁺ double positive cells) were scored both out of total cells counted on each slide and out of the total T cells counted. CD163 cells were scored out of total cells counted per slide. Overall granzyme expression was scored in relation to the total T-cell and CD163⁺ cell count. Granzyme B expression on CD8⁺ cells was scored out of the total CD8⁺ cells. Median scoring value was used for each patient per time point and two-sided Mann-Whitney U test was used for statistical analysis of the data.

The mIF and mIHC antibody panels were designed to evaluate T-cell subsets, B-cells, myeloid cells, and GZMB expression. This was conducted given 1) double positive (CD8⁺CD4⁺) T-cells with high degrees of TCR clonality have previously been described in ccRCC¹¹⁷; 2) myeloid inflammation has been associated with blunting of anti-tumour T-cell activity in metastatic ccRCC²⁰; and 3) high tumour infiltration with B-cells and plasma cells have previously been shown to correlate with favourable clinical outcomes across cancer types ^118-120.

Immunohistochemistry

Formalin-fixed paraffin-embedded tissue sections of clear-cell renal cell carcinoma and normal tonsil tissues were subjected to haematoxylin and eosin and multiplex immunostaining. The primary antibodies used for multiplex immunolabeling are as follows: CD19 (Rabbit monoclonal, SP291, 1:10000 dilution), CD138 (Mouse monoclonal, MI15, 1:100 dilution), PD-1 (Mouse Monoclonal, NAT105/E3, 1:2 dilution).

To establish optimal staining conditions each antibody was tested and optimized on 2-4 um cut tissue sections of human reactive tonsil and normal kidney by applying conventional single immunohistochemistry. In brief sections were de-waxed and re-hydrated prior to the multiplex immunolabeling whose procedure was adapted and performed according to the established protocol described elsewhere¹²¹. Total of 59 samples (40 pre-treatment and 19 on treatment samples) for the mIHC panel.

Staining assessment and data handling

Specificity of the staining was assessed by a haematopathologist with expertise in multiplex-immunostaining. Scanned slide images were obtained with the use of NanoZoomer Digital Pathology System (Hamamatsu, Japan). Total of 60 samples (41 pre-treatment and 19 on treatment samples) were used for analysis.

Flow cytometry

Renal tumour resections and normal tissue were cut into small pieces (2-3mm) by using sterile disposable scalpel plus forceps in RPMI (Sigma-Aldrich) with Collagenase I (Sigma-Aldrich) (for ADR013 tumour and normal tissue), Liberase (for ADR001 tumour tissue) and DNAse I (Roche) and was digested for 1 hour at room temperature using the gentleMACS dissociator (Miltenyi Biotec). The digest was passed through a 70-µm cell strainer by using 5-10 ml of RPMI containing 2% fetal bovine serum (FBS) to obtain a single cell suspension. Lymphocytes were obtained from the single cell suspension by using Ficoll Paque Plus (GE Healthcare) density gradient centrifugation (750g for 10 min). Isolated lymphocytes were washed with RPMI and 2%FBS and cryopreserved in 90% FBS with 10% dimethyl sulfoxide (Sigma–Aldrich). PBMCs were isolated from blood samples collected in Vacutainer EDTA blood collection tubes (BD) using Ficoll Paque Plus (GE Healthcare) density gradient centrifugation and cryopreserved in in 90% FBS with 10% dimethyl sulfoxide (Sigma–Aldrich).

Thawed lymphocytes were washed with 1x phosphate-buffered saline (PBS) and were stained with the antibodies listed below. Antibody mastermixes were prepared in Brilliant Staining Buffer (BD). eBioscience™ Foxp3 / Transcription Factor Staining Buffer Set was used for the intracellular staining. Samples were stained using the following antibodies: CD8 (RPA-T8, BUV496), CD45RA (HI100, BUV563), CD4 (SK3, BUV615), CD38 (HIT2, BUV737), CD3 (SK7, BUV705), FoxP3 (206D, BV421), CCR4 (L291H4, BV510), Viability dye (Yellow Fluorescent reactive dye, BV570), CD57 (QA17A04, BV605), Ki67(B56, BV650), CD39 (TU66, BV711), CCR7 (G043H7, BV750), CD69(FN50, BV785), CD103 (Ber-ACT8, BB515), CXR5 (RF8B2, PerCp-Cy5-5), TCF-7 (7F11A10, PE), Granzyme B (GB11, PE-CF594), CD25 (M-A251, PE-Cy5), PD-1 (EH12.2H7, PE-Cy7), TOX (REA473, APC), HLA-DR (LN3, AF700), IgG4(Biotin), 4-1BB (4B4-1, BUV661), TIM3 (7D3, BV650), KLRG1 (13F12F2, Superbright702), CD27 (O323, BV750), ICOS (C398.4A, BV785), EOMES (WD1928, Percp-eFlour710), CTLA-4 (L3D10, APC), GITR (108-17, APC/Fire-750), Streptavidin (BUV395). The samples were acquired on the BD Symphony flow cytometer. Data was analysed using the FlowJo (version 10).

PD-1 competition binding assay to evaluate anti-PD1 monoclonal antibody binding

PBMC isolated from healthy individuals were activated in vitro using plate coated anti-CD3 and soluble anti-CD28 with 100IU IL-2 per well. 50ul (5ug/ mL solution) anti-CD3 was used to coat wells of a 96 well plate which was kept at 4°C overnight. Two washes using 200ul of PBS were performed to remove unbound antibodies the next day. Subsequently, 2 × 10⁵ PBMC were added into each well with subsequent addition of soluble anti-CD28 (2ug/ mL). The plate was placed into a humidified 37°C incubator for 72 hours. Following this period, the wells containing activated PBMC were either incubated with 50ul (2.5mg) pembrolizumab or PBS control for 30 minutes. PBS washes were used to remove unbound therapeutic antibodies. Flow cytometric staining of CD3, PD-1 and anti-IgG4 was performed thereafter.

Single-cell RNA/TCR Sequencing

Tumour infiltrating lymphocytes from ADR001 and AD013 were stained with CD3 (PE, SK7 clone), IgG4 (Biotinylated) and Streptavidin (BV650) antibodies for flow cytometry. Stained cells were FACS sorted as CD3⁺IgG4^- (40,000 cells) and CD3⁺IgG4⁺ (20,000 cells) for ADR001 and CD3⁺IgG4^- (50,000 cells) and CD3⁺IgG4⁺ (90,000 cells) for ADR013. FACS sorted cells were single cell sorted using the 10X Genomic machine. The sorted cells were processed using the 10X Genomic Chromium Next GEM Single Cell 5’ Reagents Kit V2 (dual index) for 5’gene expression library construction and V(D)J library construction. The samples were sequenced on the NextSeq using the High Output Kit v2.5 (150 Cycles).

FASTQ files containing gene expression (GEX) and VDJ were demultiplexed using cellranger mkfastq (10x Genomics). GEX reads were aligned to GRCh38 and counted using cellranger count, VDJ reads were aligned to cellranger’s GRCh38 VDJ reference dataset using cellranger vdj. Expression matrices were analysed using the Seurat package¹²². To remove technical variation in the data, TCR, ribosomal and heat-shock protein genes were removed from the analysis, also cells with mitochondrial reads making up >10% total read content were removed. 8382 CD3⁺IgG4^- and 10083 CD3⁺IgG4⁺ cells in ADR013; and 4648 CD3⁺IgG4^- and 3343 CD3⁺IgG4⁺ cells in ADR001 were retained after quality control filtering. Datasets were integrated using SCTransform integration¹²³ using the recommended parameters and regressing the % mitochondrial read content. PCA and UMAP dimensional reduction (dims = 1:30) and clustering (res = 0.3) was then performed using RunPCA, RunUMAP, FindNeighbours, and FindClusters. Publicly available gene signatures for T cell states were obtained from the following publications: Schietinger et al.¹²⁴, Thommen et al.¹²⁵, Guo et al.¹²⁶, Li et al. ¹²⁷, Yost et al.¹²⁸, Miller et al.⁶⁹, Zhou et al.¹²⁹, and Litchfield et al.¹³⁰ (Supplemental Data Table 3). The proportion of reads mapping to the genes in each signature for each cell was then calculated using PercentageFeatureSet. All differential gene expression analysis were carried out on log normalised gene expression values (using NormalizeData, default parameters) using the MAST algorithm¹³¹ within FindMarkers. GOBP analysis was carried out using the XGR package¹¹⁰ using the “lea” algorithm. scTCR data was analysed using scRepertoire¹³². Cells were considered of the same clone if they contained a matching TRB sequence and CDR3 gene.