Comprehensive Analysis of circRNA Expression Pattern in Senescent EPCs

Introduction: Atherosclerosis(AS) is a major disease affecting the elderly health. CircRNAs are found to play important roles in AS. This study aims to explore the biological functions of circRNAs in senescent EPCs. Methods: Endothelial progenitor cells which can differentiate into endothelial cells and their replicative senescent cells were used. we analyzed the circRNA expression proles using RNA-seq, picked up- and down- regulated circRNAs, next done gene ontology (GO) and KEGG analysis, and last circRNA-miRNA-target interactions were predicted. Results: Many circRNAs were identied as signicantly dysregulated in senescent cells. GO analysis showed that differential expressed circRNAs origin genes were involved mainly in metabolic process, immune system process and response to stimulus. KEGG pathway enrichment is involved mainly in endocytosis, ubiquitin mediated proteolysis, mTOR signaling pathway and AMPK signaling pathway. Furthermore, we generated a dysregulated circRNA-miRNA binding and found four circRNAs (has_circ_0007710, _0067750, _0007029, _0002202) were related to atherosclerosis. Conclusions: These results show that circRNAs might play crucial roles in the development of atherosclerosis in aging.


Introduction
Cardiovascular disease (CVD) is one of the most important chronic diseases that affect aging health. Cardiovascular risk increases dramatically with age [1]. Nearly 50% of the world's population died of CVD, it brings enormous burdens to the whole world. Atherosclerosis AS is the pathological basis of CVD. It's risk factors include increased abnormal blood lipid level, smoking, diabetes and aging, of which aging is an independent risk factor [2]. Endothelial dysfunction is the initiation of AS. Endothelial progenitor cells (EPCs) are a kind of stem cells which can differentiate into endothelial cells (ECs) and repair injured ECs [3,4] . While EPCs senescence limits this repair ability, and thus promotes the progress of AS.
Circular RNA (circRNA) is a class of special non-coding RNA molecules which comes from the exon or intron sequences of a gene [5]. CircRNA differs from linear RNA structures because it has no 5' cap or polyA tail, it forms covalently closed loops, so it is more stable than linear RNA [6]. It is reported that circRNA played important roles in the occurrence and development of atherosclerosis [7,8]. CircRNA can serve as miRNA sponges to regulate translational process or directly bind to proteins to regulate protein localization and function [9]. For instance, it is con rmed that circFOXO3 in uences cell senescence, the cell cycle, and apoptosis [10]. Circular ANRIL(cANRIL), whose expression affects the human INK4a/ARF locus, is related with the risk of human atherosclerosis [11]. CircRNA is involved in the occurrence and development of atherosclerosis.
In our study, we seeked the circRNA expression pro les and identi ed differentially expressed (DE) circRNAs in the senescent EPCs by using RNA-seq analysis. Subsequently, we found the up-and downregulated circRNAs, and predicted the binding microRNAs. According to the ceRNA theory, a DEcircRNA-microRNA network was constructed and gene ontology (GO) enrichment analysis and KEGG pathway enrichment analysis were performed to explore the potential regulatory functions of circRNAs. Our ndings provide new evidence for understanding the molecular mechanisms of circRNAs in the pathogenesis of senescence and senescence related AS.

Data And Materials
Culture and identi cation of EPCs This work was approved by the Ethics Committee of Xinhua Hospital A liated to Shanghai Jiao Tong University, Shanghai, China, and consent from donors was received. Human cord blood mononuclear cells (MNCs) isolated by density gradient centrifugation with Histopaque-1077 (Sigma) were suspended in complete EGM-2 medium (Lonza Clonetics) and then seeded on bronectin (Gibco) precoated 24-well plates. Approximately 9 × 10 6 cells were collected from 20 ml of cord blood, and 1.5 × 10 6 cells were plated per well. Medium was changed every 3 days until the rst passage. For identi cation of an EPC, cells were stained with antibodies against endothelial markers CD34 and VEGFR-2 as well as stem/progenitor marker CD133 by uorescence microscopy, and by dual staining for acetylated lowdensity lipoprotein (Dil-ac-LDL; Molecular Probes) uptake and UEA-1 (Sigma) binding.
Replicative senescence model and SA-β-Gal assay Through repeated subculture to the 20th passage (P20), the EPCs gradually appeared large in size and irregular in shape, and some were tree branch-like and polygonal or long spindle-shaped. SA-β-Gal staining (Cell Signaling Technology, America) was performed according to the manufacturer's instructions, 1 × 10 5 cells were rinsed three times with PBS, xed with 0.5 ml of xative solution for 20 min at room temperature, and rinsed again with PBS. Then the cells were incubated with the staining solution overnight at 37℃, and senescent cells were identi ed as bluish green stained cells under a phase-contrast microscope.
Total RNA isolation and sequencing Total RNA was isolated from six EPCs samples, including three S4 samples and three S20, TRIzol reagent (Thermo Fisher, American) was used according to the manufacturer's protocol. After the total RNA was extracted, the mRNA was enriched with a magnetic bead with Oligo (dT) and the appropriate amount of interrupt reagent was added to the resulting mRNA under high temperature conditions to segments, a cDNA strand was synthesized as the interrupted mRNAs for a template, then two-stranded cDNA was synthesized by the two-stranded synthesis reaction system, and the kit was used to pure recycling, viscous end repairment, 3' end of cDNA plusing base "A" and connecting the connector, then selecting fragment size, and nally performing PCR amplifying. Con gured libraries were quali ed by Agilent 2100 Bioanalyzer and ABI StepOnePlus Real-Time PCR System, and quali ed for sequencing.

RNA-seq data analysis
The raw data produced by sequencing is called raw reads. First, we lter out low-quality, connector contamination, and unknown base N-content high reads, ltered data is called clean reads. The clean reads are then compared to the reference genome, via CIRI and nd_circ softwares to predict circRNA.
After combining the results of the two softwares, the circRNA is quanti ed and differentially expressed analysis, and the circRNA origin genes are commented for GO, KEGG. Finally, the genes from the DEcircRNAs are done for GO functional enrichment analysis and Pathway functional enrichment analysis.

Data ltering
The data is ltered by the software SOAPnuke (https://github.com/BGI-exlab/SOAPnuke). The steps include: remove reads of containing connectors (connector contamination), remove reads with an unknown base of concluding more than 5% N, remove low-quality reads (the proportion of a base with a value less than 15 is more than 20% of the total bases of the read is called low-quality read). Filtered clean reads is saved as FASTQ format.

Prediction and annotation of circRNA
Two softwares CIRI and nd_circ (https://sourceforge.net/projects/ciri/; https://github.com/marvin-jens/ nd_circ) are used to predict circRNA, and last the results are integrated two software results according to the start and end locations of circRNA (combining circRNAs within the front and last 10 bases in the start and end positions into one class).

CircRNA expression analysis
This study calculates the expression of circRNA based on the number of junction reads at both ends of circRNA, as two softwares CIRI and nd_circ are used, the nal number of junction reads takes the average of the results of both. It is used RPB (junction reads per billion mapped reads, compared to all reads on the upper genome standardized to a billion after crossing the backspliced number of bits of the junction reads) to homogenize for each sample.

Detection of DEcircRNA
This project uses the DEGseq algorithm for the detection of differential expressed circRNAs. DEGseq [12] method is based on the Poisson distribution.

GO function analysis of DEcircRNA origin genes
Based on the results of GO annotations and the o cial classi cation, we functionally classify genes from DEcircRNAs, use phyper function from R software to do enrichment analysis. The p-value calculation is as follows: (https://en.wikipedia.org/wiki/Hypergeometric_distribution) Pathway functional analysis Based on the results of the KEGG annotations and the o cial classi cation, we classify the genes for biological pathways, using the the phyper function from the R software. The p-value calculation is the same as the GO function analysis.

Culture and identi cation of EPCs
Under the present culture conditions, the adherent cells grew in colonies and showed a cobblestone morphology after two weeks of isolation ( Figure 1A), they were double-positive stained for Dil-Ac-LDL and FITCUEA-I ( Figure 1B). EPCs were con rmed by examining the expression of endothelial cell surface antigen VEGFR-2 and progenitor cell surface antigens CD34 and CD133, as shown in Fig. 1C. All these are considered important markers of EPCs.

Successful establishment of senescent EPCs
The senescent model was veri ed by examining the morphology and markers of senescent cells. The fourth passage (P4) was small in size and had cobble-stone like morphology. The EPCs, repeated subculture to the 20th passage (P20), appeared large in size and irregular in shape, and some were tree branch-like and polygonal or long spindle-shaped (Figure 2A). The positive staining rate of P4 was less than that of P20 by Senescence-associated β-galactosidase (SA-β-Gal) assay (P<0.05) ( Figure 2B).

Identi cation of differentially expressed RNAs in senescent EPCs
In the present study, we used CIRI, nd_circ softwares to predict circRNAs and integrated the both results based on circRNA start and end location. As is shown in gure3A, 2433 in CIRI, 6770 in nd_circ, and 9717 were in both softwares. The circRNA ID number will be provided if it is included in circBase. CircRNA is divided into introgenic circRNA and intergenic circRNA according to its' start and end location in chromosome. In total, we identi ed 18920 (8909 circBase and 10011 non-circBase) circRNAs from P20, including 18199 introgenic circRNAs (96.2%) and 721 intergenic circRNAs (3.8%) ( Figure 3B). In addition, the circRNA transcripts were broadly distributed in all chromosomes ( Figure 3C). To investigate the possible biological functions of circRNAs in senescent EPCs, we analyzed circRNA expression pro le data from senescent samples and normal controls. RNAs with a fold-change >2 and q-value ≤ 0.001 were identi ed as signi cantly differentially expressed. RNA-seq analysis identi ed 11572 and 12721 circRNAs upregulated and downregulated respectively. The top 10 dysregulated circRNAs based on q-value were summarized in Table 1. Hierarchical clustering showed that DEcircRNA expression patterns among samples were distinguishable( gure 3D,E,F). The data suggested that the expressions of circRNAs in senescent EPCs are different from those in control samples.

GO function analysis of the DEcircRNA
According to the results of the DEcircRNA analysis, Gene Ontology (GO) function analysis was performed ( Figure 4). GO can be divided into molecular function, cellular component and biological process. Functional enrichment analysis based on gene ontology biological process terms was performed.
Figure4A shows the number of genes involved in the three GO terms respectively. and Figure 4C shows the top 10 terms according to up-and down-regulated circRNAs. The top ten enriched GO-BP terms, included cellular process, cellular component organization or biogenesis, cell proliferation, biological regulation, metabolic process, immune system process, response to stimulus, regulation of biological process, signaling and pigmentation. We also visualized and clustered the enriched GO-BP terms ( Figure   4B). We found that the enriched GO-BP terms were mainly clustered into metabolic process, intracellular organelle part in cellular component, immune system process, response to stimulus. Furthermore, these enriched GO terms have been implicated in many physiological and pathophysiological activities of atherosclerosis. Therefore, the results indicated that circRNAs in the pathogenesis of atherosclerosis. Therefore, it is reasonable to predict that these circRNAs might contribute to the molecular regulation of atherosclerosis.

KEGG analysis ofthe DEcircRNA
Based on the results of the DEcircRNA analysis, KEGG pathway classi cation and enrichment analysis were done. The KEGG pathway involved mainly in cellular processes, environmental information processing, genetic information processing, human diseases, metabolism and organismal systems and so on. The top 5 enrichment pathway is: Endocytosis, protein processing in endoplasmic reticulum, ubiquitin mediated proteolysis, mTOR signaling pathway and AMPK signaling pathway, and other pathway include insulin resistance, lysin degradation and circadian rhythm. These enriched pathway results have been indicated in process of senescence. Therefore, it is con rmed circRNAs involved in senescence.

CircRNA-miRNA binding analysis
The widely known effect of circRNA is ceRNA and forms circRNA-miRNA. After analyzed the DEcircRNAs, we predicted the probable binding miRNAs. has_circ_0007710 was found upregulated in senescent EPCs, it can bind miR-486-5p, and miR-486-5p was reported an aging-associated miRNA regulator of stem cell aging [14]. miR-486-5p also was reported presumably regulated autophagy by a mTOR dependent mechanism [15], and we know aging is the result of mitochondrial damage, miR-486-5p is capable to modulate aging through the regulation autophagy of mitochondria [16]. has_circ_0067750 was found upregulated in senescent EPCs, while has_circ_0067750 can form has_circ_0067750-miR-141-3p, and miR-141-3p was con rmed altered in response to cellular senescence [17]. has_circ_0007029 was found downregulated in senescent EPCs. It can bind miR-199a-3p, and (COX)-2 of chondrocytes in ammation and degeneration can act as anti-aging effects of polyphenols exactly by targeting mir-199a-3p [18]. It is con rmed that calorie-restricted feeding upregulated the expression of miR-16-5p, and RAW264 cells transfected with a miR-16-5p mimic signi cantly decreased the mRNA expression of IL-1β, IL-6, and TNFα under LPS stimulation, so it is concluded miR-16-5p might be a critical factor involving the antiin ammatory effects of calorie-restricted feeding [19], while has_circ_0002202 can bind miR-16-5p to form ceRNA, and we found has_circ_0002202 upregulated in senescent EPCs. All these results suggest that many circRNAs play critical roles in senescent EPCs, and further studies are required to con rm the potential values of these circRNAs.

Discussion
Atherosclerosis is a kind of chronic in ammatory process, participated by a variety of cells (vascular endothelial cells, smooth muscle cells, Macrophages) in the walls of blood vessels, it begins with the vascular endothelial cell injury, then forms atherosclerotic plaques, eventually results to thrombosis and a series of adverse events [20]. Therefore, it is necessary to nd new mechanisms of atherosclerosis to prevent its occurrence and progression. Many studies have shown that epigenetics acted a key role in regulating the development of atherosclerotic disease [21,22]. During recent years, great efforts have been made to discover new therapeutic targets and biomarkers of atherosclerosis, but the results were limited.
More recently, circRNAs have received attention as new diagnostic markers for diseases including cancers and aging [23,24]. circRNA expression pro les and functions in senescence remain to be determined. circRNAs-miRNAs binding participates in large-scale ceRNA. It has exciting implications for post-transcriptional gene regulation during multiple physiological and pathophysiological processes. Endothelial dysfunction is the initiation of AS. Some studies have indicated that aging stem cell dysfunction is epigenetic [25], we con rmed the mechanism of DNA methylation, histone modi cation in senescent EPCs dysfunction in our former studies [26,27]. And now we found circRNAs also play important roles in senescent EPCs. Our nding added another epigenetic mechanism of senescent EPCs and strengthened the point aging was epigenetic.
There are some limitations, rst, this study only seeks some circRNAs, and through software prediction and connection, but it lacks deep veri cation. It concludes the probable connection with aging and atherosclerosis and the connection with MeCP2, but it doesn't clear the mechanism. These all need deep study.
In conclusion, it is the rst report systematically analyze the dysregulated circRNA-related miRNAs and regulation in senescent EPCs. We characterized a pro le of dysregulated circRNAs that might be prospective clinical markers associated with the development of aging and atherosclerosis. However, based on these results, future work is needed to uncover the underlying molecular mechanisms of circRNAs in aging and atherosclerosis.

Declarations
Ethics approval and consent to participate This work was approved by the Ethics Committee of Xinhua Hospital A liated to Shanghai Jiao Tong University, Shanghai, China, and consent from donors was received.

Consent for publication
No individual personal data included.

Availability of data and material
Not available

Competing interests
The authors declare that they have no con icts of interest.
Funding This paper is funded by National Natural Science Foundation of China(No. 81471399).

Authors' contributions
Chen SY designed the study approved the nal version. Wang F also edited the manuscript and approved the nal version. Wang CL carried out the analysis and wrote/edited the manuscript. Liu F helped interpret the study results, edited the manuscript and approved the nal version.  Figure 1 Cultivation and identi cation of EPCs derived from umbilical cord blood. A, EPCs exhibited a cobblestonelike cell monolayer at 14 days after seeding (×100). B, Uptake of Dil-Ac-LDL and binding of FITC-UEA-1 observed with a uorescence microscope (×400). C, Typical uorescence images of expression of CD34, CD133, and VEGFR-2 in EPCs (×400).

Figure 2
The characteristic of senescent EPCs. A, The morphology of two groups(X100). B, Beta-gal staining of the two groups (X100)