Proteome and phosphoproteome signatures of recurrence for HPV+ head and neck squamous cell carcinoma.

Background: Head and neck squamous cell carcinoma (HNSCC) is the sixth most common cancer worldwide and the human papillomavirus (HPV+)-driven subtype is the fastest rising cancer in North America. Although most cases of HPV+ HNSCC respond favorably to the treatment via surgery followed by radiochemotherapy, up to 20% recur with a poor prognosis. The molecular and cellular mechanisms of recurrence are not fully understood.
Methods: To gain insights into the mechanisms of recurrence and to inform patient stratification and personalized treatment, we compared the proteome and phosphoproteome of recurrent and non-recurrent tumors by quantitative mass spectrometry.
Results: We observe significant differences between the recurrent and non-recurrent tumors in cellular composition, function, and signaling. The recurrent tumors are characterized by a pro-fibrotic and immunosuppressive tumor microenvironment (TME) featuring markedly more abundant cancer-associated fibroblasts, extracellular matrix (ECM), neutrophils, and suppressive myeloid cells. Defective T cell function and increased epithelial-mesenchymal transition potential are also associated with recurrence. These cellular changes in the TME are accompanied by reprogramming of the kinome and the signaling networks that regulate the ECM, cytoskeletal reorganization, cell adhesion, neutrophil function, and coagulation. Conclusions: In addition to providing systems-level insights into the molecular basis of recurrence, our work identifies numerous mechanism-based, candidate biomarkers and therapeutic targets that may aid future endeavors to develop prognostic biomarkers and precision-targeted treatment for recurrent HPV+ HNSCC.

Introduction

HPV infection causes approximately 25% of all HNSCC cases and the rate has been rising in recent years[1]. HPV+ HNSCC is biologically and clinically distinct from non-HPV driven (HPV-) HNSCC, which is typically associated with tobacco and alcohol consumption[2-5]. Despite the availability of an HPV vaccine, modelling using current vaccination rates suggests a continued rise in cases until at least 2045[6], indicating that HPV+ HNSCC will be a clinical challenge for the foreseeable future. The most common treatment for HPV+ HNSCC is surgery and high-dose cisplatin chemotherapy given concurrently with radiation[7,8]. While this treatment strategy cures the majority of patients[9,10], significant long-term patient burdens can occur, including difficulty swallowing, kidney damage, osteoradionecrosis, and delayed death due to chronic aspiration[11]. Because many of these patients with HPV+ HNSCC will survive, they must cope with the toxicity of their treatment potentially for decades. Conversely, as the number of patients continue to rise, encountering patients with relapsed disease is becoming an increasingly common clinical scenario. Thus, there is an unmet need to develop accurate biomarkers to de-intensify treatment for HPV+ HNSCC with low risk of relapse. In the meantime, it is important to develop novel and more effective therapies for patients with recurrent HNSCC because of their generally poor response to conventional treatment regimens.

The mechanism of HNSCC recurrence is not fully understood to date. Emerging evidence suggests that the amount, type, and levels of activation of tumor-infiltrating lymphocytes are associated with response to conventional therapy and immunotherapy or survival[12-18]. For example, a recent multi-omics study carried out on HPV− HNSCC identified actin dysregulation and immunosuppression, due to widespread deletion of immune regulatory genes, as important determinants of disease pathology and response to treatment[19]. Intriguingly, these insights were obtained primarily from analysis of the tumor proteome rather than the transcriptome as the protein data substantially outperformed RNA data in co-expression-based function prediction[19]. This is reminiscent of another recent study demonstrating that the proteome, but not the transcriptome, of pancreatic ductal adenocarcinoma (PDAC), was able to differentiate different tumor groups[20]. Furthermore, large-scale proteogenomic studies estimates that only 20% of protein-coding genes have high correlations between the transcript and protein abundance[21]. The poor-to-moderate correlation between the transcriptome and proteome data underscores the importance of characterizing the tumor proteome independently or together with the transcriptome. In this regard, we note that no systematic proteomic study has been carried out to date on HPV+ HNSCC and consequently, the proteome of this cancer subtype remains unexplored.

