Increased ATG5 Expression Predicts Poor Prognosis and Promotes EMT in Cervical Carcinoma

Cervical cancer has the second-highest incidence and mortality of female malignancy. The major causes of mortality in patients with cervical cancer are invasion and metastasis. The epithelial–mesenchymal transition (EMT) process plays a major role in the acquisition of metastatic potential and motility. Autophagy-related genes (ARGs) are implicated in the EMT process, and autophagy exerts a dual function in EMT management at different phases of tumor progression. However, the role of specific ARGs during the EMT process has not yet been reported in cervical cancer. Based on the data from the Cancer Genome Atlas (TCGA) cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC) sequencing database, we performed the prognosis analysis for those ARGs obtained from the Human Autophagy database. ATG5 was identified as the only important harmful marker influencing survival of cervical cancer patients by univariate Cox regression (HR 1.7; 95% CI: 1.0–2.8, p = 0.047), and the 5-years survival rate for the high- and low-ATG5 expression groups was 0.486 (0.375–0.631) and 0.782 (0.708–0.863), respectively. TCGA CESC methylation data showed that eight methylation sites of ATG5 could also be significantly associated with the overall survival (OS) of cervical cancer patients. Single-sample gene-set enrichment and gene functional enrichment results showed that ATG5 was correlated with some cancer-related pathways, such as phagocytosis-related genes, endocytosis-related genes, immune-related genes, EMT score, and some EMT signature-related genes. Next, cell migration and invasion assay and Western blot were applied to detect the function of ATG5 in EMT of cervical cancer. In cervical cancer cells, ATG5 knockdown resulted in attenuation of migration and invasion. The functional study showed that knockdown of ATG5 could reverse EMT process by P-ERK, P-NFκBp65, P-mTOR pathways, and so on. In conclusion, the present study implies that ATG5 was a major contributor to EMT regulation and poor prognosis in cervical cancer.


INTRODUCTION
Cervical cancer secondary to breast cancer has the second-highest incidence and mortality of cancer in the female with over 604,127 new cases and 341,831 deaths worldwide per year (Sung et al., 2021). Cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC) are the main pathological types of cervical cancer. Squamous cell carcinomas account for about 75% of invasive cervical carcinoma cases. The main treatment methods for cervical cancer patients include surgery or concurrent radiotherapy and chemotherapy regimens, such as chemotherapy concurrent with radiotherapy and afterloading brachytherapy (Small et al., 2017). Clinical trials including the adoptive T-cell therapy, checkpoint inhibitors, and human papillomavirus (HPV) vaccines have shown satisfactory results. Based on these treatments, 3-years local control rate of patients with early-stage CESC is above 87%. The 5-years survival rate is 16.5%, and the median survival time is only 8-13 months for metastatic cervical cancer because of no standard treatment due to their heterogeneous manifestations (van Meir et al., 2014).
In most cases, the patient has already advanced to the moderate or late stages when first diagnosed. Therefore, invasion and metastasis are the major characteristics and causes of mortality in patients with advanced cervical cancer. Cancer cells acquire metastatic potential and gain motility mostly through the epithelial-mesenchymal transition (EMT) process. EMT is a cellular process characterized by the transformation from an epithelial phenotype to a mesenchymal phenotype (Nieto et al., 2016). This process plays a significant role in cancer initiation, invasion, metastasis, tumor immunosuppression, and immune evasion (Nieto et al., 2016;Jiang and Zhan, 2020). Targets associated with EMT regulation preventing cancer invasion and metastasis could be found to improve the prognosis of cervical carcinoma patients.
Novel increasing evidence showed that the process of EMT could be promoted or inhibited by autophagy (Catalano et al., 2015;Dash et al., 2018;Daskalaki et al., 2018;Wang et al., 2019;Wu et al., 2019;Liu et al., 2020). Autophagy is a catabolic procedure involved in maintaining homeostasis by metabolizing intracellular proteins and organelles under various conditions of cellular stress (Dikic and Elazar, 2018). There was a complicated regulatory network between the regulation of EMT and autophagy. EMT-involved oncogenic signal proteins such as SNAI1, SLUG, ZEB1/2, and NOTCH1 were functionally related to autophagy (Zada et al., 2021). Autophagy plays a tumorinhibitory effect in the initial stages of carcinogenesis and a tumor-promotion function during the cancer progression (White, 2012). By its controversial role in tumor, the acquisition of EMT phenotype required autophagy activation during cancer progression and could be inhibited by autophagy inducers in the early stage of tumor formation (Hu et al., 2018;Wang et al., 2020). ARGs, especially the core mammalian autophagy proteins including ATG2-10, 12-14, 16-18, were a series of functional proteins implicated in the regulation of autophagy. However, the role of specific ARGs during the EMT process of cervical cancer has not yet been reported. DNA methylation is an epigenetic modification process that affects the expression of many genes, contributing to genome stability and the regulation of gene transcription (Moore et al., 2013). Methylation genes such as cadherin 1 (CDH1), death-associated protein kinase 1 (DAPK1), telomerase reverse transcriptase (TERT), and cell adhesion molecule 1 (CADM1) have been found as the highest methylation frequency genes (Wentzensen et al., 2009) and identified as promising markers for prognostic prediction of cervical cancer patients (Overmeer et al., 2011;Abudukadeer et al., 2012;Sun et al., 2015). Although some ARGs have been identified to prospectively predict the prognosis in patients with squamous cell cervical cancer (Chen et al., 2020;Shi et al., 2020), the DNA methylation of ARGs has not been analyzed in the prediction of cervical cancer prognosis.
In this paper, we aimed to find the clinical prognostic values of essential ARGs in cervical cancer based on the TCGA database. We found that a high level of ATG5 expression was an important poor prognostic factor in cervical cancer. In addition, further studies were performed to verify the role and the potential mechanisms of ATG5 in promoting migration and invasion in cervical cancer cell lines. Our results suggested that ATG5 acted as a poor prognostic marker and may be applied as a potential target for reversing the invasion and metastasis in cervical cancer.

