Identification of Autophagy Related circRNA-miRNA-mRNA-Subtypes Network With Radiotherapy Responses and Tumor Immune Microenvironment in Non-small Cell Lung Cancer

Lung cancer (LC) is one of the most frequently diagnosed cancers and the leading cause of cancer death worldwide, and most LCs are non-small cell lung cancer (NSCLC). Radiotherapy is one of the most effective treatments for patients with lung cancer, either alone or in combination with other treatment methods. However, radiotherapy responses vary considerably among NSCLC patients. The efficacy of radiotherapy is influenced by several factors, among which autophagy is of importance. Autophagy is induced by radiotherapy and also influences cell responses to radiation. We explored the clinical significance of autophagy-related genes (ARGs) and gene sets (ARGSs) and the underlying mechanism in NSCLC patients treated with radiotherapy. First, differentially expressed ARGs (SNCA, SESN3, DAPL1, and ELAPOR1) and miRNAs (miR-205-5p, miR-26a-1-3p, miR-6510-3p, miR-194-3p, miR-215-5p, and miR-375-3p) were identified between radiotherapy-resistant and radiotherapy-sensitive groups. An autophagy-related radiosensitivity risk signature (ARRS) by nine ARmRNAs/miRNAs and an autophagy-related overall survival risk signature (AROS) by three ARmRNAs were then constructed with estimated AUCs of 0.8854 (95% CI: 0.8131–0.9576) and 0.7901 (95% CI: 0.7168–0.8685), respectively. The correlations between ARGSs or prognostic signatures and clinicopathological factors, short-term radiotherapy responses (radiotherapy sensitivity), long-term radiotherapy responses (overall survival), and immune characteristics were analyzed. Both ARGSs and prognostic signatures were related to immune checkpoint inhibitors (ICIs), infiltration of tumor-infiltrating immune cells (TIICs), and the activity of the cancer immune cycle. Finally, after target prediction and correlation analysis, circRNA (hsa_circ_0019709, hsa_circ_0081983, hsa_circ_0112354, hsa_circ_0040569, hsa_circ_0135500, and hsa_circ_0098966)-regulated miRNA/ARmRNA axes (miR-194-3p/SESN3, miR-205-5p/ELAPOR1, and miR-26a-1-3p/SNCA) were considered potential modulatory mechanisms by influencing the regulation of autophagy, macroautophagy, and chaperone-mediated autophagy.

Introduction

With an estimated 2.2 million new cases and 1.8 million deaths, lung cancer (LC) is one of the most frequently occurring cancer and the leading cause of cancer death according to the most recent global cancer statistics (Sung et al., 2021). In most countries, the 5-year survival rate of patients with LC is only 10 to 20% during 2010 through 2014 (Allemani et al., 2018). To increase the survival rate of patients, improving therapeutic effectiveness is as important as early screening. Radiotherapy is one of the most effective treatments for patients with LC, either alone or in combination with other treatment methods. However, because of individual heterogeneity, radiotherapy responses vary among patients, especially in those with non-small cell lung cancer (NSCLC) (Baker et al., 2016), which accounts for 80% of LC. An important focus of radiation oncology research is to predict radiotherapy responses by using molecular analysis.

Autophagy, a major type of programmed cell death, has been generally regarded as a survival or cytoprotective response under stressful conditions, for example, exposure to radiation and chemicals (Murrow and Debnath, 2013). A growing body of evidence indicates that tumor resistance to anticancer therapies, such as radiotherapy, was often associated with the regulation of autophagy (Sharma et al., 2014; Tam et al., 2017). Although no consensus has been reached about the antitumor or protumor action of autophagy induction, autophagy inhibitors or promoters are potential drug-drug or drug-radiation combinations to promote therapeutic efficacy. Thus, understanding the functional relevance of autophagy within radiotherapy is critical to evade resistance and enhance the effects for NSCLC patients. In addition, few studies have discussed the selective types of autophagy, which are highlighted in our study.

Non-coding RNAs (ncRNAs), accounting for 98% of the human genome, mediate protumorigenic/antitumorigenic responses to different cancer therapies (Zhang X. et al., 2020). MicroRNAs (miRNAs) are a family of small ncRNAs of approximately 22 nucleotides that play an important role in biological pathways by silencing mRNAs and regulating the expression of genes posttranscriptionally (Ambros, 2004). Circular RNAs (circRNAs), another type of ncRNA that can act as gene regulators or even be encoded into proteins, also play vital tumor-regulated roles in numerous cancers (Chen, 2020). Many cases have shown that circRNAs can interact with miRNAs and then form a network to regulate cellular physiological and pathological activities (Hansen et al., 2013). Moreover, due to their relatively stable structure, miRNAs and circRNAs can also be used as biomarkers of cancer therapeutic effects.

In the present study, we made full use of publicly available large-scale cancer omics data, mainly The Cancer Genome Atlas (TCGA), to investigate the clinical significance and underlying mechanisms of autophagy-related genes and gene sets in radiotherapy responses of NSCLC patients. First, patients receiving radiotherapy with complete prognostic information were retrieved, and autophagy-related genes (ARGs) and gene sets (ARGSs) were identified. Then, radiotherapy sensitivity- and overall survival (OS)-related risk signatures were generated following the differential analysis of ARGs and miRNAs. Risk signatures were then subjected to correlation analysis of clinicopathologic factors, predictive value of prognosis, and characteristic analysis of the immune microenvironment. Finally, after targeting prediction and correlation analysis of expression levels, a circRNA-miRNA-ARmRNA-ARGS network was constructed to explain the potential regulatory mechanism.

Materials and Methods

Schematic Diagram of the Study Design

