Next Article in Journal
Cigarette and Cannabis Smoking Effects on GPR15+ Helper T Cell Levels in Peripheral Blood: Relationships with Epigenetic Biomarkers
Next Article in Special Issue
Characterization of a Homozygous Deletion of Steroid Hormone Biosynthesis Genes in Horse Chromosome 29 as a Risk Factor for Disorders of Sex Development and Reproduction
Previous Article in Journal
Constitutive and Plastic Gene Expression Variation Associated with Desiccation Resistance Differences in the Drosophila americana Species Group
Previous Article in Special Issue
Androgen Receptor Gene Variants in New Cases of Equine Androgen Insensitivity Syndrome
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Genes
Volume 11
Issue 2
10.3390/genes11020147
genes-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Effects of Castration on miRNA, lncRNA, and mRNA Profiles in Mice Thymus

by
Bingxin Li
1,†,
Kaizhao Zhang
1,†,
Yaqiong Ye
2,
Jingjing Xing
1,
Yingying Wu
1,
Yongjiang Ma
1 and
Yugu Li
1,*
1
College of Veterinary Medicine, South China Agricultural University, Guangzhou 510642, China
2
School of Life Science and Engineering, Foshan University, Foshan 528000, China
*
Author to whom correspondence should be addressed.
†
These authors contributed equally to this work.
Genes 2020, 11(2), 147; https://doi.org/10.3390/genes11020147
Submission received: 29 December 2019 / Revised: 16 January 2020 / Accepted: 28 January 2020 / Published: 30 January 2020
(This article belongs to the Special Issue Genomics of Sexual Development and Reproduction in Mammals)
Browse Figures
Review Reports Versions Notes

Abstract

:
Thymic degeneration and regeneration are regulated by estrogen and androgen. Recent studies have found that long non-coding RNAs (lncRNAs) and microRNAs (miRNAs) are involved in organ development. In this study, RNA sequencing (RNA-seq) results showed that ovariectomy significantly affected 333 lncRNAs, 51 miRNAs, and 144 mRNAs levels (p < 0.05 and |log2fold change| > 1), and orchiectomy significantly affected 165 lncRNAs, 165 miRNAs, and 208 mRNA levels in the thymus. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis showed that differentially expressed genes (DEGs) were closely related to cell development and immunity. Next, we constructed two lncRNA–miRNA–mRNA networks using Cytoscape based on the targeting relationship between differentially expressed miRNAs (DEMs) and DEGs and differentially expressed lncRNAs (DELs) analyzed by TargetScan and miRanda. Besides, we screened DEGs that were significantly enriched in GO and in ceRNA networks to verify their expression in thymocytes and thymic epithelial cells (TECs). In addition, we analyzed the promoter sequences of DEGs, and identified 25 causal transcription factors. Finally, we constructed transcription factor-miRNA-joint target gene networks. In conclusion, this study reveals the effects of estrogen and androgen on the expression of miRNAs, lncRNAs, and mRNAs in mice thymus, providing new insights into the regulation of thymic development by gonadal hormones and non-coding RNAs.

1. Introduction

As a primary immune organ, the thymus regulates adaptive immunity by producing naïve T lymphocytes [1]. Thymus degeneration is closely related to age, and the main age-related manifestations are thymus size becomes smaller, reduced thymic epithelial cells, and gradually decreased naive T lymphocyte output [2,3]. When naive T lymphocytes in the peripheral circulation system decrease with aging, effector-memory T lymphocytes will proliferate to maintain a constant number of total lymphocytes [4,5,6]. However, the reactivity of old lymphocytes to new infections is weak. Therefore, the reduction in the number of naïve T lymphocytes is considered to be one of the causes of the decline in resistance to infections and cancer in older adults [7,8]. Moreover, Gui found that the number of thymocytes in female and male mice decreased sharply at 3 months of age [9].
Gonadal hormone levels are also closely related to aging [10]. Estrogen and androgen can regulate thymus development, which has been widely reported [11,12]. Studies have found that pregnancy or estrogen injection can induce thymus atrophy and inhibit the expression of chemokines, such as monocyte chemoattractant protein-1 and stromal cell-derived factor 1 [13,14]. Chemokines are essential for thymocyte migration into the thymus and induce T lymphocyte differentiation [15]. In addition, androgen therapy can also induce thymus atrophy [11]. Heng et al. found that prepubertal castration significantly reduced the rate of thymic degeneration and increased early T lineage progenitors [16]. Interestingly, the thymus rapidly regenerated when ovaries or testicles were removed [17,18]. Therefore, exploring the effects of estrogen and androgen on thymic development is very valuable.
microRNAs (miRNAs) are conserved endogenous 18–26 nucleotide (nt) RNAs that can target the 3’UTR gene region to regulate gene expression [19]. In animal genomes, miRNAs are involved in regulating more than 30% of gene expression [20]. Recent studies have found that miRNAs widely participate in the regulation of thymus aging in mice [2,21]. Our previous study found that miR-195a-5p can inhibit the proliferation of thymic epithelial cells (TECs) by targeting Smad family member 7 [22]. In addition, we also found that miRNA181a-5p may enhance TECS proliferation by regulating TGF-β signal transduction [23]. Although some studies have been performed on the effects of miRNAs on thymic aging, additional miRNAs and their regulatory mechanisms require extensive investigation.
Long non-coding RNAs (lncRNAs) are important components of non-coding RNAs [24]. lncRNAs can regulate the expression of target genes in a variety of ways, such as cis, trans, or combined with miRNAs, so it is not easy to determine their target genes [25,26]. The roles of lncRNAs in organ development, cancer treatment, and aging have been widely reported over the past decades [27,28,29]. In addition, Hu et al. found that lncRNAs play an important role in the differentiation of T helper subsets [30]; Li et al. found that the levels of 17β-estradiol (E2) and progesterone in mice ovarian cells increased after overexpression of lncRNA SRA [24]. A recent study found that androgen and estrogen influence lncRNA expression during penis development in marsupials [31]. Meanwhile, our previous study found that E2 treatment can inhibit cell activity and proliferation, and influence the lncRNA expression profile in TECS [29]. However, the effects of gonadal hormones on RNAs in the thymus have not been reported. In this study, we removed mice ovaries and testes to explore the effects of estrogen and androgen on lncRNA, miRNA, and mRNA expression profiles in the thymus. Subsequently, functional analysis and annotation of differentially expressed genes (DEGs) were performed to further explore the importance of gonadal hormones in thymic degeneration. Finally, we constructed and preliminarily validated the lncRNA–miRNA–messenger RNA (mRNA) regulatory network to determine the key factors of estrogen and androgen regulation of thymic degeneration.

2. Materials and Methods

2.1. Ethics Approval

All the methods in this experiment are based on animal science experimental guidelines. In addition, all animal experiment programs were approved by the Animal Protection Committee of South China Agricultural University (Guangzhou, China), with approval number SCAU 0014.

2.2. Animals and Sample Collection

For the experiment, 48 one-month-old CD1 mice (half male and half female) were purchased from the Laboratory Animal Center of Guangzhou University of Chinese Medicine (license key: SCXK (Yue) 2013-0034). The mice were fed in a specific pathogen-free environment (12/12 h light/dark cycle, 22–24 °C, 40–60% humidity) and randomly subdivided into ovariectomy group (F3x), female control group (F3), orchiectomy group (M3x), and male control group (M3) (12 mice per group). Mice were allowed to adapt for 7 days before the experiment. On the eighth day, the mice in the treatment groups were anesthetized for ovary or testicle removal, and the control groups received the same treatment but were not neutered. At 3 months old [9], the thymuses of each group were collected aseptically, frozen immediately with liquid nitrogen, and stored at –80 °C until use.

