Identification of circRNA-miRNA-mRNA networks to explore the molecular mechanism and immune regulation of postoperative neurocognitive disorder

Postoperative neurocognitive disorder (PND) is a common complication in older patients. However, its pathogenesis has still remained elusive. Recent studies have shown that circular RNA (circRNA) plays an important role in the development of neurodegenerative diseases, such as PND after surgery. CircRNA, as a competitive endogenous RNA (ceRNA), mainly acts as a molecular sponge for miRNA to “adsorb” microRNA (miRNA) and to reduce the inhibitory effects of miRNAs on target mRNA. The sequencing data of circRNA were obtained from the Gene Expression Omnibus (GEO) database. By bioinformatic methods, circAtlas, miRDB, miRTarBase and miRwalk databases were applied to construct circRNA-miRNA-mRNA networks and screen differentially expressed mRNAs. To improve the accuracy of the data, we randomly divided aging mice into control (non-PND group) and PND groups, and used high-throughput sequencing to analyze their brain hippocampal tissue for analysis. Three key genes were cross-detected in the data of both groups, which were Unc13c, Tbx20 and St8sia2 (as hub genes), providing new targets for PND treatment. According to the results of the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses, immune cell infiltration analysis, gene set enrichment analysis (GSEA), Connectivity Map (CMap) analysis, quantitative real-time polymerase chain reaction (qRT-PCR), the genes that were not related to the central nervous system were removed, and finally, mmu_circ_0000331/miR-1224-3p/Unc13c and mmu_circ_0000406/miR-24-3p/St8sia2 ceRNA networks were identified. In addition, the CMap method was used to select the top 4 active compounds with the largest negative correlation absolute values, including cimaterol, Rucaparib, FG-7142, and Hydrocortisone.

after surgery was almost twice that of patients without PND [4]. The risk of PND in old patients was reported to be 25-40% [5]. The global prevalence of dementia was 46.8 million in 2015, and it is expected to increase to 131.5 million by 2050. The global cost of dementia was estimated to be $818 billion in 2015 [6]. PND leads to delayed recovery, prolonged hospitalization, increased medical costs, enhanced complications, and even loss of self-care ability, resulting in a series of medical, social, and economic problems [7]. Several studies have explored the possible mechanism of PND, and neuroinflammation [8], neuronal apoptosis [9], autophagy disorders [10], and synaptic plasticity [11] have been reported. At present, the pathogenesis of PND is not clear, and there is no effective treatment and specific marker for PND.
With the rapid development of molecular biology and transcriptology, the biological functions and regulatory mechanisms of non-coding RNA can explain the occurrence and development of several complex diseases. Circular RNA (circRNA) is a newly discovered class of non-coding RNA with a length of about 100 nucleotides. The circRNA molecule has no PolyA tail with 5′ and 3′ ends, and it has a closed ring structure with a high stability. The roles of circRNAs in various diseases have been widely studied due to their unique molecular mechanism and molecular function, and have become a research hotspot. Their main function is to act as a molecular sponge for microRNAs (miRNAs) to "adsorb" miRNAs and to reduce the inhibition of miRNAs on target mRNAs. The mechanism of competitive endogenous RNAs (ceRNAs) in the occurrence and development of diverse diseases has been most thoroughly studied. Studies have reported that ceRNA network could have a great potential in the treatment of the central nervous system (CNS) diseases [12]. Recent studies have shown that Figure 1. Flowchart of the comprehensive analysis process and all methods utilized in present study. AGING circRNA_22058 and circRNA_44122/EGFR ceRNA network and circRNA_22673/Prkacb ceRNA could be the mechanism of PND [13]. However, the role of circRNA in elderly patients with PND has still remained elusive.
Bioinformatics research has revealed the molecular targets of PND mechanism in a systematic, accurate, and effective manner, and clarified the theoretical basis of PND occurrence. In the present study, prefrontal cortex data of PND in aged mice obtained from GEO database were screened, and a circRNA-miRNA-mRNA network was constructed by bioinformatics method to identify differential mRNAs. To further improve the accuracy of data, the hippocampus of splenectomized aged mice were sent to Novogene Co., Ltd. (Beijing, China) for high-throughput sequencing of mRNA. The two sets of data were intersected to detect hub genes. Surgery and trauma can trigger inflammatory responses characterized by both pro-inflammatory and antiinflammatory cytokines released by the immune system [14]. A large number of studies have shown that immune cells trigger neuroinflammatory responses, leading to postoperative cognitive dysfunction [15]. Immune regulation has become a hotspot in the study of PND mechanism. The present research aimed to identify immune-related genes. The flowchart of this study is shown in Figure 1.