As shown in Figure 1, we first mined the public data for our datasets of interest. NSCLC patient information from The Cancer Genome Atlas (TCGA) project was obtained from UCSC Xena¹. The targeted screening was performed according to the following criteria: (1) patients treated with radiotherapy without additional locoregional surgical procedure; (2) patient primary therapy outcome success and overall survival information was recorded; and (3) patient tumor samples received RNA sequencing (RNA-seq) and/or miRNA sequencing (miRNA-seq). The following four levels of primary therapy outcome were assessed: complete remission/response (CR), partial remission/response (PR), stable disease (SD), and progressive disease (PD). Patients with CR and PR were classified into the radiotherapy-sensitive group, while patients with SD and PD were classified into the radiotherapy-resistant group. Eighty-seven NSCLC patients with RNA-seq and 83 NSCLC patients with miRNA-seq met the requirements (Table 1). Moreover, autophagy-related genes (ARGs) and gene sets (ARGSs) were acquired from the Gene Ontology (GO) resource². The study was then extended to thoroughly investigate the clinical significance and regulatory mechanism of autophagy in the radiotherapy response of NSCLC patients. Differential expression of RNAs and miRNAs was analyzed, and the score of ARGSs was calculated. Clinical correlation and immune microenvironment analysis were then performed at both the gene and gene set levels. Moreover, radiotherapy sensitivity- and overall survival (OS)-related risk signatures were constructed for prognostic prediction. A circRNA-miRNA-ARmRNA-ARGS network was constructed following target prediction and correlation analysis.

FIGURE 1

Figure 1. Flowchart of the study design. CR: complete remission/response, PR: partial remission/response, SD: stable disease, PD: progressive disease, DE: differential expressed, ARmRNA: autophagy-related mRNA, NSCLC: non-small cell lung cancer.

TABLE 1

Table 1. Patients’ characteristics.

Identification and Extraction of Autophagy-Related Genes and Gene Sets

A total of 537 unduplicated autophagy-related genes (ARGs) were extracted from GO:0006914, and 9 autophagy-related gene sets (ARGSs) were identified (Supplementary Figure 1 and Supplementary Table 1). The genes were related to the following types of autophagy: 78 genes were related to autophagy of mitochondrion (GO:0000422); 8 genes were related to autophagy of peroxisome (GO:0030242); 18 genes were related to chaperone-mediated autophagy (GO:0061684); 5 genes were related to late endosomal microautophagy (GO:0061738); 308 genes were related to macroautophagy (GO:0016236); 17 genes were related to autophagy of nucleus (GO:0044804); 7 genes were related to lysosomal microautophagy (GO:0016237); 336 genes were related to regulation of autophagy (GO:0010506); and 9 genes were related to modulation by symbiont of host autophagy (GO:0075071).

Evaluation of the Immune Characteristics of Tumor Microenvironment (TME)

The immune characteristics of TME included the expression level of immune checkpoint inhibitors (ICIs), infiltration of tumor-infiltrating immune cells (TIICs), and activity of the cancer immune cycle. Overall, 20 ICIs (HAVCR2, CD274, CD86, LAG3, LAIR1, PVR, IDO1, CD80, CTLA4, SNCA, TIGIT, CD200R1, CEACAM1, CD276, CD200, KIR3DL1, BTLA, ADORA2A, LGALS3, and VTCN1) with therapeutic potential (Noam et al., 2018) were identified in our study. The infiltration levels of 28 tumor-infiltrating immune cells (activated B cells, activated CD4 T cells, activated CD8 T cells, activated dendritic cells, CD56 bright NK cells, CD56 dim NK cells, central memory CD4 T cells, central memory CD8 T cells, effector memory CD4 T cells, effector memory CD8 T cells, eosinophil cells, gamma delta T cells, immature B cells, immature dendritic cells, macrophages, mast cells, MDSCs, memory B cells, monocytes, NK cells, NK T cells, neutrophils, plasmacytoid dendritic cells, regulatory T cells, T follicular helper cells, TH1 cells, TH17 cells, and TH2 cells) (Charoentong et al., 2017) were considered. The cancer immune cycle mainly comprises the following seven steps: release of cancer cell antigens (Step 1); cancer antigen presentation (Step 2); priming and activation (Step 3); trafficking of immune cells to tumors (Step 4); infiltration of immune cells into tumors (Step 5); recognition of cancer cells by T cells (Step 6); and killing of cancer cells (Step 7) (Chen and Mellman, 2013). The activity of TIICs and the cancer immune cycle were evaluated by calculating marker gene set scores based on the gene expression of individual samples.

Screening of Differentially Expressed Genes

Available RNA sequencing (RNA-seq) and miRNA sequencing (miRNA-seq) data were downloaded. We transformed miRNA-seq names into human mature miRNA names using the miRBase version 22.0 database. We then applied DESeq2, edgeR, and limma/voom to identify differentially expressed mRNAs (DEmRNAs) and miRNAs (DEmiRNAs). The criteria for determining differential DEmRNAs and DEmiRNAs were set with an adjusted p-value < 0.05 and | log fold change (FC)| > mean ± standard deviation (sd). We determined the common DEmRNAs and DEmiRNAs by utilizing the VennDiagram R package (Chen and Boutros, 2011). Volcano plots visually displaying the distribution of DEmRNAs and DEmiRNAs were generated using ggpubr R packages, and heatmaps describing the expression of differentially expressed autophagy-related mRNAs (DEARmRNAs) and miRNAs (DEmiRNAs) were generated utilizing the pheatmap R package.

Establishment of Specific Risk Signatures

We extracted the DEmRNAs and DEmiRNAs expression profiles collected from NSCLC patients receiving radiotherapy with prognostic information. Differential analysis by Student’s t-test was conducted to compare the radiotherapy-resistant and radiotherapy-sensitive groups, while the overall survival (OS) difference was calculated by the log-rank test and described by the K-M curve. The significant variables were included in a logistic or Cox regression model. Finally, we generated an autophagy-related radiosensitivity risk signature (ARRS) and an autophagy-related OS risk signature (AROS) for each sample using the following equation: ARRS or AROS = ∑_i Cofficient(RNAi) × Expression(RNAi). The receiver operating characteristic (ROC) curves with risk score against radiosensitivity and survival status were generated using ROCit and survivalROC/timeROC R packages (Blanche, 2015), respectively. Based on the mean as a cutoff point, patients were divided into high- and low-risk groups. Student’s t-test and log-rank test were used in univariate differential analysis, while multivariate logistic and Cox regression were used in independent predictor tests.

Construction of the circRNA-miRNA-mRNA-ARGS Network

