A Novel Scoring System for Pivotal Autophagy-Related Genes Predicts Outcomes after Chemotherapy in Advanced Ovarian Cancer Patients

Ovarian cancer is one of the major threats to the health of women (1, 2). In spite of the relatively low incidence, it has the highest mortality among gynecologic cancers. The performance of current tests for ovarian cancer screening and surveillance has been poor, either used alone or in combination (3). Largely dependent on the cancer antigen 125 (CA125), as well as some CA125-based indexes such as risk of malignancy index (RMI) and risk of ovarian malignancy algorithm (ROMA), the approaches to screening and surveillance of ovarian cancer still have marked limitations (4–6). As a result, a majority of patients with ovarian cancer are at advanced stage at diagnosis, namely, stage III/IV by the criteria of the International Federation of Gynecology and Obstetrics (FIGO). Ultimately, most of these patients will develop recurrence and chemoresistance even after receiving the standard treatment of optimal debulking surgery and systemic cytotoxic platinum-based chemotherapy (7, 8). Currently, the FIGO staging system has an essential role in the clinical practice of ovarian cancer, and the correlation between ovarian cancer prognosis and the FIGO stage is well-established (7, 9). However, an apparent limitation of FIGO staging system is that it mainly focuses on clinical features whereas ovarian cancer is a disease with high degree of histologic and genetic heterogeneity (10). Previous studies have reported accurate quantitative prognosis prediction paradigms based on molecular markers or critical gene profiles that provided different information from FIGO stage and could be helpful to the optimization of therapeutic regimens (11, 12).

As one of the important intracellular stress coping systems, autophagy maintains homeostasis and provides essential energy for metabolism under the stress of chemotherapeutic agents. Autophagy activation is closely associated with chemotherapy resistance, the maintenance of tumor dormancy, and the survival of cancer stem cells responsible for tumor progression and relapse (13–17). Several clinical trials conducted in different types of cancer reported the safety and promising survival benefit of the therapeutics that target autophagy (18). Furthermore, a recent haploinsufficiency network analysis illustrated that autophagy and proteostasis pathways were more severely disrupted in ovarian cancer than 20 other types of cancer (19). The significance of autophagy in ovarian cancer suggests that it has the potential to serve as a functional biomarker in clinical management.

Several previous studies attempted to develop scoring systems by evaluating critical biological processes related to cancer (20–24), such as DNA damage-repair pathway, immunologic condition, and cell stemness. Similarly, autophagy is also considered to be a target biological process in the treatment of cancer. However, the role of autophagy in ovarian cancer seems highly context-dependent, because the dual effects of promoting tumor cell survival and suppressing oncogenesis were both validated in previous studies regarding key autophagy executors (18, 25–28). As for the treatment of cancer, therapies to inhibit and stimulate autophagy have been both reported, suggesting the importance of patient evaluation and selection (18). At the same time, little research has been conducted to utilize the autophagy-related (ATG) gene expression profile as a tool to evaluate the prognosis of the patients with ovarian cancer.

In our research, a novel scoring system that consists of screened pivotal ATG genes has been developed. Our model introduces a novel perspective and a promising tool in the prognostic assessment of advanced ovarian cancer patients who received platinum-based chemotherapy.

