Skip to main content
Download PDF

Synchronous profiling and analysis of mRNAs and ncRNAs in the dermal papilla cells from cashmere goats

BMC Genomics volume 20, Article number: 512 (2019) Cite this article

Abstract

Background

Dermal papilla cells (DPCs), the “signaling center” of hair follicle (HF), delicately master continual growth of hair in mammals including cashmere, the fine fiber annually produced by secondary HF embedded in cashmere goat skins. Such unparalleled capacity bases on their exquisite character in instructing the cellular activity of hair-forming keratinocytes via secreting numerous molecular signals. Past studies suggested microRNA (miRNAs) and long non-coding RNAs (lncRNAs) play essential roles in a wide variety of biological process, including HF cycling. However, their roles and related molecular mechanisms in modulating DPCs secretory activities are still poorly understood.

Results

Here, we separately cultivated DPCs and their functionally and morphologically distinct dermal fibroblasts (DFs) from cashmere goat skins at anagen. With the advantage of high throughput RNA-seq, we synchronously identified 2540 lncRNAs and 536 miRNAs from two types of cellular samples at 4th passages. Compared with DFs, 1286 mRNAs, 18 lncRNAs, and 42 miRNAs were upregulated, while 1254 mRNAs, 53 lncRNAs and 44 miRNAs were downregulated in DPCs. Through overlapping with mice data, we ultimately defined 25 core signatures of DPCs, including HOXC8 and RSPO1, two crucial activators for hair follicle stem cells (HFSCs). Subsequently, we emphatically investigated the impacts of miRNAs and lncRNAs (cis- and trans- acting) on the genes, indicating that ncRNAs extensively exert negative and positive effects on their expressions. Furthermore, we screened lncRNAs acting as competing endogenous RNAs (ceRNAs) to sponge miRNAs and relief their repressive effects on targeted genes, and constructed related lncRNAs-miRNAs-HOXC8/RSPO1 interactive lines using bioinformatic tools. As a result, XR_310320.3-chi-miR-144-5p-HOXC8, XR_311077.2-novel_624-RSPO1 and others lines appeared, displaying that lncRNAs might serve as ceRNAs to indirectly adjust HFSCs status in hair growth.

Conclusion

The present study provides an unprecedented inventory of lncRNAs, miRNAs and mRNAs in goat DPCs and DFs. We also exhibit some miRNAs and lncRNAs potentially participate in the modulation of HFSCs activation via delicately adjusting core signatures of DPCs. Our report shines new light on the latent roles and underlying molecular mechanisms of ncRNAs on hair growth.

Background

Cashmere goats gain worldwide reputation for yielding cashmere, the fine hair fiber with excellent quality and commercial value, from the secondary hair follicle (HF) of skins [1,2,3]. An outstanding feature of cashmere growth is annual rhythm, which synchronously occurs with yearly sequential switches of HF among fast growth stage (anagen), gradual degeneration stage (catagen), and relative quiescence stage (telogen) [4,5,6]. Although cashmere elongation tightly synchronizes with the alterations of photoperiod and endocrine status, the basic rationales are similar with other mammals [7]. Cyclic transformations of follicular keratinocytes among proliferation, differentiation and apoptosis are the cellular foundations of periodical HF regrowth [8]. However, the specialized fibroblasts locating at the bottom of HF, dermal papilla cells (DPCs), are widely accepted as the controlling center of HF growth. Via spatiotemporally releasing specific signals,

DPCs can finely instruct the activities of follicular keratinocytes to reshape HF and produce new hair shaft. DPCs also serve as a transfer station to relay the effects of locally or systematically generated hormones and other molecules on hair growth [9]. Past studies suggested that the particular ability of DPCs on HF cycling resides in their characteristic gene expression profiles, and identified some involved key regulators and signaling pathways such as Sox2, Prdm1, β-catenin, Wnt/β-catenin pathway, FGF pathways and others [10,11,12,13,14]. Whereas, the current picture of overall landscape is far from complete, especially for the versatile regulatory non-coding RNAs (ncRNAs).

Long thought as “evolutionary junk”, nowadays ncRNAs have been continually implicated to play irreplaceably regulatory roles in gene expression [15]. In contrast to functionally monotonous, microRNAs (miRNAs), long non-coding RNAs (lncRNAs), which are typically defined as RNA transcripts > 200 nucleotides in length, can actively modulate gene expression through various mechanisms that are not yet explicitly understood [16]. LncRNAs extensively exert their regulatory roles at transcriptional and posttranscriptional levels and play indispensable roles in a number of vital developmental and pathological processes, such as cell differentiation [17], organogenesis [18], X-chromosome imprinting [19], stem cells pluripotent maintenance [20], and carcinogenesis [21]. Recently, a mounting quantity of research has revealed that ncRNAs are the essential participators of HF development and functionality maintenance of DPCs. Epidermal specific deletion of Dicer, the necessary miRNAs-processing enzyme, resulted in stunned HF formation and hyperproliferative follicular cells [22]. Inducible overexpression of miR-214 drastically inhibited HF development and cycling via specifically decreasing the expression of β-catenin, the key Wnt signaling mediator [23]. Moreover, alteration of miRNAs expression profiles were frequently associated with the loss of hair-modulatory functions of DPCs on humans and mice [24, 25]. Some miRNAs such as miRNA-125b and miRNA-195-5p have been confirmed as the casual factors through restraining growth factors expressions or weakening intensity of Wnt pathway, the primary pivotal signaling inside DPCs [26, 27]. LncRNAs also make a considerable contribution to ncRNAs-mediated manipulation of hair cycling, though less explored. Several reports revealed that the expression of lncRNAs experiences obvious fluctuations during stage transitions of cashmere cyclic growth [2, 5]. A few lncRNAs such as H19, HOTAIR and lncRNA-000133 have been shown to regulate genes closely involved in cyclic HF growth in DPCs [28, 29]. Apart from functioning individually, lncRNAs can serve as competing endogenous RNAs (ceRNAs) to specifically decoy miRNAs for preventing their suppression on targeted genes [30, 31]. Such feature of lncRNAs facilitates more precise and complex regulations of gene expression, which have been uncovered to relate with a series of physiological events and cancer progressions [32,33,34], and might contribute to the exquisite modulation of DPCs functionality.

Though a few studies suggest ncRNAs are important regulators of hair cycling on cashmere goats, they mainly focused on skin tissue, a complex structure comprising a dozen of cell types [3, 5, 35]. There are still no reports on the expression profiles and the potential roles of ncRNAs in goat DPCs. Otherwise, previous researchers rarely executed simultaneous mRNAs, miRNAs and lncRNAs analysis from identical DPCs sample sets either on mice or other animals, leading to insufficient interpretation of lncRNAs performing as the newly proposed ceRNAs. In the present trial, we concurrently profiled above transcripts from DPCs and dermal fibroblasts (DFs) cultivated from anagen cashmere goat skins. Using bioinformatic tools, we subsequently screened key genes and ncRNAs, and constructed their interactive networks to highlight the presumed functions of ncRNAs in hair growth, especially for hair follicle stem cells (HFSCs) activation. Our study will provide novel clues and viewpoints to deeply explore the unrecognized functions of ncRNAs in DPCs-centered hair growth.

Results

Cultivation and morphological characterization of goat DPCs and DFs

In vitro cultured DPCs specifically supported HF neogenesis and hair regrowth, whereas their relatives DFs did not [36]. To comprehensively gain insights into the innate feature of this specialized cell type, we designed the present study to culture and molecularly characterize both cell populations from lateral backsides skins of 2-year-old female cashmere goats at anagen (Fig. 1a). Both populations outgrew from their respective explant in 5–7 days, and their passaging cultures obviously exhibited distinct cellular morphologies (Fig. 1b). DPCs are flat in appearance, with spread-out surfaces and multiple cellular projections, whereas DFs typically show a bipolar and spindle-like pattern, possessing a relatively smaller cell volume. These morphological features in goats are consistent with mice and other mammals [37, 38], thereby implying the successful cultivations of related cells and the high fidelity of subsequent analysis.

Fig. 1
figure 1

Schematic presentation of the overall workflow in the present study (a) and in vitro cultivation of DPCs and DFs (b). Cells migrated out from each explant in 5–7 days, and their subcultures showed remarkably different morphologies: DPCs are relatively flat, whereas DFs show a spindle-like cellular shape. Scale bar = 1000 μm

Full size image

Comprehensive identification of ncRNAs in goat DPCs and DFs

As previously reported, when cultured in vitro goat DPCs at early passages formed a sphere-like aggregate, an indicator of their hair-inducing capacities, whereas DFs did not possess such ability [37, 39, 40]. To uncover the ncRNAs and mRNAs underlying the morphological and functional differences between DPCs and DFs, we selected both cell types at 4th passages and synchronously sequenced six cell samples (n = 3 for each type) on Illumina HiSeq 4000 (for mRNAs and lncRNAs) and HiSeq 2000 platforms (for miRNAs). Standard bioinformatic analysis procedures were carried out as previous studies [2, 5].

As a result, we obtained an average of 114,777,434 and 119,977,164,150 bp paired-end raw reads for each DPCs and DFs sample, respectively. When invalid reads were intentionally excluded as usual, filtered clean reads (accounting for 98.5% of raw reads) were mapped to the goat reference genome [41]. After the assembly and quantification of the linear transcripts using StringTie (v1.3.1) [42] and Cuffdiff (v2.1.1) [43], respectively, we distinguished mRNAs and lncRNAs through five successive sifting steps. The detailed procedures and results are shown in Additional file 1: Figure S1a). In the final step, three tools CPC [44], Pfam [45] and CNCI [46], were jointly applied to fully eliminate transcripts with uncertain coding potential, ensuring the credible discovery of fresh lncRNAs (Additional file 1: Figure S1b). Finally, 786 annotated lncRNAs and 1754 novel ones were recognized in all samples. A total of 22,972 mRNAs were also discerned.

As for miRNAs profiling, an average of 14,555,565 and 12,621,142 50 bp single-end raw reads were acquired from each DPCs and DFs sample, respectively. Filtered as reported before [2], clean reads (accounting for 97.73% of raw reads) with certain ranges of length (18–35 nt) were mapped to reference sequences [41], and mapped tags were searched within miRBase to sort out known miRNAs. In total, 397 mature known miRNAs were found. We further used two available software programs miREvo [47] and mirdeep2 [48], to speculate novel miRNAs based on a featured hairpin structure, the Dicer cleavage site, etc.; ultimately, 139 novel transcripts were predicted. Collectively, we discovered 536 expressed mature miRNAs in all samples.

Numbers of mRNAs and ncRNAs are summarized in Table 1. Genomic information of lncRNAs is detailed in Additional file 2: Table S1.

Table 1 Information of known and novel transcripts in DPCs and DFs samples
Full size table

Genomic feature comparison between lncRNAs and mRNAs

To ascertain the results of transcripts identification and supplement the annotation information of goat genome, we compared several noticeable genomic features between lncRNAs and mRNAs. The analytical results indicated that annotated lncRNAs normally contain more exons than novel lncRNAs, although fewer than mRNAs (Additional file 3: Figure S2a). Lengths of these transcripts showed an obviously distinct pattern, in which annotated and novel lncRNAs are shorter and more narrowly distributed than mRNAs (Additional file 3: Figure S2b). Moreover, the overall expression levels of both lncRNAs were significantly lower than mRNAs, and the newly identified lncRNAs were the lowest (Additional file 3: Figure S2c). These genomic traits are highly in accordance with previous reports on goats and other species [49,50,51,52], suggesting the credible profiling of lncRNAs in present trial.

Differentially expressed mRNAs and ncRNAs

To fulfill the purpose for selecting potential candidates related to hair growth, differentially expressed mRNAs and ncRNAs were examined between two sets of samples using Ballgown suite [53]. As visualized by heatmaps in Additional file 4: Figure S3, three samples perfectly gathered inside DPCs and DFs groups. Volcano plots in Fig. 2 graphically displayed fold changes and statistical significances of entire transcripts. In total, 1286 up-regulated and 1252 down-regulated mRNAs emerged in goat DPCs compared with DFs. Then, we checked the expressions of several genes frequently utilized for in vivo labelling of DPCs on mice, and discovered HOXB6, ENPP2, CRABP1 and PRDM1 are highly expressed in goat DPCs. Among them, CRABP1 is a constant marker of DPCs, and expresses throughout entire stages of hair cycling [54]. Of note, COL15A1, one of top expressed genes, is solely abundant in DPCs. Such condition is in consistent with the discovery that DFs do not produce COL15A1 protein in culture dishes [55]. Most importantly, we defined the core signatures of DPCs through overlapping with mice data [56], and as a result a total of 25 gene emerged (Table 2). The reason for a relatively fewer number of genes is that the mice markers were screened through comparing with other four cell lineages including melanocytes, outer root sheath cells, HMCs and DFs. The gene list contains LEPR, WNT5A, PRDM1, HOXC8, FGFR2, BMP4, LTBP1, and PTGFR. They almost participate in all pivotal aspects of DPCs functions in HF cycling, such as HFSCs activation, hair matrix cells (HMCs) proliferation and differentiation, angiogenesis, and hormone-mediated hair growth regulation. Moreover, we also detected that two hormone receptors PTGER4 and ESR1 are more abundant in DPCs, as reported before [57, 58].

Fig. 2
figure 2

Differentially expressed ncRNAs and mRNAs between goat DPCs and DFs. a Differentially expressed mRNAs. b Differentially expressed lncRNAs. c Differentially expressed miRNAs. Among these mRNAs and ncRNAs, 1286 mRNAs, 17 lncRNAs, 42 miRNAs were upregulated, and 1254 mRNAs, 53 lncRNAs and 44 miRNAs were downregulated in DPCs compared with DFs. Red dots represent upregulated transcripts and green dots represent downregulated transcripts in DPCs compared with DFs

Full size image
Table 2 Information of 25 core signatures in goat DPCs and DFs
Full size table

A relatively less quantity of 71 lncRNAs exhibited varied abundances between two different functional cell populations. Of 18 elevated lncRNAs in DPCs, XR_001296062.1, LNC_001710, XR_001917771.1 and LNC_000335 showed exclusive expressions. Likewise, LNC_000726, LNC_000799, XR_310887.3, LNC_000964, LNC_000901, and LNC_000299 are uniquely occupied by DFs among the remaining 53 downregulated lncRNAs. These results indicated that lncRNAs may regulate DPCs functionality and further hair cycling via their cell-type restricted presences.

As for miRNAs analysis, we found 86 miRNAs have significantly different abundances between goat DPCs and DFs. Among them, 42 and 44 miRNAs showed upregulated and downregulated expressions in DPCs compared with DFs, respectively. Some miRNAs, such as miR-125b [24] and miR-196a [26] have been documented with regulatory roles in DPCs functions on human and mice.

We listed top 20 (10 upregulated and 10 downregulated) differently expressed lncRNAs and miRNAs in Table 3. Holistic lists of mRNAs and ncRNAs and more details are provided in Additional file 5: Table S2.

Table 3 Top 20 differentially expressed ncRNAs in goat DPCs and DFs
Full size table

To assure the accuracy of sequencing strategy, bioinformatic identifications and analysis, we randomly picked several differentially expressed mRNAs (ESR1, INHBA, INHBB, IGFBP2, LEPR, BMP4, THBS1, WNT2, and WNT5A) and lncRNAs (LOC102172689, LOC106503672, JAM3, and LOC102171315), and used qPCR to ascertain their relative levels between the two sample groups. The expression levels of transcripts determined by qPCR and RNA-seq are highly consistent (Additional file 6: Figure S4), thereby implying the prominent reliability of RNA-seq data acquisitions and following analytical procedures in present study.

KEGG pathway enrichment analysis