DEARmRNAs were the key module in the ceRNA network. MiRNAs targeting DEARmRNAs were predicted by miRWalk 3.0³ (Dweep et al., 2011). These miRNAs were intersected with DEmiRNAs to obtain the miRNA module. The circRNA sponges of the miRNA modules were obtained by circBank⁴ (Liu et al., 2019). CircRNAs sponging more than one candidate miRNAs were included in the circRNA module. Correlation analysis was then performed between the ARmRNA module and miRNA module in the total population of the TCGA NSCLC project to obtain the negatively correlated ARmRNAs and miRNAs. The targeted ARGSs of ARmRNAs and circRNA sponges of miRNAs were added to construct the final regulatory network. The R package ggalluvial and Cytoscape (Shannon et al., 2003) were used to visualize the ceRNA and circRNA-miRNA-ARmRNA-ARGS network.

Additional Bioinformatics and Statistical Analyses

R software 4.0.4⁵, GraphPad Prism 9.0 (GraphPad Software Inc., San Diego, CA, United States), and Cytoscape 3.8.2⁶ were used to analyze and visualize the data. The scores of gene sets (ARGSs, immune characteristics of the tumor microenvironment) in each sample were quantified via both single-sample gene set enrichment analysis (ssGSEA) and the gene set variation analysis (GSVA) algorithm based on the bulk RNA-seq data using the GSVA R package (Hnzelmann et al., 2013; Supplementary Table 2). We used the chi-square test for correlation analysis between categorical variables, and Pearson correlation coefficients for correlation analysis between continuous variables. A p value < 0.05 was considered statistically significant.

Results

Clinical Significance of Autophagy-Related Gene Sets in NSCLC Patients Treated With Radiotherapy

We plotted heatmaps to describe the distribution of ARGS scores by ssGSEA or GSVA (Figure 2A) and performed differential analysis between the radiotherapy-resistant and radiotherapy-sensitive groups. Late endosomal microautophagy (GO:0061738) was identified as significant by both methods (Figure 2B). The correlation between clinicopathological factors and ARGS was evaluated (Figure 2C). The consistent results of ssGSEA and GSVA score showed that autophagy of mitochondrion (GO:0000422) and macroautophagy (GO:0016236) were discriminatory for different histological types, while autophagy of peroxisome (GO:0030242), autophagy of nucleus (GO:0044804), late endosomal microautophagy (GO:0061738), and symbiont of host autophagy (GO:0075071) were discriminatory for patient gender and histological type. Multivariate logistic regression analysis demonstrated that lysosomal microautophagy (GO:0016237) and late endosomal microautophagy (GO:0061738) were independent risk factors for radiotherapy sensitivity (Figure 2D), while no ARGSs were associated with OS (Supplementary Figure 2). These findings demonstrated the clinical significance of autophagy or selective types of autophagy in NSCLC patients receiving radiotherapy.

FIGURE 2

Figure 2. Clinical correlation of autophagy-related gene sets (ARGSs) score by ssGSEA and GSVA methods. (A) The distribution of ARGSs score in radiotherapy resistant and sensitive groups. (B) The differential analysis of ARGSs score between radiotherapy resistant and sensitive groups. The left is radiotherapy resistant group while the right is radiotherapy sensitive group. (C) The correlation analysis of ARGSs and clinicopathologic factors. (D) Multivariate logistic regression analysis of ARGSs with radiotherapy responses.

Correlation Analysis of ARGSs and Immune Microenvironment Characteristics in NSCLC Patients Receiving Radiotherapy

With regard to immune microenvironment characteristics (Supplementary Figure 3A and Supplementary Table 3), autophagy (GO:0006914), regulation of autophagy (GO:0010506), macroautophagy (GO:0016236), and autophagy of peroxisome (GO:0030242) were related to infiltration of central memory CD8 T cells and gamma delta T cells (Supplementary Figure 3B). In addition, macroautophagy (GO:0016236) was related to infiltration of CD56 bright NK cells, and symbiont of host autophagy (GO:0075071) was related to central memory CD4 T cells. Regarding the immune cycle, only autophagy of peroxisome (GO:0030242) correlated with the trafficking of monocytes to tumors (Step 4) (Supplementary Figure 3C). Finally, and most importantly, we evaluated the association with ICIs (Supplementary Figure 3D). Autophagy (GO:0006914) was associated with expression levels of ADORA2A, CD200R1, CD274, CD80, CD86, HAVCR2, LAIR1, LGALS3, and TIGIT; autophagy of mitochondrion (GO:0000422) was correlated with CD276 and SNCA; regulation of autophagy (GO:0010506) was correlated with ADORA2A, CD200R1, CD274, CD80, CD86, HAVCR2, LAIR1, and TIGIT; macroautophagy (GO:0016236) was correlated with CD80, CD86, HAVCR2, LAG3, and LAIR1; autophagy of peroxisome (GO:0030242) was correlated with CD276 and SNCA; chaperone-mediated autophagy (GO:0061684) was correlated with CD200R1, CD80, CD86, HAVCR2, LAG3, and LAIR1; and late endosomal microautophagy (GO:0061738) and symbiont of host autophagy (GO:0075071) were correlated with SNCA.

Identification of Differentially Expressed Autophagy-Related mRNAs (DEARmRNAs) and miRNAs (DEmiRNAs) Associated With Radiotherapy Sensitivity in NSCLC Patients

In addition to the levels of gene sets, we explored the clinical significance of autophagy at the gene level. The DEmRNAs were identified from the radiotherapy-sensitive group compared to the radiotherapy-resistant group by DESeq2, edgeR, and limma/voom seperately (Figure 3A). In total, 461 DEmRNAs (235 upregulated and 226 downregulated) were found by intersection of these three methods (Figure 3B). Then, DEARmRNAs were recognized from DEmRNAs (Figure 3C) and one downregulated gene (ELAPOR1), three upregulated genes (SNCA, SESN3, and DAPL1) were the consistent results. The same methods were used in DEmiRNA screening (Figure 3D), which identified 3 upregulated (hsa-miR-205-5p, hsa-miR-26a-1-3p, and hsa-miR-6510-3p) and 3 downregulated (hsa-miR-194-3p, hsa-miR-215-5p, and hsa-miR-375-3p) DEmiRNAs by intersection of DESeq2, edgeR, and limma/voom (Figures 3E,F).

FIGURE 3

