Circadian Clock Genes Are Correlated with Prognosis and Immune Cell Infiltration in Colon Adenocarcinoma

He, Aoxiao; Huang, Zhihao; Zhang, Rongguiyi; Lu, Hongcheng; Wang, Jiakun; Cao, Jiaqing; Feng, Qian

doi:https://doi.org/10.1155/2022/1709918

Computational and Mathematical Methods in Medicine

On this page

Abstract Introduction Materials and Methods Results Discussion Conclusion Data Availability Conflicts of Interest Authors’ Contributions Acknowledgments Supplementary Materials References Copyright Related Articles

Research Article | Open Access

Volume 2022 | Article ID 1709918 | https://doi.org/10.1155/2022/1709918

Circadian Clock Genes Are Correlated with Prognosis and Immune Cell Infiltration in Colon Adenocarcinoma

Aoxiao He,¹Zhihao Huang,¹Rongguiyi Zhang,¹Hongcheng Lu,¹Jiakun Wang,¹Jiaqing Cao,¹and Qian Feng²

Academic Editor: Peng Li

Received04 Aug 2021

Revised27 Oct 2021

Accepted30 Dec 2021

Published25 Jan 2022

Abstract

Background. Colon adenocarcinoma (COAD) is a malignancy with a high incidence and is associated with poor quality of life. Dysfunction of circadian clock genes and disruption of normal rhythms are associated with the occurrence and progression of many cancer types. However, studies that systematically describe the prognostic value and immune-related functions of circadian clock genes in COAD are lacking. Methods. Genomic data obtained from The Cancer Genome Atlas (TCGA) database was analyzed for expression level, mutation status, potential biological functions, and prognostic performance of core circadian clock genes in COAD. Their correlations with immune infiltration and TMB/MSI score were analyzed by Spearman’s correlation analysis. Pearson’s correlation analysis was performed to analyze their associations with drug sensitivity. Lasso Cox regression analysis was performed to construct a prognosis signature. Moreover, an mRNA-miRNA-lncRNA regulatory axis was also detected by ceRNA network. Results. In COAD tissues, the mRNA levels of CLOCK, CRY1, and NR1D1 were increased, while the mRNA levels of ARNTL, CRY2, PER1, PER3, and RORA were decreased. We also summarized the relative genetic mutation variation landscape. GO and KEGG pathway analyses demonstrated that these circadian clock genes were primarily correlated with the regulation of circadian rhythms and glucocorticoid receptor signaling pathways. COAD patients with high CRY2, NR1D1, and PER2 expression had worse prognosis. A prognostic model constructed based on the 9 core circadian clock genes predicted the COAD patients’ overall survival with medium to high accuracy. A significant association between prognostic circadian clock genes and immune cell infiltration was found. Moreover, the lncRNA KCNQ1OT1/hsa-miRNA-32-5p/PER2/CRY2 regulatory axis in COAD was also detected through a mRNA-miRNA-lncRNA network. Conclusion. Our results identified CRY2, NR1D1, and PER2 as potential prognostic biomarkers for COAD patients and correlated their expression with immune cell infiltration. The lncRNA KCNQ1OT1/hsa-miRNA-32-5p/PER2/CRY2 regulatory axis was detected in COAD and might play a vital role in the occurrence and progression of COAD.

1. Introduction

Colon adenocarcinoma (COAD) has become a severe health burden that causes over half a million cancer-related deaths annually [1]. COAD is the most common type of colon tumor. In many developing countries, patients with COAD are usually diagnosed at an advanced stage due to the lack of obvious symptoms at an early stage, related screening projects, and screening instruments [2]. The quality of life of advanced COAD patients is very low due to severe digestive symptoms [3]. Although surgical treatment, targeted treatment, chemotherapy, and immunotherapy have been flexibly utilized in the treatment of COAD, the prognosis of COAD patients is still disappointing [1]. It has been reported that the overall five-year survival rate of COAD patients is less than 40%, mainly due to postoperative recurrence and metastasis [4, 5]. Recently, aberrant expression of genes in tumor tissues and their potential applications in clinical practice have been a focus of clinicians. Thus, exploring potential gene signatures for therapeutic and prognostic assessments of COAD is significant clinically.

Physiological circadian rhythms are vital for various biological functions, such as metabolism, regeneration, immunology, and endocrinology [6]. All human tissues are modulated by the circadian clock through extremely complicated processes. The molecular circadian clock consists of a core clock gene loop, which mainly includes CLOCK, ARNTL, PER1, PER2, PER3, CRY1, CRY2, NR1D1, and RORA [6, 7]. In recent years, an increasing number of studies have attempted to explore the association between the circadian clock and cancer. Previous studies have found that disruption of normal rhythms and dysfunction of core clock gene loops contribute to the occurrence and development of many cancer types and influence the functions of the tumor immunity cycle [8–11]. For cancer therapy, previous studies revealed that cancer chronotherapy, i.e., timing chemotherapy administration based on circadian rhythms, might reduce the toxicity of drugs [12]. Some studies based on animal models have reported that drugs targeting biological clock targets, such as RORγ synthetic agonists and RORα synthetic agonists, could promote the activation of anticancer immunity [13, 14]. Nevertheless, the role of the circadian clock in prognosis evaluation and its clinical significance in COAD are rarely discussed.

With the rapid development of bioinformatics databases and the enrichment of genomic data, comprehensive analysis of the role of the circadian clock in COAD has become feasible. The present study explores the expression profiles, prognostic value, and related regulatory axis of the circadian clock in COAD, which may offer novel insights into prognostic biomarkers and biological clock targets for COAD.

2. Materials and Methods

2.1. Datasets

