Chi-miR-3031 regulates beta-casein via the PI3K/AKT-mTOR signaling pathway in goat mammary epithelial cells (GMECs)

Background MicroRNAs can regulate gene expression at the posttranscriptional level through translational repression or target degradation. Our previous investigations examined the differential expression levels of chi-miR-3031 in caprine mammary gland tissues in colostrum and common milk stages. Results The present study detected the role of chi-miR-3031 in the lactation mechanisms of GMECs. High-throughput sequencing was used to analyze transcriptomic landscapes of GMECs transfected with chi-miR-3031 mimics (MC) and a mimic negative control (NC). In the MC and NC groups, we acquired 39,793,503 and 36,531,517 uniquely mapped reads, respectively, accounting for 85.85 and 81.66% of total reads. In the MC group, 180 differentially expressed unigenes were downregulated, whereas 157 unigenes were upregulated. KEGG pathway analyses showed that the prolactin, TNF and ErbB signaling pathways, including TGFα, PIK3R3, IGF2, ELF5, IGFBP5 and LHβ genes, played important roles in mammary development and milk secretion. Results from transcriptome sequencing, real-time PCR and western blotting showed that chi-miR-3031 suppressed the expression of IGFBP5 mRNA and protein. The expression levels of β-casein significantly increased in the MC and siRNA-IGFBP5 groups. We observed that the down-regulation of IGFBP5 activated mTOR at the Ser2448 site in GMECs transfected with MC and siRNA-IGFBP5. Previous findings and our results showed that chi-miR-3031 activated the PI3K-AKT-mTOR pathway and increased β-casein expression by down-regulating IGFBP5. Conclusions These findings will afford valuable information for improving milk quality and contribute the development of potential methods for amending lactation performance. Electronic supplementary material The online version of this article (10.1186/s12917-018-1695-6) contains supplementary material, which is available to authorized users.


Background
Previous studies have shown that mammary epithelial cells (MECs) play key roles in milk protein and lipid biosynthesis [1][2][3]. During lactation, MECs secrete large quantities of milk proteins, notably caseins. Caseins (αS1-, αS2-, β-, and κ-CN in most species) are the main milk-specific proteins [4]. Transcriptomic regulations affect the changes of milk composition and yield in MECs. Thus, in this study, MECs were harvested from the mammary gland.
RNA sequencing (RNA-seq) provides a holistic view of transcriptomes, allows the refinement of gene structures and also uncovers novel splice isoforms, SNPs and transcribed regions [5][6][7], in addition, it also produces complete gene expression measurements instead of relative measurements, thereby providing more significant insights and accuracy than microarrays [8,9]. Using a transcriptome analysis, Suarez-Vega et al. (2017) found that RNA transport, and fatty acid elongation Jak-STAT signaling pathways could affect milk yield, protein percentage and fat percentage [10]. Via an RNA-seq analysis of mammary epithelial cells, Wanyonyi et al. (2017) suggested that seven pathways, including lysosome, galactose metabolism, p53 signaling, cell adhesion molecule (CAM), complement, coagulation and Nod-like receptor signaling pathways, are only significantly responsive to extracellular matrices in the presence of lactogenic hormones [11].
MiRNAs are noncoding small RNA molecules (about 22 nucleotides) that regulate gene expression post-transcriptionally by binding to the 3′-UTR of their target mRNAs [12,13], which refer to proliferation [14], apoptosis [15], reproduction [16] and milk performance [17,18]. The overexpression of miR-130a notably decreases the levels of cellular triacylglycerol (TAG) and suppresses lipid droplet formation, but the inhibition of miR-130a results in increased TAG accumulation and lipid droplet formation in bovine mammary epithelial cells [19]. The miRNA target network shows that bta-miR-574 can affect lactation and mammary gland development by leptin receptor, which is involved in the adipocytokine signaling pathway. MiR-3880 regulates the development of mammary gland and lactation by serine/ threonine-protein phosphatase, which participates in the oxytocin signaling pathway [20]. The miR-15b expression decreases during lactation and oppositely correlates with lipid synthesis proteins, suggesting that miR-15b may refer to milk production and lipid synthesis [21]. The overexpression of miR-145 notably increases TAG accumulation, fat droplet formation, and proportions of unsaturated fatty acids. Oppositely, silencing miR-145 impairs fatty acid synthesis [22].
The insulin-like growth factor (IGF) system mainly includes IGF-I and IGF-II ligands, IGF-IR and IGF-IIR receptors and a family of IGF-binding proteins (IGFBPs) in vertebrates [23]. This system performs indispensable roles in adjusting growth, development and reproduction [24]. The functions of IGFs are regulated by their receptors and IGF-IGFBPs that serve as major regulators of IGF activity through acting as carriers [25]. These molecules also dominate the bioactivity and bioavailability of IGFs that interact with IGF receptors to modulate downstream signal transduction networks [26]. IGFBP5 is more conserved than other IGFBPs in mammals [27,28]. A further study showed that IGFs act as important upstream activators of PI3K/AKT-mTOR signaling pathways to regulate milk protein synthesis [29][30][31]. The effect of milking frequency on milk yield is locally controlled within mammary glands and can be a function of changes in either the number or activity of secretory mammary epithelial cells. IGF-I signaling can mediate these effects because it can be locally controlled within mammary glands [32].
Our previous study examined the differential expression levels of chi-miR-3031 in caprine mammary gland tissues in colostrum (2 days postpartum) and common (90 days postpartum) milk stages [20]. Bioinformatics analysis suggested chi-miR-3031 and IGFBP5 are key signaling factors that regulate protein synthesis and cell proliferation in GMECs. IGFBP5 can regulate IGFs downstream of signal transduction networks by interacting with IGF receptors [26][27][28][29]. Until now, there have been no reports on the effects of chi-miR-3031 and IGFBP5 on milk protein synthesis. Therefore, in the present study, we used chi-miR-3031 mimics and mimic negative controls to transfect GMECs. We performed RNA-Seq to identify DEGs and the relative signal pathways affecting protein synthesis in GMECs. Among the identified DEGs, due to their involvement in lactation mechanisms, we selected IGFBP5 based on bioinformatics for further studies. This study aimed to determine the regulatory mechanisms of chi-miR-3031 on casein synthesis in GMECs.

Summary of transcriptome sequencing
RNA-Seq was used to analyze the transcriptomic landscapes of GMECs transfected with MC and NC. Total RNAs from GMECs were extracted to construct RNA libraries for Illumina sequencing. Reads with adaptors and low-quality reads were removed prior to assembly. In total, we acquired 46,351,307 and 44,748,620 clean reads from each library (Table 1). In the MC and NC groups, uniquely mapped reads totaled 39,793,503 and 36,531,517, respectively, accounting for 85.85 and 81.66% of total reads. Only the uniquely mapped reads were considered in this analysis. Distribution of reads in different regions of the reference genome was evaluated using the following standard metrics: exon, intron and intergenic reads. Table 2 summarizes these results.

Identification of DEGs
Although a number of genes were differentially expressed between the two libraries, the present study focused on genes meeting the designated criteria of p values < 0.05 and |log 2 FoldChange| > 1 (Additional file 1: Table S1). According to changes in relative gene abundance between the two libraries, 180 differentially expressed unigenes were down-regulated, and 157 unigenes were up-regulated in the MC group compared with the NC group of GMECs (Additional file 1: Table S1). Differential expression of transglutaminase 3, thromboxane synthase, ATPase, H+ transporting, lysosomal accessory protein 1-like, fragile histidine triad and seven other genes was observed at three more times in the MC and NC groups. Compared with the NC group, IGFBP5, SRY-box containing gene 11, phosphatidylethanolamine-binding protein 4, keratin 79, parkin 2 and fructosamine-3-kinase genes were upregulated in the MC group. FGF5, sprouty RTK signaling antagonist 4, mesenchymal stem cells, PTGS2, neuregulin 1 and TGFα genes were downregulated in the MC group compared with their expression in the NC group. Other differentially expressed genes between the two libraries included POU domain, class 2, transcription factor 3; IGFBP3, TIMP3, VEGFC, G protein-coupled receptor 151 and SOCS3 (Additional file 1: Table S1).

GO and KEGG pathway analysis of DEGs
To further investigate the biological function of DEGs, GO analysis was performed through running queries for each DEG against the GO database. Additional file 2: Figure S1 presents these GO functional annotations. Overall, GO terms were assigned to DEGs and classified into three independent ontology categories: (1) 334 terms in biological processes, (2) 209 terms in molecular functions and (3) 284 terms in cellular components. The most important enriched terms are bound up with biological regulation, single-organism processes and cellular processes. These terms were followed by developmental and metabolic processes. The IGFBP5 and IGF2 genes in the GO:0005576 term are known to affect mammary development and milk secretion [33,34], which were downregulated in the MC group. In the cellular component, GO terms with high levels of significance comprised cell parts, organelles and membrane parts (Additional file 2: Figure  S1). This organelle is linked with lipid secretor mechanisms of mammary epithelial cells [35]. Additionally, cell and membrane parts have effects on mammary gland involution [36].

Effect of chi-miR-3031 on IGFBP5 expression level
To examine whether chi-miR-3031 can reduce endogenous IGFBP5 expression, GMECs were transfected with MC, siRNA and NC. After 36 h of culture, total RNAs were extracted from the MC, siRNA and NC groups to determine the transfection efficiency of chi-miR-3031 and IGFBP5 mRNA expression. RT-qPCR results showed that the expression level of chi-miR-3031 in the MC group were 500 times higher than that of the NC group ( Fig. 1; P < 0.01). Additionally, the expression level of IGFBP5 mRNA in GMECs markedly decreased in the MC group ( Fig. 2a; P < 0.01), confirming transcriptome sequencing results. Western blot results showed that the expression level of IGFBP5 protein was significantly reduced in the MC group when compared with the NC group ( Fig. 3; P < 0.05). SiRNA-IGFBP5 significantly reduced the expression levels of IGFBP5 mRNA and protein ( Fig. 2b and Fig. 3; P < 0.01). These results suggest that chi-miR-3031 plays a negative regulatory role in IGFBP5 expression.

Chi-miR-3031 and IGFBP5 promote β-casein protein expression
To detect the concentrations of β-casein and κ-casein, GMECs were transfected with NC or MC by direct lysis after 36 h of culture. ELISA results showed that the expression levels of β-casein protein significantly increased in the MC and siRNA-IGFBP5 groups when compared with the NC group ( Fig. 4b; P < 0.01 or P < 0.05). However, the expression levels of κ-casein protein showed no significant differences among the NC, MC and siRNA -IGFBP5 groups (Fig. 4a). Note: Total reads: total number of sequencing reads. Total mapped reads: the reads that aligned to reference sequence and its ratio. Multiple mapped reads: in total mapped reads, reads aligned to two or more places. Uniquely mapped reads: in total mapped reads, reads aligned to only one position. Reads mapped in proper pairs: paired-end sequencing sequence aligned to the chromosome reasonable direction and position

Chi-miR-3031 induces p-mTOR expression via downregulating IGFBP5
Previous studies have suggested that PI3K/AKT-mTOR signaling pathways regulate milk protein synthesis [45]. Thus, this study investigated the expression levels of p-mTOR in GMECs transfected with NC, siRNA or MC.
Western blot results showed that chi-miR-3031 and siRNA-IGFBP5 induced the expression of p-mTOR (Fig. 5), suggesting that chi-miRNA-3031 can promote the expression of p-mTOR by downregulating IGFBP5.

Discussion
MiRNAs play key roles in the development of milk components [46,47] and mammary glands [48]. However, only a limited number of genes are known to be directly or indirectly regulated by miRNAs in GMECs. Thus, the detection and identification of relevant genes are necessary to determine their function. Our previous studies investigated the differential expression of miRNAs in caprine mammary gland tissues during colostrum and common milk stages using high throughput sequencing and bioinformatics. Herein, we focused on chi-miR-3031 due to its prominent roles in influencing mammary gland development and lactation performance [20]. In MC and NC with transfected GMECs, we selected DEGs that were directly or indirectly regulated by chi-miR-3031. A total of 157 up-regulated and 187 downregulated genes were detected in the MC group compared with the NC group.    [52]. Additionally, miR-221/222 cluster, miR-18a and transcription factors homeobox A7, and NOTCH 3 and 4 have important roles in regulating lactose. MiR-186, miR-148a and miR-200a are remarkably associated with milk yield and are likely the most important miRNAs for component traits [53]. MiR-130a overexpression decreased FABP3, PPARG, FATP1 and PLIN2 mRNA expression levels, whereas miR-130a downregulation increased PPARG, C/ EBP beta, C/EBP alpha, PLIN2, FABP3 and FATP1 mRNA expression levels [19]. In this study, results from transcriptome sequencing, real-time PCR and western blot showed that chi-miR-3031 can suppress IGFBP5 mRNA and protein expression levels. However, β-casein expression levels significantly increased in the MC and siRNA-IGFBP5 groups. IGFBP5 is an important member of the IGFBP family and has high affinity for IGF [54]. IGFs act as key upstream activators of PI3K/AKT-mTOR signaling pathways that regulate milk protein synthesis [29,31]. We observed that the down-regulation of IGFBP5 can activate the phosphorylation of mTOR at the Ser2448 site in GMECs transfected with MC and siRNA-IGFBP5. Previous findings have suggested that  the activation of p-mTOR implies activation of the PI3K pathway [55]. We speculated that chi-miR-3031 can activate the PI3K/AKT-mTOR pathway and increase β-casein expression by down-regulating IGFBP5. However, κ-casein was not affected by chi-miR-3031.

Conclusion
In this study, 157 upregulated and 187 downregulated genes were detected by transcriptome sequencing in the MC group. KEGG pathway analyses and previous studies showed that TGFα, PIK3R3, IGF2, IGFBP5, ELF5 and LHβ play key roles in mammary development and milk secretion. Thus, the differential expression of these genes can take part in the regulation of lactation or mammary gland development. Additionally, chi-miR -3031 can activate the PI3K/AKT-mTOR pathway and improve the expression of β-casein by down-regulating IGFBP5. Therefore, we conclude that chi-miR-3031 plays key roles in regulating β-casein synthesis in GMECs.

Animals and cell culture
Three Guanzhong dairy goats (3 years old, female) at peak lactation (90 days postpartum) periods were selected from a local farm in the Northwest A&F University of China. These Guanzhong dairy goats are research animals. Half-hour before the operation, these animals were anesthetized by intramuscular injection of 150 mg phenobarbital sodium. In the middle part of mammary gland, open the wound of 1 cm with a scalpel, remove about 1 cm 3 of mammary gland tissue, then sew and sterilize the wound. The special person kept the experimental animals, and after 1 week, the surgical line was removed, and all the animals recovered. Mammary gland tissues for culturing GMECs were stored in PBS containing streptomycin (100 mg/mL) and penicillin (100 U/mL). GMECs were cultured on the basis of a previous study [56].

Transfection and total RNA extraction
In accordance with the Lipofectamine RNAiMAX Reagent (Invitrogen, Carlsbad, CA, USA) protocol, GMECs were transfected with MC, NC and siRNA-IGFBP5 at final concentrations of 100 nmol/L. The MC, siRNA and NC were chemically synthesized and purified by Shanghai Gene-Pharma (Shanghai, China). GMECs were cultured for another 36 h after transfection, and TRIzol reagent was used to extract total RNA. The quality of total RNA was detected using an Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA). The RIN values of the MC group were 9.2, 9.0 and 8.8, and those of the NC group were 9.3, 9.1 and 8.9. Additional file 5: Table S3 provides the sequences of MC, NC and siRNA.

Transcriptome sequencing Library construction and sequencing
For each independent experiment, three wells of a 6-well plate were mixed to extract the total RNAs to construct the MC group. The NC group was constructed in the same way. The total RNAs for each group were homogenized to maintain the same concentrations. Sequencing libraries were constructed with NEBNext Ultra RNA Library Prep Kit for Illumina (NEB, Beijing, China). Two MC groups were independently sequenced, respectively. The NC group was sequenced in the same way. Then, in accordance with instructions for use, libraries with different indices were loaded on an Illumina HiSeq instrument (Illumina, San Diego, CA, USA). Sequencing data is processed and analyzed by Personalbio (Shanghai, China).

Mapping reads on goat reference genome and analysis of gene expression
After removal of the 3ˊadaptor sequence, low-quality and redundant reads were removed using Cutadapt (version 1.9.1) and clean data were aligned to caprine genome via Hisat2 software (v2.0.1). Caprine genome sequences and gene model annotation files were downloaded from the NCBI genome website (GCF_001704415.1_ASM170441v1).
The level of gene expression was measured by read density, and gene expression calculation was performed with Cuffdiff (v2.2.1), which calculates fragments per kilo bases per million (FPKM) reads as FPKM = total exon reads/mapped reads in millions × exon length in kb. The results from the cuffdiff analysis were further analyzed to identify genes with significant differential expression according to the criteria of fold-changes greater than 1 and p values less than 0.05.

Gene ontology (GO) and KEGG pathway analysis of DEGs
The DEGs were mapped to GO terms in the Gene Ontology database (http://www.geneontology.org/). Remarkably enriched GO terms in differentially methylated genes were defined by hypergeometric test [57][58][59]. Pathway-based analyses aid in further understanding gene function [60]. The calculation formula of KEGG (http://en.wikipedia.org/wiki/KEGG) [61] was the same with that of GO analysis. Pathway analysis identified remarkably enriched signal transduction pathways for DEGs [62]. Table 3 shows the names of miRNA and genes for reverse transcription and real-time PCR. The reverse transcription system contained the following: (1) 800 ng of total RNA, 1 μL of gDNA Eraser, 2 μL of 5 × gDNA Eraser Buffer and RNase-Free dH 2 O mixed to a final volume of 10 μL and incubated at 42°C for 2 min. (2) 4 μL of 5 × Prime Script Buffer 2, 1 μL of RT Primer Mix, 1 μL of PrimeScript RT Enzyme Mix 1 and RNase-Free dH 2 O mixed to a final volume of 20 μL and incubated at 42°C for 15 min followed by 85°C for 5 s using a PrimeScript RT reagent Kit (TaKaRa, Dalian, China). cDNA products were stored at − 20°C.

Reverse transcription and real-time PCR
Real-time PCR was performed with a 25 μL volume containing 2 μL of template cDNA, 12.5 μL of SYBR Premix Ex Taq II (TaKaRa, Dalian, China) and 1 μM of primers using the CFX Connect Real-Time PCR Detection System (Bio-Rad, CA, USA). Thermal cycling conditions were 95°C for 10 min, followed by 40 cycles at 94°C for 15 s, 60°C for 30 s and 72°C for 30 s. Primers are shown in Table 3. The β-actin or U 6 mRNA levels were used for normalization. Each experiment was independently repeated at least three times, and fold changes were analyzed via the 2 -ΔΔCt method.

Enzyme-linked immunosorbent assay
ELISA was performed to detect the concentration of caseins in cell lysates using goat β-casein and κ-casein kit (Tongwei, Shanghai, China). Briefly, GMECs were transfected with NC or MC after 36 h of culture in six-well plates and then lysed with 200 μL RIPA Lysis Buffer (Bioteke, Beijing, China). According to manufacturer instructions, 50 μL cell lysate or standards were used, and absorbance at 450 nm was measured through an Epoch microplate reader (Biotek, Winooski, USA). Concentrations were calculated after subtracting the standard curve. Kit sensitivity was 1.0 μg/mL, and mean intra-and inter-assay variable coefficient values were less than 15 and 10%, respectively.

Western blotting
GMECs were lysed in denaturing lysis buffer containing protease inhibitors (RIPA, Bioteke, Beijing, China) for 30 min on ice and centrifuged (12,000×g) for 15 min at 4°C. Protein concentrations in the lysates were determined by

Statistical analysis
The one-way ANOVA of SPSS 16.0 software package was used to analyze data. Means of multiple comparisons was analyzed by least significant difference (LSD). Data were presented as means ± standard deviation of three independent experiments. P < 0.05 or P < 0.01 was considered significant.