Figure 3. Differentially expressed (DE) analysis of mRNAs and miRNAs between rediotherapy resistant and sensitive groups. (A) Volcano plots of DEmRNAs by DESeq2, edgeR, and limma/voom. (B) Venn diagram of common DEmRNAs. (C) Heatmaps of DEARmRNAs. (D) Volcano plots of DEmiRNAs by DESeq2, edgeR, and limma/voom. (E) Venn diagram of common DEmiRNAs. (F) Heatmaps of DEmiRNAs. Red represents upregulated genes and blue indicates downregulated genes.

Validation of the Prognostic Value of DEARmRNAs and DEmiRNAs in NSCLC Patients Receiving Radiotherapy

To establish crucial miRNAs and ARmRNAs with prognostic value in NSCLC patients receiving radiotherapy, we first verified the differential expression of mRNAs and miRNAs between radiotherapy-sensitive and radiotherapy-resistant groups. Our results showed that hsa-miR-194-3p, hsa-miR-215-5p, hsa-miR-375-3p, and ELAPOR1 were upregulated in the radiotherapy-resistant group, while hsa-miR-205-5p, hsa-miR-26a-1-3p, SESN3, SNCA, and DAPL1 were upregulated in the radiotherapy-sensitive group (Figure 4). To determine whether these DERNAs are associated with the long-term prognosis of NSCLC patients treated with radiotherapy, we generated Kaplan-Meier curves to analyze differences in OS. We found that SNCA, SESN3, and DAPL1 were related to the OS of NSCLC patients (Figure 5).

FIGURE 4

Figure 4. The distribution of differentially expressed miRNAs and mRNAs between rediotherapy resistant and sensitive groups. Blue represents radiotherapy resistant group and pink indicates radiotherapy sensitive group.

FIGURE 5

Figure 5. Overall survival analysis of differentially expressed miRNAs and mRNAs. The high- and low-expression values of four autophagy-related mRNAs (ARmRNAs) and six miRNAs were compared by Kaplan-Meier survival curve for NSCLC patients. The median survival time were indicated by dashed line.

Establishment of the ARmRNA/miRNA Signature to Predict Prognosis in NSCLC Patients Receiving Radiotherapy

Based on the above results, we first established a 9 ARmRNA/miRNA signature by multivariate logistic regression to predict the radiosensitivity of NSCLC patients, and the score for each patient was calculated as follows: ARRS = 1.543103 −0.002191^∗hsa-miR-205-5p −0.500703^∗hsa-miR-215-5p +0.776517^∗hsa-miR-26a-1-3p −0.300282^∗hsa-miR-194-3p −0.041439^∗hsa-miR-375-3p −0.394497^∗ELAPOR1 −0.053145^∗SNCA +0.331343^∗SESN3 +0.181778^∗DAPL1. The ROC curve was generated, and the estimated AUC was 0.885 with a 95% CI of 0.813–0.958 (Figure 6A). The ARRS discriminated the radiotherapy-sensitive group from the radiotherapy-resistant group (p < 0.001) by higher ARRS score (Figure 6C) and related to histology and stage (Figure 6E). Besides, ARRS could serve as an independent radiotherapy sensitivity predictor for NSCLC and high ARRS score patients are more likely to get better radiotherapy sensitivity (OR:3.13[95%CI:1.66–4.96], p < 0.001) (Figure 6G).

FIGURE 6

Figure 6. Construction and clinical correlation analysis of autophagy-related prognostic risk signature. (A) The receive operator curve (ROC) analysis for autophagy-related radiosensitivity risk signature (ARRS). AUC: Area Under Curve, FPR: false positive rate, TPR: true positive rate. (B) The ROC analysis for autophagy-related overall survival risk signature (AROS). (C) The distribution of ARRS between rediotherapy resistant and sensitive groups. Blue represents radiotherapy resistant group and pink indicates radiotherapy sensitive groups. (D) The Kaplan-Meier survival curve grouping by high- and low- AROS. The median survival time were indicated by dashed line. (E) The correlation analysis of ARRSs and clinicopathologic factors. (F) The correlation analysis of AROSs and clinicopathologic factors. (G) Multivariate logistic regression analysis of ARRSs with radiotherapy sensitivity. (H) Multivariate Cox regression analysis of AROSs with overall survival.

Moreover, a 3 ARmRNA signature was generated by multivariate Cox regression to predict the OS of NSCLC patients, and the score for each patient was calculated as follows: AROS = −0.19011^∗hsa-miR-6510-3p −0.18664^∗SNCA −0.14049^∗SESN3 +0.06797^∗DAPL1. The ROC curve was generated, and the estimated AUC was 0.790, with a 95% CI of 0.717–0.869 (Figure 6B). The K-M curves were different between the high- and low-AROS groups (p < 0.001) (Figure 6D), and multivariate Cox regression revealed that AROS served as an independent predictor of OS for NSCLC patients who scored higher AROS with shorter OS (HR:3.40[95%CI:1.44–7.99], p = 0.005) (Figure 6H). The AROS was also related to histology and stage (Figure 6F).

Landscape of Immune Microenvironment Characteristics Associated With the ARmRNA/miRNA Signature

The two prognosis-related signatures (ARRS and AROS) were then estimated for immune microenvironment characteristics. ARRS positively correlated with the infiltration of CD56 bright NK cells and central memory CD8 T cells but negatively correlated with eosinophils and type 17 T helper cells (Figure 7A), while AROS negatively correlated with the infiltration of CD56 bright NK cells in both methods (Figure 7C). In terms of the immune cycle, ARRS was positively correlated with the trafficking of eosinophil cells to tumors (Step 4) (Figure 7B), while AROS was positively correlated with the trafficking of TH1 cells to tumors (Step 4) (Figure 7D) in both methods. Finally, we also investigated the relationship between ARRS or AROS and the expression levels of immune checkpoint molecules. Low ARRS indicated low expression of SNCA, CD200R1, CD276, LGALS3, and VTCN1 but high expression of CEACAM1 (Figure 7E), while low AROS represented high expression of SNCA, CD200R1, CD276, and VTCN1 (Figure 7F).

FIGURE 7