Data Source and Preprocessing
RNA-seq reads data, methylation 450 k data, and clinical data of patients with cervical carcinoma from TCGA were acquired from Broad Institute Firehose (https://xenabrowser.net/ datapages/). ARGs were obtained from the Human Autophagy database (HADb http://www.autophagy.lu/project. html). Samples without clinical survival status follow-up data were excluded from the study. Overall, a total of 296 CESC samples were enrolled in the present study. The 31 core ARGs involved in the further analysis are described in Table 1. These ARGs, identified as important regulators of autophagy, were essential for autophagosome formation and were involved in the further result analysis. Rstudio and R programming(4.03) were the principal tools for analyzing data throughout the study.

Survival Analysis of ARGs in CESC
The optimal cutoff values of ARGs based on the expression Fragments Per kilobase per Million mapped reads (FPKM) value or methylation Beta-value, the survival time, and the survival status were identified by the "surv_cutpoint" function in "survminer" package as described before (Wang et al., 2021). The cutoff value of the OS significant genes is listed in Supplementary Figure S1. The patients with relative values above or below the optimal cutoff were considered as high or low groups for each gene expression level, methylation level, or EMT score. Kaplan-Meier (KM) survival analysis with log-rank test was then used for investigating the OS difference between the abovementioned high and low groups. The correlation between the OS characteristics of CESC patients and each ARG expression was explored by a univariate Cox regression analysis. A multivariate Cox regression analysis was performed to investigate whether each ARG is an independent OS influence factor.
The patients with a relative value above or below the optimal cutoff were considered as high or low group for each gene expression level, methylation level, or EMT score. KM survival analysis with log-rank test was then used for investigating the overall survival (OS) difference between the abovementioned high and low groups.

EMT Score Evaluation
By "GSVA" R package, the EMT score for evaluating the EMT status of each CESC patient was quantified by single-sample gene-set enrichment analysis (ssGSEA) based on previously reported EMT signature genes (Choi et al., 2010;Hanzelmann et al., 2013). The correlation between ATG5 and EMT score or between ATG5 and EMT signature related genes is calculated by Pearson's correlation coefficient and plotted in heatmap by the "circlize" R package. The EMT score density between normal and tumor tissues and the survival difference between high and low groups were also analyzed.

ATG5 Co-Expression Network Selection and Gene Functional Enrichment Analysis
The genes co-expressed with ATG5 in CESC were screened according to the Pearson correlation coefficient (|cor| > 0.3, p < 0.05) and plotted in the volcano plot by the "ggplot2" package. Heatmaps of the top 50 negative or the top 50 positive ATG5 expression correlated genes were plotted by the "pheatmap" R package. Metascape (Zhou et al., 2019) was employed to gain insights into the biological functions of those co-expressed genes of ATG5.

Identification Prognostic Methylation Sites of ATG5
The methylation site information of ATG5 in TCGA CESC methylation data based on Illumina Infinium Human Methylation450 BeadChip platform was mapped by an annotation file retrieved from the hg19 GPL16304 legacy annotation file (Price et al., 2013). The methylation level of each site in ATG5 was grouped into low and high by the survminer package as mentioned above. The OS significant associated sites were screened through KM survival analysis. Besides, MEXPRESS was used for visualizing expression, DNA methylation, and clinical parameters based on TCGA data (http://mexpress.be). The correlation between the methylation level of each site of ATG5 and the expression of ATG5 is calculated by the "corrplot" package (Wei and Simko, 2021).

Gene Mutation Analysis of ARGs in CESC
Waterfall plots of ARGs mutation were produced by the oncoplot function in the "maftools" R package (Mayakonda et al., 2018) to visualize the mutation status such as classification and frequency of mutation types, frequency of variant types, and frequency of SNV classes for CESC samples. The mutational exclusion and cooccurrence mutations among those ARGs were analyzed by Fisher's exact test (Mayakonda et al., 2018). The variant allele frequency (VAF) distribution of ARGs was also plotted.

Cell Lines
Hela, Ca Ski (human cervical cancer cell lines), and PANC-1 cells (pancreatic cancer cell line) were purchased from the American Type Culture Collection (ATCC, Manassas, VA, United States). Hela cells were incubated in DMEM supplemented with 10% FBS (Gibco, Grand Island, NY). Ca Ski and PANC-1 cells were incubated in RPMI 1640 containing 10% FBS (Gibco, Grand Island, NY). PANC-1 cells were used as a positive control for ATG5 expression . Cells are cultured in an incubator at 37°C with 5% CO 2 .

Cell Transfection
According to the manufacturer's instructions (GeneChem Corporation, Shanghai, China), Hela and Ca Ski cells were transfected with the recombinant lentiviral vector LV-ATG5-RNAi labeled with green fluorescent to establish ATG5knockdown cells and named as KD. In brief, the cells were digested at the logarithmic growth phase, then resuspended in the complete medium at a concentration of 3-5 ×10 4 cells/ml, and seeded in six-well plates. At 20% confluence, cells were transfected with the recombinant lentivirus for 16 h. Subsequently, a conditioned medium containing transfection reagent was refreshed. At 72 h post-infection, the positive expression rate of GFP was used to evaluate the transfection efficiency under a fluorescence microscope (Olympus Corporation, Tokyo, Japan). siRNAs specific for ATG5 were obtained from GenePharma (Shanghai, China), and the detailed sequences were as follows: ATG5-RNAi (98917-1) sense: CCTTTCATTCAGAAGCT GTTT, ATG5-RNAi (98918-1) sense: CCTGAACAGAA TCATCCTTAA, ATG5-RNAi (98919-1) sense: GATTCATGGAATTGAGCCAAT, NC-RNAi sense: TTCTCCGAACGTGTCACGT.

Western Blot Assay
In short, cells were harvested and washed with PBS twice, and then lysed in RIPA buffer (Beyotime, Shanghai, China) on ice for 15 min. Then, BCA Protein Assay Kit (Beyotime, Shanghai, China) was used to detect the protein concentration in the supernatant from crushed cells by ultrasonication. The final concentration was adjusted to 2 μg/μl by fresh RIPA solution, followed by protein denaturation. 10% SDS-PAGE was applied to separate the protein samples and then transferred onto PVDF membranes (Millipore, Burlington, MA). After that, membranes were incubated with TBST containing 5% skim milk to block non-specific antigens. Next, the membranes were incubated with primary antibodies and HRP-conjugated second antibodies in turn. Finally, the protein bands were visualized by the electrochemiluminescence (ECL) detection system. The primary antibodies used include mouse anti-ATG5 antibody  United States) at a density of 1 × 10 5 cells/well. In invasion assay, serum-free medium was added into the lower chamber. In the migration assay, a medium containing 30% FBS was added to the lower chamber. The noninvasive/non-migrated cells were removed by a cotton swab after 48-h incubation. Then, the cells were fixed in 4% paraformaldehyde for 0.5 h and stained with 0.5% crystal violet solution, followed by taking photographs through the microscope.  Figure (N, O) According to the forest plot for the Univariate and Multivariate Cox regression analysis results, ATG5 was a risk prognostic factor, while ATG4D, GABARAP, ATG4A, and ATG3 were favorable prognostic factors in CESC patients. And, ATG5, ATG4D, and ATG4A were independent prognostic factors.

Cell Migration and Invasion Assays
Frontiers in Cell and Developmental Biology | www.frontiersin.org November 2021 | Volume 9 | Article 757184 5

Statistical Analysis
All molecular biology assays were repeated in triplicate. Experimental data are presented as the mean ± SD. R (version 4.03) and Rstudio software were applied to perform statistical analysis. The differences between the groups were analyzed via Student's t-test and Mann-Whitney U test. Univariate and multivariate Cox regression analysis were performed to identify the factors influencing prognosis. Survival analysis was done using the KM method and the log-rank test. p < 0.05 was designated as statistically significant.

Baseline Information of the CESC Patients
A total of 309 sequencing data and 312 methylation data in TCGA CESC cohort were obtained from the TCGA data. The median age of those 306 participants was 46.57 years (min 20.95 years and max 89.04 years). The median follow-up time for CESC patients with survival information was 1.87 years (range 0-17.55 years), and among them, 225 were alive and 74 were dead during the follow-up time.

Prognosis of ARGs in CESC
High expression levels of ATG4C, ATG5, BECN1, and WIPI1 were associated with poor prognosis in CESC, while CESC patients with overexpression of ATG3, ATG4A, ATG4B, ATG4D, ATG7, ATG9B, ATG13, GABARAP, and GABARAPL2 have longer OS and higher 5-years survival rate ( Figure 1). Furthermore, the univariate Cox regression results indicated that ATG5, ATG4D, GABARAP, ATG4A, and ATG3 were significantly associated with OS ( Figure 1N). ATG5 was identified as a risk prognostic predictor in CESC patients (HR 1.7; 95% CI, 1.0-2.8, p 0.047), while another four ARGs were identified as favorable prognostic factors (HR < 1, p < 0.05). Subsequently, the multivariate Cox results showed that ATG5, ATG4D, and ATG4A were independent prognostic factors affecting OS of patients with CESC ( Figure 1O). Among them, high expression of ATG5 or low expressions of ATG4D and ATG4A were associated with unfavorable prognosis. These results indicated that ATG5 among these ARGs was the most harmful factor affecting the prognosis of CESC patients. According to the KM analysis, the 5-years survival rate for high-and low-expression groups of ATG5 is 0.486 (0.375-0.631) and 0.782 (0.708-0.863) respectively. Furthermore, KM survival analysis was also performed to show that ATG5 expression could effectively predict the prognosis of patients with CESC with different clinical stages, pathological grades, and ages at diagnosis (Figures 2A-F). For patients with cervical squamous cell carcinoma, patients with ATG5 overexpression have significantly shorter OS than the ATG5 down-expression group ( Figure 2G). For cervical adenocarcinoma, due to the limited sample size of patients included in the analysis, the p-value of the km survival curve is not significant, but patients with ATG5 down-expression FIGURE 3 | Co-expression genes correlated with ATG5 in cervical cancer (LinkedOmics) (A) Pearson test was applied to analyze the correlation between ATG5 and differentially expressed genes in Cervical cancer. Blue dots (R < −0.3, p < 0.01) represent as negatively correlated genes, while red dots (R > 0.3, p < 0.01) represent as positively correlated genes (B) Heat map showed the top 50 significant genes positively correlated with ATG5 (C) Heat map showed the top 50 significant genes negatively correlated with ATG5 (The expression of each in each heatmap is normalized with Z-score).

FIGURE 4 | GO analysis and correlation analysis (A)
GO analysis revealed the biological processes and molecular functions involved in the top 100 ATG5-related co-expressed genes by "metascape" data (B) The correlation between ATG5 and co-expressed genes involved in the regulation of the immune effector process in cervical cancer patients. Blue lines represent as negatively correlated genes, while red lines represent as positively correlated genes (C-J) Among those regulation of the immune effector process genes, ATG5 expression was significantly negatively correlated with MAP3K7, STX7 and significantly positively related to A2M, C1QA, C1QB, C1QC, RPS19, TYROBP, and IL20RB.
Frontiers in Cell and Developmental Biology | www.frontiersin.org November 2021 | Volume 9 | Article 757184 7 FIGURE 5 | The correlation between ATG5 and co-expressed genes involved in the regulation of EMT in cervical cancer patients. A. the Kaplan-Meier survival plot showed that the median OS time of low EMT score patients was significant longer than high EMT score patients. B. ATG5 showed co-expressed ship with genes involved in the regulation of EMT in cervical cancer patients. C. The EMT signature showed significant difference in cervical normal and cancer tissues. showed a trend towards a longer survival with no death during the follow-up times ( Figure 2H).

Co-Expression Genes Correlated With ATG5 in CESC
Based on the TCGA database, we analyzed the co-expressed genes of ATG5 in 304 CESC patients through the Pearson correlation coefficient. As shown in Figure 3A, red dots represented 5,666 genes positively associated with ATG5, while blue dots represented 2,114 genes negatively associated with ATG5. The heat map was used to show the top 50 genes positively or negatively correlating with ATG5 ( Figures 3B,C). Next, GO enrichment analysis was performed to investigate the biological pathway of the top ATG5-related co-expressed genes. ATG5 co-expressed genes are correlated with ribosome, cytoplasmic, C1 complex, TRBP-containing complex, glucosamine-containing compound catabolic process, and so on. There are concerns that ATG5 also participated in the regulation of the immune effector process and regulation of T-cell-mediated immunity ( Figure 4A). We further found that the ATG5 co-expressed genes involved in the regulation of the immune effector process were A2M, C1QA, C1QB, C1QC, RPS19, MAP3K7, TYROBP, STX7, and IL20RB. MAP3K7, STX7, and IL20RB were involved in the regulation of T-cellmediated immunity. Moreover, ATG5 expression was inversely associated with MAP3K7 and STX7 and positively related to A2M, C1QA, C1QB, C1QC, RPS19, TYROBP, and IL20RB ( Figures 4B-J). Furthermore, it was important to investigate the relation between ATG5 and EMT-related gene signature. The result is represented in Figure 5. We found that ATG5 was significantly correlated with EMT scores quantified by ssGSEA FIGURE 7 | Genomic methylation of ARGs in CESC (A) The methylation of ATG5 gene was significantly correlated with the OS in cervical cancer patients (B) The correlation between methylation sites of ATG5 and ATG5 expression in cervical cancer patients. Purple represents negative correlation; red represents positive correlation; "×" represents no statistical significance (C-J) Kaplan-Meier survival analysis showed the methylation sites of ATG5 significantly related to OS in cervical cancer patients (p < 0.05). cg00333154, cg15336269, cg02711268, and cg21148531 located in ATG5 were beneficial for OS, while cg24211478, cg10566763, cg26178529, and cg10850197 located in ATG5 were harmful for OS.

Genomic Alterations and Methylation of ARGs in CESC
Genetic alteration results showed that 16.15% of patients own mutations of these ARGs. Among these genes, a total of 16 genes have a mutation rate ≥1%, of which ATG2B and RB1CC1 are the most frequently mutated genes (4%). Missense and splice are the two most common types of mutations ( Figure 6A). ATG5 was altered in 3 (1.03%) of the 291 CESC patients. Figure 6B shows the mutation type and domain of ATG5. Moreover, concomitant occurrence of ATG5 mutation and ULK2 mutation was found in CESC ( Figure 6C). Compared with most mutated ARGs, ATG5 owns a relatively low VAF, except ATG10. MEXPRESS was used for visualizing expression, DNA methylation, and clinical parameters based on TCGA data, and showed the dramatic relationship between ATG5 methylation status and OS in CESC patients ( Figure 7A). According to the annotation file of Methylation450 BeadChip platform, there are 19 methylation sites of ATG5, namely, cg24211478, cg14951955, cg00333154, cg07541085, cg15336269, cg07747241, cg01786360, cg16950012, cg22876713, cg10566763, cg26178529, cg01644481, cg02711268, cg14660784, cg18206,736, cg20059537, cg21148531, cg04389950, and cg10850197. As shown in Figure 7B, there was a negative correlation between ATG5 expression and the methylation sites cg02711268, cg18206,736, cg07541085, cg24211478, cg26178529, and cg10850197. The significant sites of those methylation sites in ATG5 on OS in CESC patients were shown in KM including beneficial OS relevance methylation sites such as cg00333154, cg15336269, cg02711268, and cg21148531 and the harmful OSrelated sites such as cg24211478, cg10566763, cg26178529, and cg10850197 ( Figures 7C-J).

The Regulation of ATG5 Expression in Cervical Carcinoma Cells
ATG5 mRNA showed higher expression in cervical cancer cell lines Hela and Ca Ski compared with the cervical epithelial cell line CRL2614 ( Figure 8A). The recombinant lentiviral vector containing ATG5-RNAi labeled with green fluorescence was used to downregulate ATG5 expression of Hela and Ca Ski cells. Cells presented with green fluorescence indicated the success of transfection ( Figure 8B). RT-qPCR showed that the levels of ATG5 mRNA were significantly downregulated in Hela-KD (knockdown) and Ca Ski-KD cells compared with NC (negative control) cells ( Figure 8C, p < 0.01). In addition, the sequence with the highest knockdown efficiency was selected for use in the next study. The knockdown efficiency was further confirmed by Western blot ( Figure 8D).

Knockdown of ATG5 Impedes Cervical Cancer Cells Migration and Invasion
By Transwell migration assay, ATG5 knockdown demonstrated a markedly negative effect on the migratory capacity of Hela and Ca Ski cells compared with NC ( Figure 9). Subsequently, the role of ATG5 in invasion was assessed by a Transwell invasion assay. Likewise, the data of the Transwell invasion assay showed that downregulation of ATG5 dramatically suppressed the invasion in cervical cancer cells ( Figure 9). In brief, these results exhibited the repression of ATG5 in the migration and invasion of cervical cancer cells in vitro.
Frontiers in Cell and Developmental Biology | www.frontiersin.org November 2021 | Volume 9 | Article 757184 11 β-Catenin, Slug, P-β-Catenin, c-Myc, and Snail were not affected by ATG5 knockdown. Generally, downregulation of ATG5 could reverse the EMT process by specific pathways in cervical cancer cells and resulted in attenuation of migration and invasion.

DISCUSSION
Autophagy is an essential process induced under various conditions of cellular stress and tightly regulated to respond correctly to the ever-changing environment. It is now apparent that autophagy plays different roles in cancer progression. It plays a cancer inhibition role in the onset of tumorigenesis and cancer promotion in the cancer advancement phase (White, 2012). Autophagy inhibition may be an effective anti-cancer therapy in advanced cancer (Barbeau et al., 2020). Invasion and metastasis contribute to a high recurrence rate, which is the main cause of the poor prognosis of cancer patients (Malek et al., 2014). Despite advances in cancer management, some patients with cervical cancer undergoing early metastasis, especially lymph node metastasis, would ultimately lead to poor clinical outcomes (Nanthamongkolkul and Hanprasertpong, 2018;Xu et al., 2018). EMT, defined as a process by which cancer cells change from an epithelial phenotype to a mesenchymal phenotype, promotes invasion and metastasis by enhancing the capability of cancer cells to invade and migrate (Lamouille et al., 2014). Novel increasing evidence showed that the process of EMT could be promoted or inhibited by autophagy (Catalano et al., 2015;Dash et al., 2018;Daskalaki et al., 2018;Wang et al., 2019;Wu et al., 2019;Liu et al., 2020), and the core ARGs involved in the autophagy process plays important roles in EMT regulation (Peng et al., 2013;Liu et al., 2020). In this study, we analyzed the prognosis of ARG signature in cervical cancer patients. By univariate and multivariate analysis, ATG5, ATG4D, and ATG4A were distinguished as independent prognostic predictors from the TCGA database. It was vital that ATG5 was the only important harmful marker influencing the prognosis of cervical cancer patients. ATG5 is one of the essential regulators of autophagy. Regardless of stage, grade, and pathological characteristics, cervical cancer patients with high expression of ATG5 had shorter survival. Due to the limited sample size (29 cases), there was no statistical significance in the adenocarcinoma group. However, the KM curve of the high-ATG5 group was significantly lower than that of the low-ATG5 group. Our results suggested that ATG5, as one of ARG signatures, was an essential contributor to the poor prognosis of cervical cancer. The higher expression of ATG5 was also associated with poor clinical outcomes of other cancers (Yang et al., 2017;Cheng et al., 2019). Conversely, positive expression of ATG5 predicts a favorable prognosis in patients with breast cancer and osteosarcoma Zhao et al., 2018).
Besides the expression profile of ATG5, the hereditary genetic polymorphisms of ATG5 were also recognized as prognostic predictors of early-stage ESCC patients (Yang et al., 2017). Therefore, we investigate the genomic alterations and methylation of ARGs in CESC. Eight methylation sites of ATG5, namely, four sites with good prognosis and four sites with poor prognosis, were explored to be significantly associated with OS of cervical cancer patients. Four sites with beneficial OS (cg00333154, cg15336269, cg02711268, and cg21148531 located at chr6:106774171:106809238 FIGURE 11 | The illumination of the effect of ATG5 expression and methylation in cervical cancer. ATG5 among these 32 ARGs was the most significant harmful factor affecting the prognosis of CESC patients. ATG5 co-expressed genes involved in the regulation of the immune effector process were A2M, C1QA, C1QB, C1QC, RPS19, MAP3K7, TYROBP, STX7, and IL20RB. And, MAP3K7, STX7, and IL20RB were in-volved in the regulation of T cell-mediated immunity. Moreover, ATG5 might promote the migration and invasion of CESC by the regulation of N-cadherin, Vimentin, P-ERK, P-NFκBp65, P-mTOR, MMP-9, and Twist. ATG5 was also significantly correlated with EMT signature related genes including DPH3, CD24, KRT8, CTNNAL1, EPCAM, TSPAN13, MYB, TPD52, ECHDC1, CLSPN, USP33, TTK, TOM1L1, GCA, COMMD8, CGAS, and OSTM1. ATG5 methyltaion status also influence CESC OS. Among these 19 methylation sites in ATG5, four sites significantly increased while other four sites significantly decreased OS of cervical cancer patients. Frontiers in Cell and Developmental Biology | www.frontiersin.org November 2021 | Volume 9 | Article 757184 length ≈35 kbp 10 based on illuminaMethyl450_hg19_GPL16304_TCGAlegacy annotation data) and four sites with harmful OS (cg24211478, cg10566763, cg26178529, and cg10850197 located at chr6:106773226:106773720 length 494 bp) were located in different positions of transcription regions of ATG5. The relation between DNA methylation and gene expression is importantly subjected to the genomic location of DNA methylation, which also affects the clinicopathological characteristics of cancers (van Vlodrop et al., 2011). Traditionally, cancer research predominantly focused on the prognostic significance and effects of promoter CpG island hypermethylation located in tumor suppressor genes (Jacinto et al., 2007;Amara et al., 2008;Singh et al., 2020). It is well known that DNA methylation of promoters can silence genes while growing evidence suggests that methylation may also be associated with gene activation (Smith et al., 2020). Moreover, methylation within a gene-body or transcribed region does not block the expression of the affected gene (Salem et al., 2000;Zilberman et al., 2007;Yang et al., 2014). These might explain the different prognostic significance among the eight methylation sites of ATG5. In addition, our analysis about genomic alterations and methylation of ATG5 showed the concomitant occurrence of ATG5 mutation and ULK2 mutation and indicated that methylation sites of ATG5 could also serve as prognosis predictors in cervical cancer patients. As reported by others, the inhibition of autophagy by promoter methylation of ULK2, which leads to downregulation of transcript levels, was essential for cancer growth under the genetic background of ATG5 (Shukla et al., 2014). Moreover, ATG5 was identified as one of the immune-related genes in ESCC (Li et al., 2017). Likewise, our study found that ATG5 was significantly correlated with immune-related genes involved in the immune effector process and regulation of T-cell-mediated immunity. Research has found that ATG5 silencing could abolish the EMT promotion by miR-210-5p in osteosarcoma . Also, in hepatocellular carcinoma, ATG5 downregulation drastically dampened TGF-β2-induced EMT (Dash et al., 2018). By contrast, knockdown of ATG5 increased migration and invasion in glioblastoma cells and RAS-mutated cancer cells by promoting EMT (Catalano et al., 2015;Wang et al., 2019). However, the effect of ATG5 on EMT and the prognosis in cervical cancer is unknown to date. We analyzed the relationship between ATG5 and EMT-related gene signature in cervical cancer. In addition, cervical cancer cell lines were used to investigate the effect of ATG5 on migration and invasion.
ATG5 knockdown could inhibit migration and invasion of cervical cancer cells by reversing EMT. N-cadherin, Vimentin, P-ERK, P-NFκBp65, P-mTOR, MMP-9, Fibronectin, and Twist might be the possible mechanisms involved in ATG5-dependent EMT regulation. However, function acquisition and loss of these genes were needed to be carried out in our further epidemical investigation. MMP-9, as a Matrix metalloproteinase, can destroy the integrity of the basement membrane and extracellular matrix and promote the detachment of tumor cells from the endothelium, which will lead to metastatic process (Voura et al., 2013;Rempe et al., 2016). Fibronectin is one of the ECM components and serves as a mesenchymal marker (Paolillo and Schinelli, 2019). Elevated protein expression of N-cadherin and Vimentin promotes mesenchymal phenotypic transition (Paolillo and Schinelli, 2019). ERK-, NFκBp65-, and mTOR-signaling pathways were the important mechanisms implicated in the regulation of the EMT process (Singh et al., 2018).

CONCLUSION
In summary, ATG5 involved in the EMT process and immune regulation in cervical cancer could affect the survival of cervical cancer patients by expression and methylation level, proposing that ATG5 may be a potentially powerful therapeutic target for cervical cancer (Figure 11).

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
SZ and YX designed research; XW, JD, and HY performed experiments; SZ and YX analyzed the data; YX, SZ and XW wrote the manuscript. All authors contributed to writing and critically revising the manuscript. All authors have read and agreed to the published version of the manuscript.

FUNDING
This study was supported by Program of Taizhou Science and Technology Grant (20ywa04).