Animals
Healthy male C57BL/6 mice (age, 18-month-old; body weight, 25-30 g) were purchased from Beijing Weitong Lihua Co., Ltd. (Beijing, China). All animals were acclimated for one week before the experiment. All animal experiments were performed in accordance with the Guide for Care and Use of Laboratory Animals, and approved by the Animal Care and Use Committee of the Southern Medical University (Guangzhou, China). The aging mice were divided into two groups: control group (n = 11) and surgery group (n = 11). Mice in the surgery group received splenectomy under general anesthesia. Among them, 16 mice underwent Morris water maze (MWM) test to determine the cognitive dysfunction in aging mice after surgery, and the hippocampus of the remaining 6 mice were taken out for transcriptome detection.

Establishment of animal models
Splenectomy was performed under inhalation anesthesia (maintained by a mixture of 2.0-2.5% isoflurane and oxygen). The anal temperature was maintained between 38.8-39.8°C. After undergoing anesthesia, the left subcostal incision was made. The spleen was separated and exposed, the spleen arteries and veins in the splenic pedicle region were firmly ligated, and the spleen was removed. After confirming no active bleeding, the abdomen was closed layer-by-layer with silk thread. Mice in the control group were not anesthetized or underwent surgery. Splenectomy leads to reversible postoperative learning and memory dysfunction and subsequently induces Postoperative cognitive dysfunction (POCD) [16], and is often used as a model for establishing POCD [17].

MWM test
MWM test was employed to assess abilities of mice in learning, memory, and spatial cognition. A pool with a diameter of 120 cm and a depth of 50 cm was divided into four quadrants, and a platform with a diameter of 12 cm was set 1 cm underwater in the target quadrant. Experimental mice first received positioning navigation training for 5 days. Mice were placed on the platform for 30 s, and were put into the pool from four different quadrants (N, W, S, and E). Observe and record the distance, time and speed of mice reaching the platform for the first time. If mice could not find the platform within 60 s, they were guided to the platform for 15 s. The interval between the two training sessions was at least 30 s. The time of reaching the platform for the first time, the number of times that the distance and speed could cross the platform quadrant, and the number of times that they crossed the platform quadrant were recorded. The experiment was lasted for at least 30 s. The time to reach the platform for the first time, the distance and speed to cross the platform, and the number of times to cross the platform quadrant were recorded. On day 6 and postoperative days 3, 7, and 14, the platform was removed and the MWM assay was performed to record the distance, speed, and time spent across the platform for the first time and the percentage of time spent across the target quadrant. Finally, the time point of the most severe PND in the aging mice was selected as the time point for taking the materials. The test performers and data analyzers were blinded to the mice groups.

Sample collection and sequencing
Under isoflurane deep anesthesia, 6 mice were killed by decapital method, including 3 mice in control group and 3 mice in surgery group. The skull was separated by non-sharp method, the hippocampus were placed on ice bags, and then, the hippocampus were frozen in liquid nitrogen and stored in a −80°C refrigerator. The hippocampus were randomly selected for subsequent sequencing analysis. Total RNA was extracted from mouse hippocampus using TRIzol ® reagent (Magen, AGING Guangzhou, China). Nanodrop ND-2000 (Thermo Fisher Scientific, Waltham, MA, USA) was used to detect the A260/A280 absorbance ratio of RNA samples. An Agilent 4150 Bioanalyzer (Agilent Technologies, CA, USA) was used to measure RIN values of RNA. A PE library was prepared according to ABclonal mrna-SEQ Lib Prep Kit (ABclonal, China) specification. Sequencing was performed using the Illumina Novaseq 6000/MGISEQ-T7 sequencing platform (Illumina Inc., Chicago, IL, USA).