Figure 7. Role of autophagy-related prognostic risk signature in predicting immune phenotypes. (A) Correlations between ARRS and infiltration levels of tumor-associated immune cells. (B) Correlations between ARRS and immune cycle. (C) Correlations between AROS and infiltration levels of tumor-associated immune cells. (D) Correlations between AROS and immune cycle. (E) Correlations between ARRS and immune checkpoint inhibitors. (F) Correlations between AROS and immune checkpoint inhibitors.

Construction of the circRNA-miRNA-ARmRNA-ARGS Network to Regulate the Radiation Sensitivity of NSCLC

The targeted miRNAs of DEARmRNAs were predicted by miRWalk (Supplementary Table 4) and intersected with DEmiRNAs. We obtained three candidate ARmRNAs (ELAPOR1, SESN3, and SNCA) and 5 miRNAs (hsa-miR-26a-1-3p, hsa-miR-6510-3p, hsa-miR-205-5p, hsa-miR-375-3p, and hsa-miR-194-3p) for network construction (Figure 8A). Considering the conventionally negative correlation between mRNAs and miRNAs in regulatory relationships, we used the total population of TCGA NSCLC project for correlation analysis between these three ARmRNAs and their targeted miRNAs. After secondary screening, three miRNA/ARmRNA axes were recognized, namely, miR-205-5p/ELAPOR1, miR-26a-1-3p/SNCA, and miR-194-3p/SESN3 (Figures 8B–D). We then used circBank to identify the circRNAs targeting these three miRNAs (Supplementary Table 4). To enhance the affinity between circRNAs and the ceRNA network, we sought candidate circRNAs targeting two or more miRNAs. Finally, six circRNAs (hsa_circ_0019709, hsa_circ_0081983, hsa_circ_0112354, hsa_circ_0040569, hsa_circ_0135500, and hsa_circ_0098966) were identified. Finally, ARmRNA-associated ARGSs were added to form a complete regulatory network (Figure 9). All three ARGs (ELAPOR1, SNCA, and SESN3) participate in the regulation of autophagy. ELAPOR1 and SESN3 are involved in macroautophagy, and SNCA participates in chaperone-mediated autophagy.

FIGURE 8

Figure 8. Construction and correlation analysis of the ceRNA network. (A) The alluvial diagram of regulatory network of ceRNA. (B) Correlation analysis between has-miR-205-5p and ELAPOR1. (C) Correlation analysis between has-miR-26-1-3p and SNCA. (D) Correlation analysis between has-miR-194-3p and SESN3.

FIGURE 9

Figure 9. Diagram of the schematics of the circRNA-miRNA-ARmRNA-ARGS network. Red circles represent high-risk genes; orange circles represent their, respectively, identified regulatory miRNAs; yellow circles represent sponge circRNAs; green, blue, and purple ellipse represent corresponding ARGSs.

Discussion

Although autophagy can be non-specific, there are many selective types of autophagy (Klionsky et al., 2016). For a more detailed exploration of the role of autophagy in radiotherapy responses, we referred to the GO source to identify 9 ARGSs in the present study. The differential analysis of ARGS scores revealed that late endosomal microautophagy was distinct between the radiotherapy-sensitive and radiotherapy-resistant groups. Endosomal microautophagy requires endosomal-sorting complex systems for lysosome or endosome delivery and selectively degrades KFERQ-containing proteins recognized by HSC70 (Zheng et al., 2019). Microautophagy is the least studied form of autophagy with a largely unclear cargo delivery process (Zheng et al., 2019). The significance of endosomal microautophagy in the radiotherapy sensitivity of NSCLC patients was first proposed in our study; the intrinsic mechanism is worth pursuing in the future.

Competing endogenous RNAs (ceRNAs) are transcripts that competitively bind to shared miRNAs and act as miRNA sponges to modulate each other at the posttranscriptional level (Qi et al., 2015). With the development of high-throughput sequencing technology, abundant circRNAs have been identified and have become the focus in the ceRNA family due to the abundance of conserved miRNA response elements (MREs) (Zhong et al., 2018). Previous research has demonstrated that one of the most important mechanisms of circRNAs is their action on ceRNAs. For example, circRNA hsa_circ_100395 has been demonstrated to inhibit lung cancer progression by regulating the miR-1228/TCG21 pathway (Chen et al., 2018), while circRNA_101237 promotes NSCLC progression by regulating the miR-490-3p/MAPK1 axis (Zhang Z. Y. et al., 2020). Moreover, Jin et al. (2020) revealed potential prognostic biomarkers for radiotherapies with X-rays and carbon ions in NSCLC by integrating analysis of the circRNA-miRNA-mRNA network. Overall, the role of the circRNA-miRNA-mRNA network in the radiotherapy sensitivity of NSCLC still needs further research. In our study, after generating three miRNA-ARmRNA axes (miR-194-3p/SESN3, miR-205-5p/ELAPOR1, and miR-26a-1-3p/SNCA), we obtained six circRNAs (hsa_circ_0019709, hsa_circ_0040569, hsa_circ_0081983, hsa_circ_0098996, hsa_circ_0112354, hsa_circ_0135500) that intersected these three axes and constructed a circRNA-miRNA-ARmRNA network. It is worth noting that these ARmRNAs were contained in three ARGSs, namely, regulation of autophagy, macroautophagy, and chaperone-mediated autophagy. Macroautophagy refers to the process of autophagosomes formation and fusion with late endosomes or lysosomes to form amphisomes or autolysosomes, which are the canonical and well-known participants in the autophagy process (Zheng et al., 2019). Chaperone-mediated autophagy (CMA) is another vital type of selective autophagy which selectively degrades cytosolic proteins recognized by a specific chaperone in lysosomes (Zheng et al., 2019). CMA does not rely on vesicles or membrane invagination to deliver targeted substrates and degrades 30% of cytosolic proteins. SESN3 encodes a member of the sestrin family of stress-induced proteins, which may contribute to the positive regulation of macroautophagy (Pascale et al., 2011). ELAPOR1 is an endosome-lysosome-associated apoptosis and autophagy regulator, and it may protect cells from cell death by upregulating the autophagy pathway (Deng et al., 2010). SNCA is a member of the synuclein family and negatively regulates CMA (Alvarez-Erviti et al., 2010). In summary, our study utilized circRNA-miRNA-ARmRNA network analysis to investigate the subtypes of autophagy.

