Comprehensive analysis of circRNA expression profiles and circRNA‐associated competing endogenous RNA networks in the development of mouse thymus

Abstract The thymus plays an irreplaceable role as a primary lymphoid organ. However, the complicate processes of its development and involution are incompletely understood. Accumulating evidence indicates that non‐coding RNAs play key roles in the regulation of biological development. At present, the studies of the circRNA profiles and of circRNA‐associated competing endogenous RNAs (ceRNAs) in the thymus are still scarce. Here, deep‐RNA sequencing was used to study the biological mechanisms underlying the development process (from 2‐week‐old to 6‐week‐old) and the recession process (from 6‐week‐old to 3‐month‐old) of the mouse thymus. It was found that 196 circRNAs, 233 miRNAs and 3807 mRNAs were significantly dysregulated. The circRNA‐associated ceRNA networks were constructed in the mouse thymus, which were mainly involved in early embryonic development and the proliferation and division of T cells. Taken together, these results elucidated the regulatory roles of ceRNAs in the development and involution processes of the mouse thymus.


| INTRODUC TI ON
The thymus is defined as a primary lymphoid organ with its inimitable role in T cell maturation, education and selection. 1 It is the first formed lymphoid organ, but may also the most rapidly ageing tissue in the body. 2 In response to the demand for large numbers of mature T cells, the thymus grows fast in the early stages of post-natal life and reaches a plateau at 4-6 weeks of age in the mouse. 3 Afterwards, it starts the ageing process, of which the typical phenotype is thymic involution. The significantly decreased size and structural complexity of the thymus with age were manifested by the loss of thymocytes and thymic epithelial cells and the increase in the content of adipocytes. 4 Actually, the processes about the rapid development and gradual involution of the thymus are still unclear. 3 MicroRNAs (miRNAs) are single-stranded (average size 22 nt) non-coding RNAs that regulate gene expression at the post-transcriptional level. 5 Through their broad effects on gene expression, miRNAs participate in the regulation of many cellular processes such as organismal development, 6,7 cellular differentiation and apoptosis 8,9 as well as mitochondrial metabolism. 10 MiRNAs also take part in many aspects of T cell immunity such as T cell maturation, differentiation, activation and ageing. 11 For example, some miRNAs including miR-155, miR-181c, miR-9 and miR-31 can regulate T cell activation by modulating the IL-2 signalling pathway. 12 Moreover, miR-181a-5p and miR-195a-5p expressed in thymic epithelial cells are involved in thymus involution, which directly target transforming growth factor β receptor 1 (Tgfbr1) and Smad family member 7 (Smad7) respectively. 13,14 Circular RNAs (circRNAs) are a class of single-stranded closed RNA molecules without 3′-poly (A) and 5′-cap structures, which have been widely found in plants, animals and human beings. 15 The cell-type-specific, tissue-specific or developmental stage-specific expression profiles of circRNAs suggest their regulatory functions in biological processes. 15,16 Numerous studies have reported that circRNAs are closely related to the tumorigenesis, 17,18 neurodegenerative diseases, 19 cardiovascular diseases 20 and immune diseases. 21,22 Some circRNAs might be involved in post-transcriptional regulation by functioning as 'sponges' of miRNAs. Therefore, the circRNA-miRNA-mRNA networks may influence multiple biological pathways. 23 At present, integrative analysis of circRNA-associated competing endogenous RNAs (ceRNA) networks in the development and involution of the thymus are still scarce.
To gain further insight into the molecular events associated with thymic development and involution, RNA-seq data were systematically analysed to identify aberrantly expressed circRNAs, miRNAs and mRNAs among mice thymuses at 2, 6 weeks and 3 months, respectively. In addition, the circRNA-miRNA-mRNA networks were constructed. This is the first comprehensive high-throughput sequencing analysis of circRNA, miRNA and mRNA expression profiles in the thymus, which deepen our understanding of thymic development and involution.

| Animal tissues preparation
Male C57BL/6 mice were purchased from Beijing HFK Bioscience CO. LTD and maintained in a specific-pathogen-free environment.
The mice were provided free access to standard diet until they met age requirements (2 weeks, 6 weeks and 3 months). Two biological replicates at each time point were used for sequencing. All the animal study protocols were approved by the Committee for the Ethics

| RNA-sequencing and miRNA-sequencing
Details of the RNA-seq and miRNA-seq methods are described in Table S15.

| Integrated analysis of circRNAs-miRNAs-mRNAs
CircRNAs were blasted against the circBase for annotation. Those cannot be annotated were defined as novel circRNAs. For circRNAs that have been annotated in circBase, the target relationship with miRNAs can be predicted by StarBase (v2.0). For novel circRNAs, three softwares Mireap, Miranda (v3.3a) and TargetScan (v7.0) were used to predict targets. For the prediction of mRNAs interacting with circRNAs and miRNAs, miRTarBase (v6.1) was used to predict mRNAs targeted by miRNAs sponge. The triple network was finally built based on the ceRNA theory and the resulting correlation of circRNAs-miRNAs-mRNAs can be visualized by Cytoscape 3.01.

| KEGG enrichment analysis of differentially expressed genes
To further understand the underlying biological mechanisms and pathways of ceRNA-related genes, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was implemented using the clusterProfiler R package.

| Data access
All raw and processed sequencing data have been submitted to the NCBI Gene Expression Omnibus (GEO; http://www.ncbi.nlm.nih. gov/geo/) under accession number GSE139653.

| Changes in thymus tissues
The age-related changes were observed in the thymuses derived from the 2-week-old, 6-week-old and 3-month-old male mice ( Figure S1). Consistent with the literature, the mouse thymus grows rapidly after birth, reaches the maximal size by sexual maturity and then gradually involutes. 2 The volume of the thymus increases significantly in the 6-week-old mouse compared with that in the 2-week-old mouse. As expected, the thymus shows remarkable atrophy in 3-month-old mouse compared with that in 6-week-old mouse. Cluster analyses and heat-maps about the expression of these circRNAs, miRNAs and mRNAs were conducted respectively ( Figure 1C,F,I). Detailed information is listed in Table S1-S6.

| Construction of the ceRNA network
According to the ceRNA hypothesis, competing endogenous RNAs  (Tables S7 and S8).
The ceRNA networks included both positive and negative links ( Figure 3A,B). Figure 3A shows the decreased circRNAs, increased miRNAs and decreased mRNAs and Figure 3B shows the increased circRNAs, decreased miRNAs and increased mRNAs. In the 2W/6W (−) 6W/3M (+) group, 45 DEcircRNAs, 65 DEmiRNAs and 407 DEmRNAs were used to construct the ceRNA networks (Table S9 and S10). As mentioned above, the ceRNA networks also include two types of links in this group ( Figure 4A,B). Figure 4A shows the decreased circRNAs, increased miRNAs and decreased mRNAs and Figure 4B shows the increased circRNAs, decreased miRNAs and increased mRNAs. Information about the differentially expressed transcripts used to construct the networks can be found in Tables S7-S10. And the top15 DEcircRNAs between groups was shown in Tables 1 and 2.  (Table S13 and S14).

F I G U R E 2 Venn diagrams of differentially expressed transcripts.
A, circRNA; B, miRNA; C, mRNA

| D ISCUSS I ON
The thymus is a vital immune organ in humans and other mammals which is the primary site of T cell development. 25 The thymus provides a specialized environment to facilitate T cell maturation and generate an extremely diverse T cell repertoire. It also establishes its self-tolerance to prevent autoimmune diseases. 26 The thymus begins the involution process and apparently exhibits symptoms of age-related changes at an accelerated rate relative to other organs. Previous studies have found that some miRNAs play an important role in thymic development and involution. [27][28][29] However, the role of circRNA-associated ceRNA networks in these processes is unknown. We designed this experiment to further investigate the molecular mechanisms underlining these complex processes of the post-natal rapid development and subsequent gradual decline in the thymus.  32 and also play key roles in the differentiation and proliferation of T cells. 33,34 Moreover, mmu_circ_0000909 was found to be a ceRNA of mmu-miR-92a-1-5p, which targets KMT2D KMT2D is a major enhancer of H3K4me1/2 methyltransferase which is essential for early embryonic development in mice. 35 There is evidence that KMT2D is required for cell-type-specific gene expression and cell differentiation in model systems for adipogenesis and myogenesis. 36 Its mutation is closely related to the occurrence of T and B cell-associated lymphoma. 37,38 In the 2W/6W (−) 6W/3M (+) group, mmu_circ_0000459, novel_circ_000862, novel_circ_001902 and novel_circ_002116 were identified as ceRNAs of mmu-miR-125a-5p, which targets Erap1. Erap1 is an M1 zinc metalloprotease family member which has a significant impact on peptide processing function and the repertoire of peptides presented. In addition, the deregulation of ERAP1 can affect CD8 + T cell response. 39,40 Several researches have shown that the absence of ERAP1 expression in mice results in increased production of proinflammatory cytokines. 39 The occurrence of inflammation is closely related to thymic ageing. 4,41 The down-regulation of Erap1 in 3-monthold mice may be associated with the early thymic decline. Besides, mmu_circ_0001431, novel_circ_004960, novel_circ_005000 and novel_circ_005893 were identified as ceRNAs of mmu-miR-135a-5p, which targets Celf1. Celf1, a member of the Celf family, belongs to RNAbinding proteins. Previous research has indicated that Celf1 may affect the embryonic development process 42,43 and involve in the proliferation and division of T cells. 44 Other links in these ceRNA networks are listed in Tables S7-S10.
In order to clarify ceRNA functions in mouse thymus, KEGG

TA B L E 2 Top15
DEcircRNAs in the thymus between 6-wk-old and 3-mo-old mouse

ACK N OWLED G EM ENTS
We thank Guangzhou Genedenovo Biotechnology Co., Ltd for assisting in sequencing and bioinformatics analysis.