Here we report a deep and quantitative mass spectrometry (MS) analysis of 15 HPV+ HNSCC samples, including 7 with recurrence. Our study uncovered systematic changes in the proteome and phosphoproteome between the recurrent and non-recurrence tumors that suggests extensive remodeling of the tumor microenvironment and reprogramming of the kinome associated with recurrence. Furthermore, numerous mechanism-based biomarkers and potential therapeutic targets emerged from our deep proteomic and phosphoproteomic analysis, yielding a valuable resource for future exploration.

Methods

Patient cohort

Fresh tumor samples were prospectively collected from patients with HPV+ oropharyngeal squamous cell cancer at the Victoria Hospital and London Health Science Centre, London, Ontario, Canada between 2010 and 2016. The study was approved by the Research Ethics Board at Western University (REB 7182) and informed written consent was obtained from each patient. The samples were frozen immediately after surgical resection using the optimal cutting temperature compound as cryo embedding matrix. Patient demographics and survival outcomes were prospectively collected. Frozen section analysis was carried out to confirm tumor cellularity greater than 70%. HPV status was confirmed via PCR and Sanger sequencing. All the primary tumor samples (NR1-NR8, and R1-R3) were collected before any treatment. For the samples R4-R7, the patients received cisplatin or radiation treatment until no evidence of disease was observed. However, their tumors redeveloped after the treatment was discontinued. The R4-R7 relapse tumor samples were collected by salvage surgery. None of the 15 patients received immunotherapy. Detailed clinical information is provided in the Supplementary Table 1.

Sample processing for mass spectrometry analysis

Tissue processing

The frozen tissue was washed in ice-cold phosphate-buffered saline (PBS) once and homogenized in a pre-cooled mortar and pestle, then resuspended in 1 ml cold buffer containing 8 M urea, 50 mM Tris-HCl (pH 7.6), 2% (v/v) proteinase inhibitor (Sigma, P8340) and 1 mM NaF. The mixture was rotating-mixed at 4 oC for 20 min to dissolve the proteins. Tissue debris was removed by centrifugation, and the supernatant was transferred into a fresh tube and mixed with 5-fold volume of cold precipitation mixture (containing 50% acetone, 50% ethanol, and 0.1% acetic acid) and incubated overnight at −20 °C.

Cell culture and pervanadate treatment for the booster channel

The patient-derived HNSCC cell line 93VU147T was tested for mycoplasma contamination and kept in Dulbecco’s Modified Eagle Medium (DMEM) supplemented with 10% FBS and 1% penicillin-streptomycin. T cells were isolated from blood of healthy donors. Both 83VU137T cells and T cells were treated with sodium pervanadate as follows. Briefly, the pervanadate solution was prepared by adding 10 μl of 0.1 M sodium orthovanadate to 10 μl of 0.2 M hydrogen peroxide (diluted 50-fold from a 30% stock). The solution was then incubated at room temperature for 15 min. Excess hydrogen peroxide was inactivated by adding 2 μl of catalase in PBS (10 mg/ml).

Protein processing and digestion

The protein precipitate was collected by centrifugation, washed with 75% ethanol, and then redissolved in 200 µl 9 M urea containing 50 mM HEPES, pH 8.0, and 10 mM DTT. The mixture was incubated for 1 h at room temperature. Subsequently, the sample was alkylated with 28 mM iodoacetamide for 40 min in the dark and the reaction was quenched by adding 10 mM DTT. The sample was diluted by adding 700 µl 50 mM HEPES, pH 8.0, 1 mM sodium orthovanadate, and digested with LysC (1 mAU per 25 µg protein) for 2 h at 28 °C, followed by incubating with trypsin at 1:50 enzyme-to-substrate ratio overnight at 28 °C. The sample was acidified by adding TFA to 0.5% final concentration and was desalted on a C18 SPE column (Waters WAT054955). The peptides were eluted in 70% acetonitrile/0.1% formic acid and dried by Speedvac.