With the general success of immune checkpoint inhibitor antibodies and cell-based treatments, the age of immunotherapy has arrived, which raises the question of how autophagy interacts with the immune microenvironment and contributes to cancer treatments (Li et al., 2017; Jiang et al., 2019). It remains unclear whether autophagy inhibition impairs systematic immunity. Some evidence has shown that autophagy maintains the survival of memory T cells (Puleston et al., 2014) and promotes self-renewal of B1 cells (Clarke et al., 2018), while other evidence has shown that autophagy inhibition does not impair T cell function in preclinical models of melanoma and breast cancer, including chemotherapy-treated cells (Starobinets et al., 2016). Although a greater understanding of the role of autophagy in tumor immunity is emerging, the distinction between canonical autophagy and types of selective autophagy needs to be considered. Correlation analysis of ARGSs and ICI expression levels, immune cell infiltration, and the immune cycle was conducted in our work. We found that autophagy was related to the expression levels of many ICIs and the infiltration of central memory CD8 T cells and gamma delta T cells, while peroxisome autophagy correlated with the trafficking of monocytes to tumors. Though more extensive experiments are needed to confirm this model, these results support that autophagy levels are in tune with the immune microenvironment and have the potential to contribute to monitoring and improving immunotherapy in NSCLC patients.

Despite the tremendous development of radiation technology, tumor control and survival in NSCLC patients have not substantially improved. Individual heterogeneity partly explains this. Some patients may benefit from specific treatments while others require more aggressive treatments. To improve clinical outcomes and avoid excessive medical treatment, patients should be classified into cohorts according to differences in disease susceptibility, prognosis, and likely treatment response rates (Meehan et al., 2020). Additionally, the incorporation of molecular analysis and other patient information into the prevention, investigation, and treatment of diseases is an important aspect of precision medicine (Penet et al., 2014). Some efforts have been made to identify biomarkers that could be applied to tailor radiotherapy sensitivity to individual molecular characteristics of patient tissue. Salem et al. (2017) reported a blood biomarker panel containing interleukin (IL)-1b, neutrophil count, and cytokeratin-19 antigen to predict lung cancer radiotherapy response. Saito et al. (2014) constructed a three-microRNA signature to predict responses to platinum-based doublet chemotherapy in patients with lung adenocarcinoma. Liu et al. (2019) identified a miRNA signature by an in vitro system to assess radiosensitivity for head and neck squamous cell carcinomas and validated this signature using the TCGA database (Ning et al., 2015). These studies indicate that radiotherapy sensitivity should be considered before designing the treatment plan. Furthermore, short-term radiotherapy response does not always equate to long-term treatment outcomes. Hence, we also established an OS-related signature (AROS) beyond a radiotherapy sensitivity predictive signature (ARRS). The differential expression analysis and autophagy-related gene selection provided strong background support.

With the advent of immunotherapy, the interaction of radiotherapy and the immune system has gained widespread interest, and this interaction has been increasingly reported in NSCLC (Herter-Sprie et al., 2016). Radiotherapy has been demonstrated to promote tumor cell death and enhance antitumor immune responses by converting poorly immunogenic tumors into more highly immunogenic ones, not only through immunogenic cell death (ICD) but also through the modification of the characteristics of key immune cells within the tumor microenvironment (Keam et al., 2020). However, radiotherapy may be a double-edged sword; it induces activation and infiltration of T cells to the tumor bed, but it also triggers migration of immunosuppressive cells and upregulates inhibitory ligands and receptors (Keam et al., 2020). To improve the beneficial effects and reduce the risks, the biological responses and toxicities of radiation and drugs should be accurately modeled. However, the combination and the optimal timing, dose, or schedules of radiotherapy and immunotherapy are still controversial (Aliru et al., 2018). In addition to investigating the molecular features of patients’ responses to radiotherapy, we also described their tumor microenvironment by the bulk RNA-seq results in the present study. Autophagy-related risk scores predict not only radiotherapy sensitivity and OS but also the landscape of ICIs, immune cell infiltration, and immune cycle activation. These signature models may aid treatment decision making with consideration of concurrent radiotherapy and immunotherapy.

There were several limitations to this study. First, due to the incompleteness of primary therapy outcome success data in TCGA, fewer than 90 patients met our inclusion requirements. In addition, the lack of an index may render an inaccurate interpretation. Second, the prognostic signature was not validated because of the rarity of data sets recording radiotherapy responses in NSCLC patients. Third, although a potential regulatory mechanism has been constructed, no experimental support was provided. To ameliorate the limitations described above, single institution or multicenter clinical retrospective or prospective study should be conducted to verify the predictive value of prognostic signatures, and experiments should be rigorously designed to demonstrate the regulatory network of NSCLC.

In conclusion, we examined the role of autophagy-related genes (ARGs) and gene sets (ARGSs) in the radiotherapy response of NSCLC patients by mining public data. First, we verified the clinical significance of autophagy in the radiotherapy response of NSCLC patients by analyzing the correlation between ARGs or ARGSs and clinicopathologic factors, prognosis, and the immune microenvironment. In addition, the circRNA-miRNA-ARmRNA-ARGS network was constructed to predict the regulatory mechanisms underlying the radiation response of NSCLC. In summary, our work provided useful information to introduce potential molecular classifications and regulatory mechanisms into radiotherapy short- and long-term responses of NSCLC patients.

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 authors.

Author Contributions

LF, BL, and LS conceived and designed the whole project. LF, ZL, and LS analyzed the data and wrote the manuscript. BL helped the data discussion. All authors read and approved the final manuscript.

Funding

This study was supported by National Natural Science Foundation of China (grant number 31970636), Academic Promotion Project of Shandong First Medical University.

Conflict of Interest

ZL is employed by Shandong Yidian Gene Technology Co., Ltd.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be constructed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2021.730003/full#supplementary-material

Supplementary Figure 1 | Identification of autophagy-related gene sets (ARGSs).

Supplementary Figure 2 | Multivariate Cox regression analysis of ARGSs with overall survival.

Supplementary Figure 3 | Role of ARGS in predicting immune phenotypes. (A) The content of tumor immune microenvironment. (B) Correlations between ARGSs and infiltration levels of tumor-associated immune cells. (C) Correlations between ARGSs and immune cycle. (D) Correlations between ARGSs and immune checkpoint inhibitors. Upper left of (B–D) represent the r value while upper right represent p value by calculating ARGSs ssGSEA score; lower left represent r value while lower right represent p value by calculating ARGSs GSVA score.

Supplementary Figure 4 | Correlation analysis between other predicted miRNAs and ARmRNAs.

Footnotes

1. ^ https://xena.ucsc.edu/
2. ^ http://geneontology.org/
3. ^ http://mirwalk.umm.uni-heidelberg.de/
4. ^ http://www.circbank.cn/
5. ^ http://www.r-project.org
6. ^ https://cytoscape.org/

References

Aliru, M. L., Schoenhals, J. E., Venkatesulu, B. P., Anderson, C. C., Barsoumian, H. B., Younes, A. I., et al. (2018). Radiation therapy and immunotherapy: what is the optimal timing or sequencing? Immunotherapy 10:299. doi: 10.2217/imt-2017-0082

PubMed Abstract | CrossRef Full Text | Google Scholar

Allemani, C., Matsuda, T., Carlo, V. D., Harewood, R., Matz, M., Nikšić, M., et al. (2018). Global surveillance of trends in cancer survival: analysis of individual records for 37,513,025 patients diagnosed with one of 18 cancers during 2000–2014 from 322 population-based registries in 71 countries (CONCORD-3). Lancet 391, 1023–1075.

Google Scholar

Alvarez-Erviti, L., Rodriguez-Oroz, M. C., Cooper, J. M., Caballero, C., and Schapira, A. (2010). Chaperone-Mediated Autophagy Markers in Parkinson Disease Brains. Arch. Neurol. 67, 1464–1472.

Google Scholar

Ambros, V. (2004). The functions of animal microRNAs. Nature 431, 350–355. doi: 10.1038/nature02871

PubMed Abstract | CrossRef Full Text | Google Scholar

Baker, S., Dahele, M., Lagerwaard, F. J., and Senan, S. (2016). A critical review of recent developments in radiotherapy for non-small cell lung cancer. Radiat. Oncol. 11:115.

Google Scholar

Blanche, P. (2015). timeROC: Time-Dependent ROC Curve and AUC for Censored Survival Data.

Google Scholar

Charoentong, P., Finotello, F., Angelova, M., Mayer, C., Efremova, M., Rieder, D., et al. (2017). Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep. 18, 248–262. doi: 10.1016/j.celrep.2016.12.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, D., Ma, W., Ke, Z., and Xie, F. (2018). CircRNA hsa_circ_100395 regulates miR-1228/TCF21 pathway to inhibit lung cancer progression. Cell Cycle 17, 2080–2090. doi: 10.1080/15384101.2018.1515553

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, D. S., and Mellman, I. (2013). Oncology meets immunology: the cancer-immunity cycle. Immunity 39, 1–10. doi: 10.1016/j.immuni.2013.07.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, H., and Boutros, P. C. (2011). VennDiagram: a package for the generation of highly-customizable Venn and Euler diagrams in R. Bmc Bioinformatics 12:35. doi: 10.1186/1471-2105-12-35

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, L. L. (2020). The expanding regulatory mechanisms and cellular functions of circular RNAs. Nat. Rev. Mol. Cell Biol. 21, 475–490. doi: 10.1038/s41580-020-0243-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Clarke, A. J., Riffelmacher, T., Braas, D., Cornall, R. J., and Simon, A. K. (2018). B1a B cells require autophagy for metabolic homeostasis and self-renewal. J. Exp. Med. 215, 399–413. doi: 10.1084/jem.20170771

PubMed Abstract | CrossRef Full Text | Google Scholar

Deng, L., Feng, J., and Broaddus, R. R. (2010). The novel estrogen-induced gene EIG121 regulates autophagy and promotes cell survival under stress. Cell Death Dis. 1:e32. doi: 10.1038/cddis.2010.9

PubMed Abstract | CrossRef Full Text | Google Scholar

Dweep, H., Sticht, C., Pandey, P., and Gretz, N. (2011). MiRWalk - Database: prediction of possible miRNA binding sites by “walking” the genes of three genomes. J. Biomed. Inform. 44, 839–847. doi: 10.1016/j.jbi.2011.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Hansen, T. B., Jensen, T. I., Clausen, B. H., Bramsen, J. B., Finsen, B., Damgaard, C. K., et al. (2013). Natural RNA circles function as efficient microRNA sponges. Nature 495, 384–388. doi: 10.1038/nature11993

PubMed Abstract | CrossRef Full Text | Google Scholar

Herter-Sprie, G. S., Koyama, S., Korideck, H., Hai, J., and Wong, K. K. (2016). Synergy of radiotherapy and PD-1 blockade in Kras-mutant lung cancer. JCI Insight 1:e87415.

Google Scholar

Hnzelmann, S., Castelo, R., and Guinney, J. (2013). GSVA: gene set variation analysis for microarray and RNA-Seq data. BMC Bioinformatics 14:7. doi: 10.1186/1471-2105-14-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiang, G. M., Tan, Y., Wang, H., Peng, L., and Shan, H. (2019). The relationship between autophagy and the immune system and its applications for tumor immunotherapy. Mol. Cancer 18:17.

Google Scholar

Jin, X., Yuan, L., Liu, B., Kuang, Y., and Li, Q. (2020). Integrated analysis of circRNA-miRNA-mRNA network reveals potential prognostic biomarkers for radiotherapies with X-rays and carbon ions in non-small cell lung cancer. Ann. Transl. Med. 8, 1373–1373. doi: 10.21037/atm-20-2002

PubMed Abstract | CrossRef Full Text | Google Scholar

Keam, S., Gill, S., Ebert, M. A., Nowak, A. K., and Cook, A. M. (2020). Enhancing the efficacy of immunotherapy using radiotherapy. Clin. Transl. Immunol. 9:e1169.

Google Scholar

Klionsky, D. J., Abdelmohsen, K., Abe, A., Abedin, M. J., Abeliovich, H., Arozena, A. A. et al. (2016). Guidelines for the use and interpretation of assays for monitoring autophagy (3rd edition). Autophagy 12, 1–222.

Google Scholar

