signature=44d86896a16964053cb3a2e29a9d831f,APOBEC3A is an oral cancer prognostic biomarker in Taiwan...

Patient characteristics and clinical specimens

Two cohorts of patients who are ethnically Taiwanese were enrolled in Chang Gung Memorial Hospital at Linkuo. The first cohort comprised 50 treatment-naive OSCC patients who underwent comprehensive preoperative work-ups, and curative resection with/without post-operative adjuvant therapy. Sections of tumors that contained at least 75% tumor nuclei among the total cellular nuclei were used for DNA and RNA isolation and quantification before subjected Next Generation Sequencing experiments. Data on demography, risk exposure, clinical characteristics, and histopathological features were collected and are provided in Supplementary Table 1. For the second cohort, tumor specimens were obtained from 188 consecutively enrolled patients (172 men and 16 women) between June 2006 and November 2014, and all had regular follow-up information. Patient characteristics are presented in Supplementary Table 3. Institutional research ethical approval and written informed consent were obtained for all participants in the study (Institutional Review Board of Chang Gung Memorial Hospital, Taiwan). In general, patients underwent standard preoperative work-ups according to institutional guidelines, including a detailed medical history, a complete physical examination, computed tomography or magnetic resonance imaging scans of the head and neck, chest radiographs, a bone scan, and an abdominal ultrasound. Primary tumors were excised with adequate margins under intraoperative frozen-section control. Pathological TNM classification of all tumors was established according to the American Joint Committee on Cancer Staging Manual (2010). After discharge, all patients had regular follow-up visits every 2 months for the first year, every 3 months for the second year, and every 6 months thereafter.

The OSCC primary tumors were all excised with adequate surgical margins, and the tumor margin tissue was sent for intraoperative fresh frozen section histopathology analysis. If margins were not deemed to be tumor free, further resection was performed. Nonetheless, five of the 188 resected OSCC specimens were found to be tumor-positive in surgical margins. Various types of neck dissection were performed according to the primary tumor site and clinical lymph node status. Postoperative radiotherapy was performed on patients with pathologically identified T4 tumors and positive lymph nodes within 6 weeks following surgery. Patients with any of the following pathological features received adjuvant concurrent chemoradiotherapy: metastasis in multiple neck lymph nodes, extracapsular spread, positive surgical margins, perineural invasion or nodal dissemination at level 4 or 5. The chemotherapy was a cisplatin-based regimen, and the total radiation dose was 66 Gy. The prescribed dose was delivered in fractions of 1.8–2 Gy per day for 5 days per week.

Patients underwent standard postoperative work-ups according to institutional guidelines, which included complete physical examination at regular follow-up visits. Computed tomography or magnetic resonance imaging of head and neck and chest radiographs were performed 3 months after the treatment and every 6 months for 3 years. Additional radiological examinations, bone scans, and abdominal ultrasonography were performed for any suspicious recurrence or second primary tumors noted clinically. All patients completed regular follow-up visits every 2–3 months for the 1st year after discharge, every 3–4 months for the second and third years, and every 6 months thereafter.

Detection of human papillomavirus by E6 nested PCR

For tissue samples, genomic DNA was extracted with a DNeasy Blood & Tissue kit (Qiagen, Chatsworth, CA, USA) according to the manufacturer’s instructions. HPV detection and typing were performed as described previously

Genomic DNA extraction

Genomic DNA was extracted from the tumor tissues using a QIAamp DNA Mini kit (Qiagen), and from whole-blood samples using a Puregene Blood Core kit (Qiagen), following the manufacturer’s protocols. The quality and DNA concentrations of the samples were assessed using a Qubit fluorometer (Thermo Fisher Scientific) and 0.8% agarose gel electrophoresis.

Next generation sequencing

For whole-exome sequencing, high-quality genomic DNA (1.5 μg per sample) was subjected to capture using a SureSelect Human All Exon v5 + UTR kit according to the manufacturer’s protocol (Agilent Technologies). The qualified exome-captured libraries were sequenced using a HiSeq 2000 with the TruSeq PE Cluster kit v3 and TruSeq SBS kit v3 (all from Illumina) according to the manufacturer’s protocol. For RNA sequencing, libraries were prepared using the TruSeq RNA Access Library Prep Guide (Part # 15049525 Rev. B; Illumina) according to the manufacturer’s instructions. Equal concentrations of each library were sequenced using a NextSeq 500 (Illumina) platform. For target sequencing, libraries were prepared from each sample using an Ion AmpliSeq Comprehensive Cancer Panel (CCP; Thermo Fisher Scientific) following the manufacturer’s instructions. Equal concentrations of each library were sequenced using an Ion Proton System (Thermo Fisher Scientific).

WES data pre-processing

For all 50 donors, DNAs from tumor tissues and matched PBMCs were sequenced. The sequenced reads were mapped to the hg19 reference genome using BWA memhttp://broadinstitute.github.io/picard/). The average coverage was 244 × , the average mapping rate was 99.63%, and the average on-target rate was 69.78%. Supplementary Data 1 presents the coverage, mapping rate, and reads-on-target rate for each sample. BAM-formatted files, including read sequences, sequence qualities, and alignment status, were generated by reference genome mapping using BWA memhttps://www.broadinstitute.org/cancer/cga/indelocator), and Oncotatorhttps://wiki.nci.nih.gov/display/TCGA/Mutation+Annotation+Format+(MAF)+Specification). In total, 50 sets of tumor and matched normal samples were analyzed; from them, a total of 24,051 somatic mutations and 9883 reported germline single nucleotide variants (previously deposited in dbSNP3 lists all the somatic mutations found.

Identification of putative driver genes

The MutSigCV algorithmq-value, which was generated as the Benjamini–Hochberg false discovery rate. All parameters were set to default values, and the q-value was set to 

Validation of somatic mutations

To validate the somatic mutations detected by our exome sequencing, we used the Ion PGM System (Life Technology) to investigate mutations located in 409 genes included in the Ion AmpliSeq Comprehensive Cancer Panel (CCP, Life Technology). The average sequencing depth was >1000 ×. The Ion Reporter software (variant caller V5.0.3.5; Life Technology) was used to perform variant calling and identify mutations across 49 matched tumor and normal samples. The CCP was found to contain 816 somatic SNPs. Details for these mutations is presented in Supplementary Data 4. For validation, we used 49 matched tumor and normal sample pairs and a total of 428 mutation sites covered by the CCP target-sequencing regions.

Identification of copy CNVs

After the sequence reads were mapped to hg19 with BWA

Identification of DEGs in OSCC samples from OSCC-Taiwan

In first cohort, we have carried out RNA-Seq in normal/tumor paired sample in 39 out of the 50 patients. All RNA reads were first trimmed by Trimmomatic2. The expression levels of genes in the 39 normal/tumor pairs were estimated by RSEMA3A, A3B, and A3A_B, the genotype of each sample was taken into consideration for quantification. Genes having median transcripts per million (TPM) larger than 0.5 were used in DEG detection. A total of 3548 DEGs were selected with criteria adjusted p-value  2 with Partek Genomics Suite software (Inc. P. Partek Genomics Suite. St. Louis). The 3548 identified DEGs were further used for hierarchical clustering analysis.

Pathway analysis of the identified DEGs

The 3548 DEGs were subjected to pathway enrichment analysis using MetaCore from Thomson Reuters, which uses hypergeometric testing to examine whether genes from various pathways are found more frequently than expected in a list of DEGs.

Luciferase reporter assay

A3A 3′UTR and A3B 3′UTR were PCR amplified from genomic DNA of OC3 oral cancer cell line, using the same forward primer 5′-ACTAGTAGGATGGGCCTCAGT CTCTAAG-3′; reverse primer 5′-AAGCTTAGTGTTTGTGGAAACTCTTGCAATT C-3′ is for A3A 3′UTR, and 5′-AAGCTTAGTGTTTGTGGAAACAATTATGGAAG-3′ is for A3B 3′UTR. The amplified products were cloned into the pMIR-Report-Vector (Ambion). Precursor of miR-409 was amplified from OC3 cell genomic DNA, and cloned into pcDNA6.2-GW/EmGFP-miR (Invitrogen). 293T cells (2 × 105) were subjected to transient calcium-phosphate-mediated transfection with the following: 10 ng of pMIR-Lusiferase-A3B 3′UTR; 1 μg of the vector control or expression vectors encoding miR-409 and 10 ng of pCMV-Renilla (Promega). Luciferase and renilla activities were measured using the Dual-Luciferase Reporter Assay System (Promega). The luciferase values were normalized to those of renilla, and the results are presented as the luciferase/renilla ratio.

Comparison of DEGs in OSCC-Taiwan and OSCC-TCGA

To examine whether the identified DEGs exhibited similar expression changes in the TCGA transcriptome database (http://cancergenome.nih.gov/), we downloaded Level 3 RNASeqV2 data for HNSC. This dataset contains 315 OSCC samples that were taken from anatomic sites within the oral cavity and designated as squamous cell carcinoma. Among them, 30 samples also had transcriptomic data available from their normal tissue sample. We used the Partek Genomics Suite software to compare the expression levels of genes sharing the same Entrez IDs in all (315 tumor + 30 normal) samples.

Mutational signature analysis

The Wellcome Trust Sanger Institute (WTSI) Mutational Signature Framework

RNA-Seq-based detection of A3A_B fusion transcripts

To detect A3A_B fusion transcripts, we adapted a previously described RNA-Seq read mapping approachA3A, A3B, and A3A_B, we used the Sashimi plot of Integrated Genome Viewer 2.3.77A3A_B fusion form.

Statistical analysis

Patient characteristics were stratified using various clinicopathological factors and evaluated by the chi-square test or Wilcoxon test. Multivariate analyses were applied to define OS, DSS and DFS. Survival rates were estimated by Kaplan–Meier plotting and compared by log-rank test. Statistical analyses were performed using SAS software (v.9.3). All patients received follow-up consultations at our outpatient clinic until March 2016 or death. All p-values were two-sided, and the significance level was set at p 

Protein quantification

Each tissue specimen was homogenized with 0.1% RapiGest buffer (Waters, Bedford, MA, USA) using a bead-beating homogenizer (Precellys 24; Bertin Technologies, Ozyme, France). The tissue extracts were then digested with trypsin (Promega, Madison, WI, USA) and labeled with the iTRAQ reagent according to the manufacturer’s protocol. The iTRAQ 114 reagent was used to label the digestion product of 30-pairs of OSCC tissues, and iTRAQ 115 and 116 reagents were used to label peptides of non-tumor tissue and tumor tissue, respectively. The labeling mixtures were desalted with an in-house-built C18-microcolumn and analyzed using on-line 2D-HPLC (Dionex Ultimate 3000; Thermo Fisher) coupled with a 2D linear ion trap mass spectrometer (LTQ-Orbitrap ELITE; Thermo Fisher)

CNV analysis for OSCC-TCGA

Somatic CNVs for HNSC samples were available in TCGA database. As this project focused on OSCC, we selected 315 samples taken from the alveolar ridge, floor of mouth, buccal mucosa, hard palate, lip, tongue, or oral cavity. To assess the selected samples for CNVs, the segment means derived using Affymetrix Genome-Wide SNP Array 6.0 were downloaded from TCGA (http://cancergenome.nih.gov/) and analyzed with GISTIC 2.0

RNA extraction, quantification, and fusion gene detection

Paired OSCC tumor/pericancerous normal tissues were individually homogenized in liquid nitrogen and subjected to total RNA extraction using RNAzol B (Tel-Test, INC.), according to the manufacturer’s protocol. First-strand complementary DNA (cDNA) was synthesized from 5 μg of total RNA using oligo dT primers and SuperScript III RT and qRT-PCR analysis was carried using the QuantStudio 12K Flex Real-time PCR System (Applied Biosystems). In 10 μl of reaction volumes containing 2 μl of 10-fold diluted cDNA, 5 μl of 2 × master mix (Roche), and 0.15 μM of primers. The utilized primers included: A3A-q2F (5′-ATGGCATTGGAAGGCATAAG-3′) and A3A-q2R (5′-CAAAGAAGGAACCAGGTCCA-3′) for detection of A3A; A3B-qrtF (5′-GACCCTTTGGTCCTTCGAC-3′) and A3B-qrtR (5′-GCACAGCCCCAGGAGAAG-3′) for detection of A3B; and TBP-F849 (5′-TGCTCACCCCACCAACAATTTAG-3′) and TBP-R969 (5′-CTGGGTTTGATCATTCTGTAGATTAA-3′) for detection of TPB as an internal control. The qRT-PCR conditions were as follows: 95 °C 10 min, followed by 40 cycles of 95 °C for 30 s, 60 °C for 1 min.

A3A_B fusion gene was amplified from cDNA by PCR using primers A3A_B5′ (5′-ATGGAAGCCAGCCCAGCA-3′) and A3A_B3′ (5′-TCAAATTAAAATTAATTGACTCTGATT-3′) and the PCR conditions: 95 °C for 10 min, followed by 35 cycles of 95 °C for 30 s, 60 °C for 30 s, and 68 °C for 30 s, and a final soak at 72 °C for 10 min. The 800-bp PCR product was examined by gel electrophoresis, followed by validation with Sanger sequencing.

Genotyping of the APOBEC3B-deletion polymorphism

For OSCC-Taiwan cohort, genomic DNA was obtained from PBMCs of OSCC patients, and 20 μg of each sample was used for PCR-based A3A genotyping. Primers Del F (5′-TAGGTGCCACCCCGAT-3′) and Del 2R (5′-CTAAATATGGAGCCAATTAA-3′) were used to generate a 757-bp product as an indicator of the deletion polymorphism, while primers Ins 2F (5′-TGTCCCTTTTCAGAGTTTGAGTA-3′) and Ins 2R (5′-TCCTTAGAGACTGAGGCCCAT-3′) were used to generate a 449-bp product that represented the wild-type chromosome. The PCR conditions were as follows: 95 °C for 10 min, followed by 37 cycles of 95 °C for 15 s, 55 °C (for the deletion polymorphism) or 60 °C (for the wild-type sequence) for 30 s, and 72 °C 30 s, followed by a final soak at 72 °C for 10 min. The obtained products were analyzed by agarose gel electrophoresis.

For all the HNSC samples having WES data available in TCGA (phs000178.v9.p8) include 318 unique samples taken from anatomic sites within the oral cavity and designated as representing OSCC. We download controlled-access WES bam files for the OSCC samples of these 318 cases from the GDC websites (https://gdc.cancer.gov/, Data Release V2). To identify APOBEC3B-deletion polymorphism from the bam files, we used the expectation-maximization approach to model the copy number of each sample based on the ratio of the sequencing depths inside and outside the deletion region

Data availability

The sequencing data have been deposited in the NCBI Sequence Read Archive (SRA) database under the accession code SRP078156. The data of WES, RNASeqV2 and Affymetrix Genome-Wide SNP Array 6.0 referenced during the study are available in a public repository from the TCGA (http://cancergenome.nih.gov/) website. The authors declare that all the other data supporting the findings of this study are available within the article and its supplementary information files and from the corresponding author upon reasonable request.