It is well-established that ovarian and breast cancer have a high degree of hereditary and clinical homology. Previous studies on hereditary and risk factors revealed the close relationships between the family history and elevated risks of these 2 cancers (29, 30) and also proved the increased risk of developing breast cancer as a second primary malignancy among patients with ovarian cancer (31). In addition, pan-cancer analysis (32) and research regarding the important clinical dimensions such as morbidity (30), diagnosis (33), and prognosis (34) supported the homology between ovarian and breast cancer as well. Thus, considering the genetic and clinical similarities between ovarian cancer and breast cancer, we collected and analysed the data from these two sets (ovarian and breast cancer patients) to develop the ATG scoring system. The normalized gene expression level data (RNASeqV2 data) and clinical information of 156 ovarian serous cystadenocarcinoma samples and 1,215 breast cancer samples were downloaded from UCSC Cancer (https://genome-cancer.ucsc.edu/2016/08/21) and The Cancer Genome Atlas (TCGA; https://www.nationwidechildrens.org/research/resources-infrastructure/core-facilities/biospecimen-core-resource/the-cancer-genome-atlas), respectively. The patients who were diagnosed with FIGO stage III/IV ovarian cancer and treated with platinum drugs were included. Those patients who had received immunotherapy, targeted molecular therapy, hormone therapy, or neoadjuvant chemotherapy were excluded from this study. The 156 patients with ovarian cancer obtained from TCGA were regarded as the training dataset for the following clinical significance analysis of the newly constructed scoring system in ovarian cancer.

After a substantial literature review and database searching, a set of 135 autophagy-lysosome pathway associated genes were obtained from the study of Perera and colleagues (35). Referring to the lysosome proteomics, lysosomal diseases and autophagy interactome, such a gene set represented a comprehensive combination of pivotal autophagy-related genes. Considering the feasibility of future clinical application, we deemed these 135 genes acceptable as our preliminary candidate gene list in both quality and quantity. Combined with the TCGA data of ovarian cancer and breast cancer, 131 genes were available for subsequent analysis. Next, a univariate Cox regression analysis (UVA) of the 131 genes was performed, in which univariate P value was obtained to screen the genes that had significant correlation (P value < 0.1) with the overall survival (OS). A total of 43 genes were found significantly related to the OS of either patients with ovarian or breast cancer and constituted the final ATG scoring system. Such genes were divided into hazardous genes (HR > 1) and protective genes (HR < 1). Notably, 7 genes were found significant in both ovarian and breast cancer with the same directional effect. In contrast, we also examined the candidate genes in another major gynecologic malignancy, uterine corpus endometrial carcinoma (UCEC), and found that three genes showed dual or triple significance among UCEC and ovarian or breast cancer. However, all of the 3 genes showed opposite roles between UCEC and the other 2 kinds of cancers (data not shown). Such findings also supported the homology and similarity between ovarian and breast cancer. Next, according to a published and effective strategy of score construction (21), we established our ATG scoring system as follows. For each hazardous gene of a patient, if its expression level was above the group median, the patient got +1 point from that gene; otherwise, the patient got −1 point. As for the protective genes, the scoring rules were just the opposite. A patient would get −1 or +1 point for a high or low expression level of a protective gene, respectively. For each patient, the initial score was 0, and the 43 genes were scored separately and then added up as the final ATG score.

OS was considered as the main outcome in our study. Kaplan–Meier analysis and log-rank test were adopted to compare the prognosis of patients with different ranges of score. UVA and Multivariate Cox regression analyses (MVA) were conducted to evaluate the significance of the ATG score in prognosis prediction compared with other classic clinical characteristics. OR and 95% CI were calculated by the conditional maximum likelihood estimate (MLE). According to Wang and colleagues (36), a nomogram was then constructed and the receiver-operating characteristic (ROC) curves were adopted for the comparison of predictive significance of different models. To assess the discrimination efficacy, the concordance index (C-index) was also calculated in each model. Higher C-index represents higher discrimination ability, the range of which is 0.5 to 1. The integrated discrimination improvement index (IDI) and net reclassification index (NRI) were calculated to compare the predictive power of the nomogram and traditional clinical factors as risk prediction models. Kaplan–Meier analyses and Cox regression models were also adopted to explore whether the ATG score was linked to the recurrence-free survival (RFS). All of the data were analysed by R 3.2.2 and SPSS 19.0 (IBM SPSS Inc.). A P value < 0.05 was considered to be statistically significant.

The results were verified in a local validation dataset. We used frozen tissue specimens which were supplied by Zhejiang Cancer Hospital Biospecimen Repository. PCR array was conducted to detect the ATG expression profile. In total, 86 eligible patients with serous ovarian cancer with detectable and qualified PCR array results were included in the analyses. All of these 86 patients with ovarian cancer were diagnosed, underwent primary debulking surgery and completed the whole chemotherapy regimens at Zhejiang Cancer Hospital, Hangzhou, China. The time of first surgery ranged from May 2008 to December 2013. The chemotherapy regimens for all the patients comprised paclitaxel plus cisplatin or carboplatin. This study was conducted in accordance with the Declaration of Helsinki. The study was approved by the local ethical committees of Zhejiang Cancer Hospital and the Zhejiang Cancer Hospital Biospecimen Repository. Written informed consent was obtained from all of the patients. In addition, the gene expression profile (Affymetrix Human Genome U133 Plus 2.0 Array) and OS information of a second validation dataset was acquired from Gene Expression Omnibus (GEO) database (accession number GSE9891). The original study was designed to identify molecular subtypes of ovarian tumor by using microarrays to examine the expression profile of 285 randomly selected ovarian samples from the AOCS (Australian Ovarian Cancer Study; ref. 37). A total of 135 patients who were in accordance with the criteria of this study were selected and adopted as the second validation dataset. Because of insufficient data, we only performed the validation of the main outcomes in this dataset. Of note, all of the patients in the 3 databases that were included in this study were with the histologic types of serous ovarian cancer.

The demographics and major clinical characteristics of the patients in the TCGA training dataset and the local and GEO validation datasets are listed in Table 1. On the basis of the TCGA training dataset, we screened 43 pivotal ATG genes and constructed the ATG scoring system. The HRs, P values, and the median expression levels in the training dataset of the 43 pivotal ATG genes are listed in Table 2. The expression profiles of the selected ATG genes in the training dataset, local validation dataset, and GEO validation dataset are shown in Supplementary Fig. S1. According to the HR and relative expression level of each gene, the final scores of all the patients in the 3 datasets were calculated. The score ranges in TCGA, local, and GEO dataset were (−28) to 26, (−19) to 20, and (−21) to 19, respectively.

In the TCGA training dataset, we observed that the OS of patients with high scores was poor (low [≤0] vs. high [>0] score, median OS = 56.6 months [95% CI, 50.9-62.3 months] vs. 30.6 months [95% CI, 20.9–40.3 months]; P < 0.001) (Fig. 1A), which was verified in our local validation set (low [≤0] vs. high [>0] score, median OS = 65.0 months [95% CI, 57.0–73.0 months] vs. 48.0 months [95% CI, 34.1–61.9 months], P = 0.027) (Fig. 1B) and the GEO dataset low [<−5] vs. high [≥−5] score, median OS = 64.4 months [95% CI, 38.0–90.8 months] vs. 37.3 months [95% CI, 30.2–44.4 months], P = 0.044) (Supplementary Fig. S2). Such results were consistent with our expectation, which can also be naturally deduced according to the scoring algorithm. Thus, we identified a negative role of high ATG score in advanced ovarian cancer patients treated with chemotherapy. Since the surgical outcome plays a vital role in the survival of ovarian cancer patients, subgroup analyses were performed to identify the effects of ATG score on the patients with different sizes of residual diseases in the TCGA dataset (Supplementary Fig. S3). We noticed that the higher ATG score predicts poor survival regardless of the residual size.

