Abstract Purpose Physiologic neuronal apoptosis, which facilitates the developmental maturation of the nervous system, is regulated by neuronal activity and gene expression. Circular RNA (circRNA), a class of non-coding RNA, regulates RNA and protein expression. As the relationship between circRNA and apoptosis is unknown, we explored changes in expression patterns of circRNA during physiologic neuronal apoptosis. Methods High-throughput sequencing was used to explore changes in the expression of circRNA in the postnatal developing rat retina. Neuronal apoptosis was determined with immunohistochemistry and terminal deoxynucleotidyl transferase dUTP nick end labeling (TUNEL) in the rat retinal ganglion cell layer. Results In total, 2,654, 7,201, and 5,628 circRNA species were detected in the postnatal day (P)3, P7, and P12 rat retina, respectively. Of these circRNA species, 1,371 changed statistically significantly between P3 and P7 and 1,112 changed statistically significantly between P7 and P12. Normal developmental apoptosis, measured with the ratio of apoptotic (caspase-3- or TUNEL-positive) cells to normal cells, showed an increase from P3 to P7 and then a reduction from P7 to P12. In addition, 15 circRNAs whose host genes were associated with apoptosis were differentially expressed during the early development period. Conclusions These results associate circRNAs with neuronal apoptosis, providing potential mechanisms and treatment targets for physiologic and drug-induced apoptosis in the developing nervous system. Introduction During the period of synaptic formation and refinement, especially at the peak of synaptogenesis, many neuronal cells undergo apoptosis, or programmed cell death (PCD, also called physiologic apoptosis), to ensure the establishment of accurate neuronal connections and networks [[32]1-[33]4]. In the rodent brain, this process mainly occurs within the 2 weeks after birth, primarily during postnatal days 4 to 14 (P4–P14) [[34]5,[35]6]. The nervous system is also vulnerable to ethanol, general anesthetics, and other substances during this period [[36]7,[37]8]. The underlying mechanism regulating the neuronal apoptotic process remains unclear, but numerous studies demonstrate that neuronal activity and genetic modulation play highly significant roles. In the developing rat retina, the peak of physiologic apoptosis coincides with the transition from early cholinergic-driven, synchronized spontaneous network activity to glutamatergic-driven activity [[38]9,[39]10]. Blockade of this neural activity by ketamine during this transition period aggravates physiologic apoptosis [[40]11]. During this transition, some relevant receptors, signaling pathways, and apoptosis-related genes also experience dramatic changes [[41]10,[42]12-[43]15]. Previous studies on mRNA or transcription factors to cleave or downregulate mRNA have shown that noncoding microRNA (miRNA) species are involved in neurogenesis, proliferation, axon extension (e.g., mir-9) [[44]16,[45]17], and neuronal apoptosis (e.g., mir-21, -23, -26, -27, and -29) [[46]18-[47]21]. Circular RNA (circRNA), formed by back-splicing, has been reported for decades because of splicing errors [[48]22,[49]23]. More recently, studies have shown that circRNAs are expressed in a variety of eukaryotic organisms, including mammals, and are involved in various physiologic or pathological processes. CircRNAs are widely expressed and show temporal and spatial changes during development [[50]24-[51]27]. This large class of RNA has regulatory abilities, acting as protein or miRNA “sponges” and regulating mRNA transcription or translation [[52]28,[53]29]. Recent studies have shown that circRNAs are remarkably enriched in the nervous system, especially in synapses, during development; circRNAs were also found to regulate synaptic function and neuronal plasticity [[54]30,[55]31]. However, the relationship between circRNAs and physiologic neuronal apoptosis is largely unknown. Development-related neuronal apoptosis in the central nervous system of rodents mainly occurs within 2 weeks after birth and peaks at about P7 [[56]11]. Thus, we chose to analyze the rat retina at the P3, P7, and P12 time points, to determine any possible associations between physiologic neuronal apoptosis and changes in the expression of circRNA during postnatal development. Methods Animals and tissue dissection Sprague-Dawley rat pups (male and female) aged P3, P7, and P12 days were provided by the Experimental Animal Center at the Shanghai General Hospital in Shanghai, China. Eight rat pups were used for each group of P3, P7, and P12. Rat pups were housed with their dam under a 12 h:12 h light-dark cycle at room temperature 35–37 °C, with food and water available ad libitum until the experiments were conducted. All experimental procedures were reviewed and approved by the Animal Care Committee at the Shanghai General Hospital, Shanghai Jiao Tong University School of Medicine, and were conducted under the guidelines of the Care and Use of Laboratory Animals published by the U.S. National Institutes of Health (National Institutes of Health Publication No. 85–23, revised in 1996) and the ARVO Statement for the Use of Animals in Ophthalmic and Vision Research. Every effort was made to minimize the number of animals used, and their discomfort, during all experimental procedures. Because of the limited size of the rat retina, we combined the retinas of each group of P3, P7, and P12 rats for high-throughput sequencing of circRNA and immunohistochemistry or terminal deoxynucleotidyl transferase dUTP nick end labeling (TUNEL) assay. The basic experimental protocol was slightly modified from that in previous studies [[57]32]. Briefly, the pups were euthanized instantaneously by decapitation, and then their eyes were rapidly dissected with fine scissors and transferred to an ice-cold (0‒4 °C) bath of artificial cerebrospinal fluid (ACSF: composition in mM: NaCl 119, KCl 2.5, K[2]HPO[4] 1.0, CaCl[2] 2.5, MgCl[2] 1.3, NaHCO[3] 26.2, and D-Glucose 11). The bath was continuously bubbled with a 95% O[2]/5% CO[2] gas mixture. About one fifth of the eyeball circumference at the edge of the cornea and sclera was resected using ophthalmological scissors to facilitate the perfusion of ACSF through the retina. After 1 h of recovery in ACSF (37 °C) bubbled with a 95% O[2]/5% CO[2] gas mixture, the retinas were dissected and then either frozen quickly with liquid nitrogen for circRNA high-throughput sequencing or fixed in 4% paraformaldehyde at 4 °C for 24 h for immunohistochemistry and the TUNEL assay. Cleaved caspase-3 immunohistochemistry and TUNEL After the retinas were embedded in paraffin, retinal sections of 4 µm thickness were prepared. Three percent hydrogen peroxide was used to inactivate endogenous peroxidases, and EDTA, adjusted to pH 9.0, was used for heat-induced antigen retrieval for 8–10 min. After rinsing with PBS (1X; 137 mM NaCl, 2.7 mM KCl, 10 mM Na[2]HPO[4], 2 mM KH[2]PO[4], pH 7.4), retinal sections were incubated with rabbit anti-cleaved caspase-3 primary antibody (#9661s, dilution: 1/300, Cell Signaling Technology, Danvers, MA) was performed at 4 °C overnight, followed by incubation with goat anti-rabbit immunoglobulin G (IgG; PV-9001, ZSGB-BIO, Beijing, China) secondary antibody at 37 °C for 1 h. After washing, 3,3′-diaminobenzidine (DAB, ZLI-9017, ZSGB-BIO) was used to visualize immunoreactivity. All sections were then counterstained with hematoxylin to stain the nuclei. Finally, sections were dehydrated and mounted for microscopic examination. For the TUNEL assay, retinal sections were deparaffinized and rehydrated. After rinsing with PBS, the sections were treated with proteinase K (Roche Applied Science, Indianapolis, IN) and quenched with 3% hydrogen peroxide; then they were incubated in a terminal deoxynucleotidyl transferase (TdT) reaction mix (Roche Applied Science) for 1 h at 37 °C. Sections were then incubated for 5 min with 4',6-diamidino-2-phenylindole (DAPI) at room temperature. CircRNA enrichment and high-throughput sequencing RNA was isolated using the RNeasy system (Qiagen, Dusseldorf, Germany), including on-column DNase I digestion (Qiagen) to eliminate DNA. Then the total RNA was depleted of ribosomal and linear RNA using the RiboMinus kit (Thermo Fisher Scientific, Waltham, MA) and RNase R (Epicenter, Madison, WI), respectively. CircRNAs were validated with quantitative reverse transcription (RT)–PCR using the VAHTSTM Library Quantification Kit (Vazyme Biotech Co., Ltd, Piscataway, NJ). Specifically, samples were heated to 70 °C for denaturing and then cooled to 40 °C on a thermocycler with RNase R buffer to digest linear RNA. After RNase R digestion and circRNA reverse transcription, the cDNA library was produced by PCR (conditions: 98 °C for 30 s, followed by 15 cycles of 98 °C for 10 s, 60 °C for 30 s, and 72 °C for 30 s, then 72 °C for 5 min) following depletion of the second cDNA chain with the USER enzyme. The sequencing libraries were constructed using the VAHTSTM Library Quantification Kit (Vazyme Biotech Co) and were sequenced with an Illumina HiSeq™ 2000 flowcell sequencer (Illumina, San Diego, CA). Data analysis The RNA-sequenced (RNA-seq) data were analyzed using FastQC software to verify the quality of the data [[58]33]. The BWA-MEM alignment algorithm of the software BWA (the [59]Burrows-Wheeler Alignment tool) was used to match the RNA-seq data to the genome [[60]34]. To annotate the validated circRNAs, we used CIRI ([61]CircRNA Identifier) software to predict the circRNAs [[62]35] present in the sample and sequenced them using circBase [[63]36]. The false discovery rate (FDR) was used to correct for multiple hypothesis testing, and a twofold change and FDR≤0.001 was considered statistically significant [[64]37]. Caspase-3-positive or TUNEL-positive cells were counted in a double-blinded manner from five discontinuous views randomly selected from each sample (n = 5 retinas). Data in the figures were analyzed with GraphPad Prism 5 software (GraphPad Software Inc., San Diego, CA) and presented as mean ± standard error of the mean (SEM). A one-way ANOVA (ANOVA) followed by the Bonferroni correction, or the appropriate equivalent non-parametric test, was used for comparisons among groups, and a p-value of < 0.05 was denoted as statistically significant. Results CircRNAs change during early development in the rat retina We analyzed the raw reads or raw data of the expression of circRNA at P3, P7, and P12, using high-throughput sequencing followed by filtration and quality control (QC; [65]Figure 1A, [66]Table 1) to obtain clean reads as in a previous study [[67]35]. Because exon–exon junction reads ([68]Figure 1B-C) are important features of circRNAs that can map to the junction points of exons, the BWA software was used to match the clean reads to the genome and identify the junction reads. The results showed that the clean reads were expressed abundantly (7–9 × 10^7/group) and largely mapped to the genome (>75% of all reads) in the early developmental retina ([69]Figure 2A, [70]Table 2). Most of the clean reads mapped on the exon regions could be visualized at P7 (56.41% of all reads) but were present in only small amounts at P3 (15.21% of all reads) and P12 (17.78% of all reads; [71]Figure 2A, [72]Table 2). After the CIRI predictions and false-positive read filtration, validated circRNAs were obtained and annotated according to circBase ([73]Figure 1E, [74]Table 3, Appendix 1, Appendix 2 and Appendix 3). Figure 1. [75]Figure 1 [76]Open in a new tab The basis and pipeline of circRNAs identification and analysis. A: After high throughput sequencing, raw data were collected and subsequently filtered and analyzed for quality control (QC). Clean data were mapped on the genome, and then the circular RNAs (circRNAs) were predicted via CIRI (CircRNA Identifier) and the circBase database. Finally, analyses of the differential expression between the time points, as well as function and origin analyses, were performed. B: Two segments of junction reads of classic circRNAs match to the relevant sequence in reverse orientation, albeit separately. C: If a flanking part of the junction is a short exon, the segment can align to the flanking exon of the circRNAs. D: If the circRNAs length is short, the read will possibly contain two terminal segments that align to the area of junction within circRNA. E: To reduce the false positive rate, candidate circRNAs are filtered via criteria as in a previous study. Table 1. Expression, filtration and quality control of raw data during early development in rat retina. Sample Raw Reads Clean reads Raw data(G) Clean data(G) Error(%) Q20(%) Q30(%) P3 __________________________________________________________________ 96,172,040 __________________________________________________________________ 79,680,804 __________________________________________________________________ 14.43 __________________________________________________________________ 11.95 __________________________________________________________________ 0.03 __________________________________________________________________ 96.91 __________________________________________________________________ 92.68 __________________________________________________________________ P7 __________________________________________________________________ 96,890,346 __________________________________________________________________ 73,230,288 __________________________________________________________________ 12.84 __________________________________________________________________ 9.67 __________________________________________________________________ 0.05 __________________________________________________________________ 93.26 __________________________________________________________________ 87.04 __________________________________________________________________ P12 104,779,880 89,480,408 15.72 13.42 0.02 97.41 93.56 [77]Open in a new tab Raw Reads: original sequence data. Clean reads: data after raw data filtration. Raw data: the number of original sequence multiplied by the length of the sequence. Clean data: the number of filtrated sequence multiplied by the length of the sequence. Error: base-rate error. Q20 Q30:ratio of bases with Phred quality score >20 or 30 in total base. Figure 2. [78]Figure 2 [79]Open in a new tab The expression of circRNA reads during postnatal development in the rat retina. A: Classification of the reads during early development in the rat retina. B: Classification of the origins of the circular RNA (circRNA) species during early development. C: Distribution of the total circRNA species at time points P3, P7, and P12. D: Distribution of the circRNA species after filtration (>twofold change, false discovery rate (FDR)≤0.001 for P3 versus P7 and P7 versus P12) at P3, P7, and P12. Table 2. The origin of clean reads in different stage of developing rat retina. Sample P3 P7 P12 Total number of clean reads __________________________________________________________________ 79,680,804 __________________________________________________________________ 73,230,288 __________________________________________________________________ 89,480,408 __________________________________________________________________ Number of clean reads mapped on genome __________________________________________________________________ 61,697,391
(77.43%) __________________________________________________________________ 56,325,618
(76.92%) __________________________________________________________________ 72,183,843
(80.67%) __________________________________________________________________ Number of clean reads mapped on intron regions __________________________________________________________________ 49,562,407
(62.20%) __________________________________________________________________ 15,004,813
(20.49%) __________________________________________________________________ 56,256,016
(62.87%) __________________________________________________________________ Number of clean reads mapped on intergenic regions __________________________________________________________________ 19,376
(0.02%) __________________________________________________________________ 11,008
(0.02%) __________________________________________________________________ 15,521
(0.02%) __________________________________________________________________ Number of clean reads mapped on exon regions 12,115,608
(15.21%) 41,309,797
(56.41%) 15,912,306
(17.78%) [80]Open in a new tab The number of clean reads in P3, P7 and P12 rat retina and the clean reads origination such as exon regions, intron regions or intergenic regions. Table 3. The annotation of circRNA. circRNA_ID Chr circRNA_start circRNA_end #junction
reads #non_junction
reads CircRNAs
type gene_ID chr1:1901431|1901647 __________________________________________________________________ Chr1 __________________________________________________________________ 1,901,431 __________________________________________________________________ 1,901,647 __________________________________________________________________ 5 __________________________________________________________________ 34 __________________________________________________________________ exon __________________________________________________________________ Ppil4 __________________________________________________________________ chr1:1901431|1909369 __________________________________________________________________ Chr1 __________________________________________________________________ 1,901,431 __________________________________________________________________ 1,909,369 __________________________________________________________________ 5 __________________________________________________________________ 22 __________________________________________________________________ exon __________________________________________________________________ Ppil4 __________________________________________________________________ chr1:2375315|2377265 __________________________________________________________________ Chr1 __________________________________________________________________ 2,375,315 __________________________________________________________________ 2,377,265 __________________________________________________________________ 7 __________________________________________________________________ 7 __________________________________________________________________ exon __________________________________________________________________ Ust __________________________________________________________________ chr1:2839790|2849731 __________________________________________________________________ Chr1 __________________________________________________________________ 2,839,790 __________________________________________________________________ 2,849,731 __________________________________________________________________ 2 __________________________________________________________________ 8 __________________________________________________________________ intron __________________________________________________________________ SNORA32 __________________________________________________________________ chr1:3879152|3915233 Chr1 3,879,152 3,915,233 5 79 exon Stxbp5 [81]Open in a new tab circRNA_ID: code of circRNAs. Chr: host chromosome of circRNAs. circRNA_start and circRNA_end: the starting and ending location of circRNA in host chromosome. #junction reads: number of junction reads formed by back-splicing #non_junction_reads: number of reads which mapped on the flanking area of the head and tail end. CircRNA type: the head and tail end of circRNAs that could be both mapped on the exon area belong to exon type, otherwise belong to intron type. gene_ID: corresponding host gene code of circRNAs according to the location information. To identify circRNA changes in the developing rat retina, the total circular junction reads of the different stages were further separated with BWA. CircRNAs with the same sequence were regarded as one circRNA species. The number of circular junction reads at P12 (143,815) was greater than that at P3 (78,758) and P7 (78,385), while the number of circRNA species at P7 (7,021) was greater than that at P3 (2,654) and P12 (5,628; [82]Figure 2B, [83]Table 4). Table 4. The expression of circRNAs in different stage of developing rat retina. Sample P3 P7 P12 Number of circular junction reads __________________________________________________________________ 78,758 __________________________________________________________________ 78,385 __________________________________________________________________ 143,815 __________________________________________________________________ Number of circRNA species __________________________________________________________________ 2,654 __________________________________________________________________ 7,021 __________________________________________________________________ 5,628 __________________________________________________________________ Number of circRNA species originated from exon regions __________________________________________________________________ 2,134
(80.41%) __________________________________________________________________ 5,576
(79.42%) __________________________________________________________________ 4,509
(80.12%) __________________________________________________________________ Number of circRNA species originated from intron regions __________________________________________________________________ 517
(19.48%) __________________________________________________________________ 1,435
(20.44%) __________________________________________________________________ 1,108
(19.69%) __________________________________________________________________ Number of circRNA species originated from intergenic regions 3
(0.11%) 10
(0.14%) 11
(0.20%) [84]Open in a new tab The expression of circular junction reads and circRNAs species in P3, P7 and P12 rat retina. And the expression of circRNAs origination such as exon regions, intron regions or intergenic regions. After the data were integrated, a total of 10,890 circRNAs were obtained; additionally, the overall expression of circRNAs was different among the P3, P7 and P12 time points. For example, there were 4,123 circRNAs expressed solely at P7 ([85]Figure 2C). By comparing the expression of circRNA between the time points, we detected differential expression (greater than a twofold change, FDR≤0.001; [86]Figure 3A,B, Appendix 4 and Appendix 5) of 1,371 circRNAs (out of 8,223 circRNAs) between P3 and P7 and 1,112 circRNAs (out of 10,016 circRNAs) between P7 and P12. We found that 558 circRNAs and 813 circRNAs (out of 1,371 circRNAs) were upregulated and downregulated, respectively, from P3 to P7. The circRNAs upregulated from P3 to P7 included those originating from the genes Akt3, Ccdc6, Gnaz, and Mpped2; the circRNAs downregulated from P3 to P7 originated from the genes Ly6g6f, Melk, Casc5, Smarcc1, and Pank2 ([87]Figure 3C, Appendix 4). We found that 434 or 678 circRNAs (out of 1,112 circRNAs) were upregulated or downregulated, respectively, from P7 to P12. Upregulated circRNAs included those derived from the genes Pde6a, Unc13b, Pde4b, Pi4ka, Mtr, and Dnmbp, while downregulated circRNAs included those derived from the genes Pank2, Smarcc1, Melk, Ly6g6f, Casc5, Dst, and Pak3 ([88]Figure 3D, Appendix 5). After the circRNAs expressed at only P3 or P12 were filtered out, a total of 438 circRNAs were retained for measuring significant changes in the circRNAs across the time points (Appendix 6, [89]Figure 2D). Of these 438 circRNAs, ten were gradually downregulated, and 17 were gradually upregulated from day P3 to P12 (Appendix 7). In addition, 145 circRNAs were present in only small amounts or were undetectable (48 circRNAs) at P7 relative to P3 and P12. In addition, 266 circRNAs of the 438 circRNAs were present at P7 but were absent or diminished at P3 and P12 (172 circRNAs appeared only at P7, [90]Figure 2D). Figure 3. [91]Figure 3 [92]Open in a new tab The changes in circRNAs during postnatal development in the rat retina. A, B: Expression profile of circular RNAs (circRNAs) in the early developing rat retina between P3 and P7, and P7 and P12, respectively. C, D: Bar graphs showing the expression changes in the circRNAs from P3 to P7 and P7 to P12, respectively. Physiologic neuronal apoptosis during early development in the rat retinal ganglion cell layer Most neurons undergoing PCD exhibit the morphological and biochemical hallmarks of apoptosis, such as cell shrinkage, activation of caspases, and DNA fragmentation [[93]5]. Cleaved caspase-3 is produced early in the apoptotic process and is responsible for executing the apoptotic process, making it a strong biomarker of apoptosis [[94]38]. The TUNEL assay detects apoptosis by labeling the fragmented DNA generated during PCD. For enhanced reliability of the results, we used anti-cleaved caspase-3 immunohistochemistry and the TUNEL assay to detect apoptosis in retinas from rats at ages P3, P7, and P12. The cleaved caspase-3-positive-to-negative cell ratio increased from 1.69 ± 0.25% to 3.19 ± 0.42% between P3 and P7 (p = 0.005); it then decreased to 1.69 ± 0.25% between P7 and P12 (p = 0.008; [95]Figure 4A,C). The TUNEL-positive cell ratio increased from 6.75 ± 0.88% to 12.77 ± 1.17% between P3 and P7 (p<0.001) and then decreased to 7.16 ± 0.62% between P7 and P12 (p<0.001; [96]Figure 4B,D). The results confirmed that normal age-dependent physiologic apoptosis of neurons was occurring. Figure 4. [97]Figure 4 [98]Open in a new tab Physiologic neuronal apoptosis occurred in an age-dependent manner during the postnatal development of the rat retina. Neuronal apoptosis of rat retinas was detected with cleaved caspase-3 immunolabeling (cleaved caspase-3-positive) and with the terminal deoxynucleotidyl transferase dUTP nick end labeling (TUNEL) assay (TUNEL-positive; n = 5 retinas per group). A: Representative photomicrograph of caspase-3-positive cells (brown) in the rat retinal GCL. Scale bar = 50 μm. B: Respective photomicrograph of TUNEL-positive cells (red) in the rat retinal GCL. Scale bar = 25 μm. C: The mean percentages of caspase-3 positive cells at P3, P7, and P12. D: The percentages of TUNEL-positive cells in the retinal GCL. #p<0.05, ##p<0.01 versus P3, #p<0.05, ##p<0.01 versus P12. One-way ANOVA followed by Bonferroni’s post-hoc test. ONL, outer nuclear layer; INL, inner nuclear layer; GCL, ganglion cell layer. Candidate circRNAs from host genes associated with apoptosis We identified changes in the expression of circRNA that could be linked to normal developmental apoptosis during in the neonatal rat retina. A total of 15 circRNAs from genes associated with apoptosis were obtained after further filtration and screening ([99]Table 5). We found that five circRNAs (from the genes Optn, Thoc1, Rbm5, Ddx19a, and Bnip2) could not be detected, and two circRNAs (from the genes Pias1 and Naa35) were not expressed at P7 but were abundantly expressed at P3 and P12. Three circRNAs (from the genes Rere, Arhgap10, and Birc6) were expressed specifically at P7, and three circRNAs (from the genes Ranbp9, Epha7, and Braf) could be visualized at P7 but were absent or diminished at P3 and P12 ([100]Figure 5). Table 5. Candidate circRNAs and the host gene may link to neuronal physiologic apoptosis. circRNA_ID Gene_ID chr17:77176705|77178100 __________________________________________________________________ Optn __________________________________________________________________ chr18:1164675|1169033 __________________________________________________________________ Thoc1 __________________________________________________________________ chr8:116505048|116505353 __________________________________________________________________ Rbm5 __________________________________________________________________ chr19:43254907|43259022 __________________________________________________________________ Ddx19a __________________________________________________________________ chr8:76405782|76406571 __________________________________________________________________ Bnip2 __________________________________________________________________ chr8:67801800|67812158 __________________________________________________________________ Pias1 __________________________________________________________________ chr17:5429865|5436758 __________________________________________________________________ Naa35 __________________________________________________________________ chr5:167486007|167558256 __________________________________________________________________ Rere __________________________________________________________________ chr19:34233020|34253181 __________________________________________________________________ Arhgap10 __________________________________________________________________ chr6:22033342|22038870 __________________________________________________________________ Birc6 __________________________________________________________________ chr17:24052253|24058139 __________________________________________________________________ Ranbp9 __________________________________________________________________ chr5:43659900|43661144 __________________________________________________________________ Epha7 __________________________________________________________________ chr4:67449579|67450836 __________________________________________________________________ Braf __________________________________________________________________ chr13:95193009|95196402 __________________________________________________________________ Akt3 __________________________________________________________________ chr5:59831651|59833591 Melk [101]Open in a new tab The circRNA_ID of candidate circRNAs and the corresponding host gene. Figure 5. [102]Figure 5 [103]Open in a new tab Expression changes in candidate circRNAs derived from apoptosis-related genes during postnatal development in the rat retina. The number of circular junction reads of a circular RNA (circRNA) was normalized to get the transcripts per million (TPM) that was used to compare the expression of circRNAs between groups. The TPM is used as the unit of the y-axis. Dotted lines represent the “y=0” line, and the open circles on the dotted lines mean that the circRNA candidates were not detected. Discussion In this study, we explored physiologic neuronal apoptosis and the expression of circRNAs during early development in the rat retina. Hundreds of circRNAs were abundantly and differentially expressed during the early development of the rat retina. Expression of circRNA followed several patterns: 1) stage-specific expression or suppression, 2) gradual increases or decreases between time points, or 3) stable expression at all stages of early development. We confirmed age-dependent physiologic neuronal apoptosis was occurring in the developing rat retina. In addition, circRNAs originated from genes that were associated with apoptosis (such as Thoc1, Akt3, Rere, Rbm5, Ranbp9, Pias1, Optn, Naa35, Melk, Epha7, Ddx19a, Braf, Bnip2, Birc6, and Arhgap10) changed dramatically in the developing rat retina, suggesting a link between circRNA and physiologic neuronal apoptosis. Numerous changes in the nervous system occur during the period of peak developmental apoptosis (P4 to P14) [[104]9,[105]10]. For example, the N-Methyl-D-aspartate (NMDA) receptor intracellular signaling pathway switches from the ERK1/2 to the p38 mitogen-activated protein kinase pathway [[106]14]. CircRNA is non-coding RNA mainly composed of exons and is abundant in eukaryotic cells, with highly conserved spatiotemporal specificity [[107]25,[108]26,[109]39]. Several studies have shown that circRNAs are associated with physiologic and pathological processes, such as development of the nervous system [[110]28,[111]40], atherosclerosis [[112]27,[113]41], and cancer [[114]42-[115]44]. The present data indicate that circRNAs may play an important role in the development of the rat retina. According to these results, apoptosis increased to a peak at or near P7 in the developing rat retinal ganglion cell layer [[116]11]. Although the mechanisms underlying circRNA function are not well-known, studies by Hansen et al. and Memczak et al. showed that circRNA participates in post-transcriptional regulation by acting as a miRNA sponge, interfering with cognate linear isoforms [[117]28,[118]29,[119]42]. In addition, researchers showed that circRNA can bind with Argonaute proteins [[120]29] or RNA-binding proteins, or can base-pair with other RNAs [[121]45,[122]46], to interfere with the transcription and translation of mRNA. Some circRNAs can even produce proteins [[123]47]. Through a bioinformatics analysis of circRNAs at different ages of the developing rat retina, the expression of circRNA in 15 circRNAs originating from genes associated with apoptosis was found to change over time. Previous studies have showed large numbers of circRNAs that are endogenous to mammalian cells, and many of these circRNAs are much more stable than linear RNAs [[124]30, [125]48]; therefore, we believe the observed differences in the expression of transcripts per million (TPM) of the candidate circRNAs were mostly validated rather than simply caused by outliers. As circRNAs can interfere with their cognate linear RNAs, these circRNAs may be linked to physiologic neuronal apoptosis in the developing retina. For example, the Pias1 gene (also known as STAT1) encodes a protein involved in the G1/S transition of the mitotic cell cycle and is a negative regulator of apoptotic processes. The miRNA rno-miR-204-3p is one of the predicted targets of circRNA from Pias1; upregulation of this miRNA induces glioma cell apoptosis [[126]49]. The present study detected low expression of the circRNA from this gene at P7, which may lead to high levels of this miRNA and silencing or cleavage of its cognate mRNA, thus facilitating the apoptotic process. In another example, the Rere gene encodes a member of the atrophin family of arginine-glutamic acid dipeptide repeat-containing proteins; this protein colocalizes with a transcription factor in the nucleus that triggers apoptosis when overexpressed. Notably, overexpression of Rere-derived circRNA may absorb (sponge) rno-miR-128-3p miRNA, which may result in overexpression of the relevant protein (Mapk14, a proapoptotic protein) and subsequent apoptosis [[127]50]. Thus, we believe that the identified circRNAs are potentially important to the processes of physiologic neuronal apoptosis during maturation of the nervous system, acting as microRNA sponges or through other mechanisms mentioned above. The present study has the following limitations. We did not precisely classify the data or screen other circRNAs that may be related to apoptosis, nor did we verify the candidate circRNAs and their specific effects in physiologic neuronal apoptosis. In addition, we did not analyze the expression patterns of the host genes of the candidate circRNAs, and their effects on microRNA or mRNA, and the expression profiles of the circRNAs in the adult retina. Such studies should be performed in the future. However, we did predict possible functions or pathways through Gene Ontology (GO) enrichment analysis, pathway enrichment analysis, and miRNA-binding target prediction (data not shown), which provide fertile areas for further study. In conclusion, we confirmed the temporal dependence of developmental physiologic neuronal apoptosis and found 15 circRNAs that had significant differential expression patterns and may be involved in developmental apoptosis. As the nervous system is highly vulnerable to ethanol, general anesthetics, and other substances during this postnatal developmental period, these results may offer potential mechanisms for regulating drug-induced apoptosis in the nervous system. Acknowledgments