Data collection and identification of DEcircRNAs and DEmRNAs
In the present study, GSE174410 dataset was downloaded from the GEO database (http://www.ncbi.nlm.nih.gov/geo), as well as annotation platforms for GPL21273, including 3 normal tissue samples and 3 post-operative neurocognitive impairment tissue samples, accounting for a total of 6 samples. In the present experiment, high-throughput sequencing was performed on the mRNA of hippocampus of aging mice, including 6 samples from control group (n = 3) and PND group (n = 3). The "Limma" package was used to standardize the data of the two groups, and the circRNA/mRNA was significantly differentially expressed in the control group and the PND group. A gene with |LogFC|>1 and P < 0.05 was considered to be significantly differentially expressed, and it was displayed through a volcano plot.

Biological functional analysis
In order to obtain the biological functions and signaling pathways involved in the occurrence and development of PND, the Metascape database (http://www.metascape.org) was used to enrich and analyze 530 mRNAs with differences in the experimental data for annotation and visualization. Specific genes were analyzed by the Gene Ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses. Min overlap ≥3 and P ≤ 0.05 were considered statistically significant.

Construction of the circRNA-miRNA-mRNA network
We used circAtlas database to predict circRNA-miRNA interaction pairs. In addition, the interaction between miRNA and mRNA was predicted by combining data collected from miRDB, miRTarBase, and miRwalk databases. Targeted mRNAs identified by more than two databases were selected for further analysis. Then, the circRNA-miRNA-mRNA network was established by combining circRNA-miRNA interaction and mRNA-miRNA interaction, and Cytoscape software was used for visualization of biological networks.

Immune infiltration analysis
The CIBERSORT algorithm was used to analyze RNAseq data of mice in the control and PND groups to infer the relative proportions of 25 immune infiltrating cells. The "Pheatmap" package was used to draw the heatmap of immune cell infiltration to explore the distribution of immune cells. The "Corrplot" package was utilized to analyze the interaction between immune cells and to further analyze the influence of the interaction between immune cells. The "Vioplot" package was employed to plot the relative content of immune cells. P < 0.05 was considered statistically significant. Then, the Spearman correlation analysis was conducted between HUB gene expression level and immune cell content, and "GGploT2" package was used for visualization. Chemokines refer to small cytokines or signal proteins that have the function of making immune cells targeted to chemotaxis, which can control the migration and residence of all immune cells. Autophagy, as a defense mechanism of the collective, in the body's immunity, inflammation, neural degeneration disease, aging, etc., has shown a very important role in the pathogenesis, and it is closely associated with the immune cells. We, in the present study, further discussed the hub genes and chemotactic factors, as well as the relationship between the autophagy using the Pearson correlation analysis.

Gene set enrichment analysis (GSEA) of hub genes
Unc13c, Tbx20, and St8sia2 mRNA were analyzed by the GSEA. GSEA analyzes the expressions of a group of functionally-related genes based on gene expression profile data. The principle of GSEA is to search for sets of genes that are significantly over-represented in a given list of genes. According to the number of genes contained in the gene set and the expressions of genes, the normalized enrichment score (NES), false discovery rate (FDR), and adjusted P value of the signaling pathway were calculated. Genes with P < 0.05 were considered significantly enriched after 1000 permutations.

Identification of co-expressed genes (Unc13c, Tbx20, and St8sia2)
The PPI networks of Unc13c, Tbx20, and St8sia2 were constructed and visualized separately based on STRING database. The correlation analysis was performed between hub genes and their co-expressed gene expression levels.

RNA extraction and quantitative real-time polymerase chain reaction (qRT-PCR)
RNA was extracted and validated from the cerebral cortex of 5 control and 5 PND aging mice using qRT-PCR. Total RNA was extracted using TRIzol reagent (Cat. No. 1596-026; Thermo Fisher Scientific). The reverse transcription reaction system was prepared according to the instructions of Bestar qPCR RT Kit.
The total system was 20 μL, and the first strand of cDNA was synthesized. Follow the instructions presented for DBI Bestar SybrGreen qPCR masterMix reaction kit, annealing was carried out at 94°C for 2 min, at 94°C for 15 s, and at 60°C for 15 s, at 72°C for 15 s, with 40 cycles. GAPDH was used as an internal reference gene. A QuantStudio3 fluorimeter (Thermo Fisher Scientific) was used for fluorescence quantitative analysis. Relative gene expression was calculated by the 2 −ΔΔCT method. The primers are shown in Table 1.

Connectivity Map (CMap) analysis
CMap is a database of expression profiles based on the expression of intervening genes developed by the Broad Institute. Lamp proposed that it could be used to discover the association between drugs, genes and diseases, and constructed the CMap database when the gene changes were obtained after cell disturbance [18]. CMap was used to compare the data list of differentially expressed genes measured by our research group with the database reference data set. Finally, the database will obtain a correlation score (−100-100) according to the enrichment of differentially expressed genes in the reference gene expression spectrum. The positive value indicated that the up-regulated and down-regulated differential expressions were similar, while the negative value indicated that the up-regulated and down-regulated differential expressions were opposite, and the reference gene expression was sequenced according to the value.

Statistical analysis
All data were analyzed using Graphpad Prism 8 software, and the measures were expressed according to the mean ± SD. Data were compared between groups using analysis of variance (Two-way ANOVA) P-value < 0.05 was considered statistically significant with great importance.

Availability of data and materials
Part of datasets generated during and/or analyzed during the current study are available in the Gene Expression Omnibus (GEO) datasets (http://www.ncbi.nlm.nih.gov/geo/) and the other part of the data came from the high-sequencing results of our experiment.

Behavioral test results
During the 6 days of preoperative training, all experimental mice showed a gradual downward luminal trend in the time and distance to find the platform, indicating that the experimental mice were able to improve their spatial learning memory ability by constantly repeating finding the platform. The aging mice showed a significantly impaired memory on postoperative day 3, as the time and distance to find the platform ( Figure 2B, 2C) both significantly increased; the number of times to pass the platform within 60 s and the percentage of target quadrants ( Figure 2D, 2E) both significantly decreased, while the average speed did not significantly differ ( Figure 2F), and basically recovered on postoperative day 14. Comparing postoperative day 3 with preoperative training day 6, postoperative day 7, and postoperative day 14, postoperative day 3 was the day with the worst memory ability and the most obvious cognitive impairment, as evidenced by a significant increase in the time and distance to find the platform for the first time and a significant decrease in the number of times to pass the platform within 60 s and the percentage of target quadrants. Therefore, we selected the hippocampus of aging mice for transcriptomic assay on the 3rd postoperative day when memory impairment was most severe. Timeline diagram of the Morris Behavior Test is shown in Figure 2A.

Differential expression of circRNAs and mRNAs
A total of 18 differentially expressed circRNAs, including 17 up-regulated circRNAs and 1 downregulated circRNA, were identified from the circRNA transcriptome data of GSE174410 downloaded from the GEO database ( Figure 3A). Compared with the control group, the difference in circRNA expression in the surgical group was statistically significant. A total of 530 differentially expressed mRNAs, including 266 upregulated mRNAs and 264 down-regulated mRNAs, were screened, which are represented as volcano plots ( Figure 3B).

The GO and KEGG pathway enrichment analyses
We further conducted the GO and KEGG pathway enrichment analyses of 530 differentially expressed mRNAs through the Metascape database, and the results showed that, the main gene enrichment pathways were pheromone receptor activity, DNA-binding transcription factor activity, RNA polymerase II-specific and B cell proliferation, female sex differentiation, spinal cord development, anion transport, positive regulation of transmembrane receptor serine/threonine kinase signaling pathway, detection of chemical stimulus, and endocrine system development ( Figure 3C).

Establishment of the circRNA-miRNA-mRNA network
In order to find out the downstream genes of circRNA and to construct the ceRNA network, we first used the online circAtlas database (circatlas.biols.ac.cn) to predict miRNA targets for the 5 differentially up-regulated circRNAs with the largest change. The results showed that there were 96 circRNA-related targeted miRNAs, including 105 circRNA-miRNA pairs ( Figure 4A). A total of 4,386 mRNAs and 6,725 pairs of miRNA-mRNAs were predicted ( Figure 4B). Besides, miRwalk, miRDB, and miRTarBase databases were further used to predict miRNA-related targeted mRNAs, the list of miRNA-mRNA was shown in Appendix 1.
Finally, the intersection of 4386 predicted mRNAs and 530 differentially measured mRNAs was obtained, and a total of 13 mRNAs were detected, a Venn diagram was shown in Figure 4C. Then, 19 ceRNA networks were successfully constructed for circRNA-miRNA-mRNA pairs and visualized by the Cytoscape software ( Figure 4D). From the Figure 4E we can see that Unc13c, Tbx20, and St8sia2 are located in the largest yellow, pink, and purple regions respectively, so we selected these three genes with the most miRNA relevance as the core genes (mRNA relationship pairs ≥3). We conducted ROC curve analyses to indicate whether three hub genes could better predict PND. The results showed that AUC values of the three hub genes were 1.000 ( Figure 4F), 0.889 ( Figure 4G), and 1.000 ( Figure 4H), respectively. These three hub genes had high diagnostic potential in differentiating aging mice with PND from healthy control mice.

Immune cell infiltration analysis
Immune microenvironment is mainly composed of immune cells, extracellular matrix, various growth factors, inflammatory factors, and special physical and chemical characteristics, significantly influencing the diagnosis and clinical treatment of diverse diseases. After analyzing the relationship between differentially expressed genes and immune invasion in PND data, the potential molecular mechanism of differentially expressed genes influencing the progression of immune invasion was further discussed. The results of immune cell infiltration assay showed the contents of immune cells in each sample, of which M0 macrophage and M2 macrophage accounted for the highest proportions ( Figure 5A), and the correlation between immune cells is shown in Figure 5B. The number of mast cells was significantly lower in patients with disease than that in normal patients ( Figure 5C).
We additionally explored the relationship between core genes and immune cells, and the results showed that UNC13C and St8sia2 were significantly positively correlated with Treg cells and mast cells, respectively, and significantly negatively correlated with natural-killer (NK) cells and naive CD4 T cells ( Figure 5D, 5F). Tbx20 expression in surgery group significantly increased, it was significantly positively correlated with the number of naive CD4 T cells, and it was significantly negatively correlated with the number of mast cells ( Figure 5E). These results confirmed that key genes were closely related to the level of immune cell infiltration and play key roles in the immune microenvironment.

Correlation analysis of hub genes with chemokines and autophagy-related genes
Correlation analysis of hub genes and chemokines showed that Tbx20 expression was significantly positively correlated with Ccl1 expression (Pearson R = 0.98), while Tbx20 expression was significantly negatively correlated with CCL9 expression (Pearson R = −0.83) ( Figure 6A). Correlation analysis of hub genes and autophagy-related genes showed that Tbx20 and Epg5 were significantly negatively correlated (Pearson R = −0.91), while Unc13c and Atg4a were significantly positively correlated (Pearson R = 0.84) ( Figure 6B).

GSEA
GSEA of the three hub genes showed that high expression of Unc13c was mainly enriched in TUBULIN_BINDING, n-acetyltransferase_activity, VESICLE, and other pathways. The high expression of Tbx20 was mainly concentrated in FK506_BINDING, Protein_Peptidyl-prolyl_isomerization, MIDDLE_ EAR_MORPHOGENESIS, and other signaling pathways. The high expression of St8sia2 was mainly concentrated in gamma-tubulin_binding, VESICLR, and response_to_intracestar signaling pathways. Besides, Unc13c, Tbx20, and St8sia2 could regulate the development of PND and affect the course of disease and prognosis of elderly mice through these signaling pathways (Figure 7).

Identification of the co-expressed genes
The co-expressed genes (St8sia2, Tbx20, and Unc13c) were retrieved from the String database, and the confidence score was set to 0.4 to establish the corresponding protein interaction network ( Figure 8A-8C). In addition, the correlation analysis of the three hub genes and their corresponding co-expressed genes was performed separately, and the results showed that Unc13c and Rims2 were significantly positively correlated, and Tbx20 and Myh6 were significantly negatively correlated. The correlation coefficient between Unc13c and Rims2 was 0.965 (P = 0.002) ( Figure 8D). The correlation coefficient between Tbx20 and Myh6 was 0.823 (P = 0.044) ( Figure 8E).

Identification of four bioactive compounds by the CMap analysis
The CMap method was used to select the top 4 active compounds with the largest negative correlation absolute values, including cimaterol, Rucaparib, FG-7142, and Hydrocortisone ( Figure 10A-10D), as shown in Appendix 2. These four active compounds could be therapeutic targets of PND.

DISCUSSION
The course of delirium is acute in PND, mainly appearing 1-3 days after surgery and anesthesia, which is coincident with establishing aging mouse models with PND to screen the maximum severity of PND. However, devastating consequences were obtained, such as the increased mortality in the first year after surgery [19], the decreased quality of life, and the increased long-term risk of Alzheimer's disease (AD) [20]. In adult mice, circRNA is more abundant in the brain than that in other organs (e.g., heart, liver, and lung) [21]. Cerebral circRNA enrichment is associated with the neurotransmitter function, neuronal maturation, and synaptic activity [22]. CircRNAs have been reported to target aging-related mRNAs in the brain and to regulate the aging process by changing the expressions and availability of specific mRNAs [23]. Based on the important role of circRNAs in regulating the CNS diseases, the current experiment first extracted circRNAs from the GSE174410 dataset, which was combined with our experimental data. Four circRNAs were detected, including mmu_cirC_0000400, mmu_cirC_0000331, mmu_cirC_0000406, and mmu_cirC_0000798, as well as 3 hub genes (Unc13c, Tbx20, and St8sia2).
We first reviewed the literature to understand the three hub genes. Unc13c was widely expressed in 11 tissues, including brain, cerebellum, heart, and liver [24]. Unc13c is a protective gene against AD, and it may be involved in synaptic plasticity and synaptic transmission [25]. This was consistent with the results of the present study, in which Unc13c was down-regulated in the PND group. Unc13c was significantly down-regulated in spinal cord tissues of patients with amyotrophic lateral sclerosis [26]. It was reported as one of the hub genes in a genome-wide association study of post-traumatic stress disorder in Iraq-Afghanistan veterans [27]. Unc13c has also been found to be associated with neurodegeneration in a genetic dementia in Finland [28]. These studies suggested that Unc13c could be closely associated with the CNS, and it could be a target gene for reducing PND. Tbx20 has been detected in unique neurons and epithelial cells of mouse and human embryonic eye tissues [29], associating with the CNS development [30]. Song studied Tbx20 as a gene associated with motor neurons. However, most studies on Tbx20 have concentrated on the cardiovascular system, and Tbx20 was found to play a variety of basic roles in cardiovascular development and homeostasis in response to pathophysiological stress, as well as cardiac remodeling [31,32]. However, the relationship between Tbx20 and PND needs to be studied by more experiments. St8sia2 syntheses polysialic acid (PSA), which is crucial for cerebral development and is closely related to synaptic plasticity [33]. Lukasz et al. suggested that St8sia2 is a candidate gene that is associated with brain development and plasticity in schizophrenia. Sebastian et al. demonstrated that St8sia2 promotes oligodendrocyte differentiation and myelin repair [34]. Reported in a patient with behavioral disorders, epilepsy and autism spectrum disorders, St8sia2 was expressed in the developing brain, and it showed to play an important role in neuronal migration, axon guidance, and synaptic plasticity [35]. A study on the role of neurodevelopment and genes in psychiatric comorbidity and the regulation of inflammatory process in AD revealed that St8sia2 was only associated with the existence of clinical depression, rather than with AD, and St8sia2 could regulate inflammatory factors (IL-6, IL-1β, etc.) [36].  TANTRA [37] showed that St8sia2 could protect juvenile mice from delayed cognitive impairment several years after cannabis exposure. Taken together, St8sia2 could regulate the cognitive function of the brain, while the specific regulatory mechanism needs more experimental studies. As a result, we found that Unc13c and St8sia2 were closely associated with the CNS and cognitive function, and they could regulate PND.
The results of GO and KEGG pathway enrichment analyses revealed that, 530 differentially expressed genes were enriched in B cell proliferation and transmembrane receptor protein serine/threonine kinase signaling pathway functions. This suggested that PND could be related to immunity. Few studies have systematically screened biomarkers associated with PND immune infiltration. To further explore the role of immune cell infiltration in PND, CIBERSORT was used to comprehensively evaluate PND. The results showed that the number of mast cells in PND group was significantly reduced compared with that in the normal group. This is consistent with the results of a study on depression: isolated social stress caused depression in mice, significantly reducing the total number of mast cells in the brain by 90% on the first day [38]. There are also conflicting and controversial findings about the role of mast cells.
On the one hand, mast cells seem to mediate neuroinflammation protection. Mast cells release specific lyase McP-4 to degrade inflammatory cytokines in a mouse model of traumatic brain injury (TBI) to modulate the CNS disease [39]. Some proteases released by mast cells play homeostasis, protection, and even anti-inflammatory roles [40]. On the other hand, palmitoethanolamide (PEA) induced attenuation of the number of mast cells, amylase and trypsin, and increased cerebral edema, infarct volume, and brain injury in experimental TBI mice [41]. Mast cells release mediators into the CNS, promoting neurogenesis, such as serotonin and IL-6, providing neuroprotection (e.g., IL-1β), and maintaining the integrity of the blood-brain barrier (BBB), such as histamine. However, excessive levels of these mediators have detrimental effects on the integrity of neurons and the BBB. These results seem to be contradictory, whereas the effect of inflammation on TBI may be beneficial or detrimental depending on the time after injury and the stage of TBI. Logically speaking, the role of mast cells in these different stages is two-fold [42]. The results of the correlation between Unc13c and St8sia2 and immune cells showed that the combination of Unc13c and St8sia2 significantly decreased in the surgery group compared with that in the control group, suggesting that Unc13c and St8sia2 could protect PND by regulating immunity. In recent years, the role of T cells in neurodegenerative diseases has noticeably attracted scholars' attention [15]. In healthy cerebrospinal fluid without inflammation, 90% of cells were reported to be T cells, mainly CD4 cells [43]. In pathological conditions, T cells can penetrate the brain parenchyma, while CD4+ T cells may play a role in the pro-inflammatory process of postoperative cognitive dysfunction [44]. Regulatory CD4 cells (Treg cells) provide neuroprotection by attenuating microglial activation in the CNS diseases [45]. Zhu et al. [44] demonstrated that PND is related to the increased number of T cells and NK cells in the hippocampus of aging mice after surgery and anesthesia. NK cell involvement in neurodegenerative diseases, such as multiple sclerosis and AD, has been reported [46]. The reduction of NK cell activity was previously found in AD patients [47]. It is noteworthy that NK cell infiltration into the BBB may not be a defense response, while it could be the result of PND progression, leading to immune system activation. The results of qRT-PCR and Western blotting of Tbx20 combined with the correlation of Tbx20 immune cells suggested that Tbx20 might be a hub gene, promoting PND. Chemokines are a large family of small cytokines that control the migration and residence of all immune cells. The CXC family includes chemokines CXCL1-CXCL17, and the CC family includes CCL1-CCL28. CC chemokines not only stimulate monocytes, but also basophils, eosinophils, T lymphocytes, and NK cells [48]. Studies have shown that changes in chemokines and chemokine receptors play an important role in the inflammatory response to AD, in which CCL1 and CCR8 prevent bacterial infection from participating in the demyelination of AD. The results of the present study showed that Tbx20 was positively correlated with CCL1, suggesting that Tbx20 could regulate immunity and participate in the development of neurodegenerative diseases. Unc13c and St8sia2 were significantly correlated with CCL9, which could be closely related to PND [49,50]. These evidences suggested that the three hub genes, which are closely related to levels of immune cell infiltration, play a critical role in the immune microenvironment. It is possible to influence the occurrence of PND by regulating the whole immune process.
Autophagy is an important pathway to eliminate abnormal protein aggregation in mammalian cells, which is related to protein homeostasis and neuronal health. Several studies have recently shown a close relationship between autophagy and PND [10,51,52]. Based on the close relationship between the transmembrane receptor protein serine/threonine kinase pathway of GO function and mTOR and autophagy, we conducted correlation analysis of hub genes and autophagy-related genes, and found that Unc13c, Tbx20, and St8sia2 were significantly correlated with different autophagy-related factors. Therefore, we speculated that these hub genes might regulate PND in aging mice through autophagy.

Limitations
As of the time of our analysis of this project, only one circRNA dataset (GSE174410) on postoperative neurocognitive disorder was available in the GEO public database. Therefore, subsequent studies by our group will also be conducted based on this dataset. If more circRNA datasets can be selected in the GEO database, the bias of the experiment can be reduced. In addition, the four drugs we screened will be validated in the follow-up study; and the mechanism of action of screening PND-related ceRNAs will also be studied in more depth.

CONCLUSIONS
In summary, we detected three mRNAs and four circRNAs that were associated with immune regulation of PND using publicly available data and combined mRNA analysis of high-throughput sequencing of aging mouse hippocampus after splenectomy by our experimental data. The results of the qRT-PCR verified the differential expression of these seven genes in the control and surgical groups, and Western blotting reverified the differential expression of the three mRNAs in the two groups. After reviewing the literature to remove the genes that were not associated with the CNS, two new ceRNA hub regulatory networks were finally constructed. The data of differentially expressed genes measured by our experiments were analyzed using CMap, and four active compounds (cimaterol, rucaparib, FG-7142, and hydrocortisone) were detected to interfere with gene expression. The results of the the present study may provide new insights into the immunomodulation-related pathogenesis and potential therapeutic targets of PND.