Circadian rhythms of core circadian genes in LUAD
Twelve circadian clock genes,(ARNTL, CLOCK, CRY1, CRY2, NR1D1, NR1D2, NPAS2, PER1, PER2, PER3, RORA, and TIMELESS) were used in this study. Through the Circadb database,we used RNA sequences to analyze the expression levels of 11 biological clock genes (except TIMELESS) at different times in lung tissue through Circadb database12(Pizarro et al., 2013). The results was shown that the expression levels of these genes significantly fluctuated with time (Supplementary Figure 1).
Defining circadian clock genes in LUAD and normal tissue
In the UALCAN database, compared to normal tissues, the expression levels of CLOCK, NPAS2, NR1D1, and TIMELESS were found to be up-regulated inLUAD patients, while the expression levels of CRY1, CRY2, NR1D2, PER1, PER2, PER3, and RORA were down-regulated, with statistically significant differences (P < 0.05) (Figure 1). The expression levels of circadian clock genes in LUAD tissues and normal tissues were detected by immunohistochemistry (IHC). Expression levels of CLOCK, NPAS2, NR1D1, and TIMELESS in LUAD tissues were found to be higher than thosein normal tissues while the the expression levels of CRY1, CRY2, and NR1D2 in LUAD tissues were lower than thosein normal tissues(Figure 2). Regarding the correlation between differentially expressed core circadian clock genes in LUAD and normal tissues, TIMELESS, NR1D1, NPAS2, PER1, and PER2 were partially negatively expressed in tumor tissues while NPAS2 and PER2 were partially negatively expressed in normal tissues(Figure 3).
Prognostic values of circadian clock genes in LUAD patients
Circadian genes have a definite prognostic value. Table 1 shows the models for predicting prognosis using circadian genes.. Up-regulated genes, including CRY1(HR =0.58, P＜0.001), CRY2(HR =0.51, P＜0.001), NR1D2(HR =0.49, P＜0.001), PER2(HR =0.44, P＜0.001), PER3(HR =0.5, P＜0.001), RORA (HR =0.45, P＜0.001) were correlated with improved overall survival; while up-regulated CLOCK (HR =1.4, P＜0.001), NPAS2(HR =2.19, P＜0.001), and TIMELESS (HR =1.88, P＜0.001) were correlated poor survival outcomes (Supplementary Figure 2). The time to first progression of the following up-regulated genes was also found to be increased: CRY1(HR =0.61, P=0.002), CRY2(HR =0.5, P＜0.001), NR1D2(HR =0.63, P＜0.001), PER1(HR =0.58, P＜0.001), PER2(HR=0.51, P＜0.001), PER3(HR =0.47, P＜0.001), and RORA (HR =0.57, P＜0.001); while the time to first progression of the following up-regulated genes was decreased: NPAS2(HR =1.74, P＜0.001), NR1D1(HR =1.54, P=0.007), and TIMELESS (HR =1.8, P＜0.001) (Supplementary Figure 3). We also validated the Overall Survival (OS) of the genes of the genes by the PrognoScan database(Table 2). There was a significant correlation between the expression levels of CLOCK, CRY2, NPAS2, NR1D1, NR1D2, PER1, PER3, RORA, and TIMELESS and OS. Moreover, the OS prognostic values CRY2, NPAS2, NR1D2, PER3, RORA, and TIMELESS genes were comparable to those of the Kaplan-Meier Plotter database. LUAD patients with elevated expressions of CRY2, NR1D2, PER3, and RORA and suppressed expressions of NPAS2 and TIMELESS exhibited better survival outcomes(P<0.05).
Circadian clock genes associated with cancer hallmarks in LUAD
To reveal the potential impact of circadian dysregulation in LUAD, we evaluated the role of the circadian clock in LUAD. SNV analysis revealed that missense mutations and splicing sites were the most common types of circadian clock gene variations in LUAD (Figure 4).In addition, the methylation levels of CRY2 and NR1D1DNA in LUAD tissues were lower than those in normal tissues, while the methylation levels of NPAS2, PER1, PER3, and RORA in LUAD tissues were higher than those in normal tissues, with statistically significant differences (P < 0.05)(Figure 5).
Then, we determined the expression levels of circadian clock genes in the typical cancer-related signaling pathways (Apoptosis, Cell Cycle, DNA Damage Response, Epithelial-mesenchymal transition (EMT), Hormone androgen receptor (AR), Estrogen receptor (ER), Phosphatidylinositol 3‑kinase (PI3K)/Protein kinase B (AKT), RAS/ Mitogen-activated protein kinase (MAPK), Receptor tyrosine kinase (RTK), Tuberous sclerosis complex (TSC)/ Mammalian target of rapamycin (Mtor)). The circadian clock genes were closely associated with most of the typical cancer-related marker pathways. Among them, CRY2 activates PI3K/AKT. RORA, PER2, CRY2, and NPAS2 activate RAS/MAPK, while most of the Circadian clock genes are involved in inhibiting Apoptosis, Cell cycle, and DNA Damage Response Pathway( Figure 6A). We also established the correlation between circadian clock gene expression and drug sensitivity. Elevated expression levels of CLOCK, CRY1, and NR1D2 as well as suppressed expression levels of PER2 and CRY2 were associated with drug resistance ( Figure 6B).
Enrichment analysis of circadian clock in LUAD
We used GO and KEGG enrichment analysis to evaluate the potential roles of the circadian clock genes. These genes exert a major impact on the circadian rhythm and biological process, such as on the transcription factor complex of cellular component, and on transcription corepressor binding, transcription cofactor binding, and E-box binding of molecular function (Figure 7A). KEGG enrichment analysis showed that circadian clock genes are involved in circadian regulation of gene expression, circadian rhythms, and rhythmic process (Figure 7B). PPI network results showed that they are associated with circadian rhythm, rhythmic Process, circadian regulation of gene expression, and regulation of circadian rhythm and regulatory region DNA binding (Figure 7C). These findings imply that circadian rhythm genes can regulate the transcription factor activity and circadian rhythm in LUAD.
Immune infiltrates analysis of circadian clock in LUAD patients
Few studies have evaluated the role of the biological clock in the LUAD immune microenvironment. We analyzed the correlations between clock genes and the abundance of six types of immune cells(Table 3). The expression levels of ARNTL and CLOCK were positively correlated with the abundance of infiltration of B cells (r =0.158, P＜0.001/ r =0.287, P＜0.001), CD8+T cells (r =0.221, P＜0.001/ r =0.144, P =2.45e-03), CD4+T cells (r =0.232, P＜0.001/ r =0.272, P＜0.001), Macrophages (r =0.207, P＜0.001/ r =0.478, P＜0.001), Neutrophils (r =0.321, P＜0.001/ r =0.446, P＜0.001), and Dendritic cells (r =0.315, P＜0.001/ r =0.361, P＜0.001) (Supplementary Figure 4A/B). However, CRY1 expression was only positively correlated with the infiltration abundance of Macrophage (r = 0.134, P = 3.03e-03) and Neutrophils (r = 0.155, P＜0.001) (Supplementary Figure 4C). And CRY2 expression also was not correlated with CD8 + T cells and neutrophil, but positively correlated with B cells (r = 0.281, P＜0.001), CD4 + T cells (r = 0.297, P＜0.001), Macrophages (r = 0.171, P＜0.001) and Dendritic cells (r = 0.187, P＜0.001) (Supplementary Figure 4D). In addition, the expression levels of NPAS2 and NR1D1 were positively correlated with CD4 + T cells (r = 0.187, P＜0.001/ r= 0.251, P＜0.001) and Dendritic cells (r = 0.118, P = 9.33e-03/ r= 0.24, P＜0.001) (Supplementary Figure 4E/F). As found in ARNTL and CLOCK, NR1D2 expression was positively correlated with B cells (r = 0.204, P＜0.001/ r = 0.339, P＜0.001/ r = 0.368, P＜0.001), CD8 + T cells (r = 0.269, P＜0.001/ r = 0.134, P＜0.001/ r = 0.258, P＜0.001), CD4 + T cells (r = 0.239, P＜0.001/ r = 0.4, P＜0.001/ r =0.394, P＜0.001), Macrophages (r = 0.259, P＜0.001/ r = 0.23, P＜0.001/ r = 0.23, P＜0.001), Neutrophils (r = 0.298, P＜0.001/ r = 0.232, P＜0.001/ r = 0.325, P＜0.001) and Dendritic cells (r = 0.348, P＜0.001/ r = 0.419, P＜0.001/ r = 0.419, P＜0.001) (Supplementary Figure 4G/J/K). On the contrary, PER1 expression was not correlated with the infiltration abundance of immune cells (Supplementary Figure 4H). PER2 expression was only positively correlated with CD4 + T cells (r = 0.163, P＜0.001), but was not correlated with the other immune cells (Supplementary Figure 4I). As for TIMELESS, its expression positively correlates the abundance of Neutrophils (r = 0.114, P＜0.001) (Supplementary Figure 4L).
Besides, changes in somatic copy number (SCNA) of circadian rhythm genes inhibited immune infiltration, including CD8+T cells, neutrophils, dendritic cells, macrophages, CD4+T cells, and B cells. These results suggest a possible association between the circadian clock and immune infiltration in LUAD (Figure 8). We also assessed the association between clock genes with a different expression and immune cell infiltration. Using the Cox proportional risk model, we adjusted for confounding factors: B cells, CD4 + T cells, macrophages, neutrophils, dendritic cells, ARNTL, CLOCK, CRY1, CRY2, NR1D1, NR1D2, PER1, PER2, RORA, NPAS2, and TIMELESS. LUAD patients were correlated with elevated expression levels of B cells (P < 0.001), CD4+ T cells (P= 0.047), CRY2 (P= 0.02), PER3 (P=0.014), RORA (P= 0.013), NPAS2 (P < 0.001), and TIMELESS (P= 0.025) (Table 4).