To better understand the functioning approaches of DPCs in hair growth, KEGG analysis for upregulated mRNAs is performed, and finally 247 pathways were enriched. Of the 20 top pathways, focal adhesion and ECM-receptor interaction frequently appeared in transcriptomic studies on goat skin tissues and DPCs on mice [3, 59, 60]. Most importantly, several hormone-related pathways (Estrogen signaling pathway, Thyroid hormone signaling and Adipocytokine signaling pathway) draw our special attention, though they are not significantly enriched. There are two reasons: one is all of them were emphasized during stage transitions of cashmere growth in skins by researchers before [2, 5]; the other is corresponding hormone-binding receptors (i.e., ESR1, IGTAV and LEPR) were on the list of upregulated genes, which strongly suggests goat DPCs are targets of estrogen, thyroid hormone, and leptin in the HF. Full list of pathways and their involved genes are shown in Additional file 7: Table S3.

Prediction of lncRNAs and miRNAs on the regulation of core signatures

MiRNAs exert their physiological roles via specifically binding, repressing translation and promoting degradation of their targeted mRNAs [61]. To infer the mighty functions of miRNAs in DPCs, we constructed their interactive network with core signatures (Fig. 3; Additional file 8: Table S4). As a result, we found that two HFSCs activation-related genes are regulated by several miRNAs. HOXC8, a crucial gene for HFSCs activation [60], is the target of chi-miR-144-5p. Another gene RSPO1 with similar character in hair cycling is negatively modulated by chi-miR-874-3p, chi-miR-23b-5p and other five miRNAs. There are several ways for lncRNAs to exert their modulatory roles on gene expression, including as an enhancer, epigenetic modifier and transcriptional regulator [62]. Commonly, two modes cis-role (in which lncRNAs act on adjacent genes with 100 kb distances) and trans-role (in which the Pearson correlation of mutual expression levels is ≥0.95 or ≤ − 0.95) are widely adopted by researchers to forecast lncRNAs-gene interactive pair [2, 5]. Here, we focused on the possible relations of lncRNAs to HFSCs vitalization, and picked the two functional genes HOXC8 [60] and RSPO1 [63] as the targets. Finally, our results suggested that four lncRNAs fit the criterion of cis-role, but only LNC_000354, who situates 63,083 bp upstream of HOXC8 locus, displayed differential expression (downregulated in DPCs), potentially imposing an inhibitive effect on HOXC8. For trans-acting pattern, we detected 13 lncRNAs, comprising 6 (e.g., XR_001296062.1 and LNC_001710) positively and 7 (e.g., XR_001919776.1 and XR_001918151.1) negatively correlated, respectively (Fig. 4a). Among these lncRNAs, some such as XR_001296062.1 are uniquely expressed in DPCs, whereas LNC_000901 possesses a DFs-specific expression. As for RSPO1, none lncRNAs act as cis-regulator and 26 lncRNAs work in trans-manner (Fig. 4b). Some regulatory lncRNAs also showed a specific expressing pattern in goat DPCs and DFs. Altogether, we implied that lncRNAs might affect HFSCs status via regulating gene expressions inside DPCs.

Fig. 3
figure 3

Interactive map of miRNAs with core signatures. Triangles with grey color indicate downregulated miRNAs in DPCs compared with DFs, and circles with yellow color indicate core signature of goat DPCs. The size of each shape positively correlates with the number of connected lines

Full size image
Fig. 4
figure 4

Potential lncRNAs targeting HOXC8 (a) and RSPO1 (b). Multiple lncRNAs target the two genes. HOXC8 and RSPO1 were indicated as circle with yellow color. LncRNAs were indicated as squares, yellow and grey color represents upregulated and downregulated lncRNAs in goat DPCs compared with DFs, respectively

Full size image

We also provide entire interactive relationships among core signatures and ncRNAs in Additional file 9: Figure S5.

LncRNAs serve as ceRNAs to indirectly regulate HOXC8 and RSPO1

Extensive studies discovered that lncRNAs harbor microRNA (miRNA)-response elements (MRE), and can serve as ceRNAs to decoy miRNAs, resulting in indirect upregulation of associated mRNAs [64]. Extensive pivotal roles concerning lncRNAs as ceRNAs have been reported on goats and other animals [21, 62, 65]. Here, we identified lncRNAs that might function as ceRNAs using a combination of bioinformatic tools, and focused their regulatory relationships with HOXC8 and RSPO1. Our results suggested that some lncRNAs specifically bind to certain miRNAs targeting the genes to relief associated mRNAs from suppression, and thus modulate gene expressions as ceRNAs, such as XR_310320.3-chi-miR-144-5p-HOXC8 and XR_311077.2-novel_624-RSPO1 interactive lines (Fig. 5). XR_310320.3 and XR_311077.2 act as ceRNAs to positively regulate HOXC8 and RSPO1 expressions, respectively. We showed the possibility that lncRNAs act as ceRNAs to participate in the regulation of HFSCs activation. Whole list of interactive pairs was provided in Additional file: Table S6.

Fig. 5
figure 5

LncRNAs serve as ceRNAs. a XR_310320.3 functions as ceRNAs to upregulate HOXC8. b XR_001295577.2, XR_311077.2 and other eight lncRNAs function as ceRNAs to upregulate RSPO1. Circle represents mRNAs, square represents lncRNAs and triangle represents miRNAs. Yellow color indicates upregulated transcripts and blue color represent downregulated transcripts. LncRNAs were assumed to specifically sponge miRNAs to relief their suppressive roles on targeted mRNAs expression. The targeted relationship of miRNAs with mRNAs/lncRNAs were predicted using miRanda

Full size image

Discussion

Mammalian mature HF is a complicated miniorgan, mainly composing diverse types of epithelial cells. Although it is largely epithelium originated, HF contains a flock of specialized fibroblasts at its bottom, DPCs, who play pivotal roles in the regulation of continual regeneration of HF [66]. As hair regrowth initiates, signals secreted by DPCs are deemed to direct HFSCs locating in the bulge region of HF to proceed transit division. Then, stem cell progenies migrate downward into the base of HF, where they encircle the DPCs, forming the HMCs. Receiving further instructive signals from DPCs, HMCs serially undergo rapid proliferation and terminal differentiation to form keratinized hair shaft [7, 67]. Though unquestionable importance of DPCs in hair cycling, regulatory mechanisms regarding the expressions of such molecules inside DPCs are not well understood yet. In the present study, we executed a comprehensive ncRNAs and coding genes expression analysis between DPCs and DFs of cashmere goats, providing new views of ncRNAs in hair cycling.

We initially acquired DPCs and DFs from goat skins, and exhibited they possess obviously heterogeneous external appearances. This inconformity was identical with human [68], mice [36] and other animals [69, 70]. Previous studies found that skin implantation of cultured DPCs induced new hair growth, whereas DFs can’t [36, 71]. Above facts built a solid foundation for subsequent trial. Next, we carried out a comprehensive RNA-seq of both cell populations, and proceeded a series of bioinformatic analysis. Functional specificity of a cell is determined by its unique gene expression pattern, which is further controlled by diverse factors [56, 72]. Several reports on humans and mice have unveiled the crucial characters of ncRNAs in modulating the functionality of DPCs during hair growth. Indeed, altered profiles of miRNAs have been frequently connected to dysfunctions of human and mice DPCs, and recently some lncRNAs were also found as casual elements of such conditions [28]. However, previous studies mainly focused on humans and mice, little information was known about the functions of ncRNAs in goats, especially for less conservative lncRNAs [62]. In present study, we discovered that a plenty of miRNAs and lncRNAs were abundantly present in DPCs and DFs, and displayed remarkably differential expressions between two cell types. Unique expressions of certain miRNAs (e.g., novel_144 and chi-let-7c-3p) and lncRNAs (e.g., LNC_001710 and LNC_000799) in individual cell lineage were also discerned (Table 2), suggesting they might participate in the regulation of cellular functions in a cell type-restricted manner [15]. Our further analysis showed that a sum of 86 miRNAs, 71 lncRNAs and 2538 mRNAs transcripts were differentially expressed between DPCs and DFs. Some miRNAs (e.g., chi-miR-1271-5p, chi-miR-9-5p, and chi-miR-22-3p) and lncRNAs (e.g., LNC_000349, XR_001297172.2, and XR_310056.3) were previously identified as differentially expressed transcripts among stage transitions of cashmere growth [2, 3, 5]. However, additional cautions should be needed because these studies profiled ncRNAs from skin tissue, a complex structure comprising diverse kinds of cell including epidermal keratinocytes, fibroblasts, intradermal adipocytes, as well as HF itself. The majority of them experience drastic fluctuation in morphologies and gene expressions during hair cycling [7]. We also found that two miRNAs miR-125b (downregulated in DPCs) and miR-335 (upregulated in DPCs) deserve more attention. miR-125b was shown to decrease level of FGF7, a potent growth factor secreted by goat DPCs [27]. miR-335 was positively correlated with the hair-inducing capacity of mice DPCs [26]. Despite the truth that the functions of ncRNAs in DPCs have not been extensively explored yet, their remarkably distinct expression patterns between DPCs and DFs suggest they should be a pivotal participator of hair growth regulation. Compared with previous studies on complex skin tissue of goats, the present single cell type transcriptomes could offer reliable and valuable candidates for downstream functional verification. On the other hand, synchronous detection of ncRNAs and mRNAs enables us to faithfully construct relevant interplay network among diverse transcripts and further deduce unearthed mechanisms of ncRNAs, which have been pervasively done in other cells or tissues [30, 32, 73], however sparsely reported in mammalian DPCs.

The HF is extraordinarily sensitive to locally and systematically generated hormones that modulate hair growth, leading to an adaptive accommodation to external environment, such as climate change [74]. Hormone-participated regulation makes sure that the HF timely transits among different developmental stages [75]. We found that three hormone-related pathways Estrogen signaling pathway, Thyroid hormone signaling and Adipocytokine signaling pathway were enriched in our KEGG analysis and corresponding hormone receptors (i.e., ESR1, IGTAV and LEPR) were significantly upregulated in DPCs compared with DFs. These pathways were also frequently appeared in recent reports on cashmere periodical development [2, 3, 5], implying DPCs are the follicular targets of these hormones. These results are also in consistent with the discoveries on humans and mice [57, 76]. Leptin and thyroid hormone markedly prolonged growing phase of HF [77, 78], whereas estrogen exerted an inhibitive effect on hair growth [79]. Though crucial roles of these hormones on hair development, working approaches in DPCs are not clear yet. A series of intraceullar kinases and transcription factors are responsible for corresponding signaling transductions and gene expression modulations of the hormones [80, 81]. Some of them (e.g., JAK2, MAPK1, PIK3R1, JUN1 and SP1) were highly expressed in goat DPCs, thus providing useful clues for further researches.

We further defined 25 core signatures of DPCs through overlapping present result with a previously published data of mice [56]. Among these genes, the functions of the majority of them were well-characterized on mice models. Induced activation of HFSCs by signals emanating from DPCs drives HF reconstruction and hair regrowth [82]. Here, we found that HOXC8, RSPO1 and BMP4 are the involved molecules. RSPO1, a secreted Wnt/β-cateinin agonist, is highly enriched in mice DPCs, and its recombinant protein injection into mice skins resulted in vitalizing HFSCs and precocious anagen entry from mid-telogen [63]. The transcription factor, HOXC8, exerted similar influences on hair cycling through acting as an upstream positive regulator of Wnt signaling [60]. Of greatly noticeable, the hypermethylation status of HOXC8 exon1 is associated with shorter fleece length on cashmere goats [6], and first exon methylation tightly concerns with transcriptional silencing [83]. These observations demonstrate that artificial alteration of HFSCs status through manipulating gene expression of DPCs is a hopeful strategy for enhancing fiber productions in farm animals. Another gene, BMP4 was reportedly suppressing HFSCs in quiescent state [84]. Interestingly, melatonin, the famous hormone affects seasonal growth of cashmere, seems to intimately interplay with both Wnt/β-cateinin and BMP4 in other tissues and systems [85]. Whether it acts in the same manner in DPCs should be determined in future.

Accumulating evidences suggest miRNAs and lncRNAs impose an indispensable modulatory role on gene expression [21]. Compared with a sole approach of miRNAs, lncRNAs can perform as enhancers, transcriptional regulators, miRNAs sponges and others, which are generally categorized as cis- and trans- regulators [16]. Several paradigms have been well appreciated such as Xist acting as a cis-modulator for genomic imprinting and HOTAIR serving as trans-regulator for HOXD gene silencing [19, 86]. Past reports have identified a few ncRNAs and uncovered their vital functions in DPCs, such as miR-125b, miR-195p, HOTAIR and LncRNA-000133 [26,27,28,29]. Whereas, there are still no reports focused on the connections of ncRNAs with DPCs-mediated HFSCs activation. In this study, we established the regulatory relationships of miRNAs and lncRNAs with HOXC8 and RSPO1. Our results suggested that some miRNA and lncRNAs can specifically target the genes and impose negative or positive effects of their expressions, implying the possible characters of ncRNAs in HFSCs vitalization. Otherwise, the crucial characters of lncRNAs functioning as ceRNAs in stem cell biology have also been demonstrated, such as linc-RoR specifically decoys miR-145, upregulates Oct4, Nanog and Sox2 and contributes to the self-renewal of embryonic stem cells [87]. We also identified lncRNAs might serve as ceRNAs to sponge miRNAs and modulate the levels of two genes, and built related interactive lines such as XR_310320.3-chi-miR-144-5p-HOXC8 and XR_311077.2-novel_624-RSPO1. XR_310320.3 and XR_311077.2 may work as potent ceRNAs to particularly enhance the expressions of HOXC8 and RSPO1, resulting indirect vitalization of HFSCs. Overall, we exhibited multiple facets of lncRNAs in DPCs, and their critical importance in adjusting HFSCs status. These will greatly facilitate the understanding of ncRNAs in DPCs and HF development.

Conclusion

In present study, we reported the expression profiles of miRNAs and lncRNAs, as well as mRNAs in goat DPCs and DFs. Through definition of core signatures of DPCs and construction of their regulatory network with ncRNAs, we exhibited the specific modulation of ncRNAs on HOXC8 and RSPO1, the genes involved in HFSCs activation. Our study provided specific transcriptomic data for cashmere goat research and hinted the potential roles of ncRNAs in HF development and cycling.

Methods

Cell cultivation

The Shanbei cashmere goats were carefully managed in a farm located in the Yangling District, Shaanxi Province, China. We selected three healthy 2-year-old female Shanbei White Cashmere goats with similar body weights (~ 35 kg) and uncrossed lineage records for the current study. Subcutaneous procaine injection was performed to alleviate animal suffering. In September (anagen phase of cashmere growth), approximately 2 cm2 skin samples on lateral backsides of the goats were surgically removed under standard anesthetic and aseptic procedures. After suturing all wounds, we intentionally took good care of experimental goats to boost their recovery from surgery. The excised skin samples were washed three times with sterile phosphate-buffered solution (PBS) to clean off blood or other contaminants. Then, the tissues were cut at the interface of the subcutaneous adipose layer and dermis. An incised subcutaneous adipose layer was applied for isolation of intact HF and in vitro cultivation of dermal papilla cells, as previously stated [39]. The remaining tissues were submerged in 0.25% dispase II overnight at 4 °C to separate the epidermis and dermis, then the dermis was cut into small blocks approximately 1 mm3 in size and used as explants for DFs culture. Normally, DPCs and DFs grew from each explant about five days after initial adhesion and were passaged when they reached 100% confluency. Typically, DPCs and DFs from the 4th passages were used for all of the experiments. Primary and passaging cultures of each cell type were maintained in DMEM/F12 supplemented with 10% FBS, 100 units/ml penicillin, and 100 μg/ml streptomycin, and the cultures were incubated at 37 °C and 100% humidity in a 5% CO2/ 95% air incubator.

Total RNA extraction, library preparation, and sequencing

Six strains of cells (three for DPCs or DFs) without cross contamination were selected, and their total RNAs were extracted using the RNA isolation kit mentioned above according to the manufacturer’s instructions. To ensure the fidelity of co-expression analysis, two independent libraries were simultaneously constructed for each strain: the lncRNAs library (for mRNAs and lncRNAs) and the miRNAs library. The preparation of both libraries was performed as previously reported [2], and the resulting libraries were further sequenced on Illumina HiSeq 4000 (for mRNAs and lncRNAs) or HiSeq 2000 platforms (for miRNAs), which generated 150-bp paired-end and 50-bp single-end reads, respectively.

Sequencing data analysis

To determine the landscape of transcript expression in all samples, clean reads were acquired by removing reads containing poly-N, along with low-quality and other invalid reads among all of raw reads. At the same time, calculations of the Phred score (Q20, Q30) and GC content of the clean data were also performed before downstream analysis. For lncRNAs library mapping, Bowtie2 (v2.2.8) was employed to assess the index of the reference genome and paired-end clean reads were aligned to the reference genome using HISAT2 (v2.0.4) [42]. At the same time, Bowtie (v0.12.9) was used to map small RNA tags to the reference sequence without mismatch to analyze their expression and distribution on the reference. Assembly of mapped reads was achieved using StringTie (v1.3.1) in a reference-based approach. Cuffdiff (v2.1.1) was used to calculate the fragments per kilo-base of exons per million fragments mapped (FPKMs), and the Ballgown suite was adopted for the determination of differentially expressed lncRNAs and mRNAs transcripts [42]. LncRNAs were distinguished based on several key well-known criteria: exon number (> 2); transcripts length (> 200 nt); abundance (FPKM > 0.5); and protein-coding potential, which were excluded by the CNCI (v2) [46], CPC (cpc-0.9-r2) [44] and Pfam-scan (v1.3) [45] tools. For target gene prediction of lncRNAs, the cis-role was set as lncRNAs act on neighboring target genes, which were defined as genes located 100 k upstream and downstream of lncRNAs locus. Expression correlations between lncRNAs and coding genes were calculated to predict the trans-role of lncRNAs.

MiRBase20.0 was chosen as the reference to filter known miRNAs, and novel miRNAs were predicted using miREvo [47] and mirdeep2 [48] through the exploration of the secondary structure, the dicer cleavage sites and others. Prediction of the target gene of miRNAs was carried out by miRanda [88], and miRNA expression abundances were estimated by transcript per million (TPM) using the following criteria: normalized expression = mapped readcount/total reads*1000000. Differential expression of miRNAs among samples was analyzed using the DESeq R package (1.8.3) [89].

cDNA synthesis and quantitative real-time PCR (qPCR)

Total RNA extracted from cell samples was used for qPCR validation of sequencing data or other data. First strand cDNA was synthesized using the cDNA Synthesis Kit according to the manufacturer’s instructions and then was subjected to qPCR experiments on a Bio-Rad CFX96 Real-Time PCR Detection System. qPCR assays were executed as the manufacturer suggested. Detailed sequences of primers used for this study are listed in Additional file 10: Table S6, and GAPDH was used to normalize gene expression. Three biological and technical replicates were set for all experiments. The 2-ΔΔCt method was employed to calculate relative transcript expression.

KEGG pathway analysis

As reported in previous studies [90, 91], we used the KOBAS3.0 [92] program to arrange the differentially expressed transcripts into their involved KEGG pathways.

Construction of interactive networks among differentially expressed transcripts

To better understand the regulatory relationships of noncoding genes with mRNAs, and further reveal their functions in hair biology, we constructed interaction networks among these transcripts. As many researchers did before [2, 51], the predictions of the target transcripts of miRNAs, including lncRNAs and mRNAs, were achieved using miRanda [88], which is a program exclusively developed for animals. Based on co-location and co-expression, we constructed an interacting network of mRNAs-lncRNAs. Likewise, based on the prediction of miRNA binding sites inside mRNAs or lncRNAs, and their expression level relevance, we built networks of mRNAs-miRNAs and lncRNAs-miRNAs. The other types of interacting pairs were created with the popular ceRNAs hypothesis in which lncRNAs act as sponges to absorb miRNAs and then free its suppressive roles on mRNA transcription or translation. Finally, lncRNAs-miRNAs-mRNAs interactive lines were fabricated. All results were graphically displayed using Cytoscape3.6.

Availability of data and materials

The transcriptomic data generated for this study were deposited in the NCBI Sequence Read Archive (http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/sra). All accession numbers are listed as follows: SRR8138850, SRR8138851 and SRR8138852 for lncRNAs of DPCs; SRR8138849, SRR8138855 and SRR8138856 for lncRNAs of DFs; SRR8138853, SRR8138854 and SRR8138858 for miRNAs of DPCs; SRR8138857, SRR8138859 and SRR8138860 for miRNAs of DFs.

Abbreviations

BMP4 :

Bone morphogenetic protein 4

COL15A1 :

Collagen type XV alpha 1 chain

CRABP1 :

Cellular retinoic acid-binding protein 1

DFs:

Dermal fibroblasts

DPCs:

Dermal papilla cells

ESR1 :

Estrogen receptor alpha

FGF:

Fibroblast growth factor

FGFR2 :

Fibroblast growth factor receptor 2

FPKM:

Fragments per kilobase of transcript per million mapped reads

HF:

Hair follicle

HFSCs:

Hair follicle stem cell

HMCs:

Hair matrix cells

HOTAIR :

Transcript antisense RNA

HOX :

Homeobox protein

ITGA :

Integrin alpha

LEPR :

Leptin receptor

lncRNAs:

Long non-coding RNAs

LTBP1 :

Latent-transforming growth factor beta-binding

ncRNAs:

Non-coding RNAs

PRDM1 :

PR domain zinc finger protein 1

PTGER4 :

Prostaglandin EP4 receptor

PTGFR :

Prostaglandin F receptor

RSPO1 :

R-spondin 1

SOX2 :

SRY (sex determining region Y)-box 2

TMEM100 :

Transmembrane protein 100

WNT5A :

Wnt Family Member 5A

References

  1. Li C, Li Y, Zhou G, Gao Y, Ma S, Chen Y, et al. Whole-genome bisulfite sequencing of goat skins identifies signatures associated with hair cycling. BMC Genomics. 2018;19. https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-018-5002-5.

  2. Zhou G, Kang D, Ma S, Wang X, Gao Y, Yang Y, et al. Integrative analysis reveals ncRNA-mediated molecular regulatory network driving secondary hair follicle regression in cashmere goats. BMC Genomics. 2018;19. https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-018-4603-3.

  3. Liu Z, Yang F, Zhao M, Ma L, Li H, Xie Y, et al. The intragenic mRNA-microRNA regulatory network during telogen-anagen hair follicle transition in the cashmere goat. Sci Rep. 2018;8. https://0-doi-org.brum.beds.ac.uk/10.1038/s41598-018-31986-2.

  4. Yuan C, Wang X, Geng R, He X, Qu L, Chen Y. Discovery of cashmere goat (Capra hircus) microRNAs in skin and hair follicles by Solexa sequencing. BMC Genomics. 2013;14:511.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Wang S, Ge W, Luo Z, Guo Y, Jiao B, Qu L, et al. Integrated analysis of coding genes and non-coding RNAs during hair follicle cycle of cashmere goat (Capra hircus). BMC Genomics. 2017;18:767.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  6. Bai WL, Wang JJ, Yin RH, Dang YL, Wang ZY, Zhu YB, et al. Molecular characterization of HOXC8 gene and methylation status analysis of its exon 1 associated with the length of cashmere fiber in Liaoning cashmere goat. Genetica. 2017;145:115–26.

    Article  CAS  PubMed  Google Scholar 

  7. Paus R, Foitzik K. In search of the “hair cycle clock”: a guided tour. Differentiation. 2004;72:489–511.

    Article  CAS  PubMed  Google Scholar 

  8. Galbraith H. Fundamental hair follicle biology and fine fibre production in animals. animal. 2010;4:1490–509.

    Article  CAS  PubMed  Google Scholar 

  9. Paus R, Langan EA, Vidali S, Ramot Y, Andersen B. Neuroendocrinology of the hair follicle: principles and clinical perspectives. Trends Mol Med. 2014;20:559–70.

    Article  CAS  PubMed  Google Scholar 

  10. Yang C-C, Cotsarelis G. Review of hair follicle dermal cells. J Dermatol Sci. 2010;57:2–11.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Schmidt-Ullrich R, Paus R. Molecular principles of hair follicle induction and morphogenesis. BioEssays. 2005;27:247–61.

    Article  CAS  PubMed  Google Scholar 

  12. Clavel C, Grisanti L, Zemla R, Rezza A, Barros R, Sennett R, et al. Sox2 in the dermal papilla niche controls hair growth by fine-tuning BMP signaling in differentiating hair shaft progenitors. Dev Cell. 2012;23:981–94.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Tao Y, Yang Q, Wang L, Zhang J, Zhu X, Sun Q, et al. β-Catenin activation in hair follicle dermal stem cells induces ectopic hair outgrowth and skin fibrosis. J Mol Cell Biol. 2019;11:26–38.

    Article  PubMed  Google Scholar 

  14. Telerman SB, Rognoni E, Sequeira I, Pisco AO, Lichtenberger BM, Culley O, et al. Dermal Blimp1 acts downstream of epidermal TGFβ and Wnt/β-catenin to regulate hair follicle formation and growth. J Invest Dermatol. 2017;137:2270–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Fatica A, Bozzoni I. Long non-coding RNAs: new players in cell differentiation and development. Nat Rev Genet. 2014;15:7–21.

    Article  CAS  PubMed  Google Scholar 

  16. Quinn JJ, Chang HY. Unique features of long non-coding RNA biogenesis and function. Nat Rev Genet. 2016;17:47–62.

    Article  CAS  PubMed  Google Scholar 

  17. Kretz M, Siprashvili Z, Chu C, Webster DE, Zehnder A, Qu K, et al. Control of somatic tissue differentiation by the long non-coding RNA TINCR. Nature. 2012;493:231–5.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  18. Young TL, Matsuda T, Cepko CL. The noncoding RNA taurine upregulated gene 1 is required for differentiation of the murine retina. Curr Biol. 2005;15:501–12.

    Article  CAS  PubMed  Google Scholar 

  19. Sahakyan A, Yang Y, Plath K. The role of Xist in X-chromosome dosage compensation. Trends Cell Biol. 2018;28:999–1013.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Yang YW, Flynn RA, Chen Y, Qu K, Wan B, Wang KC, et al. Essential role of lncRNA binding for WDR5 maintenance of active chromatin and embryonic stem cell pluripotency. eLife. 2014;3:e02046.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  21. Kim J, Piao H-L, Kim B-J, Yao F, Han Z, Wang Y, et al. Long noncoding RNA MALAT1 suppresses breast cancer metastasis. Nat Genet. 2018. https://0-doi-org.brum.beds.ac.uk/10.1038/s41588-018-0252-3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Andl T, Murchison EP, Liu F, Zhang Y, Yunta-Gonzalez M, Tobias JW, et al. The miRNA-processing enzyme dicer is essential for the morphogenesis and maintenance of hair follicles. Curr Biol. 2006;16:1041–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Ahmed MI, Alam M, Emelianov VU, Poterlowicz K, Patel A, Sharov AA, et al. MicroRNA-214 controls skin and hair follicle development by modulating the activity of the Wnt pathway. J Cell Biol. 2014;207:549–67.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Goodarzi HR, Abbasi A, Saffari M, Fazelzadeh Haghighi M, Tabei MB, Noori Daloii MR. Differential expression analysis of balding and nonbalding dermal papilla microRNAs in male pattern baldness with a microRNA amplification profiling method: MicroRNA differential expression in baldness dermal papillae. Br J Dermatol. 2012;166:1010–6.

    Article  CAS  PubMed  Google Scholar 

  25. Midorikawa T, Chikazawa T, Yoshino T, Takada K, Arase S. Different gene expression profile observed in dermal papilla cells related to androgenic alopecia by DNA macroarray analysis. J Dermatol Sci. 2004;36:25–32.

    Article  CAS  PubMed  Google Scholar 

  26. Zhu N, Huang K, Liu Y, Zhang H, Lin E, Zeng Y, et al. miR-195-5p regulates hair follicle inductivity of dermal papilla cells by suppressing Wnt/β-catenin activation. Biomed Res Int. 2018;2018:1–13.

    Google Scholar 

  27. Zhou G, Yuan C, He X, Kang D, Wang X, Chen Y. Effect of miR-125b on dermal papilla cells of goat secondary hair follicle. Electron J Biotechnol. 2017;25:64–9.

    Article  CAS  Google Scholar 

  28. Lin C, Liu Y, Huang K, Chen X, Cai B, Li H, et al. Long noncoding RNA expression in dermal papilla cells contributes to hairy gene regulation. Biochem Biophys Res Commun. 2014;453:508–14.

    Article  CAS  PubMed  Google Scholar 

  29. Zheng Y, Wang Z, Zhu Y, Wang W, Bai M, Jiao Q, et al. LncRNA-000133 from secondary hair follicle of cashmere goat: identification, regulatory network and its effects on inductive property of dermal papilla cells. Anim Biotechnol. 2019:1–13. https://0-www-tandfonline-com.brum.beds.ac.uk/doi/abs/10.1080/10495398.2018.1553788.

  30. Liang W-C, Fu W-M, Wong C-W, Wang Y, Wang W-M, Hu G-X, et al. The lncRNA H19 promotes epithelial to mesenchymal transition by functioning as miRNA sponges in colorectal cancer. Oncotarget. 2015;6:22513–25.

    PubMed  PubMed Central  Google Scholar 

  31. Wan DC, Wang KC. Long noncoding RNA: significance and potential in skin biology. Cold Spring Harb Perspect Med. 2014;4:a015404.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  32. Smillie CL, Sirey T, Ponting CP. Complexities of post-transcriptional regulation and the modeling of ceRNA crosstalk. Crit Rev Biochem Mol Biol. 2018;53:231–45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Tay Y, Rinn J, Pandolfi PP. The multilayered complexity of ceRNA crosstalk and competition. Nature. 2014;505:344–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Zhu Z, Li Y, Liu W, He J, Zhang L, Li H, et al. Comprehensive circRNA expression profile and construction of circRNA-associated ceRNA network in fur skin. Exp Dermatol. 2018;27:251–7.

    Article  CAS  PubMed  Google Scholar 

  35. Zhu YB, Wang ZY, Yin RH, Jiao Q, Zhao SJ, Cong YY, et al. A lncRNA-H19 transcript from secondary hair follicle of Liaoning cashmere goat: identification, regulatory network and expression regulated potentially by its promoter methylation. Gene. 2018;641:78–85.

    Article  CAS  PubMed  Google Scholar 

  36. Jahoda CAB, Horne KA, Oliver RF. Induction of hair growth by implantation of cultured dermal papilla cells. Nature. 1984;311:560–2.

    Article  CAS  PubMed  Google Scholar 

  37. Jahoda C, Oliver RF. The growth of vibrissa dermal papilla cells in vitro. Br J Dermatol. 1981;105:623–7.

    Article  CAS  PubMed  Google Scholar 

  38. Katsuoka K, Schell H, Hornstein OP, Deinlein E, Wessel B. Comparative morphological and growth kinetics studies of human hair bulb papilla cells and root sheath fibroblasts in vitro. Arch Dermatol Res. 1986;279:20–5.

    Article  CAS  PubMed  Google Scholar 

  39. Ma S, Zhou G, Chen Y. Effects of all-trans retinoic acid on goat dermal papilla cells cultured in vitro. Electron J Biotechnol. 2018;34:43–50.

    Article  CAS  Google Scholar 

  40. Zhu B, Xu T, Yuan J, Guo X, Liu D. Transcriptome sequencing reveals differences between primary and secondary hair follicle-derived dermal papilla cells of the cashmere goat (Capra hircus). PLoS One. 2013;8:e76282.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Bickhart DM, Rosen BD, Koren S, Sayre BL, Hastie AR, Chan S, et al. Single-molecule sequencing and chromatin conformation capture enable de novo reference assembly of the domestic goat genome. Nat Genet. 2017;49:643–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Pertea M, Kim D, Pertea GM, Leek JT, Salzberg SL. Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat Protoc. 2016;11:1650–67.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28:511–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  44. Kong L, Zhang Y, Ye Z-Q, Liu X-Q, Zhao S-Q, Wei L, et al. CPC: assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007;35(suppl_2):W345–9.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Finn RD, Coggill P, Eberhardt RY, Eddy SR, Mistry J, Mitchell AL, et al. The Pfam protein families database: towards a more sustainable future. Nucleic Acids Res. 2016;44:D279–85.

    Article  CAS  PubMed  Google Scholar 

  46. Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, et al. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013;41:e166.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Wen M, Shen Y, Shi S, Tang T. miREvo: an integrative microRNA evolutionary analysis platform for next-generation sequencing experiments. BMC Bioinformatics. 2012;13:140.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Friedländer MR, Mackowiak SD, Li N, Chen W, Rajewsky N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. 2012;40:37–52.

    Article  PubMed  CAS  Google Scholar 

  49. Nie Y, Li S, Zheng X, Chen W, Li X, Liu Z, et al. Transcriptome reveals long non-coding RNAs and mRNAs involved in primary wool follicle induction in carpet sheep fetal skin. Front Physiol. 2018;9. https://0-doi-org.brum.beds.ac.uk/10.3389/fphys.2018.00446.

  50. Wang J, Hua L, Chen J, Zhang J, Bai X, Gao B, et al. Identification and characterization of long non-coding RNAs in subcutaneous adipose tissue from castrated and intact full-sib pair Huainan male pigs. BMC Genomics. 2017;18:542.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  51. Gao X, Ye J, Yang C, Zhang K, Li X, Luo L, et al. Screening and evaluating of long noncoding RNAs in the puberty of goats. BMC Genomics. 2017;18. https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-017-3578-9.

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

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  53. Frazee AC, Pertea G, Jaffe AE, Langmead B, Salzberg SL, Leek JT. Ballgown bridges the gap between transcriptome assembly and expression analysis. Nat Biotechnol. 2015;33:243–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Driskell RR, Juneja VR, Connelly JT, Kretzschmar K. Tan DW-M, watt FM. Clonal growth of dermal papilla cells in hydrogels reveals intrinsic differences between Sox2-positive and -negative cells in vitro and in vivo. J Invest Dermatol. 2012;132:1084–93.

    Article  CAS  PubMed  Google Scholar 

  55. Sorrell JM. Fibroblast heterogeneity: more than skin deep. J Cell Sci. 2004;117:667–75.

    Article  CAS  PubMed  Google Scholar 

  56. Rendl M, Lewis L, Fuchs E. Molecular dissection of mesenchymal–epithelial interactions in the hair follicle. PLoS Biol. 2005;3:e331.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  57. Movérare S, Lindberg MK, Ohlsson C, Faergemann J, Gustafsson J-Å. Estrogen receptor α, but not estrogen receptor β, is involved in the regulation of the hair follicle cycling as well as the thickness of epidermis in male mice. J Invest Dermatol. 2002;119:1053–8.

    Article  PubMed  Google Scholar 

  58. Torii E, Segi E, Sugimoto Y, Takahashi K, Kabashima K, Ikai K, et al. Expression of prostaglandin E2 receptor subtypes in mouse hair follicles. Biochem Biophys Res Commun. 2002;290:696–700.

    Article  CAS  PubMed  Google Scholar 

  59. Gao Y, Wang X, Yan H, Zeng J, Ma S, Niu Y, et al. Comparative transcriptome analysis of fetal skin reveals key genes related to hair follicle morphogenesis in cashmere goats. PLoS One. 2016;11:e0151118.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  60. Yu Z, Jiang K, Xu Z, Huang H, Qian N, Lu Z, et al. Hoxc-dependent mesenchymal niche heterogeneity drives regional hair follicle regeneration. Cell Stem Cell. 2018. https://0-doi-org.brum.beds.ac.uk/10.1016/j.stem.2018.07.016.

    PubMed  Google Scholar 

  61. Ha M, Kim VN. Regulation of microRNA biogenesis. Nat Rev Mol Cell Biol. 2014;15:509–24.

    Article  CAS  PubMed  Google Scholar 

  62. Beermann J, Piccoli M-T, Viereck J, Thum T. Non-coding RNAs in development and disease: background, mechanisms, and therapeutic approaches. Physiol Rev. 2016;96:1297–325.

    Article  CAS  PubMed  Google Scholar 

  63. Li N, Liu S, Zhang H-S, Deng Z-L, Zhao H-S, Zhao Q, et al. Exogenous R-Spondin1 induces precocious Telogen-to-Anagen transition in mouse hair follicles. Int J Mol Sci. 2016;17:582.

    Article  PubMed Central  CAS  Google Scholar 

  64. Karreth FA, Pandolfi PP. ceRNA cross-talk in Cancer: when ce-bling rivalries go awry. Cancer Discov. 2013;3:1113–21.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  65. Jiapaer Z, Li G, Ye D, Bai M, Li J, Guo X, et al. LincU preserves Na ï ve pluripotency by restricting ERK activity in embryonic stem cells. Stem Cell Rep. 2018. https://0-doi-org.brum.beds.ac.uk/10.1016/j.stemcr.2018.06.010.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  66. Sennett R, Rendl M. Mesenchymal–epithelial interactions during hair follicle morphogenesis and cycling. Semin Cell Dev Biol. 2012;23:917–27.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Stenn KS. Molecular insights into the hair follicle and its pathology: a review of recent developments [J]. Int J Dermatol. 2003;42(1):40–43.

    Article  CAS  PubMed  Google Scholar 

  68. Magerl M, Kauser S, Paus R, Tobin DJ. Simple and rapid method to isolate and culture follicular papillae from human scalp hair follicles. Exp Dermatol. 2002;11:381–5.

    Article  PubMed  Google Scholar 

  69. Thornton MJ, Kato S, Hibberts NA, Brinklow BR, Loudon ASI, Randall VA. Ability to culture dermal papilla cells from red deer (Cervus elaphus) hair follicles with differing hormonal responses in vivo offers a new model for studying the control of hair follicle biology. J Exp Zool. 1996;275:452–8.

    Article  CAS  PubMed  Google Scholar 

  70. Withers AP, Jahoda CAB, Ryder ML, Oliver RF. Culture of wool follicle dermal papilla cells from two breeds of sheep. Arch Dermatol Res. 1986;279:140–2.

    Article  CAS  PubMed  Google Scholar 

  71. Reynolds AJ, Jahoda CA. Cultured dermal papilla cells induce follicle formation and hair growth by transdifferentiation of an adult epidermis. Development. 1992;115:587–93.

    CAS  PubMed  Google Scholar 

  72. Rezza A, Wang Z, Sennett R, Qiao W, Wang D, Heitman N, et al. Signaling networks among stem cell precursors, transit-amplifying progenitors, and their niche in developing hair follicles. Cell Rep. 2016;14:3001–18.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Salmena L, Poliseno L, Tay Y, Kats L, Pandolfi PP. A ceRNA hypothesis: the Rosetta stone of a hidden RNA language? Cell. 2011;146:353–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Lynch P, Russel AJF. Hormonal manipulation of cashmere growth and shedding. Proc Br Soc Anim Prod 1972. 1990;1990:58.

    Google Scholar 

  75. Lee J, Tumbar T. Hairy tale of signaling in hair follicle development and cycling. Semin Cell Dev Biol. 2012;23:906–16.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. Iguchi M, Aiba S, Yoshino Y, Tagami H. Human follicular papilla cells carry out nonadipose tissue production of leptin. J Invest Dermatol. 2001;117:1349–56.

    Article  CAS  PubMed  Google Scholar 

  77. Watabe R, Yamaguchi T, Kabashima-Kubo R, Yoshioka M, Nishio D, Nakamura M. Leptin controls hair follicle cycling. Exp Dermatol. 2014;23:228–9.

    Article  CAS  PubMed  Google Scholar 

  78. van Beek N, Bodó E, Kromminga A, Gáspár E, Meyer K, Zmijewski MA, et al. Thyroid hormones directly Alter human hair follicle functions: Anagen prolongation and stimulation of both hair matrix keratinocyte proliferation and hair pigmentation. J Clin Endocrinol Metab. 2008;93:4381–8.

    Article  PubMed  CAS  Google Scholar 

  79. Oh HS, Smart RC. An estrogen receptor pathway regulates the telogen-anagen hair follicle transition and influences epidermal cell proliferation. Proc Natl Acad Sci. 1996;93:12525–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  80. Björnström L, Sjöberg M. Mechanisms of estrogen receptor signaling: convergence of genomic and nongenomic actions on target genes. Mol Endocrinol. 2005;19:833–42.

    Article  PubMed  CAS  Google Scholar 

  81. Jiang L, Li Z, Rui L. Leptin stimulates both JAK2-dependent and JAK2-independent signaling pathways. J Biol Chem. 2008;283:28066–73.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  82. Li L, Clevers H. Coexistence of quiescent and active adult stem cells in mammals. Science. 2010;327:542–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  83. Brenet F, Moh M, Funk P, Feierstein E, Viale AJ, Socci ND, et al. DNA methylation of the first exon is tightly linked to transcriptional silencing. PLoS One. 2011;6:e14524.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  84. Plikus MV, Mayer JA, de la Cruz D, Baker RE, Maini PK, Maxson R, et al. Cyclic dermal BMP signalling regulates stem cell activation during hair regeneration. Nature. 2008;451:340–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  85. Luchetti F, Canonico B, Bartolini D, Arcangeletti M, Ciffolilli S, Murdolo G, et al. Melatonin regulates mesenchymal stem cell differentiation: a review. J Pineal Res. 2014;56:382–97.

    Article  CAS  PubMed  Google Scholar 

  86. Woo CJ, Kingston RE. HOTAIR lifts noncoding RNAs to new levels. Cell. 2007;129:1257–9.

    Article  CAS  PubMed  Google Scholar 

  87. Wang Y, Xu Z, Jiang J, Xu C, Kang J, Xiao L, et al. Endogenous miRNA sponge lincRNA-RoR regulates Oct4, Nanog, and Sox2 in human embryonic stem cell self-renewal. Dev Cell. 2013;25:69–80.

    Article  CAS  PubMed  Google Scholar 

  88. Enright AJ, John B, Gaul U, Tuschl T, Sander C, Marks DS. MicroRNA targets in drosophila. Genome Biol. 2003;5:R1.

    Article  PubMed  PubMed Central  Google Scholar 

  89. Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11:R106.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  90. Yue Y, Guo T, Yuan C, Liu J, Guo J, Feng R, et al. Integrated analysis of the roles of long noncoding RNA and coding RNA expression in sheep (Ovis aries) skin during initiation of secondary hair follicle. PLoS One. 2016;11:e0156890.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  91. Li Y, Zhou G, Zhang R, Guo J, Li C, Martin G, et al. Comparative proteomic analyses using iTRAQ-labeling provides insights into fiber diversity in sheep and goats. J Proteome. 2018;172:82–8.

    Article  CAS  Google Scholar 

  92. Xie C, Mao X, Huang J, Ding Y, Wu J, Dong S, et al. KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic Acids Res. 2011;39(suppl_2):W316–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

Not applicable.

Funding

The present work was supported by the National Natural Science Foundation of China (Program No. 31872332, 31772571, 31572369, and 31372279), the National Key Research and Development Program of China (Program No. 2016YFD0500508), and China Agriculture Research System (Program No. CARS-39-12). The funders of this study had no role in study design, data collection, data analysis and interpretation of data and in writing the manuscript.

Author information

Authors and Affiliations

  1. College of Animal Science and Technology, Northwest A&F University, Yangling, 712100, Shaanxi, China

    Sen Ma, Ying Wang, Guangxian Zhou, Yi Ding, Yuxin Yang, Xiaolong Wang, Enping Zhang & Yulin Chen

  2. Department of Animal Science, Guangdong Ocean University, Zhanjiang, 524088, Guangdong, China

    Guangxian Zhou

Authors
  1. Sen Ma
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Ying Wang
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Guangxian Zhou
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Yi Ding
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Yuxin Yang
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Xiaolong Wang
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Enping Zhang
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Yulin Chen
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Conceived and designed the experiments: SM YY XW EZ YC. Performed the experiments: SM GZ YD YW. Wrote the manuscript: SM. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Enping Zhang or Yulin Chen.

Ethics declarations

Ethics approval and consent to participate

All experimental procedures in this study were approved by and performed according to the guidelines for the care and use of experimental animals established by the Animal Care and Use Committee of the Northwest A&F University (Approval ID: 2013–31101684).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

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

Additional files

Additional file 1:

Figure S1. Identification of lncRNAs in goat DPCs and DFs. a Five successive screening steps of valid lncRNAs transcripts from all samples: step 1: exon counts (> 2); step 2: transcript length (> 200 nt); step 3: expression abundance (FPKM > 0.5); step 4: preclusion of known lncRNAs transcripts in the database and step 5: the absence of protein-coding potential. b Coding potential prediction of lncRNAs candidates using CPC, CNCI and PFAM tools. Overall, 1754 lncRNAs transcripts that are not supposed to code protein according to any tool were thought to be novel lncRNAs (overlapped region) and proceeded to downstream analysis. (TIF 6882 kb)

Additional file 2:

Table S1. Information of lncRNAs transcripts. (XLSX 209 kb)

Additional file 3:

Figure S3. Comparison of the genomic features of lncRNAs and mRNAs. Density distribution of exon counts (a) and transcript length (b) in annotated lncRNAs, novel lncRNAs and mRNAs. Usually, mRNAs contain more exons and nucleotides than lncRNAs. c Box plot showing the abundances of transcripts, and the expression levels of mRNAs are higher than annotated lncRNAs or novel lncRNAs as indicated. (TIF 8461 kb)

Additional file 4:

Figure S4. Heatmaps of differentially expressed transcripts. a Heatmap of differentially mRNAs. b Heatmap of differentially lncRNAs. c Heatmap of differentially miRNAs. Compared with DFs, 1286 mRNAs, 18 lncRNAs, and 42 miRNAs were upregulated, while 1254 mRNAs, 53 lncRNAs and 44 miRNAs were downregulated in DPCs. (TIF 22373 kb)

Additional file 5:

Table S2. Lists of differentially expressed transcripts. (XLSX 261 kb)

Additional file 6:

Figure S5. Validation of relative abundances of transcripts in DPCs and DFs by qPCR. GAPDH was adopted as the internal reference to normalize all gene expression. Relative gene expressions were acquired through the classic comparative cycle threshold method. Fold change was calculated as the mean of three independent samples. Similar trends were found in all examined genes. (TIF 3019 kb)

Additional file 7:

Table S3. Enriched pathways of upregulated genes in DPCs. (XLS 144 kb)

Additional file 8:

Table S4. Interactive relationship of ncRNAs with core signatures. (XLSX 60 kb)

Additional file 9:

Table S5. Interactive lines of ncRNAs with core signatures. (XLSX 38 kb)

Additional file 10:

Table S6. Primer information used in present study. (XLSX 9 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

Ma, S., Wang, Y., Zhou, G. et al. Synchronous profiling and analysis of mRNAs and ncRNAs in the dermal papilla cells from cashmere goats. BMC Genomics 20, 512 (2019). https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-019-5861-4

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://0-doi-org.brum.beds.ac.uk/10.1186/s12864-019-5861-4

Keywords

BMC Genomics

ISSN: 1471-2164

Contact us