Li, C. J., Liao, W. T., Wu, M. Y., and Chu, P. Y. (2017). New Insights into the Role of Autophagy in Tumor Immune Microenvironment. Int. J. Mol. Sci. 18:1566. doi: 10.3390/ijms18071566

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, M., Wang, Q., Shen, J., Yang, B., and Ding, X. (2019). Circbank: a comprehensive database for circRNA with standard nomenclature. RNA Biol. 16, 899–905. doi: 10.1080/15476286.2019.1600395

PubMed Abstract | CrossRef Full Text | Google Scholar

Meehan, J., Gray, M., Martínez-Pérez, C., Kay, C., and Turnbull, A. K. (2020). Precision Medicine and the Role of Biomarkers of Radiotherapy Response in Breast Cancer. Front. Oncol. 10:628. doi: 10.3389/fonc.2020.00628

PubMed Abstract | CrossRef Full Text | Google Scholar

Murrow, L., and Debnath, J. (2013). Autophagy as a Stress-Response and Quality-Control Mechanism: implications for Cell Injury and Human Disease. Annu. Rev. Pathol. Mech. Dis. 8, 105–137. doi: 10.1146/annurev-pathol-020712-163918

PubMed Abstract | CrossRef Full Text | Google Scholar

Ning, L., Boohaker, R. J., Jiang, C., Boohaker, J. R., and Xu, B. (2015). A radiosensitivity MiRNA signature validated by the TCGA database for head and neck squamous cell carcinomas. Oncotarget 6, 34649–34657. doi: 10.18632/oncotarget.5299

PubMed Abstract | CrossRef Full Text | Google Scholar

Noam, A., Gao, Z., Sang, L. J., Frederick, D. T., Miao, B., Tabea, M., et al. (2018). Robust prediction of response to immune checkpoint blockade therapy in metastatic melanoma. Nat. Med. 24, 1545–1549. doi: 10.1038/s41591-018-0157-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Pascale, G., Livstone, M. S., Lewis, S. E., and Thomas, P. D. (2011). Phylogenetic-based propagation of functional annotations within the Gene Ontology consortium. Brief. Bioinform. 12, 449–462. doi: 10.1093/bib/bbr042

PubMed Abstract | CrossRef Full Text | Google Scholar

Penet, M. F., Krishnamachary, B., Chen, Z., Jin, J., and Bhujwalla, Z. M. (2014). Molecular imaging of the tumor microenvironment for precision medicine and theranostics. Adv. Cancer Res. 124:235. doi: 10.1016/b978-0-12-411638-2.00007-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Puleston, D. J., Zhang, H., Powell, T. J., Elina, L., Stuart, S., Isabel, P., et al. (2014). Autophagy is a critical regulator of memory CD8+ T cell formation. ELife 3:e03706.

Google Scholar

Qi, X., Zhang, D. H., Wu, N., Xiao, J. H., Wang, X., and Ma, W. (2015). ceRNA in cancer: possible functions and clinical implications. J. Med. Genet. 52, 710–718. doi: 10.1136/jmedgenet-2015-103334

PubMed Abstract | CrossRef Full Text | Google Scholar

Saito, M., Shiraishi, K., Matsumoto, K., Schetter, A. J., Ogata-Kawata, H., Tsuchiya, N., et al. (2014). A three-microRNA signature predicts responses to platinum-based doublet chemotherapy in patients with lung adenocarcinoma. Clin. Cancer Res. 20, 4784–4793. doi: 10.1158/1078-0432.ccr-14-1096

PubMed Abstract | CrossRef Full Text | Google Scholar

Salem, A., Mistry, H., Backen, A., Hodgson, C., Koh, P., Dean, E., et al. (2017). Cell Death, Inflammation, Tumor Burden, and Proliferation Blood Biomarkers Predict Lung Cancer Radiotherapy Response and Correlate With Tumor Volume and Proliferation Imaging. Clin. Lung Cancer 19, 239–248.e7.

Google Scholar

Shannon, P., Markiel, A., Ozier, O., Baliga, N. S., Wang, J. T., and Ramage, D. (2003). Cytoscape: a Software Environment for Integrated Models of Biomolecular Interaction Networks. Genome Res. 13, 2498–2504. doi: 10.1101/gr.1239303

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharma, K., Le, N., Alotaibi, M., and Gewirtz, D. (2014). Cytotoxic Autophagy in Cancer Therapy. Int. J. Mol. Sci. 15, 10034–10051. doi: 10.3390/ijms150610034

PubMed Abstract | CrossRef Full Text | Google Scholar

Starobinets, H., Ye, J., Broz, M., Barry, K., and De Bnath, J. (2016). Antitumor adaptive immunity remains intact following inhibition of autophagy and antimalarial treatment. J. Clin. Investig. 126:4417. doi: 10.1172/jci85705

PubMed Abstract | CrossRef Full Text | Google Scholar

Sung, H., Ferlay, J., Siegel, R. L., Laversanne, M., Soerjomataram, I., Jemal, A., et al. (2021). Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J. Clin. 71, 209–249. doi: 10.3322/caac.21660

PubMed Abstract | CrossRef Full Text | Google Scholar

Tam, S. Y., Wu, V., and Law, H. (2017). Influence of autophagy on the efficacy of radiotherapy. Radiat. Oncol. 12:57.

Google Scholar

Zhang, X., Xie, K., Zhou, H., Wu, Y., Li, C., Liu, Y., et al. (2020). Role of non-coding RNAs and RNA modifiers in cancer therapy resistance. Mol. Cancer 19:47.

Google Scholar

Zhang, Z. Y., Gao, X. H., Ma, M. Y., Zhao, C. L., and Guo, S. S. (2020). CircRNA_101237 promotes NSCLC progression via the miRNA-490-3p/MAPK1 axis. Sci. Rep. 10:9024.

Google Scholar

Zheng, K., He, Z., Kitazato, K., and Wang, Y. (2019). Selective Autophagy Regulates Cell Cycle in Cancer Therapy. Theranostics 9, 104–125. doi: 10.7150/thno.30308

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhong, Y., Du, Y., Yang, X., Mo, Y., Fan, C., Xiong, F., et al. (2018). Circular RNAs function as ceRNAs to regulate and control human cancer progression. Mol. Cancer 17:79.

Google Scholar