Abstract Propofol is a frequently used intravenous anesthetic agent. The impairment caused by propofol on the neural system, especially the hippocampus, has been widely reported. However, the molecular mechanism underlying the effects of propofol on learning and memory functions in the hippocampus is still unclear. In the present study we performed lncRNA and mRNA analysis in the hippocampi of adult mice, after propofol sedation, through RNA-Sequencing (RNA-Seq). A total of 146 differentially expressed lncRNAs and 1103 mRNAs were identified. Bioinformatics analysis, including gene ontology (GO) analysis, pathway analysis and network analysis, were done for the identified dysregulated genes. Pathway analysis indicated that the FoxO signaling pathway played an important role in the effects of propofol on the hippocampus. Finally, four lncRNAs and three proteins were selected from the FoxO-related network for further validation. The up-regulation of lncE230001N04Rik and the down-regulation of lncRP23-430H21.1 and lncB230206L02Rik showed the same fold change tendencies but changes in Gm26532 were not statistically significant in the RNA-Seq results, following propofol sedation. The FoxO pathway-related proteins, PI3K and AKT, are up-regulated in propofol-exposed group. FoxO3a is down-regulated at both mRNA and protein levels. Our study reveals that propofol sedation can influence the expression of lncRNAs and mRNAs in the hippocampus, and bioinformatics analysis have identified key biological processes and pathways associated with propofol sedation. Cumulatively, our results provide a framework for further study on the role of lncRNAs in propofol-induced or -related neurotoxicity, particularly with regards to hippocampus-related dysfunction. Keywords: propofol, long non-coding RNA, hippocampus, RNA-sequencing, neurotoxicity Introduction Propofol, the most widely used intravenous anesthetic agent in medical practice (Moseley et al., [37]1988), is noted among its analogs for its shorter onset and better recovery quality. Over the past few decades, several researchers have reported its protective effects in ischemic-reperfusion injury, including cerebral ischemic injury. Most of these reports attributed propofol's neuroprotective effect to its anti-inflammatory (Nie et al., [38]2015; Samir et al., [39]2015) and antioxidant (Ucar et al., [40]2015; Hsiao et al., [41]2016) properties. However, in the clinical setting, propofol is often administered to patients with an intact central nervous system (CNS). Thus, it is especially important to investigate the effects of propofol on an intact CNS. In the last decade, there is increasing evidence toward propofol's neurotoxic effect. Literatures describe multiple mechanisms, including the influence on dendritic development (Vutskits et al., [42]2005), induced neuronal apoptosis (Yan et al., [43]2017), increased cell death in neurons and oligodendrocytes (Krzisch et al., [44]2013), impaired maturation of neurons in newborns (Krzisch et al., [45]2013), and disturbance of the differentiation of neurons and astrocytes (Erasso et al., [46]2013). More importantly, some of these reports suggest that propofol's neurotoxicity impairs hippocampus-related learning and memory functions (Vutskits et al., [47]2005; Erasso et al., [48]2013; Krzisch et al., [49]2013; Liu et al., [50]2017; Qiao et al., [51]2017; Yan et al., [52]2017). Hippocampus-related learning and memory may be regulated by different processes, such as neurogenesis in the developmental and adult stage, long-term potential and synaptic plasticity. A recently reported study suggested that propofol induces an increase in TNF-α, short- or long-term neuronal apoptosis, neuronal loss, and synaptic loss (Chen et al., [53]2016). Another study on the same group of participants indicated that TNF-α contributes to propofol-induced neuronal apoptosis via the PI3K/AKT signaling pathway (Deng et al., [54]2017). In addition, it has been shown that maternal exposure to propofol during the late stages of pregnancy can impair learning and memory in the newborn by inactivating the BDNF-TrkB signaling pathway in the hippocampus of the fetus (Zhong et al., [55]2016). In our previous study (Fan et al., [56]2016), we demonstrated that propofol possesses the ability to influence the expression of several microRNAs in neural stem cells (NSCs). Other studies have also shown that propofol inhibits NSC neurogenesis through a mechanism involving the miR-141-3p/IGF2BP2 axis (Jiang et al., [57]2017). However, the role of long non-coding RNA (lncRNA) in propofol-induced neurotoxicity, in the hippocampus, is still unknown. LncRNA is defined as longer than 200 nucleotides in length and possessed almost no coding potential (Rinn and Chang, [58]2012). LncRNAs can regulate gene expression in four steps: epigenetic regulation, transcriptional regulation, post-transcriptional regulation, and translational regulation (Sun et al., [59]2017). Recently, several reports have provided new insights into the mechanism by which lncRNAs may regulate gene expression. Specific mechanisms include that of scaffolding and recruiting multiple regulatory proteins, genetic imprinting and chromatin remodeling, producing microRNA sponges, shaping, and utilizing three-dimensional nuclear structures and so on (Jandura and Krause, [60]2017; Sun et al., [61]2017; Bunch, [62]2018). Interestingly, in the brain, a large proportion of tissue-specific lncRNAs are preferentially expressed in specific regions or within different cell types (Derrien et al., [63]2012). These lncRNAs in the CNS participate in many aspects of brain function and in CNS development, from early neural differentiation to late-stage synaptogenesis (Briggs et al., [64]2015). Ramos et.al (Ramos et al., [65]2015) found that Pnky, a conserved lncRNA, may regulate neurogenesis in embryonic and postnatal NSC populations. Further, lncRNA-Map2k4 may regulate neuronal proliferation and apoptosis through a miR-199a/FGF1 pathway (Lv, [66]2017), and lncRNA Meg3 acts as a functional regulator in the regulation of the PTEN/PI3K/AKT signaling cascade, during the process of synaptic plasticity in neurons (Tan et al., [67]2017). Although several studies suggest that microRNA, another non-coding RNA, plays an important role in propofol-induced neurotoxicity, the relationship between propofol and an lncRNA has only been mentioned in one report, which indicates that propofol can inhibit lncRNA HOTAIR, which induces apoptosis in cervical cancer cells (Zhang et al., [68]2015). Therefore, the question remains whether propofol can affect lncRNA expression profiles in the hippocampus, and if so, how is propofol-induced lncRNA expression altered in propofol-induced neurotoxicity? In the present study, we adopted a previously reported propofol-induced neurotoxicity mouse model (Krzisch et al., [69]2013) and evaluated the extent of impairment of learning and memory functions using the Morris water maze. Then the different expression patterns of lncRNAs and mRNAs in the hippocampi of propofol-sedated mice were identified by RNA-Seq. Through bioinformatics analysis of the genes with varying expression patterns, we discovered the potential biological processes and signaling pathways that may play important roles during the process of propofol sedation, and closely connect with the hippocampus-related learning and memory functions. Materials and methods Experimental animals Eight to ten-week-old SPF (specific pathogen free) wild type C57BL6/J male mice (25–30 g) were purchased from the Laboratory Animal Center of Southern Medical University, Guangzhou, China. The mice were maintained under standard laboratory conditions [12-h light-dark cycles (lights on from 7:00 a.m. to 7:00 p.m.) with free access to food and water] unless otherwise indicated. All the protocols were approved by the Animal Use and Committee for Research and Education of Southern Medical University (Protocol number: SYXK-2015-0056). The experimental animals underwent adaptation to the environment in the course of a 2-week breeding period, and were subsequently subjected to experiments. Anesthesia procedure Adult mice were randomized into two groups (control group and propofol group; n = 8 in each group). The protocol used in this study was based on previous report (Krzisch et al., [70]2013). On day 0, mice in the Prop (propofol) group were anesthetized for 6 h; an initial intraperitoneal injection of propofol (100 mg/kg, propofol 1% medium-chain triglycerides, AstraZeneca, London, UK) was administered, followed by five subsequent injections at 50 mg/kg, at the rate of one injection per hour. The sedative effects were confirmed by noting the absence of the clip tail reflection and righting reflex. The Con (control) group was treated according to the same procedure, but instead using the same volumes of medium-chain and long-chain triglycerides (20%, Huarui, Wuxi, China). All mice were kept in a water bath (Yiheng, Shanghai, China) to maintain body temperature at 37°C during the entire period under anesthesia. Animals continued to feed till the next day, and were subjected to the behavioral test (Figure [71]1A). Figure 1. [72]Figure 1 [73]Open in a new tab Experimental design and behavioral readouts of MWM tests. (A). Schematic illustration of the experimental design; (B) The latency of mice to reach the platform, in the two groups. Two-way ANOVA showed that propofol-induced anesthesia increases the escape latency of MWM, when compared with the control group. (*P < 0.05), the Bonferroni test shows that the propofol-treated mice have a longer escape latency compared to that of the mice receiving control treatment, from day 1 to day 5 (^#P < 0.05 vs. control group at each day). Results are presented as mean ± standard error (n = 8, ^*P < 0.05); (C) Number of platform crossings in the probe trial. Results are presented as median and interquartile range (n = 8, ^*P < 0.05); (D) Representative searching swimming paths of two mice in the probe tests. (Left: Control group, Right: Propofol group). Morris water maze (MWM) The Morris water maze test (n = 8 for each group) was performed to assess the spatial learning and memory functions (Figure [74]1A). The maze consists of a circular polypropylene pool (110 cm in diameter and 20 cm in height) that was filled with tap water, up to an approximate height of 15 cm, at room temperature (23±3°C). The water was made opaque using milk powder, to ensure camouflage of the escape platform. A white escape platform (4.5 cm in diameter, 14.5 cm in height) was submerged 1 cm below the water surface and placed at the midpoint of one quadrant (the target quadrant). A colorful flag was added for visibility of the platform trial. Three pre-set extra mazes were placed on the wall of the testing room and the swimming activities were recorded via a digital video camera above the pool, which were analyzed using the Digi Behave system. During the testing period, the room was dimly lit with diffuse white light. The platform was positioned in the middle of the northwest (NW) quadrant for all mice. Each mouse was released from a predetermined, pseudo-random start location (south, east, west, or north) in the tank to receive daily training, for 4 consecutive days. Escape latency, in finding the submerged escape platform, was calculated for each trial. When the mouse failed to find the platform after 90 s, it was guided onto the platform and allowed to stay on it for ~10 s. After 4 days (day 1 to day 4) of training, mice were tested for spatial memory in a 90 s probe trial, with no platform present. The time spent in the target quadrant during the probe trial and the number of times the mice crossed the original location of the platform, were recorded. On the 5th day (day 5), mice were tested (four trials) using a visible platform that was placed in the same location within the northwest quadrant. Quadrant locations remained the same each time the mice were tested. Total RNA extraction The day after the MWM test (day 6) (Figure [75]1A), three mice from each group were randomly selected and the right hippocampus of each mice was harvested for the next experiments. Each of the hippocampus samples were washed twice, in cold PBS (Hyclone, Logan, USA), and immediately stored at −80°C. Total RNA was isolated using a Trizol reagent (Invitrogen, Carlsbad, USA), according to the manufacturer's protocol. RNA degradation and contamination were monitored on 1% agarose gels and RNA purity was analyzed using the NanoPhotometer spectrophotometer (IMPLEN, CA, USA). RNA concentration was measured using the Qubit RNA Assay Kit in Qubit 2.0 Flurometer (Life Technologies, CA, USA). RNA Nano 6000 Assay Kit of the Bioanalyzer 2100 system (Agilent Technologies, CA, USA) was used to assess the integrity of the isolated RNA. Library preparation for lncRNA sequencing A total of 3 μg of RNA, per sample, was used as input material for the RNA sample preparations. Firstly, ribosomal RNA was removed by the Epicentre Ribo-zero rRNA Removal Kit (Epicentre, Madison, USA), and the residual free rRNA was removed by ethanol precipitation. Subsequently, sequencing libraries were generated with the rRNA-depleted RNA using the NEBNext Ultra Directional RNA Library Prep Kit for Illumina (NEB, Ipswich, USA), according to the manufacturer's recommendations. Clustering and sequencing The clustering of the index-coded samples was performed on a cBot Cluster Generation System using the TruSeq PE Cluster Kit v3-cBot-HS (Illumia, San Diego, USA), according to the manufacturer's instructions. Following cluster generation, the libraries were sequenced on an Illumina Hiseq 2000 platform and 100 bp paired-end reads were generated. Quality control Clean data (clean reads) were obtained by removing reads containing adapter, reads containing ploy-N and low-quality reads from raw data. At the same time, Q20, Q30, and GC content of the clean data were calculated. The downstream analysis was done on high-quality clean data. Mapping to the reference genome The reference genome and gene model annotation files were downloaded directly from a genome website (Reference genome: [76]http://ftp.ensembl.org/pub/release-79/fasta/mus_musculus/dna/) (Annotation files: [77]http://ftp.ensembl.org/pub/release-79/gtf/mus_musculus/). The index of the reference genome was built using Bowtie v2.0.6 and the paired-end clean reads were aligned with the reference genome using TopHat v2.0.9 (Trapnell et al., [78]2012). Transcriptome assembly The mapped reads of each sample were assembled by both Scripture (beta2) (Guttman et al., [79]2010) and Cufflinks (v2.1.1) (Trapnell et al., [80]2012), following the reference. Both methods use spliced reads to determine the exon connectivity, via different approaches. Scripture uses a statistical segmentation model to distinguish expressed loci from experimental noise, and uses spliced reads to assemble the expressed segments. It reports all the statistically expressed isoforms in a given locus. On the other hand, Cufflinks uses a probabilistic model to simultaneously assemble and quantify the expression levels of a minimal set of isoforms, which provides a maximum-likelihood explanation to the expression data in a given locus (Cabili et al., [81]2011). Coding potential analysis CNCI (Coding-Non-Coding-Index) (v2) (Sun et al., [82]2013); CPC (Coding Potential Calculator) (0.9-r2) (Kong et al., [83]2007); Pfam Scan (v1.3) (Punta et al., [84]2012) and PhyloCSF (phylogenetic codon substitution frequency) (v20121028) (Lin et al., [85]2011) were used to predict the coding potential of transcripts. Once the transcripts were predicted to have coding potential by either/all of the four tools listed above, they were filtered, and those without coding potential were included as our candidate set of lncRNAs. Target gene prediction The cis and trans target mRNAs of the lncRNAs were used to predict their functions, due to a lack of adequate existing functional annotations of lncRNAs. A “cis role” indicates that the lncRNA acts on adjacent target genes. We searched the protein-coding genes 10k/100k upstream and downstream of the target lncRNA, and analyzed their functions. A “trans role” indicates that the lncRNA is used in identification, through the corresponding expression level. Although there were no more than 25 samples, we calculated the correlation of expression between lncRNAs and protein-coding genes with custom scripts; in addition, we clustered the genes from different samples with WGCNA (Langfelder and Horvath, [86]2008) to search common expression modules and then analyzed their function through functional enrichment analysis. Differential expression analysis Cuffdiff (v2.1.1) was used to calculate fragments per kilobase million (FPKMs), of both lncRNAs and coding genes, in each sample. Gene FPKMs were computed by summation of the FPKMs of transcripts in each gene group. Cuffdiff provides statistical routines for determining differential expression in digital transcript or gene expression data, using a model based on the negative binomial distribution. For biological replicates, transcripts or genes with a p-value of < 0.05 were accepted as differentially expressed. GO and KEGG enrichment analysis Gene Ontology (GO) enrichment analysis of differentially expressed genes, or lncRNA target genes, was implemented with DAVID ([87]http://david.abcc.ncifcrf.gov/). GO terms with p < 0.05 were considered significantly enriched by differentially expressed genes. Kyoto Encyclopedia of Genes, and Genomes (KEGG) is a database resource for understanding high-level functions and effects of the biological system ([88]http://www.genome.jp/kegg/). DAVID ([89]http://david.abcc.ncifcrf.gov/) was also used to test the statistical enrichment of genes, or target genes of lncRNA, with differential expression in KEGG pathways. The networks of the pathways and pathway-related genes were constructed using Cytoscape (version 3.2.1) plugin ClueGO + Cluepedia application. Construction of the co-expression network Significantly expressed mRNAs, which were involved in the FoxO signaling pathway, were superimposed onto the lncRNA-mRNA correlation network to determine their association with the lncRNAs. STRING (version 10.5) (Szklarczyk et al., [90]2017) was used to provide critical assessment and integration of protein-protein interactions. In the network, the circular nodes represent significantly expressed mRNAs and the diamond nodes represent the related lncRNAs; the black lines show connections between lncRNAs and their target mRNAs; the purple lines show integration of protein-protein interactions and the thicker lines represent a larger combined score. PCR and western blotting validation As previously mentioned, total RNA was isolated. The first strand cDNA was generated using the Reverse Transcription System Kit (Takara, China), and real-time PCR was performed using SYBR Premix Ex Taq II kit (Takara, China) on an Applied Bio systems 7500 Real-time PCR System (ABI, USA), in triplicates. The expressions of lncRNA or mRNA were normalized by U6 and GAPDH, respectively. All the sequences of primers used are listed in Table [91]S3. The fold change gene expression was calculated using the 2^−ΔΔCt method. Western blot was performed, as previously described. The hippocampus protein samples (30 μg) were separated by SDS-polyacrylamide gels and then transferred to PVDF membranes. After blocking with 5% fat free milk, the membranes were incubated with PI3K (1:750 Wanleibio, China), AKT (1:750, Wanleibio, China), FoxO3a (1:600, Wanleibio, China) and GAPDH (1:5,000, Beijing Ray Antibody, China). Membranes were then incubated in the appropriate peroxidase-conjugated secondary antibodies (1:10,000, Beijing Ray Antibody, China), at room temperature. After washing, antibody-bound proteins were detected with the ImmobilonTM Western Chemiluminescent HRP Substrate Kit (Millipore, USA) and exposed to X-ray film (Kodak, USA) for 1–2 min. The results were normalized to GAPDH and quantified using Image J (version 1.6.0_24). Statistical analysis In MWM tests, data in escape latency period was expressed as mean ± SD. The data for platform crossing times was not normally distributed and thus were expressed as median and interquartile range. There was no missing data for the variables of MWM, during the analysis. Interaction between time and group factors, in a two-way ANOVA with repeated measurements, was used to analyze the differences in learning curves (based on escape latency), between mice in the control group and mice treated with propofol in the MWM. Bonferroni analysis was used to compare the differences in escape latency between the control group and the propofol group, for each day. The Mann–Whitney test was used to determine the differences in platform crossing times, between the control and propofol conditions. An unpaired t-test was used to determine the differences in the levels of RNAs and proteins, between the control and propofol groups. Results Propofol exposure impaired spatial learning and memory in adult mice The spatial learning and memory in adult mice were tested through MWM analysis, 1 day after propofol exposure. Two-way ANOVA showed that propofol anesthesia increased the escape latency of MWM when compared with the control group (Figure [92]1B) (P = 0.0002). The Bonferroni test revealed that the mice that received propofol had a longer escape latency, compared to the mice receiving the control treatment, from day 1 to day 5 (Figure [93]1B) (P day1 = 0.001, P day2 = 0.006, P day3 = 0.007, P day4 = 0.004, P day5 = 0.001). In the probe tests, the Mann–Whitney test showed that the number of platform crossings were significantly reduced in the propofol group (Figures [94]1C,D), (P = 0.0138). The data suggests that propofol exposure impairs spatial learning and memory in adult mice. RNA-seq analysis of hippocampus transcriptome-identified mRNAs and lncRNAs Total RNA from the right hippocampi of all six treated mice was isolated and sequenced. A total of 384 million (Con group) and 417 million (Prop group) raw reads were obtained by RNA-Seq (Table [95]S1). Approximately 365 million and 397 million clean reads from the Con and Prop groups were isolated, respectively, and almost 73.96% (Con group) and 70.74% (Prop group) of the clean reads were uniquely mapped to reference genome (Table [96]S1), which were selected for subsequent experiments. The sequence reads were performed to map known gene types. Approximately 83.70% (Con group) and 82.54% (Prop group) of reads were mapped to coding genes, and around 1.50 and 1.42% of reads in the Con and Prop groups, respectively, were mapped to lncRNAs (Table [97]S2). The expression level in each sample showed no differences in terms of FKPM values (Figure [98]2B). Figure 2. [99]Figure 2 [100]Open in a new tab The comprehensive evaluation and screening of candidate lncRNAs. (A). The main workflow for screening provisional lncRNAs; (B). The box plots of FPKM distribution in different groups, showing no significant differences between the control and propofol groups; (C). The Venn diagram of the coding potential of screened transcripts using CNCI, CPC, PFAM, and phyloCSF, 9587 transcripts were predicted to have no coding potential by all of the four methods; (D) The full length of mRNAs is longer than lncRNAs; (E) The violet plots of expression levels of lncRNAs and mRNAs. The results of screening provisional lncRNAs The provisional lncRNAs were screened according to the workflow shown in Figure [101]2A. Five steps of basic filtering were applied: step 1—recurrence in ≥2 samples or by ≥2 assemblers; step2—transcript length ≥200bp, exon number ≥2; step3—minimal reads coverage ≥3; step4—filtration of known non-lncRNA annotation; step 5—classification of candidate lncRNAs. Four methods, CPC, CNCI, Pfam and PhyloCSF, were used to predict the coding potential. The transcript was eliminated if it was predicted to possess coding potential, by any or all of the four methods. The rest of the transcripts were selected as the candidate set of lncRNAs. We finally had 9587 transcripts that were predicted to have no coding potential, by all the four tools (Figure [102]2C), and all these transcripts were candidates for the following lncRNA research. There were also 80162 transcripts of mRNAs in total, with 77676 transcripts mapped to the reference genome, and the rest of them unknown novel mRNAs transcripts (data not shown). The comparison between mRNAs and unknown novel lncRNAs To distinguish between the number of unknown lncRNAs and mRNAs in introns and exons, 9587 transcripts that were predicted to have no coding potential by all the four tools, were compared with mRNAs which can be mapped to the genome, which showed that mRNAs were of longer lengths (Figure [103]2D). On the other hand, lncRNAs presented with a lower expression level than mRNAs (Figure [104]2E). The results of such a comparison between mRNAs and unknown novel lncRNAs, supports the existing information on lncRNAs. Differentially expressed genes Cuffdiff was used to detect differentially expressed genes between the Prop group and Con group, and 1249 differentially expressed transcripts were screened in total, including 433 up-regulated and 816 down-regulated ones in the Prop group (Figure [105]3A). Furthermore, there were 146 differentially expressed lncRNAs containing 77 up-regulated (20 known and 57 novel) and 69 down-regulated (19 known and 50 novel) genes, the heatmap of which is displayed in Figure [106]3B. In addition, 1103 mRNAs, containing 356 up-regulated (331 known and 25 novel mRNAs) and 747 down-regulated (678 known and 69 novel mRNAs) mRNAs, were also identified and are demonstrated as a heatmap (Figure [107]3C). Further, 6628 target genes (in 10k) and 15485 target genes (in 100k) in cis role, and 18834 target genes in trans role (data not shown) were predicted for further study. Figure 3. [108]Figure 3 [109]Open in a new tab Transcriptome profile of RNA-Seq data distinguishing control and propofol groups. (A) A volcano plot of differentially expressed lncRNAs and mRNAs between control and propofol groups; (B) Unsupervised hierarchical clustering of the expression profiles of differentially expressed lncRNAs in the propofol group, compared with the control group; (C) Unsupervised hierarchical clustering of the expression profiles of differentially expressed mRNA in the propofol group compared with the control group. The GO and KEGG enrichment analysis of target genes To illustrate the functions of the differentially expressed lncRNAs and their relationship with each other, GO and KEGG pathway enrichment analysis were conducted. In the GO analysis, we investigated the target genes of up-and down-regulated lncRNAs. All the results were ranked according to the enrichment score and the top 10 of each category are displayed in Figure [110]4. In the biological process analysis, 71 terms related to up-regulated lncRNAs were significantly enriched, while 90 processes in relation to down-regulated lncRNAs were significantly enriched. In the cellular component analysis, 12 terms associated with the up-regulated lncRNAs, and 26 terms linked to the down-regulated lncRNAs, were significantly enriched. In the molecular function analysis, 31 terms associated with the up-regulated lncRNAs, and 51 terms associated with the down-regulated lncRNAs, were significantly enriched. Results of the KEGG pathway analysis were also ranked according to the enrichment score, and the top 10 gene pathways associated with the target genes of up-regulated and down-regulated lncRNAs are listed in Figures [111]5A,B. Among these, MAPK and FoxO signaling pathways also appeared in the top 10 results of the mRNA KEGG analysis, indicating a close connection between these pathways and the effects of propofol. The network of the most enriched pathways and their related genes (Figures [112]5C,D) revealed that Ccnd2, Flt1, Crebbp, Notch3, Notch1, Ep300, and Igfbp3 were all cross-talk genes, which were associated with at least two pathways. Figure 4. [113]Figure 4 [114]Open in a new tab The top 10 enrichment scores in gene ontology (GO) enrichment analysis on target genes of differentially expressed lncRNAs. (A) Analysis of the up-regulated lncRNAs; (B) Analysis of the down-regulated lncRNAs. Red bars represent biological process terms; Green bars represent cell component terms; Blue bars represent molecular function terms. Figure 5. [115]Figure 5 [116]Open in a new tab The KEGG pathway enrichment analysis on target genes of differentially expressed lncRNAs. (A,B) The top 10 enrichment scores in the KEGG pathway analysis of the target genes of the up-regulated (A) and down-regulated (B) lncRNAs; (C,D) The network of most enriched pathways of the up-regulated (C) and down-regulated (D) lncRNAs and related genes. The GO and KEGG enrichment of differentially expressed mRNAs GO and KEGG enrichment analysis were also performed with 990 differentially expressed mRNAs. Finally, in the analysis of up-regulated mRNAs, 84, 58, and 34 GO terms were significantly enriched in the biological process, cell component and molecular function, respectively (Figure [117]6A). In the analysis of down-regulated mRNAs, 33, 7, 15 GO terms were significantly enriched in the biological process, cell component and molecular function, respectively (Figure [118]6B). The top 10 results of GO analysis on mRNA are shown in Figure [119]5. KEGG analysis showed that 29 pathways were involved in the mRNA downregulation caused by propofol (top 10 are shown in Figure [120]7A), and only 2 pathways related to propofol-regulated mRNAs were significantly enriched (top 3 are shown in Figure [121]7B). Among these significantly enriched pathways, MAPK and FoxO signaling pathways also appeared in the top 10 results of the KEGG analysis of lncRNA target genes, indicating the potential connection between these pathways and the effects of propofol (Figure [122]7). Figure 6. [123]Figure 6 [124]Open in a new tab The top 10 enrichment scores in gene ontology (GO) enrichment analysis of differentially expressed mRNAs. (A) Analysis of the up-regulated mRNAs; (B) Analysis of the down-regulated mRNAs; Red bars represent biological process terms; Green bars represent cell component terms; Blue bars represent molecular function terms. Figure 7. [125]Figure 7 [126]Open in a new tab The KEGG pathway enrichment analysis on differentially expressed mRNAs. The top 10 enrichment scores in the KEGG pathway analysis of the up-regulated mRNAs (A) and top 3 down-regulated mRNAs (B). Construction of the co-expression network A co-expression network of the dysregulated lncRNAs and their target mRNAs, which were involved in FoxO signaling pathway, was constructed. The co-expression network was composed of 13 lncRNA-mRNA predicted interactions and 13 protein-protein interactions. Sgk1, Sgk3, and Sos1 were identified as the hub nodes in the network (Figure [127]8). Moreover, lncRNA RP23-430H21.1 had three targets (Sgk1, Ccnd2, Sos1) and E230001N04Rik had two target mRNAs (Sgk1, Ccnd2), while the others had only one target in the network. Figure 8. Figure 8 [128]Open in a new tab The visualization of the lncRNA/FoxO gene co-expression network in FoxO pathways. The circular nodes represent the FoxO-genes, and the diamond nodes represent the FoxO-lncRNAs. The black lines show connections between lncRNAs and their target mRNAs; the purple line shows integration of protein-protein interactions and the thicker lines represent the larger combined score. Validation of the selected RNAs and proteins With a threshold of |FC| >1, 4 lncRNAs (E230001N04Rik, RP23-430H21.1, B230206L02Rik and Gm26532) in the network were selected for further validation by qRT-PCR. The up-regulated lncRNA, E230001N04Rik (P = 0.039), and 2 down-regulated lncRNAs, RP23-430H21.1 (P = 0.004) and B230206L02Rik (P = 0.001), showed the same fold change patterns as those in the RNA-Seq results, while down-regulated lncRNA Gm26532 (P = 0.585) did not reach statistical significance (Figure [129]9A). Quantitative analysis of FoxO pathway relative molecules showed that FoxO3a was down-regulated and PI3K/AKT were up-regulated in the Prop group (Figures [130]9A,B). Figure 9. Figure 9 [131]Open in a new tab Validation of selected lncRNAs and proteins in the FoxO pathway. (A) Three of the qRT-PCR-validated lncRNAs (E230001N04Rik, RP23-430H21.1, B230206L02Rik and Gm26532) showed the same fold change patterns as those in the RNA-Seq results. The differences in Gm26532 were not statistically significant. (B) Quantitative analysis of selected proteins in the FoxO pathway that are differentially expressed in the Prop group vs. Con group. The PI3K and AKT are significantly up-regulated, but FoxO3a is decreasing in the Prop group. Data was normalized to the house keeping gene U6 (lncRNA) or GAPDH (mRNA), ^*P < 0.05. Discussion In this study, we found that 6 h of propofol sedation impaired spatial learning and memory abilities in mice, determined by the MWM test. A total of 146 differentially expressed lncRNA and 1103 mRNAs, including known transcripts and novel transcripts, were differentially expressed in the hippocampus, identified by RNA-Seq. Bioinformatic analysis, including GO analysis, pathway analysis, and network analysis suggested that the FoxO signaling pathway played an important role in the effect of propofol on the hippocampus. Four lncRNAs were selected from the FoxO-related network for further validation through qRT-PCR and 3 of them showed the same fold change patterns as those in the RNA-Seq results. Taken together, these results suggest that lncRNAs may play a complicated role in propofol-induced hippocampal dysfunction, which may contribute to the impairment of related spatial learning and memory functions. The aim of this study was to identify the difference in the lncRNA and mRNA expression profiles in the hippocampus, in a propofol-induced neurotoxicity mice model. Propofol has been wildly used, in multiple clinical settings (Moseley et al., [132]1988), and its sedative effect may partly enhance the activity of GABA[A]R in the subsynaptic membrane (Li et al., [133]2015). Till date, the neuroprotective or neurotoxic effects of propofol on the CNS have been controversial, following multiple studies in different cerebral regions and various models (Velly et al., [134]2003; Erasso et al., [135]2013; Krzisch et al., [136]2013; Twaroski et al., [137]2015; Wang et al., [138]2015). In the present study, we adopted a propofol-induced neurotoxicity mice model which may impair adult neurogenesis in the hippocampus (Krzisch et al., [139]2013). Krzisch et al. ([140]2013) found that propofol-induced anesthesia significantly decreased the survival and maturation of adult-born hippocampal neurons, in a developmental, stage-dependent manner. The results of our Morris water maze experiments demonstrated that propofol sedation increases the escape latency of MWM, and decreases the number of platform crossings during the probe test, which suggests that propofol exposure impairs spatial learning and memory in adult mice. The results from the study by Krzisch's group and those from our behavior test, are both in support of propofol-induced neurotoxicity. The mechanism of propofol-induced neurotoxicity has gained increasing attention over the last few decades. Recently published studies suggested that propofol-induced neurotoxicity may be regulated by microRNA. For instance, propofol may induce human embryonic stem cell-derived neuronal death through a signal transducer, and activation of the transcription 3/miR-21/ Sprouty 2-dependent mechanism (Twaroski et al., [141]2014). Rno-miR-665 is involved in the neurotoxicity induced by propofol via a caspase-3, through negative regulation of BCL2L1 (Sun et al., [142]2015). However, the role of lncRNA in propofol-related neurotoxicity remains unclear. Thus, it is meaningful to investigate the differences in lncRNA expression, in mice hippocampi, following propofol exposure. LncRNAs, a class of non-coding RNA molecules, are no longer than 200 nucleotides in length and have barely discernable coding potential. They are widely involved in multiple pathophysiological processes (Guttman et al., [143]2009; Song et al., [144]2016; Wang et al., [145]2016), and also present a region-specific function (Derrien et al., [146]2012). LncRNAs in the CNS participate in many aspects of brain function and their role in CNS development range from early neural differentiation to late-stage synaptogenesis (Briggs et al., [147]2015). Research indicates that lncRNAs play an intricate part in NSCs, neuronal proliferation and apoptosis (Ramos et al., [148]2015; Lv, [149]2017). In our study, we performed RNA-Seq to screen the differential expression of lncRNA and mRNA in the hippocampus, after propofol sedation. LncRNA Malat1, one of the down-regulated lncRNAs following propofol sedation, was a highly conserved lncRNA that has been found in various cancers, such as non-small cell lung cancer (Ji et al., [150]2003) and breast cancer (Feng et al., [151]2016), and plays a role in the proliferation of myocardial cells (Zhao et al., [152]2015). More importantly, Bernard (Bernard et al., [153]2010) also discovered that Malat1 can regulate the synaptic plasticity of primary culture neurons. The results of our RNA-Seq suggested that propofol can reduce the expression of Malat1. LncRNA Kcnq1ot1 has been known to be associated with Beckwith-Wiedemann Syndrome and glioma malignancy (Gong et al., [154]2017), and depressing lung adenocarcinoma chemoresistance to paclitaxel (Ren et al., [155]2017). In addition, our results revealed that LncRNA Kcnq1ot1 significantly increased following propofol sedation in mice hippocampi. Overall, we screened 146 differentially expressed lncRNAs in the hippocampus through RNA-Seq. For the first time, through our study we have presented an lncRNA expression profile in mice hippocampus following propofol exposure. These results may also indicate the role of lncRNA and the underlying molecular mechanism of propofol-related neurotoxicity. To effectively screen the candidate lncRNAs for further validation and research, bioinformatics analysis was conducted, separately, for mRNAs and the target genes of lncRNAs. The results showed several important pathways that were significantly enriched, particularly the FoxO signaling pathway, which was significantly enriched in both mRNAs and the target genes and indicated that this pathway may play a critical role in the dysfunction of the hippocampus, after propofol exposure. The FoxO family belongs to an evolutionarily conserved group of forkhead transcription factors. Mammals have four isoforms of the FoxO family: FoxO 1 (FKHR), 3a (FKHRL1), 4 (AFX), and FoxO6 (Hwangbo et al., [156]2004). FoxO transcription factors are at the interface of diverse physiological processes, e.g., coordinating gene expression that regulates proliferation, cell-cycle, DNA repair, oxidative stress resistance as well as metabolism (Medema et al., [157]2000; Nakae et al., [158]2001). Nutrient and energy stress signaling pathways regulate FoxOs and are important for NSC maintenance (Rafalski and Brunet, [159]2011; Spéder et al., [160]2011). In the adult mouse brain, FoxO1 is most abundant in the striatum, dentate gyrus, and the ventral hippocampus (Paik et al., [161]2009), while FoxO3a is highly expressed in the cortex, hippocampus, and cerebellum (Hoekman et al., [162]2006). Previous studies have shown that the PI3K/Akt-FoxO signaling pathway plays a central role in the development of the nervous system. Lisa (Kennedy et al., [163]2013) demonstrated that the insulin/IGF-1-PI3K signaling pathway modulates the activity of the DAF-16/FoxO transcription factor to promote the anterior migrations of the hermaphrodite-specific neurons during embryogenesis of C. elegans (by signaling pathways that are conserved in humans). FoxO also restricts growth and survival of apoptosis-inhibited mushroom body neuroblasts, while activation of the growth-promoting insulin/PI3 kinase pathway sustains not only long-term survival of adult mushroom body neuroblasts, but also increases their proliferation and growth rate (Paik et al., [164]2009). However, lncRNA-mediated regulation of the FoxO signaling pathway in nervous system development is yet to be elucidated. To clarify the relationship between lncRNAs and FoxO signaling pathway, we constructed a protein-protein interaction network and validated lncRNAs (E230001N04Rik, RP23-430H21.1, B230206L02Rik, and Gm26532) and FoxO3a mRNA via qRT-PCR. Furthermore, western blotting was performed to evaluate the PI3K/AKT/FoxO3a protein. According to our bioinformatic analysis and western blotting results, we hypothesize that these lncRNAs, E230001N04Rik, RP23-430H21.1, and B230206L02Rik, may participate in the FoxO signaling pathway, to regulate propofol-induced neurotoxicity. Nevertheless, the detailed regulatory mechanism of these lncRNAs in the FoxO signaling pathway is still require further study in order to confirm these findings. In conclusion, our study reveals that propofol sedation can influence the expression of lncRNAs and mRNAs in the hippocampus, and bioinformatics analysis have identified several key biological processes and KEGG pathways associated with propofol sedation. Our results provide a framework for further study on the role of lncRNAs in propofol-induced or -related neurotoxicity, particularly with regards to hippocampus-related dysfunction. Author contributions JF and QZ: Carried out most procedures of this study, participated in the data analysis, and manuscript writing; YL and XS: Participated in design of animal study and interpreted part data of this study; JT, ZQ, and JH: Participated in the collection of data, data analysis; TT, ZQ, and JT: Designed this study, provided financial support, and wrote the manuscript. All authors read and approved the final manuscript. Conflict of interest statement The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest. Acknowledgments