Tandem Mass Tag (TMT) labelling

For mass spectrometry analysis, we labelled 15 patient and 5 healthy control samples with the 11-plex TMT isobaric labelling reagent (Thermo A37725; see Supplementary Table 2 for sample identities with TMT set/channel numbers). In addition, we employed the pervanadate BOOST channel approach[22] by including a 1:1 mixture of the pervanadate-treated 93VU147T cells and T cells as the reference channel (channel 1 of each 11-plex sample). Three sets of 11-plex reagents were used to label all samples.

The TMT labelling procedure was modified from Chua et al.[22]. Briefly, the desalted peptides were reconstituted in 0.1% formic acid to determine peptide concentration by the BCA protein assay kit (Pierce 23225). Peptides in 180 µg portions from each sample were aliquoted and vacuum-dried. A reference sample was prepared by pooling the samples C3, NR1, NR3, and NR6. Each of 0.8 mg 11-plex TMT labelling reagents was reconstituted in 41 µl acetonitrile. The peptides were reconstituted in 36 µl of 0.1 M HEPES (pH 8.5) to prepare 5 mg/ml peptide solution and were mixed with 13.5 µl of the TMT reagent to label the peptides at room temperature for two hours. After the reaction, a 1 µl aliquot was taken from each sample to check TMT labelling efficiency by mass spectrometry. The reaction was quenched by adding 2.7 µl of 5% hydroxylamine. The 11 samples in a TMT set were combined (2 mg total peptides) and desalted on SepPak C18 cartridges.

Phosphopeptide enrichment

We previously reported that engineered Src SH2 domains called the SH2 superbinder (SH2S) could be used to enrich pTyr peptides for mass spectrometry analysis[23]. SH2S-Agarose beads (Precision Proteomics, London, Canada) were used for pTyr peptide enrichment. Briefly, the TMT-labelled peptides were reconstituted in 50 mM ammonium bicarbonate and incubated with 100 µl SH2S beads for 30 min at room temperature with rotation. The flow-through fraction was saved for pSer/Thr enrichment by Ti4+-immobilized metal affinity chromatography (IMAC)[24] (a gift from Dr. Mingliang Ye, Dalian Institute of Chemical Physics). The beads were washed twice in 0.2 M ammonium bicarbonate, followed by 2x washes in 50 mM ammonium bicarbonate. The bound pTyr peptides were eluted using 0.4% trifluoroacetic acid (TFA). The eluted pTyr peptides were directly loaded onto the fractionation columns (see next section).

The flow-through fraction contains peptides not captured by the SH2S-agarose beads. Serine/threonine-phosphorylated peptides were enriched by Ti4+-IMAC resin following a published protocol[24]. Briefly, a 500 µg portion of the flow-through fraction from the SH2S enrichment step was mixed 1:1 (v/v) with 80% acetonitrile/6% TFA solution and then loaded to the IMAC resin. After incubation and wash steps (wash-1 solution: 50% acetonitrile, 6% TFA, 200 mM NaCl, wash-2 solution: 30% acetonitrile, 0.1% TFA), the peptides were eluted by 10% ammonia, and dried by Speedvac.

High-pH fractionation

The peptides were separated by Pierce high pH reversed-phase peptide fractionation kit (Thermo 84868). For proteome analysis, 10-µg portions of the peptides after TMT labelling were diluted in 300 µl 0.1% TFA for loading and 8 fractions (10–50% acetonitrile in 0.1% triethylamine) were collected. The IMAC-enriched peptides were resuspended in 300 µl 0.1% TFA and loaded onto the fractionation column. The SH2S-enriched peptides, upon elution in 0.4% TFA from the SH2S-agarose beads, were loaded directly to the fractionation column. For the IMAC and SH2S-enriched peptides, the peptides were eluted at 5, 10, 12.5, 15.0, 17.5, 20.0, 22.5 and 50% acetonitrile in 0.1% triethylamine in a step-wise manner and the resulting fractions (8) were concatenated into 4 tubes (i.e., combining fractions 1-5, 2-6, 3-7 and 4-8). The fractionated peptides were dried on Speedvac.