On May 20, 2021, COAD patient genomic data were downloaded from The Cancer Genome Atlas (TCGA) database. TCGA-COAD dataset () was selected to perform the analyses. The corresponding clinical data, mainly including sex, age, tumor grade, and survival data, were also obtained from TCGA database on May 20, 2021. Furthermore, the copy number variation (CNV) and somatic mutation data in COAD were also extracted from TCGA database using UCSC Xena (https://xena.ucsc.edu/), an easy-to-use web portal for data files derived from TCGA. Setting the value cutoff as 0.05, statistical analyses were conducted using R software V4.0.3. Before subsequent processing, the expression data were standardized to transcripts per kilobase million (TPM). To further analyze the expression of circadian clock genes, another gene expression profile of COAD was obtained from Oncomine (https://www.oncomine.org) [15], a cancer microarray database and integrated data-mining platform.

2.2. Identification of Differentially Expressed Core Circadian Clock Genes

A total of 9 genes were selected to perform further analysis, which included CLOCK, ARNTL, PER1, PER2, PER3, CRY1, CRY2, NR1D1, and RORA [6, 7]. The profiles of 24-hour oscillations of the core circadian clock genes were downloaded from CircaDB (https://circadb.hogeneschlab.org) [16]. For each tissue that was cyclic ordering by periodic structure (CYCLOPS) ordered, cosinor regression was used to test if individual genes are rhythmic or not [17]. In the 24-hour oscillation plots, the -axis represented the CYCLOPS-estimated sample phase in radians (from 0 to 2𝝿). The -axis represented the expression level (TPM). The time 𝝿 represented E-box phase (i.e., time of peak expression of E-box target genes NR1D1, NR1D2, and PER3). The datasets concerning the expression levels of the 9 core circadian clock genes were also excavated from the Oncomine (https://www.oncomine.org) database [15]. To further verify the expression level of the 9 core circadian clock genes, we also analyzed their expression in TCGA-COAD dataset. Student’s -test was utilized to compare the transcriptional levels of the core circadian clock genes in COAD tissues with those in normal colon tissues using the “reshape2” and “limma” packages in R software.

2.3. PPI Network and Mutation Analysis of Core Circadian Clock Genes

To study the correlation between each member of the 9 core circadian clock genes, we constructed a PPI network with the Search Tool for the Retrieval of Interacting Genes (STRING, https://string-db.org/), a database of predicted functional associations between genes and proteins [18]. Using R software, the overall mutation status of core circadian clock genes in COAD, including mutation categories and mutation frequencies of CLOCK, ARNTL, PER1, PER2, PER3, CRY1, CRY2, NR1D1, and RORA as well as their waterfall plot, were constructed using the “maftools” package. The locations of variants of the 9 genes on 23 chromosomes were shown using the “RCircos” package.

2.4. Functional Enrichment Analysis for Core Circadian Clock Genes

To identify the potential molecular mechanisms and biological functions of core circadian clock genes in COAD, functional enrichment analysis, including gene ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis, was performed using the “ggplot2” package in R software. In addition, GO analysis included biological process (BP), cellular component (CC), and molecular function (MF) analyses.

2.5. Construction of a Prognostic Model for Core Circadian Clock Genes

To assess the potential prognostic value of core circadian clock genes, Kaplan–Meier plots, including overall survival (OS) plots, disease-free survival (DFS) plots, disease-specific survival (DSS) plots, and progression-free survival (PFS) plots, were drawn. The hazard ratio (HR) and its 95% confidence interval (CI) are also listed in the plots. Based on the 9 core circadian clock genes, a LASSO Cox regression model was used to predict the influence of the expression of core circadian clock genes on the prognosis of COAD patients. According to the median risk score, patients were divided into the low- and high-risk groups, and the OS curves of the two groups were compared. The predictive accuracy of this prognostic model was assessed by time ROC analysis. Moreover, core circadian clock genes with potential prognostic value were selected for further analyses.

2.6. Immune Infiltration, TMB, MSI, and Drug Sensitivity Analyses

TIMER (https://cistrome.shinyapps.io/timer/) is a bioinformatics platform that is aimed at visualizing the association between immune cell infiltration and target gene expression in cancer cells [19]. The “Gene” module in TIMER was utilized to determine the correlation between the expression of prognostic circadian clock genes and the immune infiltration levels of B cells, CD8+ T cells, CD4+ T cells, neutrophils, dendritic cells, and macrophages. Based on the “ggstatsplot” package in R, tumor mutation burden (TMB) and microsatellite instability (MSI) analyses were used to assess the relationship between the expression of prognostic circadian clock genes and TMB and MSI scores, respectively. Spearman’s correlation test was selected for the above analysis with a value cutoff of 0.05. Moreover, to study the effect of prognostic circadian clock gene expression on drug resistance, drug sensitivity data in the Genomics of Drug Sensitivity in Cancer (GDSC) database were extracted.

2.7. Construction of the mRNA-miRNA-lncRNA Network

An mRNA-miRNA-lncRNA network was constructed to determine the potential regulatory axis of prognostic circadian clock genes in COAD. The miRNA targets binding to prognostic circadian clock genes were identified by the miRTarBase (http://mirtarbase.mbc.nctu.edu.tw/) [20] and starBase (https://starbase.sysu.edu.cn/index.php) [21] databases. miRTarBase and starBase are bioinformatics tools that can be used to predict miRNA targets of mRNAs. The mRNA-miRNA regulatory network was constructed using Cytoscape software. The Cytohubba plug-in was used to select the most significantly connected miRNAs. Next, according to the selected miRNAs, the lncRNA targets interacting with miRNAs were extracted from starBase (https://starbase.sysu.edu.cn/index.php) and lncBase (https://diana.e-ce.uth.gr/lncbasev3) [22], and the associations between the miRNAs and the lncRNAs were visualized. starBase and lncBase are bioinformatics tools that can be used to predict lncRNA targets of miRNAs. To further confirm the significant role of miRNAs and lncRNAs, we also evaluated the expression levels and prognostic performances of the targeted miRNAs and lncRNAs using TCGA-COAD dataset.

3. Results

3.1. Work Flow of the Current Study

Work flow of the present study is shown in Figure 1. Firstly, the core circadian clock genes in human were identified based on the previous studies [6, 7]. To reveal their rhythmicity in human colon tissues, the 24-hour oscillations of these genes were listed based on the expression profiles in CircaDB. Next, to find out the differentially expressed genes, we downloaded the RNA-sequencing data in Oncomine and TCGA-COAD dataset, respectively. The datasets in Oncomine were from the previous publications. Through comparing the expression profiles in the two databases, the expression status of core circadian clock genes in COAD was more trustworthy. And then, based on clinical data in TCGA-COAD database, the potential prognostic value of these genes in COAD was also analyzed by Kaplan-Meier methods, which included OS, DFS, PFS, and DSS. After screening out the prognostic circadian clock genes in COAD, we were also interested on the correlation between prognostic circadian clock genes expression and immune cell infiltration; these results might assist researchers to explore potential immunotherapeutic target for COAD patients. Immune infiltration analysis was conducted by TIMER. MSI and TMB could predict the efficacy of cancer immunotherapy. In theory, cancer patients with high TMB/MSI score might get better therapeutic effect. Thus, we detected the influence of prognostic circadian clock genes expression on TMB/MSI score. Meanwhile, the association between gene expression and drug sensitivity of colon tissues was also analyzed. To explore potential genes that might influence the development of COAD, we performed the subgroup analysis to show the correlation between gene expression and clinical stage. In addition, we attempted to explore a potential mRNA-miRNA-lncRNA regulatory axis concerning these prognostic circadian clock genes in COAD. To comprehensively investigate these 9 genes in COAD, the following analyses were performed. Mutation analysis helped readers to intuitively observe the overall mutation status of circadian clock genes in COAD. As the core genes in circadian clock loop, the potential biological functions and molecular mechanisms of these 9 genes were investigated by functional enrichment analysis. These results would provide novel insights on the functional studies of circadian clock genes. And in recent years, machine learning and big data technique were used to predict prognosis. Thus, LASSO Cox regression model based on the core circadian clock genes was constructed. Finally, the PPI network was constructed to visualize the correlation between each member of the 9 core circadian clock genes. Correlation analysis was performed to show the coexpression among these 9 genes in COAD, thus revealing synergistic interaction on expression among these genes.

3.2. Expression of Core Circadian Clock Genes in COAD and Normal Colon Tissues

The 24-hour oscillations of CRY1 (Figure 2(a)), CRY2 (Figure 2(b)), NR1D1 (Figure 2(c)), PER1 (Figure 2(d)), PER2 (Figure 2(e)), and PER3 (Figure 2(f)) in human colon tissues were visualized according to the expression profiles in CircaDB. However, there were no CLOCK, ARNTL, or RORA expression level data in CircaDB. Based on the Oncomine database, the number of datasets that contained the expression profiles of CLOCK, ARNTL, PER1, PER2, PER3, CRY1, CRY2, NR1D1, and RORA is shown in Figure 3. In colorectal cancer, the transcriptional levels of CLOCK, CRY1, and NR1D1 were upregulated in 1, 1, and 2 datasets, respectively. The transcriptional levels of ARNTL, CRY2, PER1, PER3, and RORA were downregulated in 2, 6, 2, 5, and 5 datasets, respectively. Using TCGA-COAD dataset, the mRNA levels of the 9 core circadian clock genes in TCGA-COAD tissues and the corresponding mRNA levels in normal colon tissues were identified and are shown in Figure 4(a). Specifically, the mRNA levels of CLOCK, CRY1, and NR1D1 were significantly higher in COAD tissues, while the mRNA levels of ARNTL, CRY2, PER1, PER3, and RORA were significantly decreased in COAD tissues (all values < 0.001). At the transcriptional level, these core circadian clock genes were greatly altered in COAD tissues. Moreover, a protein–protein interaction (PPI) network was constructed to illustrate the associations among the 9 core circadian clock genes (Supplementary Figure 1(a)). Correlation analysis of core circadian clock genes showed significant positive correlations with each other apart from NR1D1 in TCGA-COAD dataset, suggesting that these core circadian clock genes had synergistic interactions on expression (Supplementary Figure 1(b)).

(a)

(b)

(c)

(d)

(e)

(f)

(a)

(b)

(c)

(d)

(e)

Figure 4

Landscape of genetic variation and expression of core circadian clock genes in COAD. (a) The mRNA levels of core circadian clock genes in COAD and normal colon tissues. (b, c) The mutation frequency and classification of core circadian clock genes in COAD. (d) The location of CNV alteration of core circadian clock genes on 23 chromosomes in COAD. (e) The CNV variation frequency of core circadian clock genes in COAD. The height of the column represents the alteration frequency. CNV: copy number variation; COAD: colon adenocarcinoma. .

3.3. Landscape of the Genetic Variation of Core Circadian Clock Genes in COAD Tissues

Next, we summarized the somatic variations (SNVs) and copy number variations (CNVs) of the 9 circadian clock genes in COAD. As shown in Figures 4(b) and 4(c), genetic alterations were detected in 59 (88.06%) of 67 COAD samples. The categories of genetic mutations included missense mutation, nonsense mutation, frame-shift insertion, frame-shift deletion, in-frame deletion, splice-site variation, and multi-hit mutation. The most common variant classification was missense mutation. In addition, SNPs were the most common variant type, and C>T ranked as the top SNV class. Of these genes, PER3 was the gene with the highest mutation frequency. The locations of the CNV variations on chromosomes are also presented in Figure 4(d). The CNV variation frequencies of the 9 core circadian clock genes in COAD are also summarized in Figure 4(e). Copy number amplification was shown in NR1D1 and CLOCK, while copy number deletion was shown in PER3, PER2, RORA, ARNTL, CRY1, and PER1.

3.4. Functional Enrichment Analysis of Core Circadian Clock Genes in COAD Tissues

To clarify the potential biological functions and molecular mechanisms of core circadian clock genes in COAD, GO and KEGG pathway analyses were performed (Figure 5). GO analysis, including BP, CC, and MF analyses, illustrated that these 9 genes mainly participated in the regulation of circadian rhythm, circadian regulation of gene expression, rhythmic process, glucocorticoid receptor signaling pathway, corticosteroid receptor signaling pathway, transcription corepressor binding, and E-box binding (Figure 5(a)). KEGG pathway analysis revealed that the 9 core circadian clock genes were mainly involved in the regulation of circadian rhythm, circadian entrainment and dopaminergic synapses (Figure 5(b)).

(a)

(b)

3.5. Prognostic Model Construction of Core Circadian Clock Genes in COAD Tissues

Kaplan–Meier plots, including OS, DFS, PFS, and DSS plots, were used to assess the potential prognostic value of the 9 core circadian clock genes. The values and their HRs with 95% CIs are shown in Table 1. COAD patients with high CRY2 expression showed worse OS (Figure 6(a)) ( value = 0.032, HR with 95% (1.04-2.30)), DFS (Figure 6(c)) ( value = 0.040, HR with (1.04-6.29)), DSS (Figure 6(e)) ( value = 0.049, HR with (1.00-2.74)), and PFS (Figure 6(f)) ( value = 0.022, HR with (1.06-2.18)). COAD patients with high NR1D1 expression showed worse DFS (Figure 6(d)) ( value = 0.006, HR with 95% (1.50-10.90)) and PFS (Figure 6(g)) ( value = 0.003, HR with (1.22-2.54)). In addition, COAD patients with high PER2 expression showed worse OS (Figure 6(b)) ( value = 0.047, HR with (1.00-2.21)). The above results demonstrated that CRY2, NR1D1 and PER2 might be of great significance in the prognostic evaluation of COAD patients. Moreover, a prognostic model including the 9 core circadian clock genes was constructed by Lasso Cox regression analysis (Figure 7(a)–7(b)): According to the risk score, patients were divided into the high- and low-risk groups. The score distribution, survival status, and expression levels of PER2, NR1D1, and CLOCK are shown in Figure 7(c). When the risk score increased, the risk of death of COAD patients increased, while the survival time decreased (Figure 7(c)). Compared with the low-risk group, COAD patients with high risk showed a worse OS ( vs. 8.2 years, value = 0.015, Figure 7(d)). Furthermore, the AUCs of the 5- and 10-year ROC curves were 0.604 (95% CI: 0.512-0.697) and 0.696 (95% CI: 0.536-0.856), respectively (Figure 7(e)), which showed medium to high accuracy.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(a)

(b)

(c)

(d)

(e)

3.6. Immune Infiltration Analysis for Prognostic Circadian Clock Genes in COAD Tissues

The above study revealed that CRY2, NR1D1, and PER2 were potential biomarkers for COAD. We then selected CRY2, NR1D1, and PER2 for further analysis because of their prognostic performance. The association between the infiltration level of immune cells and the expression level of prognostic circadian clock genes in COAD was analyzed. As shown in Figure 8(a), the expression level of CRY2 was positively correlated with the infiltration levels of CD4+ T cells (, value = ), macrophages (, value = ), neutrophils (, value = ), and dendritic cells (, value = ). The expression level of PER2 was positively associated with the infiltration levels of B cells (, value = ), CD8+ T cells (, value = ), CD4+ T cells (, value = ), macrophages (, value = ), neutrophils (, value = ), and dendritic cells (, value = ) (Figure 8(b)). In addition, the expression level of NR1D1 was positively associated with the infiltration level of CD4+ T cells (, value = ) and negatively associated with the infiltration levels of B cells (, value = ) and CD8+ T cells (, value = ) (Figure 8(c)). We were also interested in the correlation between the common immune biomarkers and the expression levels of prognostic circadian clock genes in COAD. As shown in Table 2, CRY2 expression was positively correlated with most of the gene markers on immune cells. However, PER2 and NR1D1 expression levels were significantly correlated with only 24 and 29 kinds of gene markers, respectively. These results revealed that CRY2 might be a potential immunotherapeutic target for COAD patients.

(a)

(b)

(c)

3.7. TMB, MSI, and Drug Sensitivity Analyses of Core Circadian Clock Genes in COAD Tissues

Tumor mutation burden (TMB) and microsatellite instability (MSI) are of great significance in predicting the efficacy of cancer immunotherapy. The above results demonstrated that the prognostic circadian clock genes were associated with immune cell infiltration. Therefore, we conducted TMB and MSI analyses to assess the potential clinical applications of CRY2, NR1D1, and PER2 in COAD immunotherapy. In TMB analysis, negative correlations were detected between TMB and CRY2 expression (Figure 9(a), value = 0.002) and NR1D1 expression (Figure 9(c), value = 0.048), while no significant correlation was found between TMB and PER2 expression (Figure 9(b), value = 0.210). In MSI analysis, unexpectedly, no significant associations between MSI and CRY2 expression (Figure 9(d), value = 0.769), NR1D1 expression (Figure 9(e), value = 0.948), or PER2 expression (Figure 9(f), value = 0.588) were detected. Next, drug sensitivity analysis was performed to detect the potential functions of prognostic circadian clock genes in drug screening. Based on data from the GDSC database, high NR1D1 expression in COAD tissues was mainly correlated with drug resistance to 38 kinds of drugs and drug sensitivity to 6 kinds of drugs. High PER2 expression was primarily associated with drug sensitivity in 58 kinds of drugs and drug resistance in 13 kinds of drugs, while high CRY2 expression was mainly associated with drug sensitivity in 60 kinds of drugs and drug resistance in 18 kinds of drugs (Figure 9(g)), indicating that PER2 and CRY2 were potential biomarkers for drug screening in COAD.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

3.8. Construction of a mRNA-miRNA-LncRNA Network

The relationship between clinical stage and the expression of prognostic circadian clock genes was analyzed. As shown in Figure 10, the expression levels of CRY2 (Figure 10(a), value = 0.00016), NR1D1 (Figure 10(b), value = ), and PER2 (Figure 10(c), value = 0.0089) were associated with the clinical stage of COAD patients, which demonstrated that these 3 genes might be involved in the development of COAD. However, for potential applications in drug screening and prognostic performance, PER2 and CRY2 were selected to construct an mRNA-miRNA-lncRNA network. The miRNAs targeting PER2 and CRY2 identified by starBase and miRTarBase are shown in Figure 11(a). Among these miRNAs, hsa-miRNA-32-5p, hsa-miRNA-340-5p, hsa-miRNA-20b-5p, and hsa-let-7b-5p were identified as the top highly connected miRNAs (Figure 11(b)), which were upregulated or downregulated in COAD tissues. We then detected the expression profiles and prognostic value of these four genes. The results showed that only hsa-miRNA-32-5p was differentially expressed in tumors and was significantly associated with patient prognosis. More specifically, the expression of hsa-miRNA-32-5p was upregulated in COAD tissues compared with normal colon tissues ( value = , Figure 11(c)). COAD patients with high hsa-miRNA-32-5p levels showed worse OS ( value = 0.03453, Figure 11(d)). Thus, hsa-miRNA-32-5p was selected as the most promising miRNA target of PER2 and CRY2. Next, we further detected the upstream lncRNA targets using lncBase and starBase. According to the data from the above databases, a total of 11 lncRNA targets were identified (Figure 11(e)), including NDUFA6-AS1, SNHG14, CTBP1-AS2, LINC01128, INTS6-AS1, XIST, LINC00365, KCNQ1OT1, OIP5-AS1, DAPK1-IT1, and PITPNA-AS1. The expression levels of the above lncRNA targets in COAD tissues and normal colon tissues were analyzed. Similarly, the results showed that only lncRNA KCNQ1OT1 was differentially expressed in tumors and was significantly associated with patient prognosis. The data demonstrated upregulation of the lncRNA KCNQ1OT1 in COAD tissues compared with normal colon tissues ( value = , Figure 11(f)). COAD patients with high KCNQ1OT1 expression had a worse OS than those with low KCNQ1OT1 expression (median survival time = 5.3 years vs. 7.9 years, value = 0.039, Figure 11(g)). Thus, the lncRNA KCNQ1OT1/hsa-miRNA-32-5p/PER2/CRY2 regulatory axis might influence the occurrence and progression of COAD.

(a)

(b)

(c)

(a)

(b)

(c)

(d)

(e)

(f)

(g)

4. Discussion

In 1971, scientists first found three variations in one gene on the X chromosome of Drosophila melanogaster; one caused physiological response and behavior to become completely asynchronous with the circadian clock, and the other two caused the circadian clock to change to 19 and 28 hours, respectively [23]. By analyzing the entire X chromosome, researchers identified the location of the gene and named the gene PER. After this groundbreaking study, many studies attempted to expand the circadian clock gene family and determine the potential roles of circadian clock genes in the occurrence and development of disease and their potential applications in disease diagnosis, treatment and prevention. In the past 20 years, the circadian clock has been studied at the genetic level. The genes that participated in the construction of the circadian clock loop mainly included CLOCK, ARNTL, PER1, PER2, PER3, CRY1, CRY2, NR1D1, and RORA [6, 7]. Previous studies reported that dysfunction and abnormal expression of these genes were involved in the occurrence and progression of various human cancer types [8, 9, 24, 25]. Some studies have focused on the potential of core circadian clock genes to be cancer biomarkers [26–28]. Nevertheless, the specific functions of core circadian clock genes in the prognosis and therapy of COAD are rarely discussed. Thus, we adopted bioinformatics methods to illustrate the potential role of core circadian clock genes in COAD.

To a certain extent, the expression of core circadian clock genes reflects the status of the circadian clock and influences its functions [29]. In mammals, interlocked transcriptional-translational feedback loops regulate the circadian clock. In normal tissues, the positive limb of this loop is composed of CLOCK and ARNTL, which induce the expression of CRY1-2 and PER1-3. In turn, overexpressed PER and CRY proteins suppress their own transcription, thus forming the circadian rhythm [30, 31]. Compared to those in normal colon tissues, the transcriptional levels of CLOCK, CRY1, and NR1D1 in COAD tissues were significantly higher, while the mRNA levels of ARNTL, CRY2, PER1, PER3, and RORA in COAD tissues were downregulated in the current study. The dysexpression of these genes in COAD tissues might be involved in the disruption of circadian rhythms, thus leading to the occurrence and progression of COAD.

Researchers are also interested in the potential prognostic value of core circadian clock genes in COAD. In the present study, COAD patients with high CRY2 levels showed worse OS, DFS, DSS, and PFS than those with low CRY2 levels. COAD patients with high NR1D1 expression had worse DFS and PFS. In addition, high PER2 expression predicted a worse OS for COAD patients. A previous study based on human colorectal cancer tissues indicated that CRY2 expression was decreased and that high CRY2 expression predicted worse OS, which supported the results of the current study [26]. Na et al. indicated that NR1D1 was considered a biomarker predicting good prognosis in breast cancer patients [27]. PER2 is considered a prognostic predictor in several human cancer types, including lung cancer, colorectal cancer, and gastric cancer [28, 32, 33]. Moreover, the overall prognostic performances of CLOCK, PER1, PER3, and CRY1 were also determined in a previous meta-analysis [34].

Tumor-infiltrating immune cells play pivotal roles in eliminating cancer cells and hindering cancer progression. Exploring the association between core circadian clock genes and tumor-infiltrating immune cells is of significant importance. A previous study concluded the role of the circadian clock and its core genes in cancer immunity [35]. Xu et al. reported an M1-like proinflammatory phenotype of macrophages in mice with relatively disrupted PER1 and PER2 clock genes [36]. Cao et al. reported that CRY proteins regulated the process of autoimmunity, including B cell development, the BCR signaling pathway, and C1q expression, in CRY DKO (double knockout) mice [37]. However, more high-quality studies are needed to illustrate the mechanisms and specific processes of core circadian clock genes in cancer immunity. In addition, we detected the associations between the expression of core circadian clock genes in COAD and biomarkers in immune cells, which might assist researchers in determining potential immunotherapy targets for COAD. RORγ synthetic agonists could increase the expression levels of IL-17A, IL-17F, GM-CSF, and IL-22, thus boosting Th17 cells and Tregs to block immunosuppression [38]. The expression of BMAL1, CLOCK, Rev-Erbα, and PER2 in intraperitoneal macrophages regulated the expression levels of F4/80 and CD11c in tumor tissues [39]. However, immunotherapy based on core circadian clock gene targets is still at the initial stage. Thus, more complementary studies aimed at core circadian clock genes are urgently needed.

GO and KEGG pathway analyses could predict the biological functions and molecular mechanisms of the target gene set. In our study, the main functions and mechanisms of the 9 core circadian clock genes were involved in the regulation of circadian rhythms and some signaling pathways. The corticosteroid receptor signaling pathway includes the glucocorticoid receptor signaling pathway and mineralocorticoid receptor signaling pathway. In a mouse model, So et al. identified glucocorticoid response elements in the molecular circadian clock and found that PER2 was an integral component of the glucocorticoid regulatory pathway and was directly regulated by the glucocorticoid receptor [40]. Dickmeis et al. described multiple mechanisms by which molecular circadian clock genes regulate rhythms of glucocorticoid release and modulate glucocorticoid signaling [41]. For the mineralocorticoid receptor signaling pathway, circadian gene CLOCK signaling contributed to mineralocorticoid receptor-mediated cardiac inflammation and fibrosis [42]. The above findings revealed the significance of circadian clock genes in the corticosteroid receptor signaling pathway.

Through mRNA-miRNA-lncRNA network construction, the lncRNA KCNQ1OT1/hsa-miRNA-32-5p/PER2/CRY2 regulatory axis was identified. Previous studies have reported that miRNA-32-5p regulates the progression of gastric cancer and epithelial-mesenchymal transition and metastasis in lung adenocarcinoma [43, 44]. However, studies concerning the role of miRNA-32-5p in COAD are lacking. The upstream lncRNA KCNQ1OT1 was also detected and is considered a potential diagnostic biomarker for colon and rectal cancer and enhanced drug resistance in colon cancer cells [45, 46]. In the present study, we also studied the prognostic performance of miRNA-32-5p and lncRNA KCNQ1OT1 and found that they were correlated with COAD patient prognosis. The above findings revealed that the lncRNA KCNQ1OT1/hsa-miRNA-32-5p/PER2/CRY2 regulatory axis might be of significant importance in the occurrence and progression of COAD and may have potential applications in the diagnosis and treatment of COAD. Nevertheless, in vivo and vitro studies are urgently needed to clarify the role of this regulatory axis in COAD.

Some limitations still exist in the current study. First, our study and analysis were mainly performed at the transcriptional level, and the results may not be applicable in studies based on the protein level. Second, related fundamental and clinical studies that focus on the molecular mechanisms of circadian clock genes in colon or rectal cancer are rare. Third, the genetic background and etiology of COAD patients are influenced by many factors, such as patient race, sex, and age. Thus, more in-depth studies are necessary to validate the role of circadian clock genes in COAD.

5. Conclusion

In this study, a comprehensive bioinformatics analysis was used to identify the prognostic circadian clock gene signature in COAD, which included CRY2, NR1D1, and PER2. The lncRNA KCNQ1OT1/hsa-miRNA-32-5p/PER2/CRY2 regulatory axis was identified and might be of great significance in the occurrence and progression of COAD.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

All authors declare that no conflicts of interest in the present study.

Authors’ Contributions

Aoxiao He and Zhihao Huang drafted the manuscript and designed the research. Rongguiyi Zhang, Qian Feng, and Hongcheng Lu participated in statistical analysis. Qian Feng and Jiakun Wang reviewed the manuscript and provided the idea. Jiaqing Cao reviewed the manuscript and provided modifications in the review process. All the authors read and approved the final version. Aoxiao He and Zhihao Huang have contributed equally to this work.

Acknowledgments

The current study was supported by the Department of Education Project of Jiangxi Province (No. GJJ180083), the Science and Technology Project of Health Commission of Jiangxi Province (No. 20204266), and the Traditional Chinese Medicine Science and Technology Project of Jiangxi Province (No. 2019A015).

Supplementary Materials

Supplementary Figure 1: the PPI network and correlation analysis of core circadian clock genes in COAD. (Supplementary Materials)

References

H. Sung, J. Ferlay, R. L. Siegel et al., “Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries,” CA: a Cancer Journal for Clinicians, vol. 71, no. 3, pp. 209–249, 2021.
View at: Google Scholar
N. Dalal, R. Jalandra, M. Sharma et al., “Omics technologies for improved diagnosis and treatment of colorectal cancer: technical advancement and major perspectives,” Biomedicine & Pharmacotherapy, vol. 131, article 110648, 2020.
View at: Publisher Site | Google Scholar
O. Husson, B. H. de Rooij, J. Kieffer et al., “The EORTC QLQ-C30 summary score as prognostic factor for survival of patients with cancer in the "real-world": results from the population-based PROFILES registry,” The Oncologist, vol. 25, no. 4, pp. e722–e732, 2020.
View at: Google Scholar
K. Mody and T. Bekaii-Saab, “Clinical trials and progress in metastatic colon cancer,” Surgical Oncology Clinics of North America, vol. 27, no. 2, pp. 349–365, 2018.
View at: Publisher Site | Google Scholar
P. Favoriti, G. Carbone, M. Greco, F. Pirozzi, R. E. M. Pirozzi, and F. Corcione, “Worldwide burden of colorectal cancer: a review,” Updates in Surgery, vol. 68, no. 1, pp. 7–11, 2016.
View at: Google Scholar
C. Ramanathan, H. Xu, S. K. Khan et al., “Cell type-specific functions of period genes revealed by novel adipocyte and hepatocyte circadian clock models,” PLoS Genetics, vol. 10, no. 4, article e1004244, 2014.
View at: Google Scholar
R. Matsumura, Y. Tsuchiya, I. Tokuda et al., “The mammalian circadian clock protein period counteracts cryptochrome in phosphorylation dynamics of circadian locomotor output cycles kaput (CLOCK),” The Journal of Biological Chemistry, vol. 289, no. 46, pp. 32064–32072, 2014.
View at: Google Scholar
D. Z. Kochan and O. Kovalchuk, “Circadian disruption and breast cancer: an epigenetic link?” Oncotarget, vol. 6, no. 19, pp. 16866–16882, 2015.
View at: Google Scholar
Y. M. Lin, J. H. Chang, K. T. Yeh et al., “Disturbance of circadian gene expression in hepatocellular carcinoma,” Molecular Carcinogenesis, vol. 47, no. 12, pp. 925–933, 2008.
View at: Google Scholar
A. C. Silver, A. Arjona, M. E. Hughes, M. N. Nitabach, and E. Fikrig, “Circadian expression of clock genes in mouse macrophages, dendritic cells, and B cells,” Brain, Behavior, and Immunity, vol. 26, no. 3, pp. 407–413, 2012.
View at: Google Scholar
S. Hemmers and A. Y. Rudensky, “The cell-intrinsic circadian clock is dispensable for lymphocyte differentiation and function,” Cell Reports, vol. 11, no. 9, pp. 1339–1349, 2015.
View at: Google Scholar
G. L. Re, D. A. Santeufemia, F. L. Re et al., “Interleukin-2 chronotherapy for metastatic renal cell carcinoma: results of a phase I-II study,” Cytokine, vol. 128, article 154984, 2020.
View at: Google Scholar
X. Hu, X. Liu, J. Moisan et al., “Synthetic RORγ agonists regulate multiple pathways to enhance antitumor immunity,” Oncoimmunology, vol. 5, no. 12, article e1254854, 2016.
View at: Publisher Site | Google Scholar
I. K. Lee, H. Song, H. Kim et al., “RORα regulates cholesterol metabolism of CD8 T cells for anticancer immunity,” Cancers, vol. 12, no. 7, 2020.
View at: Publisher Site | Google Scholar
D. R. Rhodes, J. Yu, K. Shanker et al., “ONCOMINE: a cancer microarray database and integrated data-mining platform,” Neoplasia, vol. 6, no. 1, pp. 1–6, 2004.
View at: Google Scholar
A. Pizarro, K. Hayer, N. F. Lahens, and J. B. Hogenesch, “CircaDB: a database of mammalian circadian gene expression profiles,” Nucleic Acids Research, vol. 41, no. D1, pp. D1009–D1013, 2012.
View at: Publisher Site | Google Scholar
R. C. Anafi, L. J. Francey, J. B. Hogenesch, and J. Kim, “CYCLOPS reveals human transcriptional rhythms in health and disease,” Proceedings of the National Academy of Sciences of the United States of America, vol. 114, no. 20, pp. 5312–5317, 2017.
View at: Google Scholar
C. V. Mering, M. Huynen, D. Jaeggi, S. Schmidt, P. Bork, and B. Snel, “STRING: a database of predicted functional associations between proteins,” Nucleic Acids Research, vol. 31, no. 1, pp. 258–261, 2003.
View at: Google Scholar
T. Li, J. Fan, B. Wang et al., “TIMER: a web server for comprehensive analysis of tumor-infiltrating immune cells,” Cancer Research, vol. 77, no. 21, pp. e108–e110, 2017.
View at: Google Scholar
H. Y. Huang, Y. C. D. Lin, J. Li et al., “miRTarBase 2020: updates to the experimentally validated microRNA-target interaction database,” Nucleic Acids Research, vol. 48, no. D1, pp. D148–D154, 2020.
View at: Publisher Site | Google Scholar
J. Li, S. Liu, H. Zhou, L. H. Qu, and J. H. Yang, “starBase v2.0: decoding miRNA-ceRNA, miRNA-ncRNA and protein-RNA interaction networks from large-scale CLIP-Seq data,” Nucleic Acids Research, vol. 42, no. D1, pp. D92–D97, 2014.
View at: Publisher Site | Google Scholar
M. D. Paraskevopoulou, I. S. Vlachos, D. Karagkouni et al., “DIANA-LncBase v2: indexing microRNA targets on non-coding transcripts,” Nucleic Acids Research, vol. 44, no. D1, pp. D231–D238, 2016.
View at: Google Scholar
R. J. Konopka and S. Benzer, “Clock mutants of Drosophila melanogaster,” Proceedings of the National Academy of Sciences of the United States of America, vol. 68, no. 9, pp. 2112–2116, 1971.
View at: Google Scholar
N. Wang, M. Mi, X. Wei, and C. Sun, “Circadian clock gene Period2 suppresses human chronic myeloid leukemia cell proliferation,” Experimental and Therapeutic Medicine, vol. 20, no. 6, 2020.
View at: Publisher Site | Google Scholar
L. Zhou, Y. Yu, S. Sun, T. Zhang, and M. Wang, “Cry 1 regulates the clock gene network and promotes proliferation and migration via the Akt/P53/P21 pathway in human osteosarcoma cells,” Journal of Cancer, vol. 9, no. 14, pp. 2480–2491, 2018.
View at: Google Scholar
G. Mazzoccoli, T. Colangelo, A. Panza et al., “Deregulated expression of cryptochrome genes in human colorectal cancer,” Molecular Cancer, vol. 15, no. 1, 2016.
View at: Publisher Site | Google Scholar
H. Na, J. Han, N. L. Ka et al., “High expression of NR1D1 is associated with good prognosis in triple-negative breast cancer patients treated with chemotherapy,” Breast Cancer Research, vol. 21, no. 1, p. 127, 2019.
View at: Publisher Site | Google Scholar
B. Liu, K. Xu, Y. Jiang, and X. Li, “Aberrant expression of Per1, Per2 and Per3 and their prognostic relevance in non-small cell lung cancer,” International Journal of Clinical and Expermental Pathology, vol. 7, no. 11, pp. 7863–7871, 2014.
View at: Google Scholar
F. Hernández-Rosas, C. A. López-Rosas, and M. V. Saavedra-Vélez, “Disruption of the molecular circadian clock and cancer: an epigenetic link,” Biochemical Genetics, vol. 58, pp. 189–209, 2020.
View at: Google Scholar
A. M. Sangoram, L. Saez, M. P. Antoch et al., “Mammalian Circadian Autoregulatory Loop: a timeless ortholog and mPer1 interact and negatively regulate CLOCK-BMAL1-induced transcription,” Neuron, vol. 21, no. 5, pp. 1101–1113, 1998.
View at: Publisher Site | Google Scholar
M. Astiz, I. Heyde, and H. Oster, “Mechanisms of communication in the mammalian circadian timing system,” International Journal of Molecular Sciences, vol. 20, no. 2, 2019.
View at: Publisher Site | Google Scholar
T. Oshima, S. Takenoshita, M. Akaike et al., “Expression of circadian genes correlates with liver metastasis and outcomes in colorectal cancer,” Oncology Reports, vol. 25, no. 5, pp. 1439–1446, 2011.
View at: Google Scholar
H. Zhao, Z. L. Zeng, J. Yang et al., “Prognostic relevance of Period1 (Per1) and Period2 (Per2) expression in human gastric cancer,” International Journal of Clinical and Experimental Pathology, vol. 7, no. 2, pp. 619–630, 2014.
View at: Google Scholar
J. Zhang, H. Lv, M. Ji, Z. Wang, and W. Wu, “Low circadian clock genes expression in cancers: a meta-analysis of its association with clinicopathological features and prognosis,” PLoS One, vol. 15, no. 5, article e0233508, 2020.
View at: Publisher Site | Google Scholar
Z. Zhang, P. Zeng, W. Gao, Q. Zhou, T. Feng, and X. Tian, “Circadian clock: a regulator of the immunity in cancer,” Cell Communication and Signaling, vol. 19, no. 1, p. 37, 2021.
View at: Publisher Site | Google Scholar
H. Xu, H. Li, S. L. Woo et al., “Myeloid cell-specific disruption of Period1 and Period2 exacerbates diet-induced inflammation and insulin resistance,” Journal of Biological Chemistry, vol. 289, no. 23, pp. 16374–16388, 2014.
View at: Google Scholar
Q. Cao, X. Zhao, J. Bai et al., “Circadian clock cryptochrome proteins regulate autoimmunity,” Proceedings of the National Academy of Sciences of the United States of America, vol. 114, no. 47, pp. 12548–12553, 2017.
View at: Google Scholar
X. Hu, Y. Wang, L. Hao et al., “Sterol metabolism controls T_H17 differentiation by generating endogenous RORγ agonists,” Nature Chemical Biology, vol. 11, no. 2, pp. 141–147, 2015.
View at: Publisher Site | Google Scholar
A. A. S. Teixeira, L. A. Biondo, L. S. Silveira et al., “Doxorubicin modulated clock genes and cytokines in macrophages extracted from tumor-bearing mice,” Cancer Biology & Therapy, vol. 21, no. 4, pp. 344–353, 2020.
View at: Publisher Site | Google Scholar
A. Y. So, T. U. Bernal, M. L. Pillsbury, K. R. Yamamoto, and B. J. Feldman, “Glucocorticoid regulation of the circadian clock modulates glucose homeostasis,” Proceedings of the National Academy of Sciences of the United States of America, vol. 106, no. 41, pp. 17582–17587, 2009.
View at: Publisher Site | Google Scholar
T. Dickmeis, B. D. Weger, and M. Weger, “The circadian clock and glucocorticoids - interactions across many time scales,” Molecular and Cellular Endocrinology, vol. 380, no. 1-2, pp. 2–15, 2013.
View at: Publisher Site | Google Scholar
E. K. Fletcher, J. Morgan, D. R. Kennaway et al., “Deoxycorticosterone/salt-mediated cardiac inflammation and fibrosis are dependent on functional CLOCK signaling in male mice,” Endocrinology, vol. 158, no. 9, pp. 2906–2917, 2017.
View at: Google Scholar
Q. Wang, Y. He, W. Kan et al., “microRNA-32-5p targets KLF2 to promote gastric cancer by activating PI3K/AKT signaling pathway,” American Journal of Translational Research, vol. 11, no. 8, pp. 4895–4908, 2019.
View at: Google Scholar
J. Zhang, W. Yang, J. Wu et al., “MicroRNA-32-5p inhibits epithelial-mesenchymal transition and metastasis in lung adenocarcinoma by targeting SMAD family 3,” Journal of Cancer, vol. 12, no. 8, pp. 2258–2267, 2021.
View at: Publisher Site | Google Scholar
X. Qi, Y. Lin, X. Liu, J. Chen, and B. Shen, “Biomarker discovery for the carcinogenic heterogeneity between colon and rectal cancers based on lncRNA-associated ceRNA network analysis,” Frontiers in Oncology, vol. 10, article 535985, 2020.
View at: Publisher Site | Google Scholar
Y. Li, C. Li, D. Li, L. Yang, J. Jin, and B. Zhang, “lncRNA KCNQ1OT1 enhances the chemoresistance of oxaliplatin in colon cancer by targeting the miR-34a/ATG4B pathway,” Oncotargets and Therapy, vol. 12, pp. 2649–2660, 2019.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2022 Aoxiao He et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

473

Downloads

998

Citations