2.3. Total RNA Extraction and Sequencing

Referring to known research and economic factors, we mixed the tissues of each group equally for sequencing by a service provider (LC-BIO Bio-tech ltd, Hangzhou, China). According to the manufacturer’s instructions, total RNA was extracted from mice thymus using Trizol reagent (Invitrogen, Carlsbad, CA, USA). A Bioanalyzer 2100 (Agilent, Santa Clara, CA, USA) and agarose gel electrophoresis were used to assess the quantity, purity, and integrity of total RNA. Total RNA with RNA integrity number greater than 7 was selected for further experiments. Then a TruSeq Small RNA Sample Prep Kit (Illumina, San Diego, CA, USA) was used to construct miRNA libraries according to the manufacturer’s instructions. Finally, we performed single-ended sequencing on an Illumina HiSeq 2500 (Illumina, San Diego, CA, USA).
Ribosomal RNA (rRNA) was depleted from DNA-free RNA by an Epicentre Ribo-Zero Gold Kit (Illumina, San Diego, CA, USA). Then rRNA-free RNA was reverse-transcribed to create the final complementary DNA (cDNA) libraries by mRNA-Seq sample preparation kit (Illumina, San Diego, CA, USA). Finally, the libraries were paired-end sequenced on an Illumina HiSeq 4000 (Illumina, San Diego, CA, USA).

2.4. Transcript Assembly

Cutadapt [32] and FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) were used to remove reads containing adaptors and low-quality bases, and to verify sequence quality (Figure 1). Then clean data were mapped to the mice genome using Bowtie 2 [33] and Tophat 2 [34]. Finally, StringTie [35] and Ballgown [36] were used to merge all transcriptomes from this experiment and estimate their expression levels by calculating fragments per kilobase per million reads.

2.5. miRNA Identification

The raw data were processed using the ACGT101-miR internal program (LC Sciences, Houston, TX, USA) to remove adaptor dimers, garbage, low complexity, common RNA families (rRNA, tRNA, snRNA, snoRNA), and repeats. Then we used miRBase 21.0 to define known miRNAs. Finally, unmapped sequences were predicted using RNAfold software (http://rna.tbi.univie.ac.at/cgi-bin/RNAfold.cgi).

2.6. lncRNA Identification

First, we removed all transcripts that overlapped with known mRNAs and were shorter than 200 bp. Then CPC [37], CNCI [38], and Pfam [39] were used to predict transcripts with coding potential. All transcripts with CPC score <–1 and CNCI score <0 were deleted. Finally, the transcripts annotated as i, j, o, u, and x were retained. These letters represent transfrags falling entirely within a reference intron; potentially novel isoform; generic exonic overlap with a reference transcript; and unknown, intergenic, and exonic overlap with reference on the opposite strand, respectively.

2.7. RNA Expression and Functional Analysis

Differentially expressed lncRNAs (DELs), miRNAs (DEMs), and DEGs were screened by the R package Ballgown with log2 (fold change (FC)) > 1 or log 2 (FC) < –1 and statistical significance p < 0.05. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) in the DAVID online tool (https://david.ncifcrf.gov/summary.jsp) were used to analyze the functions of DEGs; p < 0.05 represents significant difference.

2.8. lncRNA–miRNA–mRNA Network Analysis

Recent studies have shown that lncRNAs and mRNAs containing the same miRNA binding site can regulate each other’s expression levels by competitively binding miRNAs. To determine the interactions between DELs, DEMs, and DEGs after ovariectomy or orchiectomy, miRWalk, miRanda, RNAhybrid, and Targetscan were used for screening. Then, the intersections of co-expressed lncRNAs, miRNAs, and mRNAs in the 4 software programs were used to construct a lncRNA–miRNA–mRNA network. Finally, the lncRNA–miRNA–mRNA network was visualized using Cytoscape v3.7 software (https://cytoscape.org/). As competing endogenous RNA (ceRNA) analysis is based on miRNAs, the relationship between lncRNAs and mRNAs needs further investigation.

2.9. Cell Isolation and Culture

Thymocytes were isolated from 1-month-old mice. According to the manufacturer’s instructions, thymocytes were obtained by centrifugation using the lymphocyte separation solution (Dakewe, Shenzhen, China). Red blood cell lysate (Solarbio, Beijing, China) was used to remove red blood cells. Then, cells were diluted with RPMI-1640 (Gibco, Waltham, MA, USA) containing 10% fetal bovine serum (Gibco, Waltham, MA, USA) and 0.2% penicillin / streptomycin (Invitrogen, Carlsbad, CA, USA) and culture for 24 h in 37 °C at 5% CO2 to remove adherent cells. TECS were isolated from mice thymus in our previous experiments [40]. TECS were cultured in DMEM (Gibco, Waltham, MA, USA) containing 10% fetal bovine serum (Gibco, Waltham, MA, USA) and 0.2% penicillin / streptomycin (Invitrogen, Carlsbad, CA, USA) at 37 °C in 5% CO2. Finally, thymocytes and TECS were collected and immediately placed in liquid nitrogen and stored at −80 °C until use.

2.10. Quantitative Real-Time PCR Analysis

Total RNA in mice thymus, thymocytes and TECS were extracted with Trizol reagent (Takara, Kusatsu, Japan). First-strand cDNA was synthesized with a ReverTra Ace quantitative real time-PCR (qRT-PCR) RT Kit (Toyobo, Osaka, Japan), with random hexamers (for mRNAs and lncRNAs) and stem-loop RT primers (for miRNAs, purchased from RiboBio) according to the manufacturer’s instructions. Primers of mRNAs and lncRNAs (Table S1) were designed by Primer Premier 5.0 software (Premier Biosoft International, USA), and miRNA primers were purchased from RiboBio (Guangzhou, China). β-actin (for mRNAs and lncRNAs) and U6 (for miRNAs) were used to normalize the data. qRT-PCR was performed by a Bio-Rad CFX96 Real-Time PCR system (Bio-Rad, Hercules, CA, USA) using SYBR Green Real-Time PCR Master Mix (Toyobo) following the manufacturer’s instructions. All qRT-PCR assays were conducted in triplicate, and the relative levels were measured in terms of threshold cycle (Ct) and calculated using the formula 2−∆∆Ct [41].

2.11. Analysis of Transcription Factor Binding Sites of DEGs

The putative promoter sequence was obtained from Browser1 of the UCSC genome. TESS software v6.0 was used to searched Position-weigh matrices in the TRANSFAC database [42]. The relative fraction line was set to 0.9. Next, the internal PERL script was used to conduct hypergeometric testing. p < 0.05 was defined as enriched transcription factor.

2.12. Identification of Transcription Factor-Related miRNAs and Their Joint Target Genes

According to the analysis of DEGs-related transcription factors, we used CircuitsDB (http://biocluster.di.unito.it/circuits/) to obtain information about miRNAs and joint target genes related to these transcription factors. Cytoscape was used to build transcription factor-miRNAs-joint target genes networks

3. Results

3.1. Expression Profiles of miRNAs, lncRNAs, and mRNAs in Thymus

To identify putative transcripts in the thymus, five mice each from the F3x group, F3 group, M3x group, and M3 group were mixed together for sequencing. In this experiment, we obtained a total of 380 million uniquely mapped lncRNA reads (Table 1). In addition, we got 56,673 unique miRNA readings in the Rfam database (Table 2, Supplementary Figure S1). In total, we identified 17,498 candidate lncRNAs, 1528 miRNAs, and 46,595 mRNAs in all chromosomes (Supplementary Tables S2–S4). Venny analysis showed that 13,343 lncRNAs, 878 miRNAs, and 37,643 mRNAs were expressed in all groups, accounting for 76.3%, 57.5%, and 80.8% of the total readings, respectively (Figure 2a–c). Chromosome analysis showed that lncRNAs were unevenly distributed on the chromosome and were basically consistent with the trends of chromosome length (Figure 2d,e). Interestingly, the minimum length, mean length, median length, and N50 of lncRNAs were very close to those of mRNAs in this experiment (Figure 2f, Table 3). In addition, the distribution of miRNAs showed that miRNA length was mainly concentrated in 20–24 nt (Figure 2g).

3.2. Identification of DELs, DEMs, and DEGs

Data analysis showed that ovariectomy significantly affected the levels of 333 lncRNAs, 51 miRNAs, and 144 mRNAs in mice thymus (p < 0.05, FC > 2 or FC < 0.5) (Figure 3a–c, Supplementary Tables S5–S7). Orchiectomy significantly affected the levels of 165 lncRNAs, 165 miRNAs, and 208 mRNAs (Figure 3d–f, Supplementary Tables S8–S10). To further explore the effects of gender on thymic degeneration, we analyzed differential RNAs in the F3x and M3x groups. The results showed that 416 lncRNAs, 68 miRNAs, and 225 mRNAs were remarkably differentially expressed after ovariectomy and orchiectomy (Supplementary Figure S2a–c, Supplementary Tables S11–S13). In addition, the levels of 136 lncRNAs, 129 miRNAs, and 124 mRNAs were remarkably different between the F3 and M3 groups (Supplementary Figure S2d–f, Supplementary Tables S14–S16). The analysis results show that most DEGs were upregulated after ovariectomy when compared with the F3 group, whereas most miRNAs and lncRNAs were downregulated (Figure 3a–c). It is worth noting that most DELs were upregulated after orchiectomy when compared with the M3 group (Figure 3d–f). Finally, 37 DELs, 10 DEMs, and 18 DEGs in the ovariectomy and female control groups and 27 DELs, 14 DEMs, and 36 DEGs in the orchiectomy and male control groups were used to construct lncRNA–miRNA–mRNA networks (Figure 3g–l).
In order to verify our RNA-seq data, we randomly selected 3 DEGs, 3 DELs, and 3 DEMs for qRT-PCR analysis (Figure 4). The results showed that the differentially expressed RNAs had the same expression trends in qRT-PCR and RNA-seq (R = 0.747, p = 0.003).

3.3. Functional and Pathway Annotation

GO terms and KEGG pathway annotations of DEGs were enriched via DAVID. The enriched GO terms were divided into 3 categories: biological process (BP), cellular component (CC), and molecular function (MF). In this study, the BP enrichment of DEGs after ovariectomy includes extracellular matrix organization (p = 0.015), peptide transport (p = 0.024), extracellular structure organization (p = 0.040), and cellular cation homeostasis (p = 0.042) (Figure 5a). KEGG enrichment of DEG after ovariectomy showed no significant enrichment pathway. Meanwhile, the BP enrichment of DEGs after orchiectomy includes immune response (p = 0.00000054), response to virus (p = 0.000025), defense response (p = 0.000036), inflammatory response (p = 0.0017), regulation of nitric oxide biosynthetic process (p = 0.0091), positive regulation of defense response (p = 0.0094), cellular alkene metabolic process (p = 0.011), positive regulation of response to stimulus (p = 0.014), production of nitric oxide during acute inflammatory response (p = 0.015). KEGG enrichment of DEG after orchiectomy includes arachidonic acid metabolism (p = 0.014), viral myocarditis (p = 0.020), and intestinal immune network for IgA production (p = 0.045).

3.4. lncRNA–miRNA–mRNA Network Construction and Visualization

As a member of ceRNAs, lncRNA can affect the expression of miRNA target genes by adsorbing miRNA. In this study, we screened co-expressed genes from differentially expressed RNAs to construct ceRNA networks. As shown in Figure 6a, the female lncRNA–miRNA–mRNA network included 37 lncRNAs, 10 miRNAs, and 18 mRNAs. The male lncRNA–miRNA–mRNA network included 27 lncRNAs, 14 miRNAs, and 36 mRNAs (Figure 6b). The results of the node degree distribution analysis in the ceRNAs network showed that the slope of the power law distribution in female was −1.115 and R2 = 0.723, and the slope in the male was -1.306 and R2 = 0.795, which indicates that these two ceRNA networks have typical bio-network scale-free feature (Figure 6c,d).

3.5. Expression of DEGs in Different Cells of Thymus

In order to further explore the function of DEGs in thymus development, we screened DEGs that were significantly enriched in GO and in the ceRNA network for verification. Then, qRT-PCR was used to detect the expression of these DEGs in thymocytes and TECS. As shown in Figure 7, H2-M2 was highly expressed in thymocytes. SLC7A2, CAV1, CACNB4, and PTGS2 were highly expressed in TECS. The expression trend of POLR3K in thymocytes and TECS was consistent.

3.6. Transcription Factors of DEGs

In order to explore the regulatory mechanism of DEGs, we used TESS software to predict enriched transcription factors of upregulated and downregulated DEGs. The results showed that the binding sites of KROX, Lmo2, PEA3, STAT1, and STAT4 were significantly overexpressed in the downregulated DEGs after ovariectomy (Figure 8a), and the binding sites of AP-1, AP-2rep, c-Myc, E2F-1 DP-2, ETF, MOVO-B, Pax, and ZF5 were significantly overexpressed in the upregulated DEGs (p < 0.05) (Figure 8b). Meanwhile, the binding sites of NF-kappaB, Pbx-1, SMAD, and SREBP-1 were significantly overexpressed in the downregulated DEGs after orchiectomy (Figure 8c), and the binding sites of Hoxa4, ICSBP, IRF, ISGF-3, MOVO-B, NF-γ, Oct-4, and TEF-1 were overexpressed in the upregulated DEGs (p < 0.05) (Figure 8d).

3.7. DEGs Transcription Factor-related DEMs and Joint Target Genes prediction

We use CircuitsDB to predict the transcription factors of miRNAs associated with DEGs, and predict the joint target genes. We found that one transcription factor (PEA3), which was related to DEGs after ovariectomy, was associated with mmu-miR-135a and 15 joint target genes (Figure 9a). In addition, we found that 5 DEGs-related transcription factors are associated with DEMs after orchiectomy. PBX-1 was associated with 3 DEMs, including mmu-miR-342, mmu-miR-135a, and mmu-miR-135b, and 5 joint target genes (Figure 9b). HOXA4 was associated with mmu-miR-135b and 2 joint target genes (Figure 9c). IRF was associated with mmu-miR-15a and 2 joint target genes (Figure 9d). SREBP-1 was associated with mmu-miR-200b and 3 joint target genes (Figure 9e). NF-γ was associated with 3 DEMs, including mmu-miR-342, mmu-miR-206, and mmu-miR-144, and 24 joint target genes (Figure 9f).

4. Discussion

The thymus, as a central immune organ, is the main site for generation of naive T lymphocytes [43]. Thymus atrophy is considered to be a biomarker of immune system aging and is influenced by many factors, such as age, hormones, and gender [44]. Studies have extensively demonstrated that ovariectomy or orchiectomy can result in transient thymic regeneration. However, the mechanisms by which gonadal hormones regulate thymic development remains unclear [17,18]. In addition, the role of non-coding RNA in regulating thymic development is revealed. In this study, 333 DELs, 51 DEMs, and 144 DEGs, and 165 DELs, 165 DEMs, and 208 DEGs were identified after ovariectomy and orchiectomy, respectively. In order to identify the accuracy of the RNA-seq data in this experiment, we randomly selected some differentially expressed RNAs for qRT-PCR identification. The results showed that DELs, DEMs, and DEGs had the same expression trends in qRT-PCR and RNA-seq (R = 0.747, p = 0.003), which parts supported the high reliability of our RNA seq data.
Next, GO and KEGG in DAVID were used to enrich the function of DEGs after ovariectomy and orchiectomy. DEGs after ovariectomy were mainly concentrated in extracellular structure organization, peptide transport, cellular cation homeostasis, and extracellular matrix organization in BP. Upregulation of peptide transport may help regulate lymphocyte proliferation [45]. Besides, cellular cation homeostasis is necessary for normal life activities of cells. For instance, K+ is used to compensate negative charges and activate protein translation processes [46]; Mg2+ is an important component of cell proliferation-related cofactors [47]. It is worth noting that the thymic extracellular matrix can increase thymocyte output in vivo and promote TEC differentiation in vitro [48]. Meanwhile, DEGs after orchiectomy were mainly concentrated in immune response and inflammatory response in BP. Castration or androgen deficiency may result in thymus enlargement, which may be related to androgen inhibiting immature thymocyte proliferation and accelerating its apoptosis [49]. In addition, androgens can affect the production of immune-related cytokines in the thymus. Upregulation of TGF-β and activation of TGF-β-induced Smad signaling pathway induces accumulation of connective fibers in the thymus [50]. As a proinflammatory cytokine, IL-6 can also stimulate the proliferation of thymocytes in vitro [51]. Besides, Oxidative stress affects the homeostasis of the immune system, causing the irreversibility of thymus degradation [52].
lncRNAs can participate in proliferation and aging physiological processes in different ways [53,54]. With the development of non-coding RNA research, the theory of ceRNA was proposed. In this theory, acting as key nodes, miRNAs can bind to lncRNAs to regulate the expressions of other non-coding RNAs or target genes with the same binding sites. In this experiment, we hypothesized that lncRNA is involved in regulating thymus development after mice castration through ceRNA. Next, we performed ceRNA analysis on the differentially expressed RNAs between castration and control groups and constructed two scale-free biological networks of ceRNAs. qRT-PCR preliminarily confirmed the expression of these genes in the ceRNA network were consistent with RNA-seq. In addition, the scale-free distribution of these two ceRNA networks indicates that there are important nodes in the regulation of gonadal hormones on thymus development.
In order to further clarify the function and expression site of DEGs in the ceRNA network. we screened DEGs (H2-M2, POLR3K, SLC7A2, CAV1, CACNB4, and PTGS2 that were significantly enriched in GO and in the ceRNA network to verify their expression in thymocytes and TECS. Cav1 is a channel for Ca2+ influx required for T cell development [55]. Cyclooxygenase 2 (COX-2) is encoded by the PTGS2 gene, which inhibits T cell proliferation by regulating NF-kappaB [56]. This experiment found that Cav1 and PTGS2 were significantly overexpressed in TECS, which may be related to the TEC providing a microenvironment for thymocyte development [23]. The qRT-PCR showed that H2-M2 and POLR3K were highly expressed in thymocytes, and SLC7A2 and CACNB4 were highly expressed in TECS. However, the roles of these genes in the thymus have not been reported. These genes may be worthy of further study.
In order to explore the regulatory mechanism of DEGs, we predicted the causal transcription factors of DEGs through enrichment experiments. We found that the binding sites of KROX, Lmo2, PEA3, STAT1, and STAT4 were significantly enriched in downregulated DEGs after ovariectomy, while AP-1, AP-2rep, c-Myc, E2F-1, DP-2, ETF, MOVO-B, Pax, and ZF5 binding sites were significantly enriched in upregulated DEGs. Meanwhile, the binding sites of NF-kappaB, Pbx-1, SMAD, and SREBP-1 were significantly enriched in downregulated DEGs after orchiectomy, while the binding sites of Hoxa4, ICSBP, IRF, ISGF-3, MOVO-B, NF-γ, Oct-4, and TEF-1 were significantly enriched in upregulated DEGs. The main function of Krox-20 is to inhibit the c-Jun NH2-terminal protein kinase-c-Jun pathway, whose activation is necessary for both proliferation and death [57]. Overexpression of LMO2 leads to delayed development of thymus progenitor cells and increases the number of T lymphocytes in lymphoid organs [58]. As members of the STAT family, TAT1 and STAT4 mainly promote the differentiation of Th1 cells [59]. AP-1 is one of the most important target members in the MAPK pathway, which can regulate cell division and apoptosis [60]. AP-2 rep is a transcriptional repressor of AP-2 and plays an important role in embryo development [61]. c-Myc makes cells sensitive to TNF stimulation and promotes apoptosis of thymocytes [62]. E2F-1 regulates apoptosis and proliferation by stabilizing p53 tumor suppressor [63]. The NF-kappaB pathway is considered to be a typical transduction pathway of proinflammatory cytokines represented by the IL-1 and TNF receptor families, which was required for T cell proliferation and emigration [64]. Smad is at the core of the TGF-β pathway [65]. HOXA4 expression is inversely related to cell cycle, metastasis, and Wnt signaling pathway, and its overexpression can inhibit tumor cell proliferation, migration and invasion [66]. Oct-4 belongs to the POU domain transcription factor family that are involved in regulating the proliferation and differentiation of various cells in tissues [67]. The function of some transcription factors in thymus development has been partially verified, but the role of PEA3, DP-2, ETF, MOVO-B, Pax, ZF5, Pbx-1, SREBP-1, ICSBP, IRF, ISGF-3, MOVO-B, NF-γ, and TEF-1 in regulating thymus development needs further investigation.
Based on the analysis of DEGs transcription factors, we used CircuitsDB to predict DEMs and associated target genes associated with them. The results showed that mmu-miR-135a was related to PEA3 after ovariectomy, and 7 DEMs (mmu-miR-342, mmu-miR-135a, mmu-miR-135b, mmu-miR-15a, mmu-miR-200b, mmu-miR-206, and mmu-miR-144) were associated with 5 DEGs-related transcription factors (PBX-1, HOXA4, IRF, SREBP, and NF-γ) after orchiectomy (Figure 9). miR-15a promotes the proliferation of mouse lymphocytes by regulating the expression of cell cycle-related genes [68]. miR-200b- 3p can be used as a biomarker for mouse thymus development and degradation [40]. These are consistent with our experimental results. However, studies on mmu-miR-342, mmu-miR-135a, mmu-miR-135b, mmu-miR-206, and mmu-miR-144 in mouse thymus have not been reported, so these genes deserve further investigation.

5. Conclusions

In summary, we identified the expression profiles of lncRNAs, miRNAs, and mRNAs in mice thymus from control groups and ovariectomy or orchiectomy group, and constructed two scale-free ceRNA networks. Bioinformatics analysis showed that DEGs were mainly enriched in extracellular matrix organization and immune response-related physiological functions. This study provides new insights into the effects of estrogen and androgen on thymic development, and the specific mechanisms need to be further explored.

Supplementary Materials

The following are available online at https://www.mdpi.com/2073-4425/11/2/147/s1, Figure S1: Reading miRNAs mapped to the reference genome; Figure S2: Volcano plots of estrogen- and androgen-related co-expression RNAs in mice thymus; Table S1: Specific primer sequences for qRT-PCR; Table S2: List of expression of lncRNAs in all groups; Table S3: List of expression of miRNAs in all groups; Table S4: List of expression of mRNAs in all groups; Table S5: List of DELs in mice thymus between ovariectomy and female control groups; Table S6: List of DEMs in mice thymus between ovariectomy and female control groups; Table S7: List of DEGs in mice thymus between ovariectomy and female control groups; Table S8: List of DELs in mice thymus between orchiectomy and male control groups; Table S9: List of DEMs in mice thymus between orchiectomy and male control groups; Table S10: List of DEGs in mice thymus between orchiectomy and male control groups; Table S11: List of DELs in mice thymus between ovariectomy and orchiectomy groups; Table S12: List of DEMs in mice thymus between ovariectomy and orchiectomy groups; Table S13: List of DEGs in mice thymus between ovariectomy and orchiectomy groups; Table S14: List of DELs in mice thymus between female control and male control groups; Table S15: List of DEMs in mice thymus between female control and male control groups; Table S16: List of DEGs in mice thymus between female control and male control group.

Author Contributions

B.L., K.Z., and Y.L. conceived and designed the experiments; B.L., Y.Y., and J.X. participated in data analysis; B.L. and Y.W. carried out experimental verification; B.L. and K.Z. wrote the original manuscript; and Y.M. edited the manuscript. All authors have read and agreed to the published version of the manuscript.

Funding

This study was funded by the National Natural Scientific Foundation of China (No. 31572475).

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Yue, S.; Zheng, X.; Zheng, Y. Cell-type-specific role of lamin-B1 in thymus development and its inflammation-driven reduction in thymus aging. Aging Cell 2019, e12952. [Google Scholar] [CrossRef] [PubMed]
  2. Guo, D.; Ye, Y.; Qi, J.; Tan, X.; Zhang, Y.; Ma, Y.; Li, Y. Age and sex differences in microRNAs expression during the process of thymus aging. Acta Bioch. Bioph. Sin. 2017, 49, 409–419. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  3. Dixit, V.D. Impact of immune-metabolic interactions on age-related thymic demise and T cell senescence. Semin. Immunol. 2012, 24, 321–330. [Google Scholar] [CrossRef]
  4. Hamazaki, Y.; Sekai, M.; Minato, N. Medullary thymic epithelial stem cells: role in thymic epithelial cell maintenance and thymic involution. Immunol. Rev. 2016, 271, 38–55. [Google Scholar] [CrossRef]
  5. Appay, V.; Sauce, D. Naive T cells: the crux of cellular immune aging? Exp. Gerontol. 2014, 54, 90–93. [Google Scholar] [CrossRef]
  6. Haynes, L.; Swain, S.L. Why aging T cells fail: implications for vaccination. Immunity 2006, 24, 663–666. [Google Scholar] [CrossRef] [Green Version]
  7. Nikolich-Žugich, J. Aging of the T cell compartment in mice and humans: from no naive expectations to foggy memories. The J. Immunol 2014, 193, 2622–2629. [Google Scholar] [CrossRef]
  8. Dorshkind, K.; Montecino-Rodriguez, E.; Signer, R.A. The ageing immune system: is it ever too old to become young again? Nat. Rev. Immunol. 2009, 9, 57–62. [Google Scholar] [CrossRef]
  9. Gui, J.; Morales, A.J.; Maxey, S.E.; Bessette, K.A.; Ratcliffe, N.R.; Kelly, J.A.; Craig, R.W. MCL1 increases primitive thymocyte viability in female mice and promotes thymic expansion into adulthood. Int. Immunol. 2011, 23, 647–659. [Google Scholar] [CrossRef] [Green Version]
  10. ThyagaRajan, S.; Hima, L.; Pratap, U.P.; Priyanka, H.P.; Vasantharekha, R. Estrogen-induced neuroimmunomodulation as facilitator of and barrier to reproductive aging in brain and lymphoid organs. J. Chem. Neuroanat. 2019, 95, 6–12. [Google Scholar] [CrossRef] [PubMed]
  11. Csaba, G. Immunity and longevity. Acat Mirobiol. Imm. H. 2019, 66, 1–17. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  12. Segner, H.; Verburg-van Kemenade, B.M.L.; Chadzinska, M. The immunomodulatory role of the hypothalamus-pituitary-gonad axis: proximate mechanism for reproduction-immune tradeoffs? Dev. Comp. Immunol. 2017, 66, 43–60. [Google Scholar] [CrossRef] [PubMed]
  13. Nordqvist, J.; Bernardi, A.; Islander, U.; Carlsten, H. Effects of a tissue-selective estrogen complex on B lymphopoiesis and B cell function. Immunobiology 2017, 222, 918–923. [Google Scholar] [CrossRef] [PubMed]
  14. Oertelt-Prigione, S. Immunology and the menstrual cycle. Autoimmun. Rev. 2012, 11, A486–A492. [Google Scholar] [CrossRef]
  15. Jin, C.; Fu, W.; Xie, L.; Qian, X.; Chen, W. SDF-1α production is negatively regulated by mouse estrogen enhanced transcript in a mouse thymus epithelial cell line. Cell Immunol. 2003, 223, 26–34. [Google Scholar] [CrossRef]
  16. Heng, T.S.; Goldberg, G.L.; Gray, D.H.; Sutherland, J.S.; Chidgey, A.P.; Boyd, R.L. Effects of castration on thymocyte development in two different models of thymic involution. J. Immunol. 2005, 175, 2982–2993. [Google Scholar] [CrossRef] [Green Version]
  17. Sutherland, J.S.; Goldberg, G.L.; Hammett, M.V.; Uldrich, A.P.; Berzins, S.P.; Heng, T.S.; Blazar, B.R.; Millar, J.L.; Malin, M.A.; Chidgey, A.P.; et al. Activation of thymic regeneration in mice and humans following androgen blockade. J. Immunol. 2005, 175, 2741–2753. [Google Scholar] [CrossRef]
  18. Windmill, K.F.; Lee, V.W.K. Effects of castration on the lymphocytes of the thymus, spleen and lymph nodes. Tissue Cell 1998, 30, 104–111. [Google Scholar] [CrossRef]
  19. Ouyang, D.; Xu, L.; Zhang, L.; Guo, D.; Tan, X.; Yu, X.; Qi, J.; Ye, Y.; Liu, Q.; Ma, Y.; et al. MiR-181a-5p regulates 3T3-L1 cell adipogenesis by targeting Smad7 and Tcf7l2. Acta Bioch. Bioph. Sin. 2016, 48, 1034–1041. [Google Scholar] [CrossRef] [Green Version]
  20. Lewis, B.P.; Burge, C.B.; Bartel, D.P. Conserved seed pairing, often flanked by adenosines, indicates that thousands of human genes are microRNA targets. Cell 2005, 120, 15–20. [Google Scholar] [CrossRef] [Green Version]
  21. Ye, Y.; Li, D.; Ouyang, D.; Deng, L.; Zhang, Y.; Ma, Y.; Li, Y. MicroRNA expression in the aging mouse thymus. Gene 2014, 547, 218–225. [Google Scholar] [CrossRef] [PubMed]
  22. Guo, D.; Ye, Y.; Qi, J.; Xu, L.; Zhang, L.; Tan, X.; Tan, Z.; Yu, X.; Zhang, Y.; Ma, Y.; et al. MicroRNA-195a-5p inhibits mouse medullary thymic epithelial cells proliferation by directly targeting Smad7. Acta Bioch. Bioph. Sin. 2016, 48, 290–297. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  23. Guo, D.; Ye, Y.; Qi, J.; Zhang, L.; Xu, L.; Tan, X.; Yu, X.; Liu, Q.; Liu, J.; Zhang, Y.; et al. MicroRNA-181a-5p enhances cell proliferation in medullary thymic epithelial cells via regulating TGF-β signaling. Acta Bioch. Bioph. Sin. 2016, 48, 840–849. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  24. Li, Y.; Wang, H.; Zhou, D.; Shuang, T.; Zhao, H.; Chen, B. Up-regulation of long non-coding RNA SRA promotes cell growth, inhibits cell apoptosis, and induces secretion of estradiol and progesterone in ovarian granular cells of mice. Med. Sci. Monitor 2018, 24, 2384–2390. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  25. Cai, B.; Li, Z.; Ma, M.; Wang, Z.; Han, P.; Abdalla, B.A.; Nie, Q.; Zhang, X. LncRNA-six1 encodes a micropeptide to activate six1 in cis and is involved in cell proliferation and muscle growth. Front. Physiol. 2017, 8, 230. [Google Scholar] [CrossRef]
  26. Ma, M.; Cai, B.; Jiang, L.; Abdalla, B.A.; Li, Z.; Nie, Q.; Zhang, X. lncRNA-six1 is a target of miR-1611 that functions as a ceRNA to regulate six1 protein expression and fiber type switching in chicken myogenesis. Cells 2018, 7, 243. [Google Scholar] [CrossRef] [Green Version]
  27. Li, Z.; Cai, B.; Abdalla, B.A.; Zhu, X.; Zheng, M.; Han, P.; Nie, Q.; Zhang, X. LncIRS1 controls muscle atrophy via sponging miR-15 family to activate IGF1-PI3K/AKT pathway. J. Cachexia Sarcopeni. 2019, 10, 391–410. [Google Scholar] [CrossRef] [Green Version]
  28. Gao, H.; Li, X.; Zhan, G.; Zhu, Y.; Yu, J.; Wang, J.; Li, L.; Wu, W.; Liu, N.; Guo, X. Long non-coding RNA MAGI1-IT1 promoted invasion and metastasis of epithelial ovarian cancer via the miR-200a/ZEB axis. Cell Cycle 2019, 18, 1393–1406. [Google Scholar] [CrossRef]
  29. Wei, C.; Guo, D.; Li, Y.; Zhang, K.; Liang, G.; Li, Y.; Ma, Y.; Liu, J.; Li, Y. Profiling analysis of 17beta-estradiol-regulated lncRNAs in mouse thymic epithelial cells. Physiol. Genomics 2018, 50, 553–562. [Google Scholar] [CrossRef]
  30. Hu, G.; Tang, Q.; Sharma, S.; Yu, F.; Escobar, T.M.; Muljo, S.A.; Zhu, J.; Zhao, K. Expression and regulation of intergenic long non-coding RNAs during T cell development and differentiation. Nat. Immunol. 2013, 14, 1190–1198. [Google Scholar] [CrossRef] [Green Version]
  31. Chen, Y.; Kuroki, Y.; Shaw, G.; Pask, A.; Yu, H.; Toyoda, A.; Fujiyama, A.; Renfree, M. Androgen and oestrogen affect the expression of long non-coding RNAs during phallus development in a marsupial. Non-Coding RNA 2019, 5, 3. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  32. Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. Embnet. J. 2011, 17, 10–12. [Google Scholar] [CrossRef]
  33. Langmead, B.; Salzberg, S.L. Fast gapped-read alignment with Bowtie 2. Nat. Methods. 2012, 9, 357–359. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  34. Kim, D.; Pertea, G.; Trapnell, C.; Pimentel, H.; Kelley, R.; Salzberg, S.L. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013, 14, R36. [Google Scholar] [CrossRef] [Green Version]
  35. Pertea, M.; Pertea, G.M.; Antonescu, C.M.; Chang, T.C.; Mendell, J.T.; Salzberg, S.L. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol. 2015, 33, 290–295. [Google Scholar] [CrossRef] [Green Version]
  36. Frazee, A.C.; Pertea, G.; Jaffe, A.E.; Langmead, B.; Salzberg, S.L.; Leek, J.T. Ballgown bridges the gap between transcriptome assembly and expression analysis. Nat. Biotechnol. 2015, 33, 243–246. [Google Scholar] [CrossRef] [Green Version]
  37. Kong, L.; Zhang, Y.; Ye, Z.Q.; Liu, X.Q.; Zhao, S.Q.; Wei, L.; Gao, G. CPC: Assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007, 35, 345–349. [Google Scholar] [CrossRef]
  38. Sun, L.; Luo, H.; Bu, D.; Zhao, G.; Yu, K.; Zhang, C.; Liu, Y.; Chen, R.; Zhao, Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013, 41, e166. [Google Scholar] [CrossRef]
  39. Punta, M.; Coggill, P.C.; Eberhardt, R.Y.; Mistry, J.; Tate, J.; Boursnell, C.; Pang, N.; Forslund, K.; Ceric, G.; Clements, J.; et al. The Pfam protein families database. Nucleic Acids Res. 2012, 40, 290–301. [Google Scholar] [CrossRef]
  40. Jia, H.L.; Zeng, X.Q.; Huang, F.; Liu, Y.M.; Gong, B.S.; Zhang, K.Z.; Zeng, J.H.; Guo, D.G.; Wang, Z.Y.; Li, Y.G. Integrated microRNA and mRNA sequencing analysis of age-related changes to mouse thymic epithelial cells. Iubmb Life 2018, 70, 678–690. [Google Scholar] [CrossRef] [Green Version]
  41. Li, B.; Li, W.; Tian, Y.; Guo, S.; Qian, L.; Xu, D.; Cao, N. Selenium-alleviated hepatocyte necrosis and DNA damage in cyclophosphamide-treated geese by mitigating oxidative stress. Biol. Trace Elem. Res. 2019. [Google Scholar] [CrossRef] [PubMed]
  42. Wingender, E.; Dietze, P.; Karas, H.; Knuppel, R. TRANSFAC: A database on transcription factors and their DNA binding sites. Nucleic Acids Res. 1996, 24, 238–241. [Google Scholar] [CrossRef] [Green Version]
  43. Cepeda, S.; Griffith, A.V. Thymic stromal cells: Roles in atrophy and age-associated dysfunction of the thymus. Exp. Gerontol. 2018, 105, 113–117. [Google Scholar] [CrossRef] [PubMed]
  44. Muller, L.; Di Benedetto, S.; Pawelec, G. The immune system and its dysregulation with aging. Subcell Biochem 2019, 91, 21–43. [Google Scholar] [PubMed]
  45. Warsi, J.; Hosseinzadeh, Z.; Dong, L.; Pakladok, T.; Umbach, A.T.; Bhavsar, S.K.; Shumilina, E.; Lang, F. Effect of Janus kinase 3 on the peptide transporters PEPT1 and PEPT2. J. Membr. Biol. 2013, 246, 885–892. [Google Scholar] [CrossRef] [PubMed]
  46. Page, M.J.; Di Cera, E. Role of Na+ and K+ in enzyme function. Physiol. Rev. 2006, 86, 1049–1092. [Google Scholar] [CrossRef] [Green Version]
  47. Pisat, N.P.; Pandey, A.; Macdiarmid, C.W. MNR2 regulates intracellular magnesium storage in Saccharomyces cerevisiae. Genetics 2009, 183, 873–884. [Google Scholar] [CrossRef] [Green Version]
  48. Hun, M.; Barsanti, M.; Wong, K.; Ramshaw, J.; Werkmeister, J.; Chidgey, A.P. Native thymic extracellular matrix improves in vivo thymic organoid T cell output, and drives in vitro thymic epithelial cell differentiation. Biomaterials 2017, 118, 1–15. [Google Scholar] [CrossRef]
  49. Gubbels, B.M.; Jorgensen, T.N. Androgen-induced immunosuppression. Front. Immunol. 2018, 9, 794. [Google Scholar] [CrossRef]
  50. Liu, R.; Zeng, Y.; Lei, Z.; Wang, L.; Yang, H.; Liu, Z.; Zhao, J.; Zhang, H.T. JAK/STAT3 signaling is required for TGF-β-induced epithelial-mesenchymal transition in lung cancer cells. Int. J. oncol. 2014, 44, 1643–1651. [Google Scholar] [CrossRef] [Green Version]
  51. Kim, K.J.; Abrams, J.; Alphonso, M.; Pearce, M.; Thorbecke, G.J.; Palladino, M.A. Role of endogenously produced interleukin-6 as a second signal in murine thymocyte proliferation induced by multiple cytokines: regulatory effects of transforming growth factor-β. Cell Immunol. 1990, 131, 261–271. [Google Scholar] [CrossRef]
  52. De la Fuente, M.; Miquel, J. An update of the oxidation-inflammation theory of aging: the involvement of the immune system in oxi-inflamm-aging. Curr. Pharm. Des. 2009, 15, 3003–3026. [Google Scholar] [CrossRef] [PubMed]
  53. Yao, J.; Shi, Z.; Ma, X.; Xu, D.; Ming, G. lncRNA GAS5/miR-223/NAMPT axis modulates the cell proliferation and senescence of endothelial progenitor cells through PI3K/AKT signaling. J. Cell Biochem. 2019, 120, 14518–14530. [Google Scholar] [CrossRef] [PubMed]
  54. Puvvula, P.K. LncRNAs regulatory networks in cellular senescence. Int. J. Mol. Sci. 2019, 20, 2615. [Google Scholar] [CrossRef] [Green Version]
  55. Jha, A.; Singh, A.K.; Weissgerber, P.; Freichel, M.; Flockerzi, V.; Flavell, R.A.; Jha, M.K. Essential roles for Cavbeta2 and Cav1 channels in thymocyte development and T cell homeostasis. Sci. Signal. 2015, 8, ra103. [Google Scholar] [CrossRef]
  56. Iniguez, M.A.; Punzon, C.; Fresno, M. Induction of cyclooxygenase-2 on activated T lymphocytes: regulation of T cell activation by cyclooxygenase-2 inhibitors. J. Immunol. 1999, 163, 111–119. [Google Scholar]
  57. Parkinson, D.B.; Bhaskaran, A.; Droggiti, A.; Dickinson, S.; D’Antonio, M.; Mirsky, R.; Jessen, K.R. Krox-20 inhibits Jun-NH2-terminal kinase/c-Jun to control Schwann cell proliferation and death. J. Cell Biol. 2004, 164, 385–394. [Google Scholar] [CrossRef] [Green Version]
  58. Wiekmeijer, A.S.; Pike-Overzet, K.; Brugman, M.H.; van Eggermond, M.; Cordes, M.; de Haas, E.; Li, Y.; Oole, E.; van IJcken, W.; Egeler, R.M.; et al. Overexpression of LMO2 causes aberrant human T-Cell development in vivo by three potentially distinct cellular mechanisms. Exp. Hematol. 2016, 44, 838–849. [Google Scholar] [CrossRef]
  59. Wu, W.X.; Zuo, L.; Dine, K.E.; Shindler, K.S. Decreased signal transducers and activators of transcription (STAT) protein expression in lymphatic organs during EAE development in mice. Immunol. Innov. 2013, 1. [Google Scholar] [CrossRef]
  60. Papoudou-Bai, A.; Hatzimichael, E.; Barbouti, A.; Kanavaros, P. Expression patterns of the activator protein-1 (AP-1) family members in lymphoid neoplasms. Clin. Exp. Med. 2017, 17, 291–304. [Google Scholar] [CrossRef]
  61. Saito, Y.; Gotoh, M.; Ujiie, Y.; Izutsu, Y.; Maeno, M. Involvement of AP-2rep in morphogenesis of the axial mesoderm in Xenopus embryo. Cell Tissue Res. 2009, 335, 357–369. [Google Scholar] [CrossRef] [PubMed]
  62. Rudolph, B.; Hueber, A.O.; Evan, G.I. Reversible activation of c-Myc in thymocytes enhances positive selection and induces proliferation and apoptosis in vitro. Oncogene 2000, 19, 1891–1900. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  63. Phillips, A.C.; Vousden, K.H. E2F-1 induced apoptosis. Apoptosis 2001, 6, 173–182. [Google Scholar] [CrossRef] [PubMed]
  64. Lawrence, T. The nuclear factor NF-kappaB pathway in inflammation. Cold Spring Harb. Perspect. Biol. 2009, 1, a001651. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  65. Massague, J.; Seoane, J.; Wotton, D. Smad transcription factors. Genes Dev. 2005, 19, 2783–2810. [Google Scholar] [CrossRef] [Green Version]
  66. Cheng, S.; Qian, F.; Huang, Q.; Wei, L.; Fu, Y.; Du, Y. HOXA4, down-regulated in lung cancer, inhibits the growth, motility and invasion of lung cancer cells. Cell Death Dis. 2018, 9, 465. [Google Scholar] [CrossRef] [Green Version]
  67. Hochedlinger, K.; Yamada, Y.; Beard, C.; Jaenisch, R. Ectopic expression of Oct-4 blocks progenitor-cell differentiation and causes dysplasia in epithelial tissues. Cell 2005, 121, 465–477. [Google Scholar] [CrossRef] [Green Version]
  68. Klein, U.; Lia, M.; Crespo, M.; Siegel, R.; Shen, Q.; Mo, T.; Ambesi-Impiombato, A.; Califano, A.; Migliazza, A.; Bhagat, G.; et al. The DLEU2/miR-15a/16-1 cluster controls B cell proliferation and its deletion leads to chronic lymphocytic leukemia. Cancer Cell 2010, 17, 28–40. [Google Scholar] [CrossRef] [Green Version]
Genes 11 00147 g001 550
Figure 1. RNA library preparation and analysis process. CPC, Coding Potential Calculator; CNCI, Coding-Non-Coding Index.
Figure 1. RNA library preparation and analysis process. CPC, Coding Potential Calculator; CNCI, Coding-Non-Coding Index.
Genes 11 00147 g001
Genes 11 00147 g002 550
Figure 2. Identification of miRNAs, lncRNAs, and mRNAs in mice thymus. (a–c) Venny maps showing specific miRNAs, lncRNAs, and mRNAs, respectively, in the thymus after ovariectomy or orchiectomy. (d) Circos plot depicting the distribution of lncRNAs in different chromosomes. (e) Columnar graph showing the distribution of lncRNAs and mRNAs in different chromosomes. (f) The size and distribution of lncRNAs and mRNAs. (g) The distribution of miRNAs.
Figure 2. Identification of miRNAs, lncRNAs, and mRNAs in mice thymus. (a–c) Venny maps showing specific miRNAs, lncRNAs, and mRNAs, respectively, in the thymus after ovariectomy or orchiectomy. (d) Circos plot depicting the distribution of lncRNAs in different chromosomes. (e) Columnar graph showing the distribution of lncRNAs and mRNAs in different chromosomes. (f) The size and distribution of lncRNAs and mRNAs. (g) The distribution of miRNAs.
Genes 11 00147 g002
Genes 11 00147 g003 550
Figure 3. Comparative analysis of estrogen- and androgen-related RNAs in mice thymus. (a–c) Volcano plots of differentially expressed lncRNAs (DELs), differentially expressed miRNAs (DEMs), and differentially expressed genes (DEGs), respectively, in ovariectomy and female control groups. (d–f) Volcano plots of DELs, DEMs, and DEGs, respectively, in orchiectomy and male control groups. (g–i) Heat maps of specific DELs, DEMs, and DEGs, respectively, in ovariectomy and female control groups. (j–l) Heat maps of specific DELs, DEMs, and DEGs, respectively, in orchiectomy and male control groups. Each row represents one RNA. Green and red represent relatively lower and higher gene expression in castration group, respectively (p < 0.05 and |log2FC| > 1).
Figure 3. Comparative analysis of estrogen- and androgen-related RNAs in mice thymus. (a–c) Volcano plots of differentially expressed lncRNAs (DELs), differentially expressed miRNAs (DEMs), and differentially expressed genes (DEGs), respectively, in ovariectomy and female control groups. (d–f) Volcano plots of DELs, DEMs, and DEGs, respectively, in orchiectomy and male control groups. (g–i) Heat maps of specific DELs, DEMs, and DEGs, respectively, in ovariectomy and female control groups. (j–l) Heat maps of specific DELs, DEMs, and DEGs, respectively, in orchiectomy and male control groups. Each row represents one RNA. Green and red represent relatively lower and higher gene expression in castration group, respectively (p < 0.05 and |log2FC| > 1).
Genes 11 00147 g003
Genes 11 00147 g004 550
Figure 4. Validation of DELs, DEMs, and DEGs by qRT-PCR. Log2 fold change (FC) were expressed as mean ± SD. n = 3. The statistical significance of all genes reached P < 0.05.
Figure 4. Validation of DELs, DEMs, and DEGs by qRT-PCR. Log2 fold change (FC) were expressed as mean ± SD. n = 3. The statistical significance of all genes reached P < 0.05.
Genes 11 00147 g004
Genes 11 00147 g005 550
Figure 5. Functional analysis of significant DEGs in thymus of mice after ovariectomy or orchiectomy. Functional enrichment of DEGs after (a) ovariectomy and (b) orchiectomy. The significance cutoff for P value was set at 0.05. The font size of enriched terms is proportional to –log10 of P value.
Figure 5. Functional analysis of significant DEGs in thymus of mice after ovariectomy or orchiectomy. Functional enrichment of DEGs after (a) ovariectomy and (b) orchiectomy. The significance cutoff for P value was set at 0.05. The font size of enriched terms is proportional to –log10 of P value.
Genes 11 00147 g005
Genes 11 00147 g006 550
Figure 6. lncRNA–miRNA–mRNA network. (a) Global view of network after ovariectomy. (b) global view of network after orchiectomy. (c) Degree distribution of the female ceRNA network. (d) Degree distribution of the male ceRNA network. Rectangle, circle, and hexagon represents lncRNA, miRNA, and mRNA, respectively. Red represents upregulation of RNAs expression after castration, green represents downregulation.
Figure 6. lncRNA–miRNA–mRNA network. (a) Global view of network after ovariectomy. (b) global view of network after orchiectomy. (c) Degree distribution of the female ceRNA network. (d) Degree distribution of the male ceRNA network. Rectangle, circle, and hexagon represents lncRNA, miRNA, and mRNA, respectively. Red represents upregulation of RNAs expression after castration, green represents downregulation.
Genes 11 00147 g006
Genes 11 00147 g007 550
Figure 7. The expression levels of DEGs at thymic epithelial cells (TECS) and thymocytes. * represents p < 0.05, *** represents p < 0.001.
Figure 7. The expression levels of DEGs at thymic epithelial cells (TECS) and thymocytes. * represents p < 0.05, *** represents p < 0.001.
Genes 11 00147 g007
Genes 11 00147 g008 550
Figure 8. Analysis of transcription factor binding sites of DEGs. (a) Sequence logos of transcription factor binding sites in downregulated DEGs after ovariectomy. (b) Upregulated DEGs after ovariectomy. (c) Downregulated DEGs after orchiectomy. (d) Upregulated DEGs after orchiectomy.
Figure 8. Analysis of transcription factor binding sites of DEGs. (a) Sequence logos of transcription factor binding sites in downregulated DEGs after ovariectomy. (b) Upregulated DEGs after ovariectomy. (c) Downregulated DEGs after orchiectomy. (d) Upregulated DEGs after orchiectomy.
Genes 11 00147 g008
Genes 11 00147 g009 550
Figure 9. Transcription factor-miRNAs-joint target genes networks. (a) Transcription factor-miRNAs-joint target genes network after ovariectomy. (b–f) Transcription factor-miRNAs-joint target genes network after orchiectomy. The purple triangles, brown circles, and green diamonds represents transcription factors, miRNAs, and joint target genes, respectively.
Figure 9. Transcription factor-miRNAs-joint target genes networks. (a) Transcription factor-miRNAs-joint target genes network after ovariectomy. (b–f) Transcription factor-miRNAs-joint target genes network after orchiectomy. The purple triangles, brown circles, and green diamonds represents transcription factors, miRNAs, and joint target genes, respectively.
Genes 11 00147 g009
Table
Table 1. Reading long non-coding RNAs (lncRNAs) mapped to reference genome.
Table 1. Reading long non-coding RNAs (lncRNAs) mapped to reference genome.
SampleM3M3xF3F3x
Raw reads147,021,844148,097,166146,835,736154,631,520
Clean reads138,324,292137,140,672138,697,620144,375,220
Clean ratio (%)94.0892.6094.4693.37
Mapped reads115,678,789112,610,952114,852,768123,422,077
Mapping ratio (%)91.3591.3391.5491.26
Uniquely mapped reads95,107,62792,239,77594,004,60398,579,424
Unique mapping ratio (%)75.1074.8174.9372.89
Table
Table 2. Reading microRNAs (miRNAs) mapped to reference genome.
Table 2. Reading microRNAs (miRNAs) mapped to reference genome.
GroupType Raw Reads3 ADT and Length FilterEfam
M3Total (%)24,408,417 (100)2,742,851 (11.24)1,487,253 (6.09)
Unique (%)1,482,996 (100)847,814 (57.17)21,384 (1.44)
M3xTotal (%)26,998,643 (100)3,258,039 (12.07)860,395 (3.19)
Unique (%)1,459,597 (100)850,698 (58.28)13,133 (0.90)
F3Total (%)22,244,392 (100)2,661,777 (11.97)703,157 (3.16)
Unique (%)1,283,499 (100)779,232 (60.71)12,022 (0.94)
F3xTotal (%)20,753,214 (100)2,332,868 (11.24)590,578 (2.85)
Unique (%)1,131,578 (100)635,854 (56.19)10,134 (0.90)
Table
Table 3. Assembly results of lncRNAs and messenger RNAs (mRNAs).
Table 3. Assembly results of lncRNAs and messenger RNAs (mRNAs).
ItemMin LengthMax LengthMean LengthMedian LengthN50
lncRNA20185,600197211013320
mRNA201124,938198412233246

Share and Cite

MDPI and ACS Style

Li, B.; Zhang, K.; Ye, Y.; Xing, J.; Wu, Y.; Ma, Y.; Li, Y. Effects of Castration on miRNA, lncRNA, and mRNA Profiles in Mice Thymus. Genes 2020, 11, 147. https://doi.org/10.3390/genes11020147

AMA Style

Li B, Zhang K, Ye Y, Xing J, Wu Y, Ma Y, Li Y. Effects of Castration on miRNA, lncRNA, and mRNA Profiles in Mice Thymus. Genes. 2020; 11(2):147. https://doi.org/10.3390/genes11020147

Chicago/Turabian Style

Li, Bingxin, Kaizhao Zhang, Yaqiong Ye, Jingjing Xing, Yingying Wu, Yongjiang Ma, and Yugu Li. 2020. "Effects of Castration on miRNA, lncRNA, and mRNA Profiles in Mice Thymus" Genes 11, no. 2: 147. https://doi.org/10.3390/genes11020147

APA Style

Li, B., Zhang, K., Ye, Y., Xing, J., Wu, Y., Ma, Y., & Li, Y. (2020). Effects of Castration on miRNA, lncRNA, and mRNA Profiles in Mice Thymus. Genes, 11(2), 147. https://doi.org/10.3390/genes11020147

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop