Prognostic Implications of Six Altered Genes in Asian Non-Surgical Esophageal Carcinoma Patients Treated with Chemoradiotherapy

Introduction

Esophageal cancer (EC) is the one of most aggressive tumors in China, with a higher ratio of age-standardized mortality rate in Chinese males compared to Chinese females as well as the UK, the USA and even worldwide.¹ Esophageal squamous cell carcinoma (ESCC) is more often observed in the Asian population, while esophageal adenocarcinoma (EAC), on the other hand, occurs mainly in North America and Western Europe.^1,2 Smoking and alcohol consumption are two high-risk factors for ESCC,³ whereas obesity is a strong risk factor for EAC.⁴ The 5-year survival rates of EC depend on several factors, including the stage of cancer at the time of diagnosis. In China, EC patients achieved a 5-year survival rate of around 40%.^5,6 In the US, the 5-year survival rate is 47% for EC patients with only localized tumors and 20% for all EC patients combined according to SEER database (cancer.org).

According to the National Comprehensive Cancer Network (NCCN) guideline, surgery is a major treatment used for locally advanced resectable EC patients with additional preoperative chemoradiation or perioperative chemotherapy to improve survival.⁷ Targeted treatment strategies have also been explored in EC patients including HER2-targeted therapy, anti-angiogenesis therapy, and immunotherapy. HER2 TKI trastuzumab has been approved by the FDA for HER2-positive advanced ECs.⁸ Ramucirumab, a VEGFR-2 antibody, has been approved for pre-treated patients with advanced or metastatic EAC initially as monotherapy and subsequently as combination therapy with paclitaxel.⁹ Pembrolizumab is approved in 2017 by the FDA for EC patients with high microsatellite instability and/or PD-L1 expression.¹⁰

With the aid of next-generation sequencing, molecular characteristics are becoming an emerging aspect in multidisciplinary treatment decision-making to create an overall therapeutic plan for cancer patients. The genomic landscape of EC genomes has been studied in both Asian and Western populations.^11,12 Highly mutated genes in EC included TP53, NOTCH1, PIK3CA, RB1, CDKN2A have been identified in different ethnic groups. The association between gene alterations and prognosis has been indicated.^13,14 Here, we performed targeted panel sequencing of tumor tissues from 69 Asian EC patients. A TCGA cohort of 143 EC patients was also included in this study for comparison. Six gene alterations were shown to be potential prognostic biomarkers for overall survival in EC patients using univariate and multivariate Cox regression analyses. The Cox proportional hazard model built on the six potential prognostic biomarkers demonstrated the ability to select patients with worse disease outcome in the EC patients.

Method

Patient Cohort

A total of 69 patients clinically diagnosed as non-surgical EC was retrospectively recruited from Shandong Provincial Hospital Affiliated to Shandong University according to the NCCN guidelines.⁷ The inclusion criteria were as follows: 1) all patients had histologically proven primary non-surgical EC and 2) the patients were treated with chemoradiotherapy. The study was approved by the Ethical Review Board of Shandong Provincial Hospital, and informed written consent was obtained from each participant.

DNA Extraction and Library Preparation

As previously described,¹⁵ the genomic DNA was extracted from formalin‐fixed and paraffin‐embedded (FFPE) using the QIAamp DNA FFPE Tissue Kit (Qiagen). The extracted DNA quantity was evaluated using a Qubit 3.0 fluorometer and extracted DNA quality was measured using Nanodrop 2000 (Thermo Fisher Scientific). Sequencing libraries were prepared using the KAPA Hyper Prep Kit (KAPA Biosystems) and sequenced with a pan-cancer panel of 422 genes.

Sequencing Data Analysis

Sequencing data analysis was performed as previously described.¹⁵ In brief, FASTQ file quality control was performed with trimmomatic¹⁶ (below 15 or N bases were removed). Reads were mapped to the reference Human Genome (hg19) using Burrows-Wheeler Aligner (BWA-mem, v0.7.12; https://github.com/lh3/bwa/tree/master/bwakit). VarScan2¹⁷ was used for somatic mutation detection. Genome Analysis Toolkit was applied to local realignment around the indels and base quality score recalibration (GATK 3.4.0; https://software.broadinstitute.org/gatk/), which was also used for detecting germline mutations. Somatic variants were called if mutant allele frequency (MAF) was at least 0.2%, and at least three supporting-reads from both directions. Common SNPs were filtered out according to dbSNP (v137) and the 1000 Genomes database, and annotated using ANNOVAR.¹⁸ Genomic fusions were identified by FACTERA¹⁹ with default parameters. Copy-number variations (CNVs) were detected using ADTEx (http://adtex.sourceforge.net) with default parameters. Somatic CNVs were identified with the cut-off of 0.65 for copy-number loss and 1.50 for copy-number gain using paired normal/tumor samples for each exon.

Data Analysis

The Kaplan–Meier survival curves were performed to estimate OS in different genomic groups. The log‐rank test was performed to analyze differences between groups. The univariate and multivariate analyses were performed to evaluate the prognostic value of clinicopathological characteristics and gene alterations on OS. Overall survival (OS) was measured from the date of pathological diagnosis of EC to the date of death or last follow-up.

Backward stepwise selection with the Akaike information criterion (AIC) was used to identify variables for the multivariate Cox proportional hazards model. Model performance was evaluated by assessing discrimination against the index of concordance (C-index) and plotting Kaplan–Meier curves over the quartiles of prediction by nomogram.²⁰ The optimal cutoff for risk score was selected using X-tile based on the best model as shown in Figure S4. X-tile was a published method for biomarker assessment and out-based cut-point optimization. X-tile can assess the robustness of the relationship between a biomarker and outcome by the construction of a two-dimensional projection of every possible subpopulation, which was also employed in other studies.^21,22 Using the optimal cut-off, Kaplan–Meier curve was generated in the test and complete sets to validate the used cutoff and risk score model.

Result

Description of Analytical Cohort

We obtained 69 tumor tissue samples from a total of 69 Asian patients with non-surgical esophageal cancer (AEC cohort, Table 1) was enrolled in this study with a median age of 64 years old (yrs), ranging from 41 to 83 yrs. More than 80% of the patients (81.60%) were male and the rest (17.39%) was female. Almost all patients were squamous cell carcinomas (SCC, 98.55%) except one patient was adenocarcinomas (ADC, 1.45%). There were 17 (24,64%) stage II, 40 (57.97%) stage III, and 12 (17.29%) stage IV patients. Among them, 65.22% of the patients were smokers and more than half of patients (52.17%) had a history of alcohol consumption. All AEC patients were treated with chemoradiotherapy.

Table 1 Demographic Characteristics of Patients in AEC Cohort

Meanwhile, a TCGA cohort consisting of 143 patients with esophageal cancer was employed in this study (Table S1). The TCGA cohort included 46 (32.17%) Asian cases, 77 (53.84%) Caucasian cases, and 20 (13.99%) cases of other races. The median age was 61 yrs, ranging from 36 to 90 yrs. Similar to our cohort, the majority of the TCGA cohort was male (87.41%). Histology subtypes included 39.86% ADC and 60.13% SCC. The tumor stage included 78 (54.55%) stage II, 56 (39.16%) stage III, and 9 (6/29%) stage IV. In the TCGA cohort, 58.74% of the patients were non-smokers and 55.94% of the patients never had alcohol.

Genomic Landscape of Asian Patients with Esophageal Cancer

The genomic landscape of the AEC cohort and TCGA cohort was shown in Figures 1 and S1, respectively. In the AEC cohort, nearly 95% of patients harbored TP53 alterations and more than half of patients were identified with NOTCH1 alterations. Co-amplifications of FGF19 and CCND1 were found in 36.2% of cases. Other altered genes of high frequencies in the AEC cohort included MCL1(39.1%), MYC (31.9%), PIK3CA (21.7%), and EP300(18.8%). In the TCGA cohort, TP53(85.3%) and PIK3CA (19.6%) were the top two frequently mutated genes. Compared to the AEC cohort, co-amplifications of FGF19 and CCND1 were found at a slightly lower ratio (33.6%). Alterations in multiple DNA damage repair genes were identified in both cohort including ATM (10.1% vs 13.3%), ATR (11.6% vs 4.9%), SMARCA4(7.2% vs 6.3%). Interestingly, considering only targetable mutations, EGFR showed the highest ratio in the AEC cohort and ERBB2 mostly amplification showed the highest ratio in the TCGA cohort.

Figure 1 Genomic landscape of AEC cohort. The type of alterations was indicated by color. Each column represented one patient.

Univariate and Multivariate Cox Regression Analyses of Prognostic Parameters

Clinicopathological features and genetic alterations are all potential predictors of prognosis in cancer treatment. Next, we examined the association of these possible prognostic features with patients’ overall survival (OS) using the univariate Cox regression model in the AEC cohort. As shown in Table 2, clinicopathological features including gender, age, smoking status, and alcohol consumption were not predictors of the OS in the AEC cohort. TNM stage showed some association with OS, with poorer outcome in stage III–IV patients compared to stage II (Figure S2). The most frequently observed alterations TP53 mutation was not able to predict OS (Figure S3), which might due to limited TP53 wild-type patients.

Table 2 Univariate Cox Regression Analyses of Prognostic Parameters

Interestingly, seven gene alterations were found to be independent markers of OS (Table 2) in the univariate model including BAP1 mutation, BRIP1 mutation, KDR mutation, MYC amplification, RB1 variant (mutation and deletion), WRN mutation, YAP1 amplification. In the AEC cohort, the frequencies of these seven gene alterations varied from 5.8% to 31.88% (Table 2). Next, we performed a multivariate cox regression analysis to predict OS using seven potential prognosis-related genes identified in univariate analysis (Table 3). Except for KDR with an insignificant p-value of 0.91, the rest six genes all showed as independent factors in multivariate analysis. Compared to wild type, patients with alterations in these genes were associated with poorer OS (Figure 2). Meanwhile, we also look into the association of POLE gene alterations and DDR pathway gene alterations with OS, however, the results were insignificant (Figure S4).

Table 3 Multivariate Cox Regression Analyses of Prognostic Parameters.

Figure 2 Survival analysis in AEC patients with different gene alterations. Kaplan–Meier survival curves for overall survival for the 69 EC patients from AEC cohort. The overall survival of patients with YAP1 amplification (A), RB1 alteration (B), BAP1 mutation (C), MYC amplification (D), BRIP1 mutation (E) and WRN mutation (F) was compared to that of patients with wild-type genes, respectively.

The Construction of Proportional Hazards Model with Six Prognosis-Related Genes and Evaluation

To distinguish EC patients with poor disease outcomes, a Cox proportional hazards model has been constructed using identified prognosis-related genes. The model was evaluated using a stepwise selection approach and the best model was then chosen with the combination of six prognosis-related genes, which minimized the Akaike Information Criterion (AIC) and maximized index of concordance (C-index). The comparison of the C-index and AIC was shown in Table 4. The six-gene model achieved the lowest AIC of 264.09 and the highest C-index of 0.75 compared to models using single-gene alterations or other gene alteration combinations. The optimal cutoff point for risk score was set at 18 according to X-tile²² (maximum χ², p<0.0001, Figure S5). Based on this cutoff, the Kaplan–Meier survival curves showed AEC patients with a risk score >18 (15 patients, median OS:10.40 months) displayed a poorer OS than a risk score˂18 group (54 patients, median OS:41.86 months) (Figure 3A). We further examine the performance of this model using the TCGA cohort. The frequency of six prognosis-related genes in TCGA was shown in Table S2. With a cutoff at 18, this model was able to distinguish six EC patients with poor OS (median OS 14.31 vs 28.09 months, p=0.0008) (Figure 3B). We also validated our model using an independent Asian cohort of 47 ESCC patients treated with dCRT.²³ As shown in Figure S6, the model identified 22 ESCC patients with a risk score higher than 18, which displayed significantly worse OS (p=0.0022) than patients with a risk score lower than 18.

Table 4 Comparison of C-Index and Akaike Information Criterion (AIC) Between Models Using Different Variables

Figure 3 Kaplan–Meier curves for overall survival according to optimal cut-off point of six gene alterations in AEC (A) and TCGA (B) cohorts.

Discussion

Here, we employed a pan-cancer NGS panel of 422 genes to study the association between clinical characteristics, gene alterations, and overall survival in the Asian EC population. The prognostic value of six potential biomarkers was further evaluated in an EC cohort from TCGA. Compared to previous studies, the genomic landscape of this AEC displayed a higher ratio of NOTCH1, MCL1, PIK3CA alterations, which is likely due to the difference in the tumor stage as well as the sequencing method used.¹¹ In this AEC cohort, clinical characteristics including gender, age, smoking status, alcohol consumption were not associated with OS. Stage II EC patients showed slightly better OS than stage III–IV patients. Meanwhile, despite their high prevalence in EC, alterations in TP53, PIK3CA, NOTCH1 as well as FGF19 /CCND1 co-amplification were not associated with OS, which is consistent with previous studies.^13,23

Several gene alterations were identified as adverse prognostic factors in the AEC cohort, some of which belonged to different oncogenic signaling pathways in cancers such as YAP1 in the Hippo pathway, MYC in the MYC pathway, RB1 in the cell cycle pathway.²⁴ YAP1 amplification has been reported to be a prognostic factor of chemoradiotherapy in nonsurgical esophageal squamous cell carcinoma and is associated with shorter local recurrent‐free survival and OS.²³ MYC amplification was also reported to be associated with poorer OS in ESCC patients.²⁵ The associations between the rest alterations and OS in EC were First reported in this study with some indications of their prognostic value in other cancer types. BAP1, WRN, and RBIP1 are homologous recombination pathway-related genes, which can be targeted for DNA repair-targeted therapy.^26,27 BAP1 is deubiquitylase associated with many cancer pathways and BAP1 mutation has been reported as prognostic factors in predicting metastatic risk.^28,29 BRIP1 and WRN encode RecQ DNA helicases and are important in the normal double-strand break repair.^30,31 Altered BRIP1 is a targetable alteration in ovarian cancers and breast cancer patients with overexpression of BRIP1 displayed a poor survival rate.^32,33 KDR, also known as VEGFR-2, plays an important role in tumor angiogenesis and potential therapeutic target for esophageal carcinoma.³⁴

The use of genetic biomarkers in clinical settings is increasing, and many studies have been carried out on the characterization of prognostic biomarkers in EC patients to predict disease outcomes.^35–37 The identification of these markers provides a basis for detecting potential therapeutic strategies for specific molecular subtypes in clinical trials and will ultimately contribute to the personalized treatment plan for EC patients. For instance, some patients in AEC and TCGA cohort were identified as high-risk for poor disease outcomes using our model (a risk score higher than 18), different strategies could be considered for these patients: firstly, targeted therapy including trastuzumab (HER2 positive), ramucirumab (VEGFR); secondly, immunotherapy can be applied if the patient had MSI-H or PD-L1 expression according to NCCN guidelines;⁷ Thirdly, as shown in Figure 1, altered DDR genes especially homologous recombination pathway genes such as BRCA2, CHEK2, ATM/ATR were also identified. Preliminary studies showed that PARP inhibitors may enhance the radiosensitivity in ESCC patients.^38,39 Therefore, PARP inhibitor together with radiotherapy may also be an option for high-risk populations. More studies are needed to validate our observation.

Major limitations of this study were the small cohort size and the lack of a large-scale Asian EC cohort to validate the six candidates for OS prediction. Most Asian EC studies were either RNA and gene expression level or the survival information of the cohort was unavailable. Here, we used a TCGA cohort, which consisted of 32.17% of Asian cases and 53.84% of Caucasian cases, whereas the AEC cohort was all Asian. Meanwhile, the distribution of histology subtypes varied in both cohorts. ESCC accounted for 60.13% of the TCGA and 97.1% of the AEC. Most of the treatment history from the TCGA cohort is unavailable, which may have some impact on the interpretation of the result. In conclusion, our study identified potential prognostic biomarkers for Asian EC patients. Further studies and validations of the prognostic value of these biomarkers in larger Asian clinical cohorts are warranted.