LC-MS/MS experiments

The fractionated peptides were reconstituted in 2% acetonitrile/0.1% formic acid (FA). The peptides were analyzed by the data-dependent acquisition (DDA) method on an EASY-nLC 1000 system coupled to a Q-Exactive Plus mass spectrometer (Thermo Fisher Scientific). The peptides were loaded on an Acclaim PepMap 100 C18 column (20 mm with 75 µm in diameter, Thermo 164946), and separated on an EASY-Spray ES803A analytical column (Thermo Fisher Scientific) at the flow rate of 300 nl/min and a linear gradient from 3 to 40% acetonitrile in 0.1% formic acid. The gradient length was 2 hours for the pTyr phosphoproteome fractions and 4 h for the proteome or the IMAC phosphoproteome fractions. See Supplementary Table 3 for parameters of mass spectrometry data acquisition.

MS data processing

FragPipe version 16.0 was used for data processing[25]. For the proteome and IMAC datasets, the pool channel was used as the reference channel for batch correction. Due to the reporter ion interference from the phosphotyrosine boost channel (channel 1) to the pool channel (channel 3) caused by isotopic impurities of the TMT reagents, the pool channel was not usable for TMT batch correction for the SH2S-enriched pTyr datasets. Instead, the technical triplicate sample (the tonsil control C2 in all three TMT batches, see Supplementary Table 2) was used as the reference channels for batch correction. The human protein sequences obtained from UniProt (20367 reviewed entries, February 2020) were used in the database search. For the proteome data processing, the TMT10-bridge workflow was used for the proteome search (by changing to the plex to TMT11). For the phosphoproteome data, the TMT10-phospho-bridge workflow was loaded, followed by changing the plex to TMT11 and the minimal peptide length to 6. The median centering normalization by FragPipe was applied to the log2 intensities.

Proteome and phosphoproteome data analysis

For data analysis, only the proteins (for proteome) or phosphosites (for phosphoproteome) observed in at least three samples in each group (control, non-recurrent (NR) or recurrent (R)) were retained. Perseus version 1.6.14.0 was used to analyze the data[26]. The VolcaNoseR server was used for drawing volcano plots. The list of the human kinases was based on Manning et al.[27]. The list of meta signatures that reflect common expression programs are taken from Puram et al.[28]. The Metascape server was used for functional gene annotation analysis[29]. Pathview Web was used for KEGG pathway analysis[30]. The heatmaps were prepared with the Morpheus server. The master protein-protein interaction network was constructed by extracting proteins and phosphosites significantly more abundant in recurrence samples (>1.5-fold more abundant in the recurrence group, and p < 0.1 between NR and R groups) in the mass spectrometry data. The STRING database[31] was used to connect the nodes on Cytoscape (confidence cutoff 0.9)[32]. For the Kaplan–Meier survival analysis, the data for 71 HPV+ HNSCC cases were retrieved from The Cancer Genome Atlas (TCGA). For each gene, the 71 patients were divided into high expression (n = 36) or low expression group (n = 35), and the disease-specific survival month was plotted using GraphPad Prism 9.

Statistical analysis

The log2 intensity values were used for statistical analysis in volcano plots and box plots. Unpaired 2-sided t-tests were conducted between the non-recurrence and recurrence groups using Perseus. GraphPad Prism 9 was used for the survival analysis, with the log-rank test. The box in the box plot extends from the 25th to 75th percentiles. The center line is plotted at the median. The whiskers go down to the smallest value and up to the largest.