To investigate the profile of miRNAs induced during sheep fat tail development, tail fat tissues of four key developmental stages of fat-rumped Altay sheep were collected. Tail fat tissues of adult thin-tailed XFW sheep were collected simultaneously to detect the DE miRNAs between adult Altay and XFW sheep.

Six pregnant adult Altay ewes with similar due dates were selected. Three of them were slaughtered on the 120th day of gestation to collect the tail fat tissues of fetuses. The rest of the three ewes were fed until they gave birth, and the tail fat of the 1-month-old lambs were collected. Meanwhile, three tail fat tissues of 1- and 2-year-old Altay ewes as well as 2-year-old XFW ewes were collected, respectively. All sheep were supplied and uniformly fed by the farm of Xinjiang Academy of Agricultural and Reclamation Sciences. Prior to sampling, the sheep were anesthetized by injecting 30 mg/kg body weight of pentobarbital sodium (Ningbo, Zhejiang, China) via the ear vein (20). Pentobarbital sodium is a type of medium-efficiency barbital hypnotics that can inhibit the uplink activation system of the brain stem reticular structure. Approximately 50 g of tail fat tissue from each sheep was rapidly collected and immersed in liquid nitrogen for transportation and then stored at −80°C in the laboratory. Tail fat tissues of three sheep were collected at each stage, with a total of 15 samples collected, including fetuses, 1-month-old lambs, 1- and 2-year-old Altay ewes, and 2-year-old XFW ewes, respectively.

Total RNA was isolated from tissues using TRIzol reagent (Invitrogen, United States) according to the manufacturer’s protocol. The concentration and integrity of the total RNA were assessed using 1% agarose gel electrophoresis and 2,100 Bioanalyzer Instrument (Agilent Technologies, CA, United States) by measuring the absorbance at 260 and 280 nm.

The concentration of three RNA samples from each developmental stage was standardized and then mixed into one sample to prepare for the miRNA sequencing library using TruSeq Small RNA Sample Prep Kits (Illumina, San Diego, CA, United States). Subsequently, five distinct libraries were constructed and designated as fetal fat rump, lamb fat rump, hogget fat rump, adult fat rump, and adult thin tail, respectively. All miRNA libraries underwent single-end sequencing on an Illumina HiSeq 2,500 at BGI (Shenzhen, China).

The raw data were processed using software developed by BGI (Shenzhen, China) to get clean reads following the criteria outlined below: (1) low-quality reads (more than 10% bases could not be determined); (2) reads with 5′ primer contaminants; (3) reads without 3′ primer; (4) reads without the insert tag; (5) reads with poly A; and (6) reads shorter than 18 nt. Subsequently, the reads with a sequence length ranging from 18 to 26 nt were retained, and the reads including common RNA families (rRNA, tRNA, snRNA, and snoRNA) and repeats were eliminated using Ensembl, RFam, and Repbase databases. Thereafter, the length distribution of these filtered valid reads was summarized, including the summary of unique tags and total tags. The small RNAs with unique sequences ranging from 18 to 26 nt were aligned to the miRNA precursor of Ovis aries and Bos taurus in miRbase 21.0 (accessed on 20 February 2021)¹ to obtain the miRNA count as well as the base bias on the first position of the identified miRNAs with certain length and at each position of all identified miRNAs, respectively. The locations of mapped miRNAs were further determined by aligning against the sheep genome Oar_v3.1 (accessed on 20 February 2021).² The unmapped unique sequences were utilized to predict novel miRNAs applying Mireap software (BGI, Shenzhen, China) by analysis of the secondary structure, identification of the Dicer cleavage site, and determination of the minimum free energy.

The miRNA expression levels in two libraries (control and treatment) were normalized to obtain the expression of transcript per million (TPM), calculate fold change and p-value from the normalized expression using DESeq2 (21). The miRNAs satisfied the following conditions: log2 (fold change) ≥ 1 or ≤ −1 and p < 0.05 were identified as significant DE miRNAs.

The reliability of RNA-seq data was validated by randomly selecting six upregulated miRNAs in fetal fat rump and six downregulated ones in hogget fat rump using stem-loop real-time quantitative PCR. The expression levels of each miRNA were calculated using the 2^-∆∆CT method, with the expressions of 5S serving as a reference for relative quantification. The primer information is shown in Supplementary Table S1.

TargetScan³ and miRanda⁴ were utilized to predict the target genes of DE miRNAs. The predicted genes, for which TargetScan Total context score < −0.2, P_CT > 0.5, and miRanda energy < −20, were considered as candidate target genes. Then, the target genes were annotated by Gene ontology (GO)⁵ and enriched using the Kyoto Encyclopedia of Genes and Genomes (KEGG) biological pathway database,⁶ respectively. The p-value was corrected using the Bonferroni method, and p ≤ 0.05 was considered as significantly enriched terms (22).

Four DE miRNAs, namely oar-miR-433-3p, oar-miR-485-3p, oar-miR-409-3p, and oar-miR-495-3p, which were significantly expressed in the tail fat tissues at different developmental stages of Altay sheep, and their candidate target genes, PPARγ, PDCD4, and PIK3R3, which were enriched in pathways related to fat deposition, were selected for further verification of their target relationship using dual-luciferase reporter assays. The prediction results of target genes indicated that there is one binding site of oar-miR-433-3p in the 3’UTR of PPARγ mRNA, one binding site of oar-miR-485-3p in the 3’UTR of PDCD4 mRNA, and three binding sites of oar-miR-485-3p, oar-miR-409-3p, and oar-miR-495-3p in the 3’UTR of PIK3R3 mRNA, respectively.

The wild-type 3’UTR of the PPARγ, PDCD4, and PIK3R3 mRNA was amplified and inserted into the psiCHECH2 vector (Promega, Madison, WI, United States) between the XhoI and NotI restriction enzyme cutting sites and named psiCHECH2-PPARγ-3’UTR-WT, psiCHECH2-PDCD4-3’UTR-WT, psiCHECH2-PIK3R3-3’UTR-WTI, psiCHECH2-PIK3R3-3’UTR-WTII, and psiCHECH2-PIK3R3-3’UTR-WTIII, respectively. Their mutant type 3’UTR was generated using the Site-Directed Mutagenesis Kit (Thermo Fisher Scientific, MA, United States) and named psiCHECH2-PPARγ-3’UTR-Mut, psiCHECH2-PDCD4-3’UTR-Mut, psiCHECH2-PIK3R3-3’UTR-MutI, psiCHECH2-PIK3R3-3’UTR-MutII, and psiCHECH2-PIK3R3-3’UTR-MutIII, respectively. The primers used in plasmid construction are shown in Supplementary Table S1.

The vectors and their corresponding miRNA mimics (GenePharma, Shanghai, China), which were fluorescently labeled with 6-FAM, along with psiCHECK2 pure vectors with negative control (NC), and psiCHECK2 pure vectors with miRNA mimics were co-transfected into 293 T cells using Lipofectamine 3,000 Transfection Kit (Thermo Fisher Scientific, MA, United States). After 6 h, the culture medium was replaced, and the transfection efficiency was assessed by detecting the fluorescence of 6-FAM under an inverted fluorescence microscope (Olympus, Japan). Then, the cells were cultured for an additional 48 h, and the relative luciferase activity in the cell was measured by Dual-Luciferase Reporter Assay System (Promega, Madison, WI, United States) according to the manufacturer’s protocol. Each treatment group underwent three replicates.

The sheep preadipocytes obtained through primary culture in our laboratory were seeded in 6-well plates at a density of 1 × 10⁵ cells/well. When the cells reached a confluence of 70%, the following molecules, the 6-FAM fluorescently labeled mimics and negative controls for oar-miR-433-3p, oar-miR-485-3p, oar-miR-409-3p, and oar-miR-495-3p, as well as inhibitors and negative controls for the same miRNAs (GenePharma, Shanghai, China), were transfected into the cells, respectively, using Lipofectamine 3,000 Transfection Kit (Thermo Fisher Scientific, MA, United States) according to the manufacturer’s protocol. The concentrations of miRNA mimics, inhibitors, and NC were maintained at 50 nM within the recommended range provided by the manufacturer. Six hours later, the culture medium was replaced, and total RNA and protein were extracted to detect the expression levels of candidate target genes PPARγ, PDCD4, and PIK3R3 when the cells were further cultured for 48 and 72 h.

The mRNA and protein levels of candidate target genes were detected by qRT-PCR and Western blot, and the primers are listed in Supplementary Table S1. The qRT-PCR analysis was carried out on a Roche 480 instrument (Roche, Mannheim, Germany) using the SYBR Green PCR Master Mix Kit (QIAGEN, Germany) in accordance with the manufacturer′s instructions. The reaction mixture consisted of 10 μL 2 × Quanti Fast SYBR Green PCR Master Mix, 1 μL cDNA (<100 ng), 0.5 μL forward and reverse primers (10 μM), respectively, and ddH₂O to a total volume of 20 μL. The following conditions were used for amplification: 95°C for 5 min, followed by 45 cycles of 95°C for 10 s, 60°C for 30 s, and 72°C for 7 min. qRT-PCR analysis was conducted in triplicate for each sample, and relative expression for each gene was estimated using the 2^-△△Ct method.

The MTT cell proliferation detection kit (Solarbio, Beijing, China) was used to assess the influence of miRNAs on the proliferative capacity of sheep preadipocytes in accordance with the manufacturer’s protocol. The sheep preadipocytes were seeded in 96-well plates at a density of 1 × 10⁴ cells/well. When the cells reached a confluence of 50%, the mimics, inhibitors, and NC of oar-miR-433-3p, oar-miR-485-3p, oar-miR-409-3p, and oar-miR-495-3p were transfected into the cells according to the methods mentioned above. After the cells were further cultured for 72 h, the culture medium was removed carefully, 90 μL of complete culture medium and 10 μL of MTT solution were added per well, cultured in a cell incubator at 37°C and 5% CO₂ for 4 h, then the culture medium was removed, 110 μL of formazan solution per well was added, and were shaken on a shaker for 10 min at low speed to fully dissolve the crystal. Then, the OD values of each well were measured at a wavelength of 490 nm. To eliminate background noise, the blank wells were zeroed by adding 110 μL of formazan solution. Subsequently, the average OD values for each treatment were calculated based on six replicates.

The sheep preadipocytes were seeded in six-well plates at a density of 1 × 10⁵ cells per well. When the cells reached a confluence of 90% ~ 100%, an adipogenic differentiation induction medium was introduced, consisting of 89% DMEM/F12 (Gibco, NY, United States), 10% FBS (Gibco, NY, United States), and 1% penicillin–streptomycin (Solarbio, Beijing, China), supplemented with 80 μmol/L oleic acid (Sigma, MO, United States), 0.25 mmol/L IBMX (Sigma, MO, United States), 5 μmol/L DEX (Sigma, MO, United States), and 8 μmol/L INS (Sigma, MO, United States). Following a 4-day induction period in this medium, the induction medium was replaced with a maintenance medium comprising 89% DMEM/F12, 10% FBS, 1% penicillin–streptomycin, and 5 μmol/L INS for an additional 2 days of culture. Thereafter, the induction and maintenance medium were changed every 2 days until the 12th day (23, 24). Meanwhile, on the 3rd, 6th, and 9th days of induced differentiation, the miRNA mimics, inhibitors, and NC of oar-miR-433-3p, oar-miR-485-3p, oar-miR-409-3p, and oar-miR-495-3p were transfected into the cells to ensure that their levels remained as elevated as possible. Subsequently, the expression levels of four lipid metabolism-related genes, FABP3, PLIN1, ADIPOQ, and LPL, were detected using Western blot.

To evaluate the influence of cellular miRNA levels on their induced lipogenic differentiation ability, the Oil Red O Staining Kit (Solarbio, Beijing, China) was used to compare the extent of lipid accumulation in cells subjected to different treatments on the 12th day of induced differentiation. The staining of Oil Red O was conducted according to the instructions of the kit. After visual observation and photography of stained cells, 1 mL isopropyl alcohol was added to each well for Oil Red O extraction, and the OD value at 510 nm was determined. Subsequently, DNA was extracted from the cells after extraction of Oil Red O using phenol–chloroform–isoamyl alcohol extraction method; the amount of DNA was calculated by measuring the OD value at 260 nm; and then, the amount of fat deposition in the induced differentiated cells was determined through calculation of Oil Red O OD510/DNA OD260.

The cells were extracted, and 1 mL of pre-cooling RIPA lysis buffer containing 1 mM PMSF was added to obtain the total protein. The total protein concentration was determined using the TaKaRa BCA Protein Assay Kit (TaKaRa, Beijing, China) following the manufacturer’s protocol. The proteins were separated using a 12% SDS-PAGE gel (Solarbio, Beijing, China) at 110 V for 90 min and subsequently transferred onto a PVDF membrane (Millipore, United States) at 400 mA for 50 min. The membrane was sealed at 4°C overnight with blocking buffer (TaKaRa, Beijing, China) and washed three times using TBST (Solarbio, Beijing, China). It was then probed with rabbit polyclonal antibody against PPARγ (1:2,000, Abcam, Cambridge, United Kingdom), PDCD4 (1:2,000, Abcam, Cambridge, United Kingdom), PIK3R3 (1:2000, Abcam, Cambridge, United Kingdom), FABP3 (1:2,000, Proteintech, Wuhan, China), PLIN1 (1:5,000, Proteintech, Wuhan, China), ADIPOQ (1:200, Proteintech, Wuhan, China), LPL (1:1,000, Proteintech, Wuhan, China), and beta-actin (1:2,000, Abcam, Cambridge, United Kingdom) at 4°C overnight. Following this incubation period, the membrane was washed three times for 10 min each. Then, it was incubated for 1 h at room temperature with a goat anti-rabbit IgG secondary antibody (1:5,000, Proteintech, Wuhan, China) labeled with HRP, followed by washing three times, 10 min each time. Finally, the reaction bands were developed using the SuperSignal West Femto Trial Kit (Thermo Fisher Scientific, MA, United States) and detected and analyzed using a chemiluminescence imaging system (Vilber, Paris, France).

Statistical analyses were conducted on a minimum of three independent experiments for each treatment, and all data were presented as the means ± standard deviation. All statistical analyses were performed using SPSS version 20 (IBM, NY, United States). Comparisons between two groups were performed by the t-test, while comparisons among multiple groups carried out using a one-way analysis of variance followed by Tukey’s multiple comparisons test; p < 0.05 and p < 0.01 were considered as statistically significant differences.

To explore the miRNA expression profiles in tail fat tissues of sheep during different development stages, five sRNA libraries, namely, fetal fat rump, lamb fat rump, hogget fat rump, adult fat rump, and adult thin tail respectively, were constructed and sequenced. The raw data were filtered to obtain the total sRNAs and unique sRNAs (Supplementary Table S2). More than 98% of the sRNAs ranged from 18 nt to 26 nt in length (Figure 2A). There were 441,215 unique sRNAs expressed exclusively during the fetal stage, 299,425 during the lamb stage, 235,156 during the hogget stage, and 190,663 during the adult stage. Additionally, there were also 3,766 unique sRNAs expressed simultaneously across these four stages (Figure 2B). Then, the unique sRNAs were classified and annotated, and 5,299, 3,330, 3,656, 3,024, and 3,446 miRNAs were identified from tail fat tissue of fetal, lamb, hogget, adult Altay sheep, and adult XFW sheep (Figure 2C).

The miRNA expression of adult Altay sheep was used as the control, while fetal, lamb, hogget Altay sheep, and adult XFW sheep were considered as treatment groups to normalize their expression. The fold change and p-value were calculated to identify the significant DE miRNAs. A total of 350 DE miRNAs, which log2 (fold change) ≥ 1 or ≤ −1 and p < 0.05, were identified (Figure 3A; Supplementary Table S3). The cluster analysis and heatmap drawing were performed using the Heatmap Plot tool in Hiplot Pro⁷ (Figure 3B). Among these DE miRNAs, 191, 60, 26, and 21 miRNAs were extremely upregulated in tail fat tissues of fetal, lamb, hogget Altay sheep, and adult XFW sheep, respectively, but were significantly downregulated in other development stages. Then, their expression profiles across different developmental stages of Altay sheep were further analyzed. From a pool of 245 DE miRNAs that exhibited significant differential expression among these four development stages, the top 20 DE miRNAs were selected for further investigation (Figures 3C,D).

To validate the sequencing data, a total of six upregulated miRNAs in fetal fat rump and six downregulated ones in hogget fat rump were randomly selected, and their expression levels were detected in tail fat tissues of corresponding developmental stages using stem-loop real-time quantitative PCR. By analyzing the qRT-PCR results and miRNA sequencing data, we found that the fold change of selected miRNAs showed similar trends (Figure 4), indicating that the miRNA sequencing data were reliable and efficient.

A total of 4,476 candidate target genes of the top 20 DE miRNAs were predicted (Supplementary Table S4), and their biological functions were analyzed by using GO and KEGG enrichment analysis (Supplementary Table S5). Ninety-nine GO terms were significantly enriched, including regulation of cellular metabolic process (GO:0031323), regulation of primary metabolic process (GO:0080090), phospholipid transporter activity (GO:0005548), and lipid binding (GO:0008289) that are associated with fat metabolism. The top 20 GO terms are shown in Figures 5A–C.

The top 20 KEGG terms are presented in Figure 5D, while the regulatory network of miRNAs, candidate target genes, and associated signaling pathways are depicted in Figure 5E. Several pathways, such as the PPAR signaling pathway (ko03320), AMPK signaling pathway (ko04152), insulin signaling pathway (ko04910), mTOR signaling pathway (ko04150), and PI3K-AKT signaling pathway (ko04151), have already been proved playing critical roles in adipogenesis (25). Meanwhile, the interaction analysis revealed that these miRNAs intricately interconnected these signaling pathways through the regulation of their respective target genes (Figure 5E).

Among these top DE miRNAs and candidate genes, the expression of miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p was significantly different in the tail fat tissues at different developmental stages of Altay sheep. Their candidate target genes, PPARγ, PDCD4, and PIK3R3, have been shown to play key roles in adipocyte proliferation, differentiation, and lipid deposition (25). Therefore, the regulatory relationship between these four miRNAs and their three candidate genes was further validated, and their functions in adipocyte development and lipid deposition were further investigated.

The miRNA mimics were labeled with 6-FAM for easy evaluation of transfection efficiency by fluorescence detection using an inverted fluorescence microscope. Fluorescence detection results showed that the mimics were efficiently transfected into 293 T cells (Figure 6A). The binding sites of miRNAs on the 3’UTR of candidate target genes are shown in Figure 6B. The luciferase assay results revealed that the relative luciferase activity was significantly lower in WT + mimics groups than Mut + mimics, psiCHECH2 + mimics, and psiCHECH2 + NC groups (p < 0.01), indicating that the mimics of miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p could efficiently bind to their target sites and suppress the expression of luciferase (Figure 6C).

Then, sheep preadipocytes were transfected with mimics, mimics negative control, inhibitors, and inhibitors negative control of miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p. The mRNA and protein levels of candidate target genes such as PPARγ, PDCD4, and PIK3R3 were detected using qRT-PCR and Western blot after the cells were further cultured for 48 and 72 h. The results revealed no significant differences in the mRNA levels of PPARγ, PDCD4, and PIK3R3 among the blank control group or any treatment groups involving these miRNAs (p > 0.05; Figure 7A). However, the protein levels of PPARγ and PIK3R3 were significantly downregulated when mimics were transfected into cells comparing to other treatment groups (p < 0.01; Figures 7B,C). The results indicated that miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p could suppress the transcription of mRNA of PPARγ and PIK3R3, but did not induce mRNA degradation because the mRNA levels did not show significant downregulation (p > 0.05). In terms of inhibitor groups, it was observed that although the single endogenous miRNA combined with its inhibitor failed to suppress the expression of its target gene, but there was no significant up-regulation in the protein levels of PPARγ and PIK3R3 (p > 0.05; Figure 7B). This suggests that PPARγ and PIK3R3 might be targeted by other multiple miRNAs.

It is evident that miR-485-3p did not suppress the expression of PDCD4 in sheep preadipocytes (Figures 7B,C), although it has target site on the 3’UTR fragment of PDCD4 (Figure 6B), and could effectively bind to it to suppress the expression of luciferase in dual-luciferase assay (Figure 6C). Thus, the luciferase assay may have a certain rate of false negatives, and the targeted regulatory relationship between miRNA and its candidate target gene should be finally validated by Western blot.

When the expression levels of miRNAs were upregulated by transfecting mimics of miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p, respectively, the preadipocytes exhibited a relatively rapid proliferation until reaching 100% confluent, and the adipogenic differentiation of cells was almost completely inhibited (Figures 8Ab,B). On the contrary, the cells in the blank control and rest treatment groups displayed slower proliferation than the mimics groups and gradually shifted to adipogenic differentiation with numerous tiny lipid droplets appearing within the cells (Figures 8Aa–e). Additionally, the number of cells was also significantly lower than that of mimics groups (p < 0.01; Figure 8B). In contrast, the inhibitors groups showed stronger adipogenic differentiation ability (Figure 8Ad) and slower rate of cell proliferation than the blank control, mimics negative control, inhibitors, and inhibitors negative control groups, but there was no significant difference among themselves (p > 0.05; Figure 8B).

When sheep preadipocytes, transfected with mimics, inhibitors, and NC of miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p, respectively, were induced to undergo adipogenic differentiation, the cells in blank control, mimics negative control, inhibitors, and inhibitors negative control groups shifted to adipogenic differentiation states quickly. The tiny lipid droplets rapidly accumulated and then fused together into larger lipid droplets, but the adipogenic differentiation of the cells transfected with mimics was almost inhibited and kept stasis (Figure 9A). Moreover, these cells exhibited significantly lower levels of lipid deposition than those in the blank control group and other treatment groups (p < 0.01; Figure 9B).

The expression levels of four lipid metabolism-related genes, FABP3, PLIN1, ADIPOQ, and LPL, were detected to further investigate the effect of up- or downregulation of these four miRNAs. Following the transfection of miRNA mimics into sheep preadipocytes, the protein levels of FABP3, PLIN1, ADIPOQ, and LPL were significantly downregulated compared to both the blank control and other treatment groups (p < 0.01; Figures 9C,D). These data indicated that miR-433-3p, miR-485-3p, miR-409-3p, and miR-495-3p could inhibit the expression of their target genes PPARγ and PIK3R3, then regulate the expression of key genes related to lipid metabolism, finally inhibit the adipogenic differentiation of cells, and promote their proliferation in vitro.