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

Screening and evaluating of long noncoding RNAs in the puberty of goats



Long noncoding RNAs (lncRNAs) are involved in regulating animal development, however, their function in the onset of puberty in goats remain largely unexplored. To identify the genes controlling the regulation of puberty in goats, we measured lncRNA and mRNA expression levels from the hypothalamus.


We applied RNA sequencing analysis to examine the hypothalamus of pubertal (case; n = 3) and prepubertal (control; n = 3) goats. Our results showed 2943 predicted lncRNAs, including 2012 differentially expressed lncRNAs, which corresponded to 5412 target genes. We also investigated the role of lncRNAs that act cis and trans to the target genes and found a number of lncRNAs involved in the regulation of puberty and reproduction, as well as several pathways related to these processes. For example, oxytocin signaling pathway, sterol biosynthetic process, and pheromone receptor activity signaling pathway were enriched as Kyoto Encyclopedia of Genes and Genomes (KEGG) or gene ontology (GO) analyses showed.


Our results clearly demonstrate that lncRNAs play an important role in regulating puberty in goats. However, further research is needed to explore the functions of lncRNAs and their predicted targets to provide a detailed expression profile of lncRNAs on goat puberty.


Puberty is a pivotal stage in female goat development. It marks the first occurrence of ovulation and the onset of reproductive capability [1]. The mechanism of puberty onset is complex and thought to be associated with environmental factors, neuroendocrine factors, genetic factors and their interactions. In general, the secretion of gonadotropin-releasing hormone (GnRH) is considered a crucial factor in puberty onset for goats [2]. A popular view is that during the prepubertal period, secreting neurons suffer persistent trans-synaptic inhibition. This means that GnRH secretions increase as long as this inhibition is eliminated, which leads to puberty [2]. However, these influences are based on substantial genetic control [3].

It was previously reported that the initiation of puberty in female rats is regulated by epigenetic mechanism of transcriptional repression [4], whereby epigenetic control was composed of several mechanisms. Two well established mechanisms include: modification of chromatin and chemical modification of the DNA (including DNA methylation and hydroxymethylation). Non-coding RNA is the most recently unveiled mechanism of epigenetic control, which affords epigenetic information by lncRNAs or microRNAs [5]. Broadly, lncRNAs are known as transcripts greater than 200 nt in length that do not appear to code proteins [6].

During the past decades of transcriptome studies, multiple lncRNAs have been discovered, such as Xist and H19. The advent of RNA-seq has been a powerful tool in exploring and quantifying lncRNAs [7], which has led to the identification of many more lncRNAs that await functional validation. Most identified lncRNAs have primarily originated from human and mouse studies [8, 9]. Recent studies in bovine [1012] and porcine species [13, 14] have enriched the mammalian lncRNAs databases, providing a promising future for further lncRNAs studies.

LncRNAs have been shown to participate in the regulation of transcriptional and post-transcriptional control [15]. In recent years, lncRNAs have proven to play roles in lactation, ovary development, and embryo and sperm maturation. Therefore, we inferred that the onset of goat puberty is also regulated by lncRNAs. In this study, we applied RNA-seq and investigated the expression profiles of mRNA and lncRNAs in pubertal and prepubertal goats to explore the association of lncRNAs with the onset of puberty.


Preparation of animals and tissues

This study was authorized and endorsed by the Animal Care and Use Committee of Anhui Agricultural University. We housed three, prepubertal (aged 2.5-3.0 months) and three, pubertal (aged 4.5-5.0 months) female Anhuai goats under the same conditions on a farm in Anhui Province, China. We determined puberty goats in studied femal by male goats detecting estrus and the appearance changes of vulva [16]. The average weight of pubertal goats was 16.17 kg compared with the prepubertal goats 8.30 kg, and the average weight of the pubertal goats’ ovary was 0.76 g compared with the prepubertal goats 0.30 g. The animals were deeply anesthetized by intravenous administration of 3% pentobarbital sodium (30 mg/kg; Solarbio, P8410, China) and sacrificed by exsanguination in a healthy physiological stage when pubertal goats were in the late follicular phase. Hypothalamic tissues were surgically removed, and frozen in liquid nitrogen immediately. These tissues were stored at −80 °C until the RNA extraction [17].

RNA sequencing and quality control

We isolated total RNA from goat hypothalamus using TRIzol Reagent (Invitrogen, Carlsbad, CA, USA), according to the standard extraction protocol. The contamination and degradation of RNA was detected by 1% agarose gels. The purity of RNA was monitored using the NanoPhotometer® spectrophotometer (Implen, Los Angeles, CA, USA). We measured the concentration of RNA using Qubit® RNA Assay Kit in Qubit® 2.0 Flurometer (Life Technologies, Carlsbad, CA, USA). The integrity of RNA was monitored as previous reported [18]. We used 3 μg RNA per sample for the RNA sample preparations. Firstly, ribosomal RNA in total RNA was removed [19], and then the residue was cleaned up by using ethanol precipitation. Then the libraries whith high strand-specificity for sequencing was generated [19], following manufacturer’s recommendations. Then the process was followed as previously described [20]. Illumina Hiseq 4000 platform was adopted on sequencing and 150 bp paired-end reads were generated. Raw reads were dealt with in-house perl scripts. The reads with more than 10% unknown bases, reads containing adapter and reads with more than 50% of low-quality bases (whose Phred scores were < 5%) were removed, yielding only the clean reads. Meanwhile, the quality of clean reads (Q20, Q30, and GC content) were detected. All the following analyses were based on high quality clean reads.

Transcriptome assembly

We used a GTF file ( with the annotation of the goat genome. Index of the reference genome was created by Bowtie v2.0.6 [21, 22] and then we aligned paired-end clean reads to the reference genome using TopHat v2.0.9 [23]. The mapped reads of each sample were assembled by both Scripture (beta2) [24] and Cufflinks (v2.1.1) [25, 26] in a reference-based approach. Both methods determined exons connectivity by spliced reads. Scripture ran using default parameters, while Cufflinks ran with min-frags-per-transfrag = 0’ and–library-type fr-firststrand’. Other parameters were set as default.

Expression and coding potential analysis of transcripts

Gene expression was calculated using FPKMs of transcripts in each sample [27]. We confirmed differential expression in gene expression data using Cuffdiff as it based on the negative binomial distribution provides statistical routines [25]. Transcripts with a P < 0.05 were assigned as significantly differentially expressed between two groups. We used three analytic tools, including Coding-Non-Coding-Index (CNCI; v2) [28], Coding Potential Calculator (CPC; 0.9-r2) [29], Pfam Scan (v1.3) [30] to screen out candidate lncRNAs. CNCI (v2) profiles distinguished protein-coding and non-coding sequences effectively by adjoining nucleotide triplets, which was independent of known annotations. CPC (0.9-r2) was mainly used to detect the extent and quality of the Open Reading Frames (ORF) in a transcript and discover the sequences in known protein database, clarifying the coding and non-coding transcripts. Each transcript was translated in all three possible frames, then any of the known protein family was identified by Pfam Scan (v1.3) in the Pfam database (release 27; adopted both Pfam A and Pfam B). The coding potential of transcripts predicted by any of the three tools above were filtered out (non-annotated transcriptional activity by identifying novel transcripts), and those without coding potential were our candidate lncRNAs for further analysis.

Target gene prediction and functional enrichment analysis

The cis role refers to the lncRNA acting on neighboring target genes [31, 32]. To predict the cis-regulated target genes of lncRNAs, we screened protein-coding genes as potential targets 10 K/100 K upstream and downstream of lncRNAs andanalyzed their function. The trans role refers to the coexpression relationship between lncRNAs and mRNA. Expression levels of lncRNAs and mRNAs were calculated for Pearson’s correlation coefficients by custom scripts (r > 0.95 or r < −0.95). The target genes of lncRNAs were performed functional enrichment analysis by clustering the genes from various samples using the DAVID platform [33]. The significance was described as a P-value, measured by the EASE score (P < 0.05 was considered significant).

Quantitative real-time PCR

We have validated the RNA-seq data by selected eight lncRNAs, two novel transcripts, and two target genes to investigate the expression patterns in the samples using qRT-PCR. Werepeated the qRT-PCR experiments three times per sample on expression from the same hypothalamic tissues of three pre vs. three pubertal goats. We designed primers online using Primer5 software and evaluated using BLAST at NCBI. A list of the primer sequences is shown in Table 1. We performed qRT-PCR using SYBR green (Vazyme, China) method. Expression levels of genes were quantified through the cycle threshold (Ct) values and evaluated as 2-ΔΔCT. The data of expression was normalized to β-action.

Table 1 qRT-PCR primer and size of the amplification products of the target and housekeeping genes

GO and pathway analysis

In this study, Gene Ontology (GO) enrichment [34] analysis of targets was performed by the GOseq R package and corrected by P (P < 0.05 were considered significantly enriched). Pathway analysis is a functional analysis in KEGG ( pathways [35]. We evaluated the statistical enrichment of lncRNAs target genes or differential expression genes in KEGG pathways using KOBAS [36] software.

Statistical analysis

We performed further analysis of RNA-seq data and graphical representations using the statistical R package (R, Auckland, NZL), adopting multiple testing and P corrections. We applied SPSS 17.0 software package (SPSS, Chicago, IL, USA) to analyze the qRT-PCR data. Differential expression levels of genes were calculated by independent-samples t-test between prepubertal and pubertal goats. Significance of data was defined as P < 0.05.


Identification of lncRNAs

We used pubertal and prepubertal goats to perform RNA-seq analysis from the hypothalamus of six female Anhuai goats. In total, 774,998,560 raw reads were produced under the Illumina HiSeq 4000 platform. We obtained 636,544,196 reads maped to goat reference genome after discarding low-quality sequences and adaptor sequences. The percentage of mapped reads among clean reads in each library ranged from 80.45% - 84.32% (Additional file 1). After the analysis of coding potential using the CNCI, CPC and Pfam-scan software, we identified 2943 lncRNAs (Fig. 1), including 2426 large intergenic noncoding RNAs, 217 anti-sense_lncRNAs, and 300 intronic_lncRNAs.

Fig. 1
figure 1

Screening of candidate lncRNAs in hypothalamus transcriptome. The coding potential of lncRNAs were analyzed by three tools (CPC, CNCI and PFAM)

Genomic features of lncRNAs

Overall, we observed a lower expression of lncRNAs compared with mRNA (Fig. 2a). The mean length of lncRNAs in our dataset was 1180 nt, and the mean mRNA length was 2869 nt (Fig. 2b). Furthermore, we detected an ORF mean length of 105 nt for our lncRNAs, which tended to be shorter than protein coding genes (Fig. 2c). We also found lncRNAs contained fewer exons than mRNA (Fig. 2d).

Fig. 2
figure 2

The comparison of features between predicted lncRNAs and mRNA. a Expression of lncRNAs and mRNA. b Length distribution of 2943 predicted lncRNAs and 30162 coding transcripts. c ORF length distribution of lncRNAs and coding transcripts. d Exon number distribution of lncRNAs and coding transcripts

Differential expression cluster analysis

Further analysis identified 1165 significant differential expression transcripts (including lncRNAs, mRNAs, and novel transcripts) (P < 0.05), 770 up-regulated and 395 down-regulated transcripts (Fig. 3), including 57 novel transcripts. Furthermore, we detected 59 lncRNAs transcripts from 58 lncRNAs gene loci significant differentially expressed lncRNAs (P < 0.05), including 29 up-regulated and 30 down-regulated lncRNAs transcripts in pubertal samples compared with prepubertal samples (Additional file 2). We validated sequencing results using qRT-PCR analysis (Fig. 4).

Fig. 3
figure 3

Volcano plots of differential expression transcripts. X-axis is fold change (log 2) and Y-axis is P (−log 10). Red points indicate up-regulated (X axis > 0) transcripts; green points indicate down-regulated (X axis < 0) transcripts

Fig. 4
figure 4

Validation of RNA-seq results by using quantitative qRT-PCR. Some lncRNAs and target genes were examined using quantitative qRT-PCR. The data are expressed as the mean ± 1 SD (n = 3). *p < 0.05, **p < 0.01

Prediction of target genes of lncRNAs in cis and trans

LncRNAs can act on target genes, either in cis (neighbor the site of lncRNA production) or in trans to coexpression whith target genes [37]. To explore whether differences in lncRNAs affects functional regulation of goat puberty, we predicted the target genes of lncRNA using the cis and trans model. To analyze the cis role of lncRNA, we screened protein-coding genes as potential targets 10 K/100 K upstream and downstream of the lncRNAs. The results indicated that there were 2012 lncRNAs that corresponded to 5412 target genes (Additional file 3). Interestingly, we observed several genes related puberty such as PRLHR, EMC3, IGF2BP1, ZNF175, and ZNF444 [3841], which were respectively a target of XLOC_1486935, XLOC_1284300, XLOC_957527, XLOC_080674, and XLOC_910648, indicating that the onset of puberty is probably regulated by the lncRNA-tatget genes. Regarding the trans role of lncRNA, our results showed that the lncRNA, XLOC_957527, acted on GnRH1 (Table 2; Additional file 4).

Table 2 LncRNAs and its potential target genes associated with puberty

GO and KEGG analysis

Our GO analysis of predicted targets demonstrated 73 significantly enriched terms (P < 0.05). The top eight terms were as follows: pheromone receptor activity, hyalurononglucosaminidase activity, hexosaminidase activity, sensory perception of taste, viral genome packaging, helicase activity, sensory perception of chemical stimulus, and sensory perception (Additional file 5). Interestingly, the signaling pathway of the pheromone receptor activity was significantly enriched, which relates to goat estrus. In addition, the sterol biosynthetic process signaling pathway was significantly enriched. DNAJB2 was the differentially expressed target gene on the pathway, which suggests that it may be a new gene involved in the regulation of puberty onset via the sterol biosynthetic process signaling pathway.

KEGG pathway analysis of lncRNAs targets showed 90 terms were enriched (Additional file 6), in which oxytocin signaling pathway was related to puberty [42]. These results suggested that lncRNAs may be cis-acting on its target genes to regulate onset of puberty (Fig. 5).

Fig. 5
figure 5

KEGG annotation for target gene functions of predicated lncRNAs. Red indicates higher expression and green indicates lower expression. The number of differentially expressed genes is shown in parentheses

We also evaluated the trans role of 2943 lncRNAs in protein-coding genes by its correlation coefficient of gene expression (Pearson correlation ≥ 0.95 or ≤ −0.95). The results showed that 2551 lncRNAs had interactions with target genes in trans of the goat genome. Functional analysis illustrated that target genes in trans were enriched (P < 0.05) in 158 GO terms including a variety of processes (Additional file 7), such as G-protein coupled receptor activity, transmembrane signaling receptor activity, receptor activity, and so on.

We identified 273 KEGG pathways (Additional file 8), several of which were associated with puberty, including ovarian steroidogenesis, GnRH signaling pathway, steroid biosynthesis, steroid hormone biosynthesis, oxytocin signaling pathway, GABAergic synapse, estrogen signaling pathway, oocyte meiosis, glutamatergic synapse, and others. These findings indicate that lncRNAs may act on the target genes associated with puberty of goat in trans.

Specific expression of lncRNAs

There were 187 specific expressions of lncRNAs in the pubertal samples, especially, XLOC_2409732, which has a lower P than other specific expression of lncRNAs. The targets of XLOC_2409732 were detected as ASB5, WDR17, SPATA4 and SPCS3 according to RNA-seq analysis. We found 243 specific expressions of lncRNAs in prepubertal samples, particularly XLOC_1498149, which has a lower P than other specific expression of lncRNAs, and CDR1 was targeted to XLOC_1498149 (Additional file 9). These two specifically expressed lncRNAs may play a pivotal role in goat puberty and, with further studies, provide crucial information regarding the regulation of puberty.


We initially performed RNA-seq to analyze lncRNAs of hypothalamus from pubertal and prepubertal female Anhuai goats. Through sequencing, we acquired 2943 predicted lncRNAs and 30162 coding transcripts. Many studies have indicated that lncRNAs have unique features compared with mRNA; for example, lncRNAs are shorter in length than protein-coding transcripts [27]. Furthermore, we found that lncRNAs in hypothalamus were shorter than in skin (1809 bp on average); however, the number of exons is similar [20]. Interestingly, the length of lncRNAs in goat hypothalamus are longer than that in human (1 kb on average) and mouse (550 nt on average), containing fewer exons than human (2.9 exons on average) and mouse (3.7 exons on average) [43].

In this study, we screened out significant differentially expressed transcripts, including 59 lncRNA transcripts from 58 lncRNA gene loci. As previous research, the functions of lncRNAs were reflected by acting on the protein-coding genes. For example, in a recent study, a muscle-specific lncRNA, linc-MD1, influenced muscle development by targeting to MAML1 [44]. Moreover, the lncRNA, Neat1, could make a difference in pregnancy by acting on corpus luteum formation in mice [45]. Therefore, we could predict the role of mammalian lncRNAs by the relevant protein-coding genes.

Here, we predicted the potential functions of lncRNAs through the protein-coding genes in cis and trans. Several genes have been confirmed to be associated with puberty onset, including Kiss1/GPR54 [4649], IGFs [50], GABA [51] and FSHR. We discovered several differentially expressed targets in cis and trans for lncRNAs in pubertal and prepubertal hypothalamus. PRLHR, EMC3, IGF2BP1, ZNF175, ZNF444 have been reported involved in the regulation of puberty and reproduction of animal [40, 52]. For example, puberty of female rats is significantly advanced by GnRH release under the stimulation of IGF-1; IGF-1 can affect the puberty-related events by hypothalamic GnRH release [53].

Moreover, previous research has demonstrated that puberty is delayed after over expression of ZNFs in the arcuate nucleus (ARC) of female rats; subsequent oestrous cyclicity is also disrupted [40]. Our results also showed the lncRNA, XLOC_957527, acted on GnRH1 through trans interactions. Consequently, we confirmed that relevant lncRNAs might play a crucial role on regulation of puberty via the above targets. However, these predicted functions of lncRNAs need further experimental verification.

In the present study, oxytocin signaling pathways were enriched in KEGG pathways. MEF2C, as the target gene of lncRNA XLOC_2123676, is an important gene in oxytocin signaling pathway associated with puberty regulation. The age that vaginal opening occurs in female rats is delayed by treatment with an oxytocin antagonist, indicating that oxytocin enhances sexual maturation [42]. We also observed that DNAJB2, the target of lncRNA XLOC_1101518, has a crucial role in sterol biosynthetic process, which is involved in regulation puberty. ESR1 is essential for multiple estrogen feedback loops and required for puberty onset in female mouse [54].

Our GO analysis of the predicted targets indicates that pheromone receptor activity signaling pathway, which relates to goat estrus, is significantly enriched (Fig. 6) [55].

Fig. 6
figure 6

GO enrichment analysis for target gene functions of predicated lncRNAs. (MF: molecular function)

In our study, the enriched KEGG pathways and GO pathways associated with reproduction and puberty clearly suggest that these lncRNAs play a vital role in regulation of puberty in goats. However, the functions of lncRNAs and their predicted targets analyses should be carefully evaluated by further experiments.


We performed RNA-seq analysis, and screened out differentially expressed lncRNAs of pubertal and prepubertal goats. We elucidated genomic differences between lncRNAs compared with mRNA. Then, we observed several target genes of lncRNAs related to puberty. Our results clearly demonstrate that lncRNAs play an important role in regulating puberty in goats.



ankyrin repeat and SOCS box containing 5


DnaJ (Hsp40) homolog, subfamily B, member 2


ER membrane protein complex subunit 3


follicle stimulating hormone receptor


gonadotropin-releasing hormone


insulin-like growth factor 2 mRNA binding protein 1


long noncoding RNA


myocyte enhancer factor 2C


prolactin releasing hormone receptor


WD repeat domain 17


zinc finger protein 175


zinc finger protein 444


  1. Cao GL, Feng T, Chu MX, Di R, Zhang YL, Huang DW, Liu QY, Hu WP, Wang XY. Subtraction suppressive hybridisation analysis of differentially expressed genes associated with puberty in the goat hypothalamus. Reprod Fertil Dev. 2015;28(11):1781–1787.

  2. Smith MJ, Jennes L. Neural signals that regulate GnRH neurones directly during the oestrous cycle. Reproduction. 2001;122:1–10.

    Article  CAS  PubMed  Google Scholar 

  3. Abreu AP, Macedo DB, Brito VN, Kaiser UB, Latronico AC. A new pathway in the control of the initiation of puberty: the MKRN3 gene. J Mol Endocrinol. 2015;54(3):R131–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  4. Lomniczi A, Loche A, Castellano JM, Ronnekleiv OK, Bosch M, Kaidar G, Knoll JG, Wright H, Pfeifer GP, Ojeda SR. Epigenetic control of female puberty. Nat Neurosci. 2013;16(3):281–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Lomniczi A, Wright H, Ojeda SR. Epigenetic regulation of female puberty. Front Neuroendocrinol. 2015;36:90–107.

    Article  CAS  PubMed  Google Scholar 

  6. Mattick JS, Rinn JL. Discovery and annotation of long noncoding RNAs. Nat Struct Mol Biol. 2015;22(1):5–7.

    Article  CAS  PubMed  Google Scholar 

  7. Atkinson SR, Marguerata S, Bähler J. Exploring long non-coding RNAs through sequencing. Semin Cell Dev Biol. 2012;23(2):200–5.

    Article  CAS  PubMed  Google Scholar 

  8. Volders PJ, Verheggen K, Menschaert G, Vandepoele K, Martens L, Vandesompele J, Mestdagh P. An update on LNCipedia: a database for annotated human lncRNA sequences. Nucleic Acids Res. 2015;43(Database issue):D174–180.

    Article  PubMed  Google Scholar 

  9. Quek XC, Thomson DW, Maag JL, Bartonicek N, Signal B, Clark MB, Gloss BS, Dinger ME. lncRNAdb v2.0: expanding the reference database for functional long noncoding RNAs. Nucleic Acids Res. 2015;43(Database issue):D168–173.

    Article  PubMed  Google Scholar 

  10. Huang W, Long N, Khatib H. Genome-wide identification and initial characterization of bovine long non-coding RNAs from EST data. Anim Genet. 2012;43(6):674–82.

    Article  CAS  PubMed  Google Scholar 

  11. Billerey C, Boussaha M, Esquerré D, Rebours E, Djari A, Meersseman C, Klopp C, Gautheret D, Rocha D. Identification of large intergenic non-coding RNAs in bovine muscle using next-generation transcriptomic sequencing. BMC Genomics. 2014;15:499.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Weikard R, Hadlich F, Kuehn C. Identification of novel transcripts and noncoding RNAs in bovine skin by deep next generation sequencing. BMC Genomics. 2013;14:789.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Zhou ZY, Li AM, Adeola AC, Liu YH, Irwin DM, Xie HB, Zhang YP. Genome-wide identification of long intergenic noncoding RNA genes and their potential association with domestication in pigs. Genome Biol Evol. 2014;6(6):1387–92.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Zhao W, Mu Y, Ma L, Wang C, Tang Z, Yang S, Zhou R, Hu X, Li MH, Li K. Systematic identification and characterization of long intergenic non-coding RNAs in fetal porcine skeletal muscle development. Sci Rep. 2015;5:8957.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Yang G, Lu X, Yuan L. LncRNA: A link between RNA and cancer. Biochim Biophys Acta. Biochim Biophys Acta. 2014;1839(11):1097–109.

    Article  CAS  PubMed  Google Scholar 

  16. Dantas A, Siqueira ER, Fernandes S, Oba E, Castilho AM, Meirelles PRL, Sartori MMP, Santos PTR. Influence of feeding differentiation on the age at onset of puberty in Brazilian Bergamasca dairy ewe lambs. Arq Bras Med Vet Zootec. 2016;68:22–8.

    Article  Google Scholar 

  17. Pan L, Ma J, Pan F, Zhao D. Gao J: long non-coding RNA expression profiling in aging rats with erectile dysfunction. Cell Physiol Biochem. 2015;37(4):1513–26.

    Article  CAS  PubMed  Google Scholar 

  18. Kiewe P, Gueller S, Komor M, Stroux A, Thiel E, Hofmann WK. Prediction of qualitative outcome of oligonucleotide microarray hybridization by measurement of RNA integrity using the 2100 Bioanalyzer capillary electrophoresis system. Ann Hematol. 2009;88(12):1177–83.

    Article  CAS  PubMed  Google Scholar 

  19. Li P, Conley A, Zhang H, Kim HL. Whole-Transcriptome profiling of formalin-fixed, paraffin-embedded renal cell carcinoma by RNA-seq. BMC Genomics. 2014;15:1087.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Ren H, Wang G, Chen L, Jiang J, Liu L, Li N, Zhao J, Sun X, Zhou P. Genome-wide analysis of long non-coding RNAs at early stage of skin pigmentation in goats (Capra hircus). BMC Genomics. 2016;17(1):67.

    Article  PubMed  PubMed Central  Google Scholar 

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

    Article  PubMed  PubMed Central  Google Scholar 

  22. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9(4):357–U354.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

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

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Guttman M, Garber M, Levin JZ, Donaghey J, Robinson J, Adiconis X, Fan L, Koziol MJ, Gnirke A, Nusbaum C, et al. Ab initio reconstruction of cell type-specific transcriptomes in mouse reveals the conserved multi-exonic structure of lincRNAs (vol 28, pg 503, 2010). Nat Biotechnol. 2010;28(7):756.

    Article  CAS  Google Scholar 

  25. 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(5):511–515.

  26. 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  CAS  PubMed  PubMed Central  Google Scholar 

  27. Flegel C, Manteniotis S, Osthold S, Hatt H, Gisselmann G. Expression profile of ectopic olfactory receptors determined by deep sequencing. PLoS One. 2013;8(2):e55368.

  28. Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, Liu Y, Chen R, Zhao Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013;41(17):e166.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. 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 

  30. Bateman A, Birney E, Durbin R, Eddy SR, Howe KL, Sonnhammer ELL. The Pfam Protein Families Database. Nucleic Acids Res. 2000;28(1):263–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Gomez JA, Wapinski OL, Yang YW, Bureau JF, Gopinath S, Monack DM, Chang HY, Brahic M, Kirkegaard K. The NeST Long ncRNA Controls Microbial Susceptibility and Epigenetic Activation of the Interferon-gamma Locus. Cell. 2013;152(4):743–54.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Lai F, Orom UA, Cesaroni M, Beringer M, Taatjes DJ, Blobel GA, Shiekhattar R. Activating RNAs associate with Mediator to enhance chromatin architecture and transcription. Nature. 2013;494(7438):497–501.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Huang DW, Sherman BT, Lempicki RA. Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists. Nucleic Acids Res. 2009;37(1):1–13.

    Article  Google Scholar 

  34. Young MD, Wakefield MJ, Smyth GK, Oshlack A. Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010;11(2):R14.

    Article  PubMed  PubMed Central  Google Scholar 

  35. Kanehisa M, Araki M, Goto S, Hattori M, Hirakawa M, Itoh M, Katayama T, Kawashima S, Okuda S, Tokimatsu T, et al. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 2008;36(Database issue):D480–484.

    CAS  PubMed  Google Scholar 

  36. Mao XZ, Cai T, Olyarchuk JG, Wei LP. Automated genome annotation and pathway identification using the KEGG Orthology (KO) as a controlled vocabulary. Bioinformatics. 2005;21(19):3787–93.

    Article  CAS  PubMed  Google Scholar 

  37. Wang Kevin C, Chang Howard Y. Molecular Mechanisms of Long Noncoding RNAs. Mol Cell. 2011;43(6):904–14.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Bachelot A, Carre N, Mialon O, Matelot M, Servel N, Monget P, Ahtiainen P, Huhtaniemi I, Binart N. The permissive role of prolactin as a regulator of luteinizing hormone action in the female mouse ovary and extragonadal tumorigenesis. Am J Physiol Endocrinol Metab. 2013;305(7):E845–52.

    Article  CAS  PubMed  Google Scholar 

  39. Greenwald-Yarnell ML, Marsh C, Allison MB, Patterson CM, Kasper C, MacKenzie A, Cravo R, Elias CF, Moenter SM, Myers Jr MG. ERalpha in Tac2 Neurons Regulates Puberty Onset in Female Mice. Endocrinology. 2016;157(4):1555–65.

    Article  CAS  PubMed  Google Scholar 

  40. Lomniczi A, Wright H, Castellano JM, Matagne V, Toro CA, Ramaswamy S, Plant TM, Ojeda SR. Epigenetic regulation of puberty via Zinc finger protein-mediated transcriptional repression. Nat Commun. 2015;6:10195.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Zhu H, Shah S, Shyh-Chang N, Shinoda G, Einhorn WS, Viswanathan SR, Takeuchi A, Grasemann C, Rinn JL, Lopez MF, et al. Lin28a transgenic mice manifest size and puberty phenotypes identified in human genetic association studies. Nat Genet. 2010;42(7):626–U106.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Parent A-S, Rasier G, Matagne V, Lomniczi A, Lebrethon M-C, Gérard A, Ojeda SR, Bourguignon J-P. Oxytocin Facilitates Female Sexual Maturation through a Glia-to-Neuron Signaling Pathway. Endocrinology. 2008;149(3):1358–65.

    Article  CAS  PubMed  Google Scholar 

  43. Cabili MN, Trapnell C, Goff L, Koziol M, Tazon-Vega B, Regev A, Rinn JL. Integrative annotation of human large intergenic noncoding RNAs reveals global properties and specific subclasses. Genes Dev. 2015;25(18):1915–27.

    Article  Google Scholar 

  44. 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(2):358–69.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Nakagawa S, Shimada M, Yanaka K, Mito M, Arai T, Takahashi E, Fujita Y, Fujimori T, Standaert L, Marine J-C, et al. The lncRNA Neat1 is required for corpus luteum formation and the establishment of pregnancy in a subpopulation of mice. Development. 2014;141(23):4618–27.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Funes S, Hedrick JA, Vassileva G, Markowitz L, Abbondanzo S, Golovko A, Yang S, Monsma FJ, Gustafson EL. The KiSS-1 receptor GPR54 is essential for the development of the murine reproductive system. Biochem Biophys Res Commun. 2003;312(4):1357–63.

    Article  CAS  PubMed  Google Scholar 

  47. Kaiser UB, Kuohung W. KiSS-1 and GPR54 as new play-ers in gonadotropin regulation and puberty. Endocrine. 2005;26(3):277–84.

    Article  CAS  PubMed  Google Scholar 

  48. d’Anglemont de Tassigny XD, Fagg LA, Dixon JP, Day K, Leitch HG, Hendrick AG, Zahn D, Franceschini I, Caraty A, Carlton MB, et al. Hypogonadotropic hypogonadism in mice lacking a functional Kiss1 gene. Proc Natl Acad Sci U S A. 2007;104(25):10714–9.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Lapatto R, Pallais JC, Zhang D, Chan YM, Mahan A, Cerrato F, Le WW, Hoffman GE, Seminara SB. Kiss1−/− mice exhibit more variable hypogonadism than Gpr54−/− mice. Endocrinology. 2007;148(10):4927–36.

    Article  CAS  PubMed  Google Scholar 

  50. Wolfe A, Divall S, Wu S. The regulation of reproductive neuroendocrine function by insulin and insulin-like growth factor-1 (IGF-1). Front Neuroendocrinol. 2014;35(4):558–72.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  51. Camille Melon L, Maguire J. GABAergic regulation of the HPA and HPG axes and the impact of stress on reproductive function. J Steroid Biochem Mol Biol. 2016;160:196–203.

    Article  CAS  PubMed  Google Scholar 

  52. Allen MP, Xu M, Zeng C, Tobet SA, Wierman ME. Myocyte enhancer factors-2B and -2C are required for adhesion related kinase repression of neuronal gonadotropin releasing hormone gene expression. J Biol Chem. 2000;275(50):39662–70.

    Article  CAS  PubMed  Google Scholar 

  53. Hiney JK, Srivastava VK, Volz CE, Dees WL. Alcohol alters insulin-like growth factor-1-induced transforming growth factor beta1 synthesis in the medial basal hypothalamus of the prepubertal female rat. Alcohol Clin Exp Res. 2014;38(10):2572–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Cheong RY, Czieselsky K, Porteous R, Herbison AE. Expression of ESR1 in Glutamatergic and GABAergic Neurons Is Essential for Normal Puberty Onset, Estrogen Feedback, and Fertility in Female Mice. J Neurosci. 2015;35(43):14533–43.

    Article  CAS  PubMed  Google Scholar 

  55. Sakamoto K, Wakabayashi Y, Yamamura T, Tanaka T, Takeuchi Y, Mori Y, Okamura H. A population of kisspeptin/neurokinin B neurons in the arcuate nucleus may be the central target of the male effect phenomenon in goats. PLoS One. 2013;8(11):e81017.

    Article  PubMed  PubMed Central  Google Scholar 

Download references


We thank all members of the Anhui Provincial Laboratory of Animal Genetic Resources Protection and Breeding for their stimulating discussions.


This work was supported by the National Natural Science Foundation of China (Grant number 31472096), the National Transgenenic New Species Breeding Program of China (No. 2014ZX08008-005-004), the National Natural Science Foundation of China (Grant number 31301934), and the Anhui Provincial Natural Science Foundation (Grant number 1408085MKL40).

Availability of data and materials

The sequencing data were submitted to the Genome Expression Omnibus (Accession Numbers GSE84301) in NCBI.

Authors’ contributions

GXX and YJ conceived of the study, participated in its design and coordination and drafted the manuscript; YC, ZKF, LXM and LL conducted qRT-PCR validation and statistical analysis; DJP and LYS performed the statistical analysis; CHG, LYH, ZXR and LY carried out the analysis of data, animal experiments, and surgical processes; FFG and ZYH participated in the design and coordination and helped to draft the manuscript. All authors read and approved the final manuscript.

Competing interests

The authors declare that they have no competing interests.

Consent for publication

Not applicable.

Ethics approval

The study was approved by the Animal Care and Use Committee of Anhui Agricultural University. The methods were carried out in accordance with the approved guidelines. All experimental procedures involving goats were performed according to the Regulations for the Administration of Affairs Concerning Experimental Animals (Ministry of Science and Technology, China; revised in June 2004).

Author information

Authors and Affiliations


Corresponding authors

Correspondence to Fugui Fang or Yunhai Zhang.

Additional files

Additional file 1:

The production of reads from the Illumina HiSeq 4000 platform. (XLSX 9 kb)

Additional file 2:

The FPKM of differential expression transcripts. (XLSX 5267 kb)

Additional file 3:

The protein-coding genes as potential targets 10K/100K upstream and downstream of the lncRNAs. (XLSX 372 kb)

Additional file 4:

The expression of protein-coding genes related puberty. (XLSX 10 kb)

Additional file 5:

GO analysis of predicted targets of lncRNAs in cis. (XLS 119 kb)

Additional file 6:

KEGG pathway analysis of predicted targets of lncRNAs in cis. (XLS 21 kb)

Additional file 7:

GO analysis of predicted targets of lncRNAs in trans. (XLS 1551 kb)

Additional file 8:

KEGG pathway analysis of predicted targets of lncRNAs in trans. (XLS 367 kb)

Additional file 9:

The specific expression of lncRNAs in pubertal/prepubertal samples. (XLSX 36 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, 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 ( 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

Gao, X., Ye, J., Yang, C. et al. Screening and evaluating of long noncoding RNAs in the puberty of goats. BMC Genomics 18, 164 (2017).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: