Skip to main content
  • Research article
  • Open access
  • Published:

Transcriptome sequencing to detect the potential role of long non-coding RNAs in bovine mammary gland during the dry and lactation period

Abstract

Background

It is known that long non-coding RNAs (lncRNAs) play an important role in various biological processes, including cell proliferation, differentiation and apoptosis. However, their functions and profiles in lactation cycle of dairy cows are largely unknown. In this study, lncRNA-seq technique was employed to compare the expression profiles of lncRNAs and mRNAs from Chinese Holstein mammary gland in dry and lactation period.

Result

Totally 3746 differentially expressed lncRNAs (DELs) and 2890 differentially expressed genes (DEGs) were identified from the dry and lactation mammary glands of Holstein cows. Functional enrichment analysis on target genes of lncRNAs indicated that these genes were involved in lactation-related signaling pathways, including cell cycle, JAK-STAT, cell adhesion, and PI3K-Akt signaling pathways. Additionally, the interaction between lncRNAs and their potential miRNAs was predicted and partly verified. The result indicated that the lactation-associated miR-221 might interact with lncRNAs TCONS_00040268, TCONS_00137654, TCONS_00071659 and TCONS_00000352, which revealed that these lncRNAs might be important regulators for lactation cycle.

Conclusion

This study provides a resource for lncRNA research on lactation cycle of bovine mammary gland. Besides, the interaction between lncRNAs and the specific miRNA is revealed. It expands our knowledge about lncRNA and miRNA biology as well as contributes to clarify the regulation of lactation cycle of bovine mammary gland.

Background

Mammary gland is an important organ for the synthesis and secretion of milk, which provides essential nutrients for human and other animal offspring. The development and regression cycles of mammary gland include pregnancy, lactation, and involution, which is regulated by growth factors, hormones, and coding genes [1,2,3]. Janus kinases and signal transducers and activators of transcription ((JAK-STAT) have been shown to function downstream of several peptide hormones and cytokines that are required for postnatal development and secretory function of mammary gland [4]. In mammary epithelial cells, the phosphorylated STAT5A and STAT5B form homodimers and heterodimers to modulate differentiation, survival and proliferation through alterations in cellular gene expression [5]. The phosphatidylinositol 3-kinase-proteinkinase B/mammalian target of the rapamycin (PI3K-Akt/mTOR) signaling pathway regulates a wide range of cellular processes, such as cell proliferation, growth, survival and metastasis [6], and it is essential for mammary gland development [7]. A conditional knockout of Akt1 averts the extended survival of mammary epithelial cells that express hyperactive STAT5, which indicates that the PI3K-Akt/mTOR pathway is a crucial downstream effector of JAK-STAT signaling [8].

In the past few years, many studies focused on the functions of protein-coding genes and microRNAs (miRNAs) for the development and regression cycles of mammary gland [9, 10]. ErbB3 played a crucial role in mammary epithelial survival and differentiation during pregnancy and lactation [11]. Forty milk lipid synthesis- and secretion-associated genes, including FADS1, AGPAT6, GPAM, LPIN1, BTN1A1, LPL, CD36, FABP3, ACSL1, ACSS2, ACACA, FASN, SCD, XDH, BDH1, INSIG1, PPARG, and PPARGC1A were verified from dry period to the end of subsequent lactation period [12]. During lactation, MAPK14, FRAP1, EIF4EBP2, GSK3A and TSC1 had an increased expression which revealed that these genes had important roles in milk protein synthesis and secretion [13]. MiRNAs are a kind of non-coding RNA, which can silence or degrade gene expression by targeting the 3’UTR region of coding gene. An increasing number of studies had demonstrated that miRNAs were involved in lactation of mammary gland by regulating their target genes [14,15,16,17,18]. It had been found that miR-27a could regulate cellular triacylglycerol synthesis by targeting PPARG gene in bovine mammary epithelial cells (BMECs) [15]. The overexpression of miR-206 changed the expression of Wnt, Tbx3 and Lef1 genes which were essential for mammary gland development, indicating that miR-206 might be a novel candidate for morphogenesis during the initiation of mammary gland formation [16]. As a downstream regulator of PTEN, miR-486 expressed higher during bovine high quality lactation period and could regulate the secretion of β-casein, lactose and triglyceride in BMECs, which indicated that miR-486 was required for the development of cow mammary gland [17]. In cattle, the expression of miR-221 was found to be higher in peak lactation than in early lactation, suggesting its role in the control of endothelial cell proliferation or angiogenesis [18]. In mice, miR-221 regulated lipid metabolism in mammary epithelial cells and was expressed differentially at various stages during mammary gland development [19]. MiR-212 and miR-132 were necessary for mice mammary gland development and growth by targeting MMP9 gene, especially for mammary epithelial ducts [20].

Long non-coding RNAs (lncRNAs) are transcripts longer than 200 nucleotides, which have uncovered new layers in the control of various biological processes, including cell proliferation, differentiation and apoptosis [21]. H19 and SRA, two of the earliest identified regulatory lncRNAs, may play a role in the developing mammary gland [22]. Initially, H19 was found to function to restrict embryonic growth, but later evidence showed that H19 had a role in long-term maintenance of adult hematopoietic stem cells [23]. Long noncoding RNA mPINC (mouse pregnancy-induced non-coding RNA) and Zfas1 (Znfx1 antisense 1) had growth-suppressive roles in mammary epithelial cells [24, 25]. Previous observation demonstrated that lncRNA Neat1 could regulate mammary gland morphogenesis and lactation by investigating the proliferation of Neat1-mutant cells in mice [26]. lncRNAs also could act as competing endogenous RNAs (ceRNAs) to control muscle differentiation and involve in goat lactation process [27, 28]. Pregnancy-induced noncoding RNA (PINC) could inhibit terminal differentiation of alveolar cells during pregnancy to prevent abundant milk production and secretion until parturition [29]. The above-mentioned studies showed that lncRNAs had crucial roles in mammary gland cell proliferation and differentiation, which would be an important issue in lactation biology. Mammary gland development and regression was directly related with cow lactation. However, the functions and profiles of lncRNAs in bovine mammary gland during dry and lactation period were largely unknown. The objective of this study was to screen the lncRNAs associated with lactation by lncRNA-seq analysis in bovine mammary gland, which would provide a new perspective for lncRNAs in lactation biology and lay the foundation for their further function study in milk composition synthesis and secretion.

Methods

Animals and mammary gland tissue collection

Eight healthy Chinese Holstein dairy cows at the third or fourth parities used in this study were housed in free stall and had access to water and feed ad libitum at the Experimental Farm of Northwest A&F University (Yangling, Shaanxi, China). After intravenous injection of lidocaine hydrochloride, approximately 4 g of mammary gland tissues were harvested via repeated biopsies from four cows at the dry period and four cows at approximately 180 days during lactation period. The tissues were dissected, frozen in RNA later (TaKaRa, Dalian, China) and stored at − 80 °C for further analysis. All experimental and surgery procedures involved in this study were approved by the Experimental Animal Manage Committee of Northwest A&F University (2011–31101684).

Total RNA isolation and quality control for library construction

Total RNA of mammary gland tissues was isolated using Trizol reagent following the manufacturer’s instructions (Invitrogen, CA, USA). The integrity of RNA was detected using RNA Nano 6000 Assay Kit on the Bioanalyzer 2100 system (Agilent Technologies, Santa Clara, USA). RNA purity and concentration was measured using Nanodrop 2000 photometer spectrophotometer (Implen, Los Angeles, USA). The 260/280 ratio for all the samples from mammary gland tissues was about 2.0, and the RNA integrity number (RIN) was ≥8.0. After the determination of RNA purity and quality, 3.0 μg RNA per sample was used and ribosomal RNA was removed using Epicentre Ribo-zero™ rRNA Removal Kit (Epicentre, Madsion, WI, USA) for library construction. The RNA from three individuals in dry and three in lactation stage were pooled, respectively. Subsequently, the libraries were generated from the rRNA-depleted RNA pools using the NEBNext® Ultra™ Directional RNA Library Prep Kit for Illumina® (NEB, USA) following the manufacturer’s recommendations. In order to select cDNA fragments of preferentially 250~ 300 bp in length, the library fragments were purified with AMPure XP system (Beckman Coulter, Beverly, USA).

RNA sequencing, transcriptome assembly, and quantification of gene expression level

The coded samples were clustered using TruSeq PE Cluster reagent (Illumina, CA, USA) according to the manufacturer’s instructions on a cBot Cluster generation system. The libraries were sequenced on an Illumina Hiseq3000 platform and 100-bp paired-end reads were generated after cluster generation. The schematic of lncRNA-seq analysis was shown in Fig. 1. Raw data were first processed using in-house Perl scripts. In this step, clean data were obtained by trimming reads containing adapter, reads containing over 10% of ploy-N, and low-quality reads (> 20% of bases whose Phred scores were < 20) from the raw data. Phred = −10log10 (e), e is defined as the error probability of sequencing for every single base. Q20, Q30 and GC content of the clean data were calculated. Subsequently, Bowtie (v2.0.6) [30] and Tophat2 [31] (v2.0.9) was used to align paired-end clean reads to the reference genome (version GCA_000003055.5_Bos_taurus_UMD_3.1.1). The default parameters for Tophat2 were set as ‘-read-mismatches=2 (≤2 mismatches are allowed) and -read-gap-length=2’ (≤2 gaps are allowed). The mapped reads were assembled using Cufflinks (v2.1.1) in a reference-based approach [32]. Cufflinks was run with ‘min-frags-per-transfrag = 0’ and ‘-library-type’, other parameters were set as default. Fragments per kb for a million reads (FPKM) of both coding genes and lncRNAs were analyzed using Cuffdiff (v2.1.1) [33]. Differential expression analysis between two groups was performed using the DESeq R package (1.8.3). The P-values were adjusted using the Benjamini and Hochberg method [34]. P-adjust (q-value) < 0.05 and the absolute value of fold change≥2 were set as the threshold for significantly differential expression.

Fig. 1
figure 1

The schematic of RNA-seq analysis

Identification the annotated and novel lncRNAs

NONCODE database was used to characterize the annotated lncRNAs in bovine from the assembled transcripts. To identify bovine novel lncRNAs, the steps were followed as Wang et al. (2017) described [35] with little modification. (1) transcripts with length < 200 bp were removed; (2) transcripts with predicted ORF > 300 were removed; (3) transcripts were compared with mRNA, rRNA, tRNA, snRNA, snoRNA and pre-miRNA (https://www.ncbi.nlm.nih.gov/) using Cuffcompare v2.1.1 to remove the same or similar transcripts [32]. (4) transcripts with FPKM < 1 were removed; (5) transcripts that did not pass the protein-coding-score test were removed using the Coding Potential Calculator (CPC) [36] and PFAM database [37].

Prediction and functional enrichment analysis of lncRNA target genes

To reveal the potential function of lncRNAs, their target genes were predicted in trans and cis. For cis role, it refers to lncRNA’s action on neighboring target genes. In this study, the coding genes from 100 kb upstream and downstream of an lncRNA were searched. The trans role refers to the influence of lncRNAs on other genes at expression level. RNAplex bioinformatics software (http://www.bioinf.uni-leipzig.de/~htafer/index.html) was used to predict lncRNA target genes in trans. Genome distribution of differentially expressed lncRNAs and mRNAs was illustrated with Circos (http://circos.ca/). GO enrichment analysis was performed with DAVID database (https://david.ncifcrf.gov/). As to KEGG analysis, differentially expressed genes were analyzed with KEGG online website (http://www.genome.jp/kegg/). Protein-protein interaction network between differentially expressed genes were analyzed by STRING database (https://string-db.org/), and further visualized with Cytoscape (http://www.cytoscape.org/).

Prediction of potential miRNAs interacted with lncRNAs

To obtain potential miRNAs interacted with lncRNAs, PicTar (https://pictar.mdc-berlin.de/), PITA (https://genie.weizmann.ac.il/pubs/mir07/mir07_prediction.html), and RNAhybrid (https://bibiserv.cebitec.uni-bielefeld.de/rnahybrid/) were used to predict the candidate miRNAs interacted with lncRNAs. The miRNAs shared by the above three tools were selected as the candidate miRNAs to assume the mechanism of lncRNAs interacted with lactation. The potential target genes of miR-221 were predicted by Targetscan (http://www.targetscan.org/vert_71/), PITA and miRanda. The genes shared by the above three tools were selected as the candidate target genes.

Verification of sequencing data using qRT-PCR

To confirm the sequencing results, quantitative real time PCR (qRT-PCR) method was used to measure the relative expression of DEGs and DELs. The total RNAs from the cows used for lncRNA-seq were also used for qRT-PCR. The first-strand cDNA was obtained using a PrimeScript RT reagent Kit (TaKaRa, Dalian, China) following the manufacturer’s instructions. The qRT-PCR was performed in triplicate using SYBR® Premix Ex Taq™ II (TaKaRa) on the Bio-Rad CFX96 Touch™ Real Time PCR Detection System (Bio-Rad, Hercules, CA, USA). The qRT-PCR protocol was initiated with an incubation of 3 min at 95 °C, followed by 40 cycles of 95 °C for 12 s and optimized annealing temperature for 40 s. The relative expression of DEGs and DELs were analyzed using the 2−ΔΔCt method and normalized by GAPDH gene [38]. The qRT-PCR primers for DELs and DEGs were designed with Primer Premier 6.0 (Premier, British Columbia, Canada) spanning at least one intron and shown in Additional file 1.

BMECs transfection and expressions of four lncRNAs and ErBb3

BMECs were cultured in Dubecco’s modified Eagle medium (DMEM) /F12 containing 10% fetal bovine serum and 1% penicillin streptomycin (All from Gibco, Grand Island, NY, USA) at 37 °C with 5% CO2. The cells were seeded in 24-well plates, then transfected with miR-221 mimic, mimic-NC, miR-221 inhibitor and inhibitor-NC (GenePharma, Shanghai, China) at approximately 50% confluence using Lipofectamine 2000 Transfection Reagent (Invitrogen, USA) according to the manufacturer’s instructions, respectively. Briefly, Lipofectamine 2000 Transfection Reagent (1 μL) was diluted in 25 μL of Opti-MEM (Invitrogen) with miR-221 mimic, mimic-NC, inhibitor, or inhibitor-NC to yield final concentrations of 20, 20, 40, and 40 nM, respectively. The mixture was incubated at room temperature for 20 min and added to the BMECs. The transfection efficiency was assessed by fluorescence microscopy after 48 h.

Total RNA was extracted from the transfected BMECs using Trizol reagent (Invitrogen) 72 h after transfection. Then, the cDNA was obtained as above mentioned method. qRT-PCR protocol was also performed as above mentioned. Primers for the four lncRNAs and ErBb3 are shown in Additional file 1.

Results

Overview of sequencing data in cow mammary gland

A total of 58,411,766 and 69,114,038 raw reads were produced from cow mammary glands using the Illumina Hiseq3000 platform in dry and lactation periods, respectively. After discarding adaptor sequences and low-quality sequences, 54,805,136 and 65,069,892 corresponding clean reads were obtained, and the percentage of clean reads was 93.83 and 94.15%, respectively (Table 1). The whole expression feature of transcripts was shown in Fig. 2a. The expression level of the transcripts in lactation was slightly higher than that in dry period. Similarly, 64,316 and 61,791 known transcripts, 44,651 and 43,094 known mRNAs were also obtained in dry and lactation period of cow mammary gland, respectively. Of those transcripts, 928 were novel in dry stage and 841 were in lactation, respectively.

Table 1 The information of sequencing data
Fig. 2
figure 2

The genome distribution and relative expression of clean reads and novel transcripts. a The relative expression of novel transcripts and b the genome distribution of clean reads in cow mammary glands in dry and lactation periods, respectively

The genome distribution of clean reads was shown in Fig. 2b. The results revealed that most clean reads in the two periods were located in exonic region, only few reads were in intergenic region. And about a quarter of the total clean reads were in the intronic region of the bovine genome.

GO and KEGG enrichment analysis of DEGs

Totally 2890 differentially expressed genes (DEGs) were obtained according the fold change≥2 and q < 0.01 (Additional file 2). And 2300 genes were down-regulated, whereas 590 were up-regulated in lactation compared with dry period. Then functional enrichment analysis of the 590 significantly up-regulated genes was performed based on their relative expression and fold changes. As illustrated in Fig. 3a, 69 significantly enriched GO terms (P < 0.05) were identified, such as negative regulation of cell growth (GO: 0030308), growth factor activity (GO: 0008083), positive regulations of ERK1 and ERK2 cascade (GO: 0070374), and sodium channel complex (GO: 0034706), etc. ERK1 and ERK2 had been suggested to play an important role in regulating cell invasion, cell proliferation, and colony formation in triple-negative breast cancer cell lines [39]. Growth factors, especially the epidermal growth factors (EGFs) and insulin-like growth factors (IGFs), were involved in development of normal mammary gland and pathogenesis of breast cancer [40].

Fig. 3
figure 3

Functional enrichment analysis of significantly up- and down-regulated DEGs. a Gene Ontology (GO) analysis of significantly up-regulated DEGs. GO terms revealed that significantly up-regulated genes might relate to lactation-associated processes, including calcium ion binding, negative regulation of cell growth and proliferation, gene expression, etc. b GO analysis of down-regulated DEGs revealed that those genes might play an important role in cell division, intracellular signal transduction chemokine-mediated, and cell surface receptor signaling pathways, etc. c KEGG pathway enrichment analysis of significantly up- and down-regulated DEGs. KEGG terms, including cell cycle, cell adhesion molecules (CAMs), cytokine-cytokine receptor interaction, cAMP, calcium, and Wnt, PPAR, PI3K-Akt, TNF, cytokine-cytokine pathways, etc. d Cluster analysis of some lactation-related genes by using heat map. The genes in red box located in PI3K-Akt signaling pathway. Similarly, the genes in green box were in PPAR signaling pathway

Two thousand two hundred eight significantly down-regulated genes were also selected to carry out the functional enrichment analysis. Forty five terms, including cell division (GO: 0051301), cell adhesion (GO: 0007155), cell communication (GO: 0007154), cation transmembrane transport (GO: 0098655), and cell surface receptor signaling pathway (GO: 0035556), were shown in Fig. 3b. It was worthy to notice that the process of cell adhesion and division were closely associated with mammary gland architecture construction, maintenance, development, and lactation [41, 42].

Meanwhile, KEGG results from the significantly up- and down-regulated genes indicated that 53 significantly signaling pathways were enriched (Fig. 3c, P < 0.05), such as cell adhesion molecules (CAMs), PI3K-Akt, PPAR, TNF, cytokine-cytokine, cell cycle, and Wnt signaling pathways. Previous studies had determined that Wnt, PPAR, CAMs, PI3K-Akt, and TNF signaling pathways could regulate mammary gland development, milk fat synthesis, mastitis, BMECs proliferation and apoptosis [43,44,45,46,47]. The expression pattern of genes involved in PI3K-Akt and PPAR signaling was displayed in Fig. 3d. The heat map showed that PIK3R3, CSF3, TNC, TLR2, GNG11, DOIT43, NOS3, THBS3, IL3RA, FN1, SORBS1, FAPPS, SLC27A, ACP7, GK, ACSL4 and ANGPTL4 genes expressed higher in lactation than that in dry period at mRNA level.

DELs and their potential interacted miRNAs involved in lactation

Totally 23,495 expressed lncRNAs were found in the two different periods, of which 5893 were novel (Additional files 3 and 4) and 17,602 were annotated in NONCODE (Additional file 5). A total of 3746 significantly differentially expressed lncRNA transcripts were found in lactation compared with that in dry period, including 2732 down- and 1014 up-regulated lncRNA transcripts (fold change≥2, P < 0.05) (Additional files 6 and 7). The genome wide distribution of DELs and DEGs was shown in Fig. 4a. It could be seen that the DELs and DEGs was harmoniously located on autosomes and X chromosome. The stage-specific expressed lncRNAs and genes were shown in Additional files 8 and 9. One hundred forty five lncRNAs were identified and subjected to further analysis according their expression fold changes and genomic locations. The potential interacted miRNAs with the selected 145 lncRNAs were predicted by miRanda, PITA and RNAhybrid, and partial result was shown in Table 2. It could be seen that miR-103, miR-21-3p, miR-27a-5p, miR-107, and miR-24-3p were predicted to interact with lncRNA TCONS_00120917. The interaction network between candidate lncRNAs and their potential miRNAs was constructed (Fig. 4b). It showed that several miRNAs might interact with multiple lncRNAs, such as miR-221 interacting with lncRNAs TCONS_00032404, TCONS_00032444, TCONS_00040268, TCONS_00137654, TCONS_00071659, TCONS_00143274 and TCONS_00000352. Coincidentally, it had demonstrated that the expression of miR-221 was found to be higher in peak lactation than in early lactation, suggesting a role in the control of endothelial cell proliferation or angiogenesis which was closely related with lactation [19], so those above mentioned lncRNAs could be considered as important candidates for lactation.

Fig. 4
figure 4

DEGs, DELs and its potential interacted miRNAs. a Chromosome distribution of DELs and DEGs by using Circos. Red columns represent DEGs, and blue bars represent DELs. b The interaction network between DELs and their potential interacted miRNAs. Those lncRNAs were depicted as pink and purple, and miRNAs as blue and yellow

Table 2 Differently expressed lncRNAs and their potential lactation-related miRNAs

The interaction between lncRNAs and miR-221 as well as miR-221 and ErBb3

To verify the interaction between lncRNAs and miR-221 and the interaction between miR-221 and ErBb3, overexpression and inhibition of miR-221 was employed in BMECs. Our data showed that the expressions of the four lncRNAs (TCONS_00040268, TCONS_00137654, TCONS_00071659, and TCONS_00000352) decreased with the overexpression of miR-221. On the contrary, their expressions were increased with the inhibition of miR-221 in BMECs (Fig. 5a), which indicated that miR-221 could interact with these four lncRNAs; Furthermore, the overexpression of miR-221 could reduce the expression of ErBb3 gene in miR-221 mimic group compared with mimic-NC group. While the inhibition of miR-221 increased ErBb3 gene expression, which indicated that ErBb3 gene was the target of miR-221 (Fig. 5b).

Fig. 5
figure 5

The interaction between lncRNAs and miR-221 as well as miR-221 and ErBb3 gene. a The relative expression of 4 lncRNAs (TCONS_00040268, TCONS_001371654, TCONS_00071659, and TCONS_00000352) changed with the overexpression and inhibition of miR-221 b miR-221 could target ErBb3 gene by transfecting miR-221 mimic and inhibitor in BMECs

Interaction network between the candidate lncRNAs and their potential target genes

To reveal the potential function of these selected lncRNAs, their target genes were predicted in cis and trans roles, and we found that the mammary gland biology related genes, such as PRLR, FAS, and MAP3K7 genes, could be targeted by several lncRNAs in cis or trans (Table 3). In addition, the network between lncRNAs and target genes was constructed. The potential target genes of those lncRNAs included S100A1, UPK1B, BARD1, FGF2, IGFBP1, ALOX15, KLH34, SYT12, WNT10A, SULF1, SLC10A6 and CLE7A (Fig. 6). Moreover, the expression level of those genes in dry period were significantly higher than that in lactation period in our study while the corresponding lncRNAs were significantly lower, suggesting an opposite expression trend exists between those genes and corresponding lncRNAs. Previous data had demonstrated that IL1B might function as a retinoic acid induced gene and could inhibit the proliferation of normal human mammary epithelial cells [48]. Therefore, it could be deduced that TCONS_00075230 might involve in regulating the proliferation of mammary epithelial cells by interacting with IL1B gene in cow mammary gland.

Table 3 DELs and their potential target genes involved in lactation
Fig. 6
figure 6

The interaction network between candidate lncRNAs and their potential target genes. The red boxes indicated up-regulated candidate lncRNAs. And the blue boxes were down-regulated lncRNAs. The up- regulated genes in the dry period were depicted as green and down-regulated genes as,yellow, respectively

GO and KEGG analysis of lncRNAs potential target genes

A total of 47 potential target genes in cis of lncRNAs (Additional file 10) were selected to carry out functional enrichment analysis. Our results showed that nine significantly enriched GO terms focused on the negative regulation of JAK-STAT cascade (GO: 0046426), positive regulation of transcription from RNA polymerase II promoter (GO: 0045944), cytokine-mediated signaling pathway (GO: 0019221), and positive regulation of phosphorylation (GO: 0042327), etc. (Fig. 7a, P < 0.05). Previous study had suggested that JAK-STAT pathway made a marked contribution to dairy milk production traits [49]. In addition, phosphorylation regulation of mammary S6 K1 and eIF2 genes could control milk protein yield [50].

Fig. 7
figure 7

Functional enrichment analysis of lncRNAs target genes in cis and trans. a Gene Ontology (GO) analysis of lncRNAs target genes in cis and trans. The results indicated that those genes seem to play an essential role in lactation-related pathways, including phosphorylation, protein glycosylation, cytokine-mediated, and positive regulation of cell division signaling pathway, etc. b KEGG pathway enrichment analysis of lncRNAs target genes in cis and trans. These KEGG terms were lactation-associated pathways, including NF-kappa B, MAPK, PI3K-Akt, prolactin, Toll-like receptor, and CAMs signaling pathway, etc.

Similarly, a total of 459 potential target genes in trans of lncRNAs (Additional file 11) were selected to perform functional enrichment analysis (Fig. 7a, P < 0.05). These terms were protein glycosylation (GO: 0006486), positive regulation of cell division (GO: 0051781), DNA-directed RNA polymerase III complex (GO: 0005666), and growth factor binding signaling pathway (GO: 0019838), etc. It was worth mentioning that cell division signaling pathway had been demonstrated to be closely associated with mammary gland architecture construction and maintenance [42]. Additionally, protein glycosylation in mammal membrane played important roles in milk quality and biomodulator properties [51].

Meanwhile, the KEGG results from the target genes of candidate lncRNAs in trans and cis illustrated 15 pathways (Fig. 6b). Among these pathways, PI3K-Akt and MAPK pathways had been confirmed to be associated with cow lactation [43, 52].

Interaction network of protein-protein

Combine the bioinformatics analysis of DEGs and target genes of DELs, the interaction network between proteins was produced using String software (Fig. 8). The results revealed that protein-protein interaction focused on CCN (CCND1, CCNA1, CCNB2, CCNA2, CCNB1, CCNE1, and CCNE2) and CDC protein family (CDC6, CDC20, CDC25B, CDC25C, and CDCA8). CCN and CDC protein families had been demonstrated to be involved in cyclin and cell division cycle, respectively [53, 54].

Fig. 8
figure 8

Protein-protein network was analyzed using the String. Cyclin-associated proteins were depicted as light blue frame, cell division cycle-associated protein as red oval symbol, cyclin dependent kinase-associated protein as yellow frame. Their associated proteins were depicted as dark blue

qRT-PCR

Seven lncRNAs and five genes, which were significantly differentially expressed between the two periods, were selected to verify the transcriptome sequencing data using qRT-PCR. The seven lncRNAs were selected based on maximal fold change (and P < 0.05) and location of the genome. For the five genes, they were selected because of their fold change (and P < 0.05) and potential function in lactation. The results showed that the relative expression of the selected lncRNAs and genes was consistent with their sequencing data (Fig. 9).

Fig. 9
figure 9

Validation of DEGs and DELs in dry and lactation periods by qRT-PCR. Relative expression of DEGs (FGFBP1, IGFBP5, LPO, CXCL10, and SAA3) and DELs (TCONS_00040268, TCONS_00137654, TCONS_00071659, TCONS_00000352, TCONS_01118313, TCONS_ 00093001 and NONBTAG0105046.2) were verified by qRT-PCR, and the results indicated that the relative expressions of TCONS_00040268, TCONS_00137654, TCONS_00071659, TCONS_00000352 and TCONS_ 00093001 in lactation were significantly higher than that in dry stage (P < 0.05). On the contrary, CXCL10 and TCONS_01118313 in lactation were lower than that in dry stage (P < 0.05). The qRT-PCR results were consistent with the sequencing data

Discussion

Increasing data had shown that lncRNAs could regulate gene expression both at transcriptional and post-transcriptional levels, including the regulation of splicing, mRNA processing, and translation [55]. It had been demonstrated that some lncRNAs implicated in breast normal development and cancer based on their expression pattern, function and localization in human and mouse [56, 57]. However, limited researches had been reported the regulation mechanism of lncRNAs on bovine lactation. Koufariotis et al. (2015) identified and annotated novel lncRNA transcripts in the bovine genome across 18 tissues (including mammary gland) via RNA sequencing. In addition, they found most lncRNA in bovine mammary gland were downregulated compared with other tissues [58]. Tong et al. (2017) identified 36 lincRNAs located in milk related quantitative trait loci (QTL) from five RNA-seq datasets of bovine mammary glands whereas one lincRNA was within clinical mastitis QTL region, which indicated that these lncRNAs were involved in many biological functions including susceptibility to clinical mastitis as well as milk quality and production [59]. In our study, lactation-related lncRNAs were assessed by directly lncRNA-seq from cow mammary glands at the dry and lactation period, respectively. Given the functional enrichment analysis of lncRNAs target genes, 24 lncRNAs were identified to have potential role in lactation (Table 3). Previous study had suggested that lncRNAs could regulate biological processes by sponging miRNAs [28], so 15 lncRNAs interacted with lactation-associated miRNAs were also considered as potential regulators for lactation (Table 2). Among those lncRNAs, five lncRNAs were shared by the above two prediction methods.

LncRNAs might regulate expression of lactation-associated genes in trans or cis

It was well recognized that a series of genes involved in lactation initiation, maintenance as well as mammary gland growth, development and breast cancer by direct or indirect regulation [8, 9, 11, 60]. In this study, the DEGs functional enrichment results showed that these genes were related to some biological processes, including cell division, adhesion, cycle, ERK1 and ERK2 cascade, PI3K-Akt, PPAR, and TNF pathway, which were closely associated with lactation [43,44,45,46].

LncRNAs can regulate neighboring gene expression, therefore high expression correlations exert between lncRNAs and mRNAs (known as in cis) [61]. LncRNAs also can change the expression of distant mRNAs through the pairing of lncRNAs-mRNA (known as in trans). Bioinformatics analysis of lncRNAs target genes in trans showed that these genes played an important role in some pathways, such as MAPK, PI3K-Akt, prolactin, NF-kappa B, and Toll-like receptor signaling pathways, which played important roles during mammary gland development and lactation [4, 8, 62]. Consequently, many lncRNAs might function through targeting mRNA which played important roles in mammary gland from non-lactation to lactation period. For example, lncRNA TCONS_00075230 and IL1B had an opposite expression trend between dry and lactation period, and IL1B had reported to inhibit the proliferation of normal human mammary epithelial cells function as a retinoic acid induced gene [48]. So TCONS_00075230 might involve in lactation process through altering the expression of IL1B.

Regulating role of lncRNA-miRNA-mRNA network

It was known that the regulatory networks of lncRNAs, miRNA, and ceRNAs communicated with each other to regulate gene expression [28]. LncRNA SNHG7 could accelerate prostate cancer proliferation and cell cycle progression through cyclin D1 by sponging miR-503 [63]. The axis of lncRNA H19-miR-675-TGFBI had possible diagnostic and therapeutic potential for advanced prostate cancer [64]. LncRNA APF could control the expression of miR-188-3p, which suppressed myocardial infarction and autophagy by targeting ATG7 [65]. A series of miRNAs had been demonstrated to regulate lactation-associated processes, including lactation cycle, milk fat accumulation, hormone receptor activity, mammary gland involution and development, and lactation activity of mammary epithelial cells [14, 66, 67]. In this study, some candidate lncRNAs, including TCONS_00040268, TCONS_00137654, TCONS_00071659, and TCONS_00000352, could interact with lactation-related miR-221, and miR-221 could target ErBb3 gene. It had known that ErBb3 gene could drive mammary epithelial survival and differentiation during pregnancy and lactation [11]. Therefore, these lncRNAs might play an important role in regulating lactation, which would provide a new insight into lactation process of cattle.

Conclusion

In this study, 3746 significantly dysregulated lncRNA transcripts were found in lactation period compared with that in dry period, including 2732 down- and 1014 up-regulated lncRNA transcripts (fold change≥2, P < 0.05). Functional enrichment analysis of target genes and interacted miRNAs prediction of 34 lncRNAs indicated these lncRNAs might be important regulators for lactation in dairy cattle. This study would provide a resource for bovine lncRNA study involving in lactation biology.

Abbreviations

BMECs:

Bovine mammary epithelial cells

ceRNAs:

Competing endogenous RNAs

DEGs:

Differentially expressed genes

DELs:

Differentially expressed lncRNAs

DMEM:

Dubecco’s modified Eagle medium

lncRNAs:

Long non-coding RNAs

miRNAs:

MicroRNAs

qRT-PCR:

Quantitative real time PCR

References

  1. Watson CJ, Khaled WT. Mammary development in the embryo and adult: a journey of morphogenesis and commitment. Development. 2008;135(6):995–1003.

    Article  PubMed  CAS  Google Scholar 

  2. Topper YJ, Freeman CS. Multiple hormone interactions in the developmental biology of the mammary gland. Physiol Rev. 1980;60(4):1049–106.

    Article  PubMed  CAS  Google Scholar 

  3. Liu X, Robinson GW, Wagner KU, Garrett L, Wynshawboris A, Hennighausen L. Stat5a is mandatory for adult mammary gland development and lactogenesis. Genes Dev. 1997;11(2):179–86.

    Article  PubMed  CAS  Google Scholar 

  4. Rädler PD, Wehde BL, Wagner KU. Crosstalk between STAT5 activation and PI3K/AKT functions in normal and transformed mammary epithelial cells. Mol Cell Endocrinol. 2017;451(8):31–9.

    Article  PubMed  CAS  PubMed Central  Google Scholar 

  5. Furth PA, Nakles RE, Millman S, Diaz-Cruz ES, Cabrera MC. Signal transducer and activator of transcription 5 as a key signalling pathway in normal mammary gland developmental biology and breast cancer. Breast Cancer Res. 2011;13(5):220.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  6. Ersahin T, Tuncbag N, Cetin-Atalay R. The PI3K/AKT/mTOR interactive pathway. Mol BioSyst. 2015;11(7):1946.

    Article  PubMed  CAS  Google Scholar 

  7. Sobolewska A, Gajewska M, Zarzyńska J, Gajkowska B, Motyl T. IGF-I, EGF, and sex steroids regulate autophagy in bovine mammary epithelial cells via the mTOR pathway. Eur J Cell Biol. 2009;88(2):117–30.

    Article  PubMed  CAS  Google Scholar 

  8. Schmidt JW, Wehde BL, Sakamoto K, Triplett AA, Anderson SM, Tsichilis PN, Leone G, Wagner KU. Stat5 regulates the phosphatidylinositol 3-kinase/Akt1 pathway during mammary gland development and tumorigenesis. Mol Cell Biol. 2004;34(7):1363–77.

    Article  CAS  Google Scholar 

  9. Dai WT, Zou YX, White RR, Liu JX, Liu HY. Transcriptomic profiles of the bovine mammary gland during lactation and the dry period. Funct Integr Genomics. 2018;18(2):125–40.

    Article  PubMed  CAS  Google Scholar 

  10. Li Z, Liu HY, Jin XL, Lo LJ, Liu JX. Expression profiles of microRNAs from lactating and non-lactating bovine mammary glands and identification of miRNA related to lactation. BMC Genomics. 2012;13(1):731.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  11. Williams MM, Vaught DB, Joly MM, Hicks DJ, Sanchez V, Owens P, Rahman B, Elion DL, Balko JM, Cook RS. ErbB3 drives mammary epithelial survival and differentiation during pregnancy and lactation. Breast Cancer Res. 2017;19(1):105.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Bionaz M, Loor JJ. Gene networks driving bovine milk fat synthesis during the lactation cycle. BMC Genomics. 2008;9(1):366.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  13. Bionaz M, Loor JJ. Gene networks driving bovine milk protein synthesis during the lactation cycle. Bioinform Biol Insights. 2011;5:83–98.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  14. Do DN, Li R, Dudemaine PL, Ibeagha-Awemu EM. MicroRNA roles in signaling during lactation: an insight form differential expression, time course and pathway analyses of deep sequence data. Sci Rep. 2017;7:44605.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  15. Tang KQ, Wang YN, Zan LS, Yang WC. miR-27a controls triacylglycerol synthesis in bovine mammary epithelial cells by targeting peroxisome proliferator-activated receptor gamma. J Dairy Sci. 2017;100(5):4102–12.

    Article  PubMed  CAS  Google Scholar 

  16. Lee MJ, Yoon KS, Cho KW, Kim KS, Jung HS. Expression of miR-206 during the initiation of mammary gland development. Cell Tissue Res. 2013;353(3):425–33.

    Article  PubMed  CAS  Google Scholar 

  17. Li D, Xie XJ, Wang J, Bian YJ, Li QZ, Gao XJ, Wang CM. miR-486 regulates lactation and targets the PTEN gene in cow mammary glands. PLoS One. 2015;10(3):e0118284.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  18. Wang M, Moisá S, Khan MJ, Wang J, Bu D, Loor JJ. MicroRNA expression patterns in the bovine mammary gland are affected by stage of lactation. J Dairy Sci. 2012;95(11):6529–35.

    Article  PubMed  CAS  Google Scholar 

  19. Chu MQ, Zhao Y, Yu S, Hao YN, Zhang PF, Feng YN, Zhang HF, Ma DX, Liu J, Cheng M, et al. MicroRNA-221 may be imvolved in lipid metabolism in mammary epithelial cell. Intl J Biochem Cell Biol. 2018;97:118–27.

    Article  CAS  Google Scholar 

  20. Ucar A, Vafaizadeh V, Jarry H, Fiedler J, Klemmt PAB, Thum T, Groner B, Chowdhury K. miR-212 and miR-132 are required for epithelial stromal interactions necessary for mouse mammary gland development. Nat Genet. 2010;42(12):1101–8.

    Article  PubMed  CAS  Google Scholar 

  21. Guttman M, Rinn JL. Modular regulatory principles of large non-coding RNAs. Nature. 2012;482(7385):339–46.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  22. Shore AN, Rosen JM. Regulation of mammary epithelial cell homeostasis by lncRNAs. Int J Biochem Cell Biol. 2014;54:318–30.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  23. Venkatraman A, He XC, Thorvaldsen JL, Sugimura R, Perry JM, Tao F, Zhao M, Christenson MK, Sanchez R, Yu JY, et al. Maternal imprinting at the H19-Igf2 locus maintains adult haematopoietic stem cell quiescence. Nature. 2013;500(7462):345–9.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  24. Ginger MR, Shore AN, Contreras A, Rijnkels M, Miller J, Gonzalez-Rimbau MF, Rosen JM. A noncoding RNA is a potential marker of cell fate during mammary gland development. Proc Natl Acad Sci U S A. 2006;103(15):5781–6.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  25. Askarian-Amiri ME, Crawford J, French JD, Smart CE, Smith MA, Clark MB, Ru K, Mercer TR, Thompson ER, Lakhani SR, et al. SNORD-host RNA Zfas1 is a regulator of mammary development and a potential marker for breast cancer. RNA. 2011;17(5):878–91.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  26. Standaert L, Adriaens C, Radaelli E, Van Keymeulen A, Blanpain C, Hirose T, Nakagawa S, Marine JC. The long noncoding RNA Neat1 is required for mammary gland development and lactation. RNA. 2014;20(12):1844–9.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  27. Cesana M, Cacchiarelli D, Legnini I, Santini T, Sthandier O, Chinappi M, Tramontano A, Bozzoni I. A long noncoding RNA controls muscle differentiation by functioning as a competing endogenous RNA. Cell. 2011;147:358–69.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  28. Yu S, Zhao Y, Lai F, Chu M, Hao Y, Feng Y, Zhang H, Liu J, Cheng M, Li L, Shen W, Min L. LncRNA as ceRNAs may be involved in lactation process. Oncotarget. 2017;8(58):98014–28.

    PubMed  PubMed Central  Google Scholar 

  29. Shore AN, Kabotyanski EB, Roarty K, Smith MA, Zhang Y, Creighton CJ, Dinger ME, Rosen JM. Pregnancy-induced noncoding RNA (PINC) associates with polycomb repressive complex 2 and regulates mammary epithelial differentiation. PLoS Genet. 2012;8(7):e1002840.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  30. Langmead RBB, Trapnell C, Pop M, Salzberg LS. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10(3):R25.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  31. Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009;25(9):1105–11.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  32. Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL, Pachter L. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc. 2012;7(3):562–78.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  33. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28:511–5.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  34. Storey JD. The positive false discovery rate: a Bayesian interpretation and the q-value. Ann Stat. 2003;31(6):2013–35.

    Article  Google Scholar 

  35. Wang SH, Ge W, Luo ZX, Guo Y, Jiao BL, Qu L, Zhang ZY, Wang X. Integrated analysis of coding genes and non-coding RNAs during hair follicle cycle of cashmere goat (Capra hircus). BMC Genomics. 2017;18(1):767.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  36. Kong L, Zhang Y, Ye ZQ, Liu XQ, Zhao SQ, Wei L, Gao G. CPC: assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007;35:W345–9.

    Article  PubMed  PubMed Central  Google Scholar 

  37. Punta M, Coggill PC, Eberhardt RY, 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(D1):D290–301.

    Article  PubMed  CAS  Google Scholar 

  38. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(T)(−Delta Delta C) method. Methods. 2001;25(4):402–8.

    Article  PubMed  CAS  Google Scholar 

  39. Park S, Sung Y, Jeong J, Choi M, Lee J, Kwon W, Jang S, Park SJ, Kim HS, Lee MH, et al. hMAGEA2 promotes progression of breast cancer by regulating Akt and Erk1/2 pathways. Oncotarget. 2017;8(23):37115–27.

    PubMed  PubMed Central  Google Scholar 

  40. Hynes NE, Watson CJ. Mammary gland growth factors: roles in normal development and in cancer. Cold Spring Harb Perspect Biol. 2010;2(8):a003186.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  41. Streuli CH. Cell adhesion in mammary gland biology and neoplasia - preface. J Mammary Gland Biol. 2003;8(4):375–81.

    Article  Google Scholar 

  42. Knight CH. The importance of cell division in udder development and lactation. Livest Prod Sci. 2000;66(2):169–76.

    Article  Google Scholar 

  43. Samant GV, Sylvester PW. Gamma-Tocotrienol inhibits ErbB3-dependent PI3K/Akt mitogenic signalling in neoplastic mammary epithelial cells. Cell Proli. 2006;39(6):563–30.

    Article  CAS  Google Scholar 

  44. Yee LD, Guo Y, Bradbury J, Suster S, Clinton SK, Seewaldt VL. The antiproliferative effects of PPAR gamma ligands in normal human mammary epithelial cells. Breast Cancer Res Treat. 2003;78:179–92.

    Article  PubMed  CAS  Google Scholar 

  45. Frey RS, Singletary KW. Genistein activates p38 mitogen-activated protein kinase, inactivates ERK1/ERK2 and decreases Cdc25C expression in immortalized human mammary epithelial cells. J Nutr. 2003;133(1):226–31.

    Article  PubMed  CAS  Google Scholar 

  46. Hojilla CV, Jackson HW, Khokha R. TIMP3 regulates mammary epithelial apoptosis with immune cell recruitment through differential TNF dependence. PLoS One. 2011;6(10):e26718.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  47. Kessenbrock K, Smith P, Steenbeek SC, Pervolarakis N, Kumar R, Minami Y, Goga A, Hinck L, Werb Z. Diverse regulation of mammary epithelial growth and branching morphogenesis through noncanonical Wnt signaling. Proc Natl Acad Sci U S A. 2017;114(12):3121–6.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  48. Liu LM, Gudas LJ. Retinoic acid induces expression of the interleukin-1b gene in cultured normal human mammary epithelial cells and in human breast carcinoma lines. J Cell Physiol. 2002;193:244–52.

    Article  PubMed  CAS  Google Scholar 

  49. Arun SJ, Thomson PC, Sheehy PA, Khatkar MS, Raadsma HW, Williamson P. Targeted analysis reveals an important role of JAK-STAT-SOCS genes for milk production traits in Australian dairy cattle. Front Genet. 2015;6(187):342.

    PubMed  PubMed Central  Google Scholar 

  50. Toerien CA, Trout DR, Cant JP. Nutritional stimulation of milk protein yield of cows is associated with changes in phosphorylation of mammary eukaryotic initiation factor 2 and ribosomal s6 kinase 1. J Nutr. 2010;140(2):285–92.

    Article  PubMed  CAS  Google Scholar 

  51. O'Riordan N, Kane M, Joshi L, Hickey RM. Structural and functional characteristics of bovine milk protein glycosylation. Glycobiology. 2014;24(3):220–36.

    Article  PubMed  CAS  Google Scholar 

  52. Jurek B, Slattery DA, Maloumby R, Hillerer K, Koszinowski S, Neumann ID, van den Burg EH. Differential contribution of hypothalamic MAPK activity to anxiety-like behaviour in virgin and lactating rats. PLoS One. 2012;7(5):e37060.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  53. Yang Y, Gu CY, Luo C, Li F, Wang M. BUB1B promotes multiple myeloma cell proliferation through CDC20/CCNB axis. Med Oncol. 2015;32(3):81.

    Article  PubMed  CAS  Google Scholar 

  54. Cooper CR, Szaniszlo PJ. Evidence for 2 cell-division cycle (CDC) genes that govern yeast bud emergence in the pathogenic fungus wangiella-dermatitidis. Infect Immun. 1993;61(5):2069–81.

    PubMed  PubMed Central  CAS  Google Scholar 

  55. Hansji H, Leung EY, Baguley BC, Finlay GJ, Askarian-Amiri ME. Keeping abreast with long non-coding RNAs in mammary gland development and breast cancer. Front Genet. 2014;5:379.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  56. Huang J, Zhou N, Watabe K, Lu Z, Wu F, Xu M, Mo YY. Long non-coding RNA UCA1 promotes breast tumor growth by suppression of p27 (Kip1). Cell Death Dis. 2014;5:e1008.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  57. Zhang YF, Wagner EK, Guo XY, May I, Cai QY, Zheng W, He CY, Long JR. Long intergenic non-coding RNA expression signature in human breast cancer. Sci Rep. 2016;6:37821.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  58. Koufariotis LT, Chen YPP, Chamberlain A, Jagt CV, Hayes BJ. A catalogue of novel bovine long noncoding RNA across 18 tissues. PLoS One. 2015;10(10):e0141225.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  59. Tong C, Chen QL, Zhao LL, Ma JF, Ibeagha-Awemu EM, Zhao X. Identification and characterization of long intergenic noncoding RNAs in bovine mammary glands. BMC Genomics. 2017;18:468.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  60. Bach K, Pensa S, Grzelak M, Hadfield J, Adams DJ, Marioni JC, Khaled WT. Differentiation dynamic of mammary epithelial cells revealed by single-cell RNA sequencing. Nat Commun. 2017;8(1):2128.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  61. Guttman M, Amit I, Garber M, French C, Lin MF, Feldser D, Huarte M, Zuk O, Carey BW, Cassady JP, et al. Chromatin signature reveals over a thousand highly conserved large non-coding RNAs in mammals. Nature. 2009;458(7235):223–7.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  62. Garcia MC, Lopez M, Gualillo O, Seoane LM, Dieguez C, Senaris RM. Hypothalamic levels of NPY, MCH, and prepro-orexin mRNA during pregnancy and lactation in the rat: role of prolactin. FASEB J. 2003;17(11):1392–400.

    Article  PubMed  CAS  Google Scholar 

  63. Qi H, Wen B, Wu Q, Cheng W, Lou J, Wei J, Huang J, Yao X, Weng G. Long noncoding RNA SNHG7 accelerates prostate cancer proliferation and cycle progression through cyclin D1 by sponging miR-503. Biomed Phamacother. 2018;102:326–32.

    Article  CAS  Google Scholar 

  64. Zhu MJ, Chen Q, Liu X, Sun Q, Zhao X, Deng R, Wang YL, Huang J, Xu M, Yan JS, et al. lncRNA H19/miR-675 axis represses prostate cancer metastasis by targeting TGFBI. FEBS J. 2014;281(16):3766–75.

    Article  PubMed  CAS  Google Scholar 

  65. Wang K, Liu CY, Zhou LY, Wang JX, Wang M, Zhao B, Zhao WK, Xu SJ, Fan LH, Zhang XJ, et al. APF lncRNA regulates autophagy and myocardial infarction by targeting miR-188-3p. Nat Commun. 2015;6:6779.

    Article  PubMed  CAS  Google Scholar 

  66. Lin XZ, Luo J, Zhang LP, Wang W, Gou DM. MiR-103 controls milk fat accumulation in goat (Capra hircus) mammary gland during lactation. PLoS One. 2013;8(11):e79258.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  67. Feuermann Y, Kang K, Shamay A, Robinson GW, Hennighausen L. MiR-21 is under control of STAT5 but is dispensable for mammary development and lactation. PLoS One. 2014;9(1):e85123.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

Download references

Funding

This study was supported by the Program for New Excellent Talents in University of China (NCET-13-0485), National natural Science Foundation of China (31772573) and Agricultural Scientific and Technological Innovation Project of Shandong Academy of Agricultural Science (CXGC2016B03).

Availability of data and materials

The RNA-seq datasets supporting the conclusions of this article are available in the NCBI SRA repository (https://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/Traces/sra_sub/). BioProject asscession: PRJNA482783. BioSamples: SAMN09714469 and SAMN09714468.

The datasets supporting the conclusions of this article are included within the article and its additional files.

Author information

Authors and Affiliations

Authors

Contributions

WX was responsible for the concept and design of this study and manuscript review. YB and JBL contributed equally to this work, and they performed statistical analyses and writing the manuscript. GW and ZXL performed data visualization. WSH and ZHB conducted surgical processes for tissue collection. All authors approved the final manuscript.

Corresponding author

Correspondence to Xin Wang.

Ethics declarations

Ethics approval

All the experimental procedures and animal care with cows involved in this study were approved by the Experimental Animal Manage Committee of Northwest A&F University (China).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

qRT-PCR primers for lncRNAs and coding genes. (DOCX 21 kb)

Additional file 2:

The expression profile of DEGs in the dry and lactation periods. (XLSX 387 kb)

Additional file 3:

The expression profile of novel lncRNAs in the dry and lactation periods. (XLSX 980 kb)

Additional file 4:

The sequence information of novel lncRNAs. (XLSX 2011 kb)

Additional file 5:

The expression profile of annotated lncRNAs in the dry and lactation periods. (XLSX 1732 kb)

Additional file 6:

The expression file of novel DELs in the dry and lactation periods. (XLSX 230 kb)

Additional file 7:

The expression file of annotated DELs in the dry and lactation periods. (XLSX 224 kb)

Additional file 8:

Stage specific expressed lncRNAs in the dry or lactation periods. (XLSX 41 kb)

Additional file 9:

Stage specific expressed genes in the dry or lactation periods. (XLSX 21 kb)

Additional file 10:

Predicted target genes and mRNAs of DELs (in cis). (XLSX 24 kb)

Additional file 11:

Predicted target genes and mRNAs of DELs (in trans). (XLSX 84 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Yang, B., Jiao, B., Ge, W. et al. Transcriptome sequencing to detect the potential role of long non-coding RNAs in bovine mammary gland during the dry and lactation period. BMC Genomics 19, 605 (2018). https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-018-4974-5

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-018-4974-5

Keywords