MiRNA Profiling and Its Potential Roles in Rapid Growth of Velvet Antler in Gansu Red Deer (Cervus elaphus kansuensis)

A significant variety of cell growth factors are involved in the regulation of antler growth, and the fast proliferation and differentiation of various tissue cells occur during the yearly regeneration of deer antlers. The unique development process of velvet antlers has potential application value in many fields of biomedical research. Among them, the nature of cartilage tissue and the rapid growth and development process make deer antler a model for studying cartilage tissue development or rapid repair of damage. However, the molecular mechanisms underlying the rapid growth of antlers are still not well studied. MicroRNAs are ubiquitous in animals and have a wide range of biological functions. In this study, we used high-throughput sequencing technology to analyze the miRNA expression patterns of antler growth centers at three distinct growth phases, 30, 60, and 90 days following the abscission of the antler base, in order to determine the regulatory function of miRNA on the rapid growth of antlers. Then, we identified the miRNAs that were differentially expressed at various growth stages and annotated the functions of their target genes. The results showed that 4319, 4640, and 4520 miRNAs were found in antler growth centers during the three growth periods. To further identify the essential miRNAs that could regulate fast antler development, five differentially expressed miRNAs (DEMs) were screened, and the functions of their target genes were annotated. The results of KEGG pathway annotation revealed that the target genes of the five DEMs were significantly annotated to the “Wnt signaling pathway”, “PI3K-Akt signaling pathway”, “MAPK signaling pathway”, and “TGF-β signaling pathway”, which were associated with the rapid growth of velvet antlers. Therefore, the five chosen miRNAs, particularly ppy-miR-1, mmu-miR-200b-3p, and novel miR-94, may play crucial roles in rapid antler growth in summer.


Introduction
The appendages of male deer, velvet antlers, are rich in blood vessels and nerves and covered with velvet skin [1]. It is the only appendage organ in mammals which can regenerate every year [2]. Androgen-regulated antlers undergo annual cyclic shedding and regrowth, which is distinct from cellular dedifferentiation in lower organisms but based on the differentiation process in stem cells [3]. Antler's growth is a very complex process. During rapid growth, the cells proliferate rapidly but remain organized and not neoplastic. The antlers eventually undergo complete ossification, resulting in the death of velvet antler tissue and effectively preventing the uncontrolled growth of the antler [4]. However, as an ideal growth system model for cartilage damage repair and cell proliferation and differentiation without carcinogenesis, the mechanism of antler rapid growth remain unclear.
The growing season of antlers varies from being slow in spring to being exponentially accelerated in summer to being slow in autumn with a typical S curve [5]. In summer,

RNA Extraction, Library Preparation, and Sequencing
Total RNA was extracted using the TRIzol reagent (Thermo Fisher Scientific, Waltham, MA, USA) in accordance with the manufacturer's instructions. DNase I (Promega, Madison, WI, USA) was used to clean up extracted RNA and remove contaminants. The concentration and integrity of the RNA samples were assessed using a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, USA), RNA Nano 6000 analysis kit, and Bioanalyzer 2100 system (Agilent Technologies, Santa Clara, CA, USA) to confirm that they complied with sequencing requirements. It was found that the extracted total RNA satisfied the quality standards and could therefore be utilized for the construction of small RNA libraries if 1.8 < OD260/OD280 < 2.0, RIN > 7.0, and rRNA 28S/18S ≥ 0.7. Utilizing the NEBNext ® Ultra™ RNA Library Prep Kit for Illumina ® (New England Biolabs, Ipswich, MA, USA), sequencing libraries were created. A Bioanalyzer 2100 (Agilent Technologies) was used to test the quality of the libraries. The miRNA libraries were sequenced using the HiSeq™ 2500 platform (Illumina, USA).

Data Processing and Analysis
By eliminating spliced sequences, reads with poly(A) tails, low-quality bases, and sequences with less than 18 or more than 30 nucleotides, raw data were initially processed to obtain clean reads. Additionally, the cleaned dataset's Q20 and Q30 scores, GC content (%), and levels of sequence repetitions were calculated. The Trinity software [25] was used to carry out de novo assembly. Secondly, clean reads were aligned with Silva (https://www.arb-silva.de, accessed on 15 September 2022), genomic tRNA (GtRNAdb; http://lowelab.ucsc.edu/GtRNAdb/, accessed on 15 September 2022), Rfam (http:// rfam.xfam.org/, accessed on 16 September 2022), and Repbase (www.girinst.org/Repbase, accessed on 16 September 2022) databases using Bowtie to remove sequences matching rRNA, tRNA, SNRNA, SnRNA, other NCRNA, and repetitive sequences. Thirdly, the secondary structures of novel miRNAs were predicted using the Randfold program (http://www.aquafold.com, accessed on 19 September 2022). Fourth, the known miRNAs were determined via the miRBase v21.0 database (http://www.mirbase.org/, accessed on 20 September 2022) using either 0 or 1 mismatches. The novel miRNAs were predicted using RNAfold software [26] for sequences that could not be aligned to the miRBase database but could be mapped to the red deer genome. Expression profiles of miRNA were computed using the normalized approach [27]. DEseq software [28] was used to identify miRNAs that were differentially expressed across the three libraries based on normalized counts. The differentially expressed miRNAs (DEMs) were detected using |log2(Fold Change, FC)| ≥ 1.00 and False Discovery Rate (FDR) ≤ 0.01 as the screening criteria. Then, data from the three sample groups were compared in pairs (d30 versus d60, d30 versus d90, and d60 versus d90), and the DEMs of each comparison group were screened, respectively.
We used the default settings of the software programs Miranda v33 (http://www. microrna.org, accessed on 25 September 2022) and TargetScan 7.2 (http://www.targetscan. org, accessed on 25 September 2022) to estimate possible target genes of miRNAs. Both approaches' predictions of genes were taken into consideration when choosing possible targets for certain miRNAs. The Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis and GO enrichment (Gene ontology) were used to annotate these candidate genes with default parameters functionally.

Real-Time PCR to Verify Differentially Expressed microRNAs
The DEMs in the antler growth centers of Gansu red deer during the three growth periods were detected via qPCR. The U6 gene was used as the internal reference, and the primers of the DEMs and U6 gene are listed in Table 1. RNA was extracted from antler growth centers using the TaKaRa MiniBEST Universal RNA Extraction Kit (TaKaRa, Code No. 9767), and 1 ug of total RNA from each sample was reverse transcribed to synthesize cDNA using PrimeScript™ RT Master Mix (Perfect Real Time) (TaKaRa, Code No. RR036Q), which used as the template for qPCR. qPCR was performed with TB Green ® Premix Ex Taq™ II (Tli RNaseH Plus) (TaKaRa, Code No. RR820A). qPCR reaction system (20 µL): 2 µL of cDNA, 1 µL of forward primer (10 µM), 1 µL of reverse primer (10 µM), 10 µL of TB Green Premix Ex Taq II (Tli RNaseH Plus) (2×), and 6 µL of ddH 2 O. qPCR reaction conditions: pre-denaturation at 95 • C for 10 min; denaturation at 95 • C for 15 s, and annealing at 60 • C for 30 s for 40 cycles; extension at 72 • C for 5 min; storage at 4 • C. The relative expression results were calculated and analyzed via the 2 −∆∆Ct method.

Sequencing Data Processing
The sequencing results showed that 14,971,388,16,985,192, and 15,732,412 raw reads were obtained from the 30-, 60-, and 90-day growth centers of red deer velvet antlers, respectively. Further processing of the raw data resulted in 9,203,007, 12,661,360, and 10,346,245 clean reads in three groups (d30, d60, and d90), with Q30 percentages higher than 98.00% for each sample. After removing rRNA, tRNA, SNRNA, SnRNA, other NCRNAs, and repeat sequences, 8,083,271, 11,313,839, and 9,089,417 unannotated reads from the three groups were compared to the red deer reference genome; then 5,047,174, 6,834,924, and 5,384,420 mapped reads were obtained in the three groups, which accounted for 52.09%, 45.12%, and 46.46% of their unannotated reads, respectively. Finally, the mapped reads were aligned with the miRBase database for miRNA identification. Based on these results, the clean reads were reliable, and data quality was consistent and good, which meant that it could be used for subsequent analysis.

MiRNA Identification
The mapped reads were aligned with the miRBase database for the identification of conserved miRNAs, and miRDeep2 software was used for the prediction of novel miRNAs. In total, 4103, 4429, and 4307 known miRNAs were identified, and 216, 211, and 213 novel miRNAs were predicted in each of the three groups. The length distribution of the identified known miRNAs and new miRNAs showed that the miRNAs detected in this study were 18-28 bp in length, and most of them were concentrated at 20-24 bp, with 22 bp being predominant ( Figure 1). The expression abundance data were generated by homogenizing miRNA expression across different samples; Tables 2 and 3 show the top 20 known and novel miRNAs with the highest levels of expression in all three groups. Among all the identified miRNAs, 4136 were expressed in d30, d60, and d90, and 41, 161, and 91 miRNAs were specifically expressed in groups d30, d60, and d90, respectively ( Figure 2).

Screening for Differentially Expressed miRNAs
The miRNAs that were differentially expressed in the different comparison groups (d30 vs. d60, d90 vs. d60, and d30 vs. d90) were filtered out. A total of 292, 102, and 60 miRNAs were found to be differentially expressed in the three comparison groups, respectively.. In total, 122 DEMs were significantly up-regulated in d30 compared to the

Screening for Differentially Expressed miRNAs
The miRNAs that were differentially expressed in the different comparison groups (d30 vs. d60, d90 vs. d60, and d30 vs. d90) were filtered out. A total of 292, 102, and 60 miRNAs were found to be differentially expressed in the three comparison groups, respectively.. In total, 122 DEMs were significantly up-regulated in d30 compared to the

Screening for Differentially Expressed miRNAs
The miRNAs that were differentially expressed in the different comparison groups (d30 vs. d60, d90 vs. d60, and d30 vs. d90) were filtered out. A total of 292, 102, and 60 miRNAs were found to be differentially expressed in the three comparison groups, respectively. In total, 122 DEMs were significantly up-regulated in d30 compared to the expression of miRNAs in d60, and 170 DEMs were significantly down-regulated. Compared with the expression of miRNA in d90, 38 DEMs were significantly up-regulated in d30, 64 DEMs were significantly down-regulated in d30; the expression of 40 DEMs was significantly up-regulated in d60, and 20 DEMs were significantly down-regulated in d60 ( Figure 3, Supplementary Table S1).

Target Genes' Prediction of Differentially Expressed miRNAs
The target genes of DEMs in the three comparison groups were predicted using two algorithms, TargetScan and miRanda. A total of 14,050 target genes were identified after the intersection, of which 13,086 genes were predicted for conserved miRNAs, and 9584 genes were predicted for novel miRNAs. Specifically, 13,162, 13,807, and 13,725 target genes were obtained in antler growth centers at days 30, 60, and 90, respectively.
In terms of differentially expressed miRNAs, there were 3419, 3719, and 3665 in the d30 vs. d60, d30 vs. d90, and d60 vs. d90 comparison groups, respectively. Specifically, in the 30d vs. 60d comparison group, up-regulated DEMs predicted 1189 target genes, and down-regulated DEMs predicted 2663 target genes. In the d30 vs. d90 comparison group, up-regulated DEMs predicted 3308 target genes, and down-regulated DEMs predicted 730 target genes. In the d60 vs. d90 comparison group, up-regulated DEMs predicted 3603 target genes, and down-regulated DEMs predicted 115 target genes.

Target Genes' Functional Analysis of Differentially Expressed miRNAs
We examined the roles played by the target genes of the DEMs to better understand the roles of miRNA target genes in antler growth centers during various growth phases. The results of the GO enrichment analysis showed that the target genes of differentially expressed miRNAs were significantly enriched in 30 GO categories (p < 0.05) in all three comparison groups, including "cellular process", "biological regulation", and "singleorganism process" in Biological Process (BP); "cell part", "cell", and "organelle" in Cellular Component (CC); and "binding" and "catalytic activity" in Molecular Function (MF) (Figure 4).

Target Genes' Prediction of Differentially Expressed miRNAs
The target genes of DEMs in the three comparison groups were predicted using two algorithms, TargetScan and miRanda. A total of 14,050 target genes were identified after the intersection, of which 13,086 genes were predicted for conserved miRNAs, and 9584 genes were predicted for novel miRNAs. Specifically, 13,162, 13,807, and 13,725 target genes were obtained in antler growth centers at days 30, 60, and 90, respectively.
In terms of differentially expressed miRNAs, there were 3419, 3719, and 3665 in the d30 vs. d60, d30 vs. d90, and d60 vs. d90 comparison groups, respectively. Specifically, in the 30d vs. 60d comparison group, up-regulated DEMs predicted 1189 target genes, and down-regulated DEMs predicted 2663 target genes. In the d30 vs. d90 comparison group, up-regulated DEMs predicted 3308 target genes, and down-regulated DEMs predicted 730 target genes. In the d60 vs. d90 comparison group, up-regulated DEMs predicted 3603 target genes, and down-regulated DEMs predicted 115 target genes.

Target Genes' Functional Analysis of Differentially Expressed miRNAs
We examined the roles played by the target genes of the DEMs to better understand the roles of miRNA target genes in antler growth centers during various growth phases. The results of the GO enrichment analysis showed that the target genes of differentially expressed miRNAs were significantly enriched in 30 GO categories (p < 0.05) in all three comparison groups, including "cellular process," "biological regulation," and "single-organism process" in Biological Process (BP); "cell part," "cell," and "organelle" in Cellular Component (CC); and "binding" and "catalytic activity" in Molecular Function (MF) (Figure 4).
(a)  The KEGG pathway enrichment analysis revealed that in the d30 vs. d60 comparison group, the functions of the target genes of DEMs were significantly enriched in 232 KEGG The KEGG pathway enrichment analysis revealed that in the d30 vs. d60 comparison group, the functions of the target genes of DEMs were significantly enriched in 232 KEGG pathways, including "Type II diabetes mellitus", "Axon guidance", and "Small cell lung cancer", among others. The target genes of DEMs were enriched in 230 KEGG pathways in the comparison group between d30 and d90, with target genes significantly enriched in "Basal cell carcinoma", "Axon guidance", "Endocytosis", etc. In the d60 vs. d90 comparison group, differential target genes were significantly enriched in "Basal cell carcinoma", "Endocytosis", "Axon guidance", etc. (p < 0.05). Figure 5 shows the top 20 signaling pathways with the highest confidence of enrichment significance (i.e., the smallest Q value) across the three comparison groups. pathways, including "Type II diabetes mellitus," "Axon guidance," and "Small cell lung cancer," among others. The target genes of DEMs were enriched in 230 KEGG pathways in the comparison group between d30 and d90, with target genes significantly enriched in "Basal cell carcinoma," "Axon guidance," "Endocytosis," etc. In the d60 vs. d90 comparison group, differential target genes were significantly enriched in "Basal cell carcinoma," "Endocytosis," "Axon guidance," etc. (p < 0.05). Figure 5 shows the top 20 signaling pathways with the highest confidence of enrichment significance (i.e., the smallest Q value) across the three comparison groups. Additionally, we were intrigued by the molecular processes that underlining the quick development of velvet antlers in the rapid growth stage (about 60 days). As miR-NAs usually negatively regulate the expression of their target genes, in this study, we screened four DEMs that were lowly expressed in antler growth centers at 60 days and highly expressed in antler growth center tissues at 30 and 90 days, including ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94. In addition, there was another miRNA (novel miR-10) with expression characteristics that were opposite to those of the above four miRNAs. Target gene prediction results showed that ppy-miR-1, mmu-miR-200b-3p, novel miR-6, novel miR-94, and novel miR-10 predicted 821, 1154, 121, 458, and 737 target genes, respectively. The results of the KEGG pathway annotation showed that the target genes of ppy-miR-1, mmu-miR-200b-3p, novel miR-6, novel miR-94, and novel miR-10 were significantly annotated to "Ras signaling pathway," "Neurotrophin signaling pathway," "Glioma," "Pathways in cancer," and "Thyroid hormone signaling pathway" (p < 0.05), respectively. Furthermore, the target genes of these miRNAs were also annotated to the "Wnt signaling pathway," "PI3K-Akt signaling pathway", "MAPK signaling pathway," "TGF-β signaling pathway," and "mTOR signaling pathway," and other classical signaling pathways (Supplementary Table S2). It is possible that these miRNAs may play important roles in the rapid growth of velvet antlers, which needs to be further verified.

Validation of miRNAs Expression
The Venn Diagram package of R language was used to analyze the differential miR-NAs, and there were 13 common DEMs in the three comparison groups. After removing the miRNAs with consistent, mature sequences, three DEMs (bta-miR-1298, bta-miR-486, and novel miR-166) were retained for verification via stem-loop RT-qPCR. In addition, to further investigate the function of the four DEMs (ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94), their relative expression in the growth centers of velvet antlers in three growth stages was also detected. Additionally, we were intrigued by the molecular processes that underlining the quick development of velvet antlers in the rapid growth stage (about 60 days). As miRNAs usually negatively regulate the expression of their target genes, in this study, we screened four DEMs that were lowly expressed in antler growth centers at 60 days and highly expressed in antler growth center tissues at 30 and 90 days, including ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94. In addition, there was another miRNA (novel miR-10) with expression characteristics that were opposite to those of the above four miRNAs. Target gene prediction results showed that ppy-miR-1, mmu-miR-200b-3p, novel miR-6, novel miR-94, and novel miR-10 predicted 821, 1154, 121, 458, and 737 target genes, respectively. The results of the KEGG pathway annotation showed that the target genes of ppy-miR-1, mmu-miR-200b-3p, novel miR-6, novel miR-94, and novel miR-10 were significantly annotated to "Ras signaling pathway", "Neurotrophin signaling pathway", "Glioma", "Pathways in cancer", and "Thyroid hormone signaling pathway" (p < 0.05), respectively. Furthermore, the target genes of these miRNAs were also annotated to the "Wnt signaling pathway", "PI3K-Akt signaling pathway", "MAPK signaling pathway", "TGF-β signaling pathway", and "mTOR signaling pathway", and other classical signaling pathways (Supplementary Table S2). It is possible that these miRNAs may play important roles in the rapid growth of velvet antlers, which needs to be further verified.

Validation of miRNAs Expression
The Venn Diagram package of R language was used to analyze the differential miRNAs, and there were 13 common DEMs in the three comparison groups. After removing the miRNAs with consistent, mature sequences, three DEMs (bta-miR-1298, bta-miR-486, and novel miR-166) were retained for verification via stem-loop RT-qPCR. In addition, to further investigate the function of the four DEMs (ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94), their relative expression in the growth centers of velvet antlers in three growth stages was also detected.
The expression levels of bta-miR-1298 and novel miR-166 increased with antler growth, and the expression levels of bta-miR-486 showed a rising then falling trend, while the expression levels of ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94 showed a falling then rising trend, which was in line with the expression trend in the sequencing data, showing that the small RNA sequencing results were accurate (Figure 6). , x FOR PEER REVIEW 13 of 22 The expression levels of bta-miR-1298 and novel miR-166 increased with antler growth, and the expression levels of bta-miR-486 showed a rising then falling trend, while the expression levels of ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94 showed a falling then rising trend, which was in line with the expression trend in the sequencing data, showing that the small RNA sequencing results were accurate ( Figure 6).

Discussion
Antler growth is primarily driven and controlled by the antler apical growth centerwhich is made up of velvet skin, mesenchyme, and cartilage tissue-by the process o endochondral ossification [1,9]. The antler skin and cartilage tissues in the antler growth center are rich in vascular tissue and differ greatly from normal skin and cartilage tissues which is inevitably associated with the astonishing growth rate of the antler and its capac ity for repair and regeneration [29]. Therefore, the investigation of changes in gene expres sion in antler growth centers during different growth periods is highly significant in orde to reveal the mechanism of rapid antler growth. The majority of organ growth is driven by the combinatorial effects of cell proliferation and cell expansion [30].
The purpose of this study was to deep sequence for miRNAs in velvet antler growth centers at different growth stages. In total, 292, 102, and 60 miRNAs were found to be differentially expressed in the d30 vs. d60, d30 vs. d90, and d60 vs. d90 comparison groups, respectively. Furthermore, 10 miRNAs were differentially expressed in the d30 vs. d60 and d60 vs. d90 comparison groups but not in the d30 vs. d90 group, among which four differentially expressed miRNAs (ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94) with high expression levels in the antler growth center of 30 and 90 day and low expression in the antler growth center of 60 days. Meanwhile, the reverse wa true for another miRNA (novel miR-10). Given that miRNAs typically function by down regulating the expression of their target genes [31], five differentiated miRNAs were screened, which could play important roles in the rapid growth of velvet antlers. Furthe validation is required to determine whether miRNAs that are highly or lowly expressed

Discussion
Antler growth is primarily driven and controlled by the antler apical growth center-which is made up of velvet skin, mesenchyme, and cartilage tissue-by the process of endochondral ossification [1,9]. The antler skin and cartilage tissues in the antler growth center are rich in vascular tissue and differ greatly from normal skin and cartilage tissues, which is inevitably associated with the astonishing growth rate of the antler and its capacity for repair and regeneration [29]. Therefore, the investigation of changes in gene expression in antler growth centers during different growth periods is highly significant in order to reveal the mechanism of rapid antler growth. The majority of organ growth is driven by the combinatorial effects of cell proliferation and cell expansion [30].
The purpose of this study was to deep sequence for miRNAs in velvet antler growth centers at different growth stages. In total, 292, 102, and 60 miRNAs were found to be differentially expressed in the d30 vs. d60, d30 vs. d90, and d60 vs. d90 comparison groups, respectively. Furthermore, 10 miRNAs were differentially expressed in the d30 vs. d60 and d60 vs. d90 comparison groups but not in the d30 vs. d90 group, among which four differentially expressed miRNAs (ppy-miR-1, mmu-miR-200b-3p, novel miR-6, and novel miR-94) with high expression levels in the antler growth center of 30 and 90 days and low expression in the antler growth center of 60 days. Meanwhile, the reverse was true for another miRNA (novel miR-10). Given that miRNAs typically function by down-regulating the expression of their target genes [31], five differentiated miRNAs were screened, which could play important roles in the rapid growth of velvet antlers. Further validation is required to determine whether miRNAs that are highly or lowly expressed in antler growth centers at 60 days of growth promote antler cell proliferation by suppressing the expression of their target genes.
It has been reported that miR-1 has a broad range of biological functions, e.g., in modulating skeletal muscle proliferation and differentiation [32,33], playing roles in heart disease [34][35][36][37], sex determination [38], and diabetes [39]. MiR-1 also plays critical roles in chondrogenesis [40]; the proliferation, differentiation, and hypertrophy of chondrocytes [41,42]; and the development of osteoarthritis [43,44]. The growth center of the antler tip drives the longitudinal endochondral ossification process of velvet antlers [45,46]. It is probable that miR-1 has a significant impact on this procedure. In one study, it was discovered that miR-1 was capable of suppressing the proliferation of chondrocytes in Chinese sika deer by specifically targeting the 3 UTR of IGF-1 [47]. In the present study, ppy-miR-1 was lowly expressed in the growth center of antlers at 60 days of growth and highly expressed at 30 and 90 days of growth, which may be closely related to the rapid growth and ossification of antlers.
MiR-200b is a miRNA with a wide range of functions, including functions in multiple cancers or tumors [48][49][50][51], endotheliogenesis [52][53][54], angiogenesis [55], etc. It has been reported that a decreased level of miR-200b-3p may result in angiogenesis in HCC [56]. The growing tips of velvet antlers are rich in blood vessels [4,19,[57][58][59]. The antler grows quickly, including its interior components, such as the blood vessels, nerves, and skin that join it [59]. Mmu-miR-200b-3p may promote angiogenesis in antlers during the rapid growth phase, which would provide nutritional and other types of support for the rapid growth of antlers.
The novel miR-94, the target genes of which are strongly annotated to the "Pathways in cancer" and "Wnt signaling pathway", is the miRNA in which we are most interested. It was shown that the proliferation of velvet antler cells was faster than that of cancer cells during the rapid growth stage of antlers [23,60,61]. Academics refer to this phenomenon as cancer-like growth [62,63]. "Pathways in cancer" is a complicated signaling pathway, forming a complex network with the "ERK signaling pathway", "PI3K-AKT signaling pathway", "Wnt signaling pathway", "Notch signaling pathway", "TGF-β signaling pathway", "HIF-1 signaling pathway", and other signaling pathways. Studies showed that the "ERK signaling pathway" [64], "PI3K-AKT signaling pathway" [65][66][67], "Notch signaling pathway" [68,69], "TGF-β signaling pathway" [70,71], "Wnt signaling pathway" [72][73][74], and "HIF-1 signaling pathway" [66] were involved in the development and regeneration of velvet antlers. Thus, more research is necessary to determine if novel miR-94 and its target genes have a major impact on the quick development and even the regeneration of velvet antlers.
The target genes of novel miR-6 annotated fewer signaling pathways than ppy-miR-1, mmu-miR-200b-3p, and novel miR-94, and the target genes are significantly enriched in "mTOR signaling pathway", "Pathways in cancer", "Rap1 signaling pathway" and "signaling pathways regulating pluripotency of stem cells". mTOR is a vital regulator of multiple cellular functions, including cell proliferation, differentiation, and protein synthesis [75]. The mTOR signaling pathway can stimulate the protein synthesis of chondrocytes to regulate mammalian limb skeletal growth [76]. Studies have indicated that mTOR is an important regulator of cartilage development [77][78][79][80]. The Rap1 signaling pathway is a cancer-related signaling pathway with a wide range of biological functions [81][82][83][84]. The longitudinal section of the antler was divided into velvet skin, mesenchyme, pre-cartilage, and cartilage from top to bottom, which clearly reflected the dynamic process of antler ossification, indicating that the growth and development of antler is essentially the process of bone tissue growth and development [85]. All factors that theoretically play a relevant role in bone tissue development may be involved in the regulation of antler growth and development. Studies on the regulation of antler growth by the mTOR signaling pathway and Rap1 signaling pathway have not been reported yet. Whether novel miR-6 regulates the proliferation of antler cells through these two pathways remains to be further verified. The antler mesenchyme is filled with mesenchymal cells (reserve mesenchyme cells, RMCs) that have the ability to divide vigorously and differentiate towards chondrocytes, which ultimately secrete the bone matrix to complete the ossification [7]. RMCs are stem cells that drive the rapid growth of velvet antlers and are stromal cells capable of self-replication and have the potential to differentiate into chondrocytes, etc. [86]. Target genes of novel miR-6 annotated to the "Signaling pathways regulating pluripotency of stem cells", this miRNA and its target genes may be involved in the differentiation of RMCs into chondrocytes. Taken together, the screened DEMs have a certain significance for the study of the rapid growth of velvet antlers, and the regulatory mechanism remains to be further investigated.

Conclusions
In summary, we detected miRNAs that were differentially expressed in antler growth centers at different stages of growth using an RNA-seq approach. Following the comparison and analysis, we screened five differentially expressed miRNAs with target genes that were significantly enriched in the "Wnt signaling pathway", "PI3K-Akt signaling pathway", "MAPK signaling pathway", "TGF-β signaling pathway", "mTOR signaling pathway", and other classical signaling pathways. These miRNAs, especially ppy-miR-1, mmu-miR-200b-3p, and novel miR-94, could be extremely important for the rapid development of antlers.