For the purpose of evaluating the independent prognostic significance of the ATG score in ovarian cancer, univariate Cox regression analyses and multivariate Cox regression analyses were then conducted, in which other classic clinical characteristics were included and all the indicators were treated as category variables. In the univariate analysis (Fig. 1C), the result that high scores were associated with poor prognosis was verified in both datasets (TCGA dataset, low [≤0] vs. high [>0] score, HR of death = 0.41 [95% CI, 0.26 to 0.65]; P < 0.001; local dataset, low [≤0] vs. high [>0] score, HR of death = 0.47 [95% CI, 0.24 to 0.93]; P = 0.031). We found that age, pathologic grade and FIGO stage were all outperformed by the ATG score. In the multivariate analysis (Fig. 1D), ATG score was the only variate that significantly correlated to OS in TCGA dataset (low [≤0] vs. high [>0] score, HR of death = 0.27 [95% CI, 0.16 to 0.47]; P < 0.001). As for our local validation dataset, although age showed prognostic significance as well, the ATG score remained to be the most powerful variate with a smallest P value (low [≤0] vs. high [>0] score, HR of death = 0.41 [95% CI, 0.21 to 0.84]; P = 0.014). Therefore, both UVA and MVA Cox regression model supported the results of Kaplan-Meier analyses and demonstrated a strong prognostic performance of the ATG score. Furthermore, we performed UVA and MVA to determine the significance of ATG score when treated as a continuous variable (Supplementary Table S1 and Supplementary Table S2). Consistently, high ATG score was a hazardous factor for prognosis in the univariate Cox regression analyses of both datasets (TCGA dataset, HR of death = 1.050 [95% CI, 1.028 to 1.074], P < 0.001; local dataset, HR of death = 1.064 [95% CI, 1.023 to 1.107], P = 0.002). In MVA, age and ATG score were treated as continuous variables and FIGO stage and pathologic grade were categorized into detailed subtypes. The ATG score was still the only indicator that showed statistical significance in OS assessment and the negative role of high ATG score was also confirmed [TCGA dataset, HR of death = 1.063 (95% CI, 1.036 to 1.090), P < 0.001; local dataset, HR of death = 1.076 (95% CI, 1.029 to 1.126), P = 0.001].

Considering that the breast cancer data were also used in the ATG score construction, we examined the prognostic value of the ATG score in breast cancer patients from TCGA as well. The breast cancer patients with high ATG scores (25.8%) had significantly poor OS (P = 0.014), which suggested the potential of the ATG score to be used in breast cancer prognosis prediction (Supplementary Fig. S4).

After assessing the performance of the ATG score in OS prediction, we examined its performance in the prediction of the RFS. Similar to the results in OS analyses, we observed that poor RFS was associated with high ATG score in TCGA dataset [low (≤0) vs. high (>0) score, median RFS = 26.6 months (95% CI, 13.6–39.6 months) vs. 18.2 months (95% CI, 12.3–24.2 months); P = 0.006; Fig. 2A] and in local dataset [low (≤0) vs. high (>0) score, median RFS = 26.0 months (95% CI, 3.3–48.7 months) vs. 20.0 months (95% CI, 16.1–23.9 months); P = 0.038; Fig. 2B). After response to primary chemotherapy, namely complete or partial remission, the rate of recurrence was 47.2% in the low-score group and 68.3% in the high-score group (odds ratio, 0.42; 95% CI, 0.18–0.92; P = 0.03) in the TCGA dataset, as shown in Table 1. In the local dataset, the rate of recurrence after response was 55% in the low-score group and 82.9% in the high-score group (odds ratio, 0.25; 95% CI, 0.07–0.82; P = 0.01).

Furthermore, ATG score was the only variate that significantly associated with RFS in the UVA [TCGA dataset, low (≤0) vs. high (>0) score, HR of recurrence = 0.47 (95% CI, 0.27–0.82); P = 0.007; local dataset, low (≤0) vs. high (>0) score, HR of recurrence = 0.56 (95% CI, 0.32–0.99); P = 0.045; Fig. 2C] and MVA [TCGA dataset, low (≤0) vs. high (>0) score, HR of recurrence = 0.29 (95% CI, 0.15–0.57); P < 0.001; local dataset, low (≤0) vs. high (>0) score, HR of recurrence = 0.55 (95% CI, 0.31–0.97); P = 0.038; Fig. 2D). Therefore, such findings revealed a critical role of the ATG score in prediction of RFS in the patients with advanced ovarian cancer. We next examined whether the ATG score was associated with the primary chemotherapy outcomes. In TCGA dataset, we found that the proportions of complete remission, partial remission, stable disease, and progressive disease were quite similar in the patients with high scores and low scores (Supplementary Table S3). Hence, it was a weakness of the ATG score that lack of the ability to predict the primary chemotherapy remission rate.

Based on our local dataset, we next constructed a prognostic nomogram combining our ATG score and several key clinical indicators, namely ascites, FIGO stage and CA125 level (Fig. 3A). The C-index of the nomogram was 0.79 (95% CI, 0.72–0.86; P < 0.001). To test the significance of ATG score, we compared the ROC curves and C-index of different prognostic models in local dataset (Fig. 3B). We found the C-index of ATG score (0.66; 95% CI, 0.57–0.75; P < 0.001) was significantly higher than that of FIGO stage (0.56; 95% CI, 0.49–0.63; P = 0.093). Of note, the C-index of the model only including ascites, FIGO stage, and CA125 level was 0.68 (95% CI, 0.60–0.76; P < 0.001), which was much lower than that of the nomogram. Because such a model was actually the result of excluding ATG score from the nomogram, it suggested that the ATG score dramatically improved the discrimination efficacy of the traditional clinical indicators. Similar results were acquired when the IDI and NRI were calculated to compare the predictive ability of the risk prediction models before and after the ATG score was included. Comparing with the model of ascites, FIGO stage and CA125 level, we noticed that for addition of ATG score, the IDI was 0.152 (95% CI, 0.041–0.279; P = 0.007) and the continuous NRI was 0.476 (95% CI, 0.013–0.697; P = 0.027), both of which showed a significant improvement of the risk prediction model after the ATG score was included. For the lack of clinical data, the ROC curves and C-index were only determined for ATG score (C-index = 0.69; 95% CI, 0.62–0.76; P < 0.001) and FIGO stage (C-index = 0.52; 95% CI, 0.46–0.58; P = 0.41) in TCGA dataset (Fig. 3C), which corroborated that the ATG score has superiority over the traditional clinical indicators. In the GEO dataset, ROC curve of the ATG score is shown in Supplementary Fig. S5 (C-index = 0.56; 95% CI, 0.51–0.61; P = 0.028).

With the difficulty of early diagnosis and the limitation in treatment, the clinical management of ovarian cancer is still not satisfactory. Because of the possibility of chemotherapy resistance, an early and accurate prediction of prognosis may inform clinical decisions. Therefore, we screened a set of pivotal ATG genes and established a novel scoring system for ovarian cancer prognosis prediction. Such a scoring system was proven to be capable as a predictive tool for OS and RFS in multidatasets. Notably, the ATG score showed superiority compared with several classic clinical indicators. Thus, we constructed a robust prognostic nomogram that incorporated the ATG score and clinical indicators which are commonly referenced in ovarian cancer. In the comparison of different prognostic models, we found that in addition to serving as a promising independent prognosis predictor, the ATG score could dramatically elevate the C-index of the model consisting of traditional indicators CA125, ascites, and FIGO stage, which indicates that the ATG scoring system and testing of CA125 level could work independently and synergistically with each other in prediction of ovarian cancer prognosis. Thus, the introduction of ATG score will contribute to establishing a diversified prognosis assessment system to optimize chemotherapeutic regimens. Moreover, the results of TCGA training dataset which mostly comprised American and European population (38) were well verified in our local Asian validation dataset. Validation in multiple nationalities strengthened the reliability of our work and suggested the potential of applying the ATG score to patients with different ethnic backgrounds (39–41).

Immunoscores and DNA damage-repair pathway scores have been commonly used for outcome assessment in cancer (20–23). However, there are few studies that have constructed a scoring system based on the ATG genes (42, 43). Recently, although An and colleagues (42) reported an autophagy-related signature consisting of eight genes (BLOC1S1, IL24, NRG4, PDK4, PEX3, PRKG1, SIRT2, and WDR45L) could act as an independent prognostic indicator in serous ovarian cancer, the genes were not strictly autophagy-related and the key executors of autophagy were not included. In contrast, we used a comprehensive but more refined gene set (35) for the primary candidate selection, and thus created a novel scoring system on the basis of screened pivotal genes in autophagy-related pathways. For instance, MAP1LC3B acted as a hazardous gene whereas SQSTM1 played a protecting role in our ATG score. Both of these 2 genes are well known as the core ATG genes in autophagy and the corresponding roles of them are supported by previous studies as well (44, 45). Of note, no overlap was observed between the component genes of our ATG score and the study by An and colleagues (42). Gu and colleagues (43) reported an autophagy-related signature to predict prognosis in breast cancer. Although we used data of breast cancer patients to construct the ATG score, there was no overlap with an autophagy-related signature previously reported by Gu and colleagues. Given our focus on core ATG genes, our scoring system likely has a closer relationship with autophagy and may manifest novel, autophagy-specific significance in ovarian cancer. Such traits might also be helpful to the identification of the patients that could benefit from autophagy targeting therapeutics. Although in this study the evidence is not enough to draw a conclusion on the relationships between the ATG score and potential effects of autophagy targeting therapy, we consider this issue intriguing for future exploration. Furthermore, validation in the local patient data greatly increased the credibility of our ATG score. We expect large-scale clinical studies to further investigate the significance of ATG score in the future.

With 43 pivotal ATG genes, we constructed a novel scoring system, which is highly autophagy-related and also powerful in the prognosis prediction of patients with advanced ovarian cancer who received chemotherapy. However, we admit that this study has some limitations and weaknesses. The strategies of gene selection and signature construction are the most important and fundamental parts for a scoring system. We screened 43 pivotal ATG genes and compiled all of them into the ATG score. Such a strategy was proven thorough and powerful in this study, whereas it is also possible to build a more refined signature with further screening within the 43 genes, which may improve the score in a more feasible and economical way in clinical practice. Another limitation is the relatively small population size. Future study with larger scale may be able to improve the score structure and find more potential significance of the ATG score, such as the primary chemotherapy remission, which unfortunately cannot be predicted by the present ATG score. We assume that the ATG score has a closer relationship with the long-term survival of the tumor cells which survived primary chemotherapy than the sensitivity to chemotherapeutic agents. In addition, it is a pity that although the ATG score showed significance in predicting OS in TCGA patients with breast cancer, local breast cancer samples are currently unavailable for us to validate such an intriguing finding. We look forward to future explorations which may discover more application scenarios of the ATG score.

Previous studies have demonstrated the association existed in autophagy, immunotherapy, and DNA damage (46–49). As for ovarian cancer, BRCA1/2 mutational status has a robust influence on the patients' outcomes, especially for the sensitivity of platin therapy (34, 50). As a major player in the homologous recombination repair pathway, BRCA1 and BRCA2 proteins are critical molecules for high-fidelity DNA double-strand breaks repair (51). It is intriguing to combine the BRCA mutational status with our ATG score, so that the interactive effects of DNA repair and autophagy on ovarian cancer might be further explored. However, such information was temporarily not accessible in this study. Thus, further research into combining the vital players in multiple biological processes may have the potential to provide more powerful tools in prognosis prediction and clinical assessment.

This study successfully established a novel scoring system based on the expression profiles of pivotal ATG genes to accurately assess the outcomes after chemotherapy in patients with advanced ovarian cancer. Our results were well confirmed in the datasets including different races and showed strong capability of prediction. Therefore, autophagy-related pathways may provide critical biomarkers which help for clinical evaluation of patients with ovarian cancer and our novel scoring system may have a satisfactory performance in clinical practice and make contributions to the individualized therapeutic regimen design.

No potential conflicts of interest were disclosed.

Conception and design: Y. Niu, K. Chen, Z. Fu, H.-M. Shen, D. Xia, Y. Wu

Development of methodology: Y. Niu, W. Sun, K. Chen, H.-M. Shen, D. Xia, Y. Wu

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): Z. Fu, Y. Chen, J. Zhu

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): Y. Niu, W. Sun, K. Chen, H. Zhang, Y. Wu

Writing, review, and/or revision of the manuscript: Y. Niu, K. Chen, H. Chen, Y. Shi, L. Wang, H.-M. Shen, D. Xia, Y. Wu

Study supervision: D. Xia, Y. Wu

The authors thank the Zhejiang Cancer Hospital Biospecimen Repository for providing the frozen tissue specimens. This work was supported by Zhejiang Provincial Natural Science Foundation of China (grant nos. LY18C060001, LY15H160007); National Natural Science Foundation of China (grant nos. 81773016, 31471297, and 81302455); WeiJian Special Foundation, Zhejiang University School of Public Health (to Y. Wu); and the Fundamental Research Funds for the Central Universities (to Y. Wu).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.