This article is mentioned in:

Introduction

The domestic goat (Capra hircus) is reared throughout the world, particularly in China, India and other developing countries (1), and serves as an important source of fibers and pelts (2). The Laiwu black goat, of which there are 230,000, is indigenous to Shandong and supplies cashmere, pelts and meat. The advantages of this breed include the high-quality cashmere produced, resistance to disease and a high reproductive rate. With regards to cashmere, the fleece fiber diameter is <12 μm, cashmere length is up to 56 mm, and the maximum cashmere yield is 300 g (3). The cashmere growth in Laiwu black goats begins in late July, continuing to March (3).

The identification of genes involved in cashmere production may provide an opportunity to improve production efficiency and product quality and diversity. This may be accomplished through breeding programs, the development of transgenic lines or via the manufacture of therapeutic agents that improve fibers by altering gene expression (4). Previous studies have investigated the effect of genetic polymorphisms (single nucleotide polymorphisms and quantitative trait loci) on cashmere growth and regulation (5,6). Goats have two types of hair follicles: The primary hair follicle, which produces coarse coat hair and the secondary hair follicle, which produces cashmere (1). Comparative transcriptomic analysis identified 51 genes that were differentially expressed (DE) in primary and secondary hair follicles (1). In addition, microRNAs and protein coding genes with potential roles in goat and sheep hair growth have been identified (7–13). The present study used the Agilent sheep gene expression array to identify genes and proteins associated with wool follicle growth and cycling (14,15). The Agilent sheep gene expression array has previously been performed on goat samples to investigate the contribution of mammary epithelial cells to the immune response during the early stages of Staphylococcus aureus infection (16).

To date, to the best of our knowledge, no microarray or other transcriptomic analysis has been conducted to investigate hair and cashmere growth in various skin areas of adult goats. In addition, only a limited number of studies have described the molecular characteristics of cashmere growth in the Laiwu black goat (11,12,17). The aim of the present study was to analyze and compare gene expression levels in body (hair and cashmere rich) and groin skin (no hair or cashmere) during the cashmere growth period using microarray technology. This may enable the identification of genes potentially associated with hair and cashmere growth regulation in the Laiwu black goat.

Materials and methods

Sampling

The protocol of the present study was approved by the Ethics Committee of Qingdao Agricultural University (Qingdao, China). The goats were obtained from Qingdao Aote Sheep Farm (Qingdao, China), and housed in semi-closed sheep houses. Sampling was performed as described previously (14). Sampling was performed in August, during the cashmere growth period, on one ram and two ewes, aged 2 years. Full thickness skin from the body (cashmere/wool-bearing) and groin (no cashmere/wool) were sampled from the three animals under local anesthesia (2 ml 0.5% procaine hydrochloride; Harbin Pharmaceutical Group Co., Ltd., Harbin, China) for complementary DNA microarray experiments.

RNA preparation, microarray hybridization and data analysis

Total RNA was isolated from the skin samples using TRIzol® reagent (Invitrogen; Thermo Fisher Scientific, Inc., Waltham, MA, USA) according to the manufacturer's instructions. RNA integrity was confirmed by denaturing agarose gel electrophoresis and RNA was quantified using a NanoDrop® ND-1000 Spectrophotometer (NanoDrop Technologies; Thermo Fisher Scientific, Inc.). Total RNA (2.5 μg) from each sample was doubled and transcribed into fluorescent complementary RNA using a Quick Amp Labeling kit (Agilent Technologies, Inc., Santa Clara, CA, USA) according to the manufacturer's instructions. Subsequently, RNA samples were incubated with the Agilent Sheep Gene Expression Microarray (Agilent Technologies, Inc.). The microarray signals were scanned and analyzed as previously described (14). Clustering analysis of DE genes was performed using Cluster software version 3.0 (bonsai.hgc.jp/~mdehoon/software/cluster/software.htm) (18,19) to analyze the similarity in the expression patterns among different skin areas. The functional annotation of DE genes was performed using the Database for Annotation, Visualization and Integrated Discovery (DAVID) gene annotation tool version 6.6 (david.abcc.ncifcrf.gov/) (20).

Reverse transcription-quantitative polymerase chain reaction (RT-qPCR) verification

Equal amounts of RNA were reverse-transcribed using Superscript III (Invitrogen; Thermo Fisher Scientific, Inc.). qPCR was performed with SYBR Green Master Mix (Roche Applied Science, Penzberg, Germany) and gene-specific primers. The cycling conditions consisted of an initial, single cycle for 5 min at 95°C, followed by 40 cycles of 20 sec at 95°C, 30 sec at 60°C and 30 sec at 72°C. The primer sequences, melting temperatures and product sizes are listed in Table I. GAPDH served as the housekeeping gene, and results were normalized using the ΔΔCq method (21).

Table I Primers used for qPCR validation.

Table I

Primers used for qPCR validation.

Gene	Primer sequence (5′–3′)		Tm, °Ca	Target size, bp
	Forward	Reverse
GAPDHb	ACAGTCAAGGCAGAGAACGG	CCAGCATCACCCCACTTGAT	60	98
AMP18 (GKN1)	TCAAGCCCTTGGTATGCTGG	TGAAGTCCGGCTTCTTGGTC	60	198
CYP1A1	TTTCACCCTCGCTCTGAAGG	AAGTTCTGTGGCCGAGATGG	60	190
Connexin 43	TTGTACCCGGGAGGAGACAT	CTGAGCCCCTCCAAAGACTG	60	101
C-jun	GGATCAAGGCGGAGAGGAAG	CTGCGTTAGCATGAGTTGGC	60	224
ITGB1	AGCACGGATGAGGTGAACAG	ATCTCACAGGTTGGCCCTTG	60	389
C3	AACAAACGGGATCCCCTGAC	GAGTTCCCCTGCGTGTTGTA	60	113

a Tm represents the optimal temperature during PCR

b RNA levels of GAPDH were assayed for normalization. qPCR, quantitative polymerase chain reaction; Tm, annealing temperature; AMP18, antrum muscle protein 18; GKN1, gastrokine 1; CYP1A1, cytochrome P450 family 1 subfamily A member 1; ITGB1, integrin subunit β1; C3, complement component 3.

Network analysis

Groups of genes were analyzed using the Meme Software Suite version 4.10.2 (meme-suite.org/) to identify common motifs and cis-regulatory elements (22). Common motifs and cis-regulatory elements were searched by Patch in the TRANSFAC® Public version 7.0 database (gene-regulation.com/pub/databases.html) to identify transcription factor binding sites (23).

Statistical analysis

Comparisons between groups were performed using Student's t-tests in Microsoft Excel 2007 (Microsoft Corporation, Redmond, WA, USA). Data are expressed as the mean ± standard deviation. P≤0.05 was considered to indicate a statistically significant difference.

Results

Comparative transcriptome analysis between body and groin skin

Skin samples collected from the two skin areas were analyzed. A total of 15,008 transcripts (98.68% of all probe sets) were detected in the two skin areas (data not shown). Following normalization and statistical analyses, 655 genes were identified as DE, with 217 upregulated and 438 downregulated in body compared with groin skin, over the threshold of fold-change (>2.0; data not shown).

To further investigate the similarity in expression patterns of gene transcripts in the two skin areas, cluster analysis was performed using the Cluster 3.0 tool. As presented in Fig. 1, cluster analysis revealed differences between body and groin skin.

Figure 1 Hierarchical cluster analysis of genes differentially expressed in two skin areas of the Laiwu black goat. Each column represents one goat, and each row refers to a gene. Red indicates genes with a greater expression relative to the geometrical means and green indicates genes with a reduced expression relative to the geometrical means. LF1, LF2, LF3 represent the groin skin from three goats and IJ1, IJ2 and IJ3 represent the body skin from three goats.

Using DAVID, 22 of the DE genes were classified into 16 categories according to their functional correlation (data not shown). A number of these genes could be classified into more than one category. The majority of the genes associated with hair growth regulation could be assigned into the intracellular, intracellular organelle, membrane-bound vesicle, cytoplasmic vesicle, pattern binding, heparin binding, polysaccharide binding, glycosaminoglycan binding and cytoplasmic membrane-bound vesicle categories. Other DE genes, not classified in DAVID, are discussed below.

Confirmation of DE genes by RT-qPCR

To confirm the microarray data, six DE genes were selected for assessment of their relative expression levels by RT-qPCR: Antrum muscle protein 18 (gastrokine 1; GKN1; Fig. 2A), cytochrome P450 family 1 subfamily A member 1 (CYP1A1; Fig. 2B), C-jun (Fig. 2C), Connexin 43 (Fig. 2D), integrin subunit β1 (ITGB1; Fig. 2E) and complement component 3 (C3; Fig. 2F). GKN1, CYP1A1 and C-jun were significantly upregulated in body compared with groin skin (P=0.013, P=0.025 and P=0.005, respectively). Connexin 43, ITGB1 and C3 were significantly downregulated in body compared with groin skin (P=0.041, P=0.022 and P=0.038, respectively). RT-qPCR results for all six genes were consistent with the microarray data, therefore validating its reliability.

Figure 2 RT-qPCR verification of the microarray data. RT-qPCR was performed on skin samples to assess the relative expression of six differentially expressed genes. Student's t-tests were conducted to compare relative gene expression in body and groin skin. (A) GKN1 (P=0.013). (B) CYP1A1 (P=0.025). (C) C-jun (P=0.005). (D) Connexin 43 (P=0.041). (E) ITGB1 (P=0.022). (F) C3 (P=0.038). RT-qPCR, reverse transcription polymerase chain reaction; GKN1, gastrokine 1; CYP1A1, cytochrome P450 family 1 subfamily A member 1; ITGB1, integrin subunit β1; C3, complement component 3; B, body skin; G, groin skin.

GKN1 (8.14-fold) was one of the most upregulated genes in body compared with groin skin in the present study (data not shown). The other two genes upregulated in body compared with groin skin, CYP1A1 (2.04-fold) and C-jun (2.04-fold), are involved in skin disease development and epidermal keratinocyte survival and differentiation (24,25). Connexin 43 (−19.35-fold) and ITGB1 (−2.45-fold) are associated with differentiation of hair follicle stem cells (26,27). The downregulated C3 gene (−8.92-fold) is a core member of the complement system (28).

Analysis of the region 1,800 bp upstream of the translation start site of DE genes

As presented in Fig. 3, three common motifs [41 bp (Fig. 3A), 50 bp (Fig. 3B) and 50 bp (Fig. 3C)] were identified in the upstream sequences of genes upregulated in body compared with groin skin. Notably, nuclear factor 1A (NF-1A), Yi, E2 factor (E2F) and cyclic adenosine monophosphate response element binding (CREB)/CREBβ binding sites are present in the promoters of a number of upregulated genes. The promoter regions of genes downregulated in body compared with groin skin demonstrate 3 common motif sequences of 29 bp (Fig. 4A), 39 bp (Fig. 4B) and 50 bp (Fig. 4C). In these regions, binding sites were identified for the transcription factors LF-A1, Yi, E2F, Collier/Olfactory-1/early B-cell factor 1 (COE1) and peroxisome proliferator-activated receptor α (PPARα). In addition, a binding site (the U site) was identified from the promoter regions of 12 downregulated genes, which has been previously demonstrated to be present in the olfactory cyclic nucleotide gated channel and type III cyclase promoters (29).

Figure 3 Demonstrates 3 common motifs in 1,800 base pairs upstream of the translation start site in genes upregulated in body compared with groin skin. Strands + and − denote the sense and antisense DNA sequence. Start denotes the number of base pairs upstream of the transcriptional start site. P-value denotes the probability of the motif occurring by random chance. At the bottom of each panel the cis-regulatory elements for various transcription factors are presented, which may bind to the common motif of the genes listed. Common motifs containing (A) NF-1A, (B) Yi and E2F, and (C) CREB/CREBβ are represented. NF-1A, nuclear factor 1A; E2F, E2 factor; CREB, cyclic adenosine monophosphate response element binding.

Figure 4 Demonstrates 3 common motifs in 1,800 base pairs upstream of the translation start site in genes downregulated in body compared with groin skin. Strands + and − denote the sense and antisense DNA sequence. Start denotes the number of base pairs upstream of the transcriptional start site. P-value denotes the probability of the motif occurring by random chance. At the bottom of each panel the cis-regulatory elements for various transcription factors are presented, which may bind to the common motif of the genes listed. Common motifs containing (A) LF-A1, E2F and Yi, (B) COE1, and (C) PPARα are represented. E2F, E2 factor; COE1, Collier/Olfactory-1/early B-cell factor 1; PPARα, peroxisome proliferator-activated receptor α.

Discussion

In the present study, the molecular events associated with wool and cashmere growth control in the Laiwu black goat were investigated using microarray technology. Statistical analysis revealed hundreds of genes DE (>2-fold-change) between the two sampled skin areas. RT-qPCR confirmed the reliability of the microarray results.

Heat shock protein 70.1 (Hsp70.1) is necessary and sufficient to accelerate depigmentation in vitiligo-prone Pmel-1 mice (30), while at least three different mutations of Agouti signaling peptide (ASIP) were hypothesized to cause the recessive black coat color pattern in sheep (31–33). A gene duplication affecting expression of the ovine ASIP gene is responsible for white and black sheep (32). In addition, Hsp70.1 was downregulated in body compared with groin skin of the Aohan fine-wool sheep (34). In the present study, Hsp70.1 in body skin may facilitate the regulation of the depigmentation of wool fiber. Upregulation of ASIP may have a similar role.

Overexpression of Hsp70.1 alone significantly inhibits aminoglycoside-induced hair cell death in mice (35). However, in the present study Hsp70.1 was downregulated in body skin of the Laiwu black goat. The reason for this contradiction remains to be elucidated.

Matrix metalloproteinase 2 (MMP2) activity is associated with the disappearance of collagen VII during the invasion of epithelial cords of hair follicles and sweat glands in human skin (36). MMP2 is involved in hair growth-associated extracellular matrix remodeling and cell migration, and may be a downstream effector, through which thymosin β4 exerts its effect on hair growth (37). However, in the present study MMP2 was significantly downregulated in body compared with groin skin. A similar result was observed in the Aohan fine wool sheep (38), suggesting that MMP2 functions may differ in goat and sheep.

Integrins, including α3, α11 and β3, have been demonstrated to be DE between primary and secondary hair follicles (11). In keratinocytes there was revealed to be a log-linear association between the relative expression levels of β1 integrins (ITGB1) on the cell surface and proliferative capacity (39). In addition, ITGB1-mediated signaling is important in human hair growth control (40). Skin and hair follicle integrity is crucially dependent on ITGB1 expression in keratinocytes (41). However, ITGB1 expression in body skin of the Aohan fine wool sheep was downregulated >10 fold (38). The results of the present study were consistent with this, although the fold change was reduced at −2.45. The underlying mechanism requires further investigation.

Rowe et al (24) identified the primary CYP1A1 location as the sebaceous gland surrounding the hair shaft. Consistent with this, in the present study CYP1A1 expression was upregulated in body skin according to the microarray and RT-qPCR data.

Connexin 43 is a gap junction protein expressed in the dermal papilla (26). The G60S connexin 43 mutant regulates hair growth and hair fiber morphology in a mouse model of human oculodentodigital dysplasia (27). However, connexin 43 was downregulated in body skin in the microarray and RT-qPCR data of the present study. A similar pattern was observed in the Aohan fine wool sheep (38).

Fibroblast growth factor (FGF) 2 and −10 were upregulated and downregulated respectively in body compared with groin skin. FGF signaling is crucial for hair growth regulation (42,43). FGF5, expressed in the outer root sheath, has been demonstrated to control anagen-catagen transition (44) via the inhibition of hair elongation (45). FGF21 and −22, which promote the transition of catagen (46), were reportedly down-regulated in primary compared with secondary hair follicles of the Yunnan black goats (1). FGF expression varied throughout the hair follicle cycle, reflecting the role of FGFs in regulatory and developmental processes, which include patterning, morphogenesis, proliferation, differentiation and migration of cells (42). Expression analysis suggested an important role for FGF7 and −10, as well as their cognate receptor FGFR2-IIIb, in the proliferation and differentiation of the mature hair follicle (46).

Homeobox (Hox) A7 and HoxC6 were downregulated in body compared with groin skin in the present study. The first Hox gene demonstrated to have a universal role in hair follicle development was Hoxc13, as Hoxc13-deficient and -overexpressing mice exhibit severe hair growth and patterning defects (47). Members of the Hox family may perform essential and functionally diverse roles in hair that require complex transcriptional control mechanisms to ensure correct spatio-temporal patterns of Hox gene expression at homeostatic levels (47,48). HOXA7 silences differentiation-specific genes during keratinocyte proliferation (49), which is consistent with the downregulation of HOXA7 in body skin compared with groin skin in the present study.

Growth hormone (GH) and GH receptor (GHR) were upregulated in body compared with groin skin, consistent with results in Aohan fine wool sheep (38). A number of growth factors have been demonstrated to be critical for hair growth, including FGFs, transforming growth factors, insulin-like growth factors (IGFs), epidermal growth factors and hepatocyte growth factor (50). The hair cycle in the dorsal skin of male GH-deficient rats enters a long-lasting telogen phase at 8 weeks of age; however, depilation induces a transient hair cycle (51). Therefore, GH and GHR may be involved in wool and cashmere growth regulation.

IGF2 is DE in anagen, catagen and telogen stages in the Shaanbei White cashmere goat (9). IGF binding protein (IGFBP) 5-mediated FGFR2-IIIb signaling is a critical regulator of the structure of the hair shaft medulla (46). The expression of IGFBP3 and IGFBP5 was decreased in early anagen and anagen phases (IGFBP3) or in catagen and telogen phases (IGFBP5) during a depilation-induced hair cycle (51). Downregulation of FGFR2-IIIb, an effect associated with increased expression levels of IGFBP5, was demonstrated to decrease the thickness of the hair shaft and to reduce the cellularity of the hair shaft medulla (46). In accordance with these results, IGFBP2 was downregulated in body compared with groin skin in the microarray data of the present study.

In addition, the present study identified NF-1A, Yi, E2F and CREB/CREBβ binding sites in the promoter region of a number of upregulated genes. NF-1A negatively regulates target gene expression via binding to a silencer element (52); therefore, the reduced expression level of the gene subset may be caused by NF-1A upregulation or activation. CREB is a transcription factor important for keratinocyte proliferation (53–56). Adenosine stimulates growth of hair follicles by triggering phosphorylation of CREB (25). Combined recruitment of CREB, CCAAT-enhancer-binding protein β and c-jun determines the activation of promoters upon keratinocyte differentiation (25). In the present study, no differential expression was observed in CREB between the two skin regions. Therefore, CREB may have important additional roles in wool and cashmere growth, potentially via phosphorylation or other modifications.

The promoter regions of downregulated genes demonstrated three common motif sequences of 29, 39 and 50 bp. In these regions, binding sites were identified for the transcription factors LF-A1, Yi, E2F, COE1 and PPARα. The transcription factor LF-A1 is required for the cell-specific expression of the human α1-antitrypsin gene in hepatocytes (57). Binding sites of LF-A1 are present in the promoter regions of various genes expressed in the liver (α1-antitrypsin, apolipoproteins A1, B1 and A4, and pyruvate kinase). In the present study, thioredoxin-related transmembrane protein 2, adiponectin, C1Q and collagen domain containing and perilipin 3 may associate with apolipoproteins. COE1 is essential for B-cell differentiation (58). Therefore, COE1 may be involved in the immune privilege process as described previously (14). Furthermore, PPARα has been reported to contribute to hair growth and epidermal healing (59).

Notably, putative Yi and E2F sites were side-by-side in the upstream regions of upregulated and downregulated genes. Yi is an inducible DNA binding activity (60), consistent with its involvement in upregulated and downregulated genes. E2F has important roles in hair follicle growth, differentiation and survival (61,62). E2F1-deficient mice have a high incidence of spontaneous epidermal tumors of hair follicle origin (63). However, the reason for the existence of putative E2F sites in upregulated and downregulated genes requires further investigation. The U site was present in the promoter regions of 12 downregulated genes, of which 11 contained COE1 binding sites, consistent with a previous study (29).

A previous study in sheep revealed that the number and fold-change of DE genes in December were markedly reduced compared with August (21), which is consistent with the decreased activity of hair follicles in winter. Di et al (64) investigated DE genes of skin tissue in fine wool sheep with various fiber diameters using the Agilent sheep gene expression microarray. Therefore, the present study validates the use of Agilent Sheep Gene Expression Microarray in wool follicle investigations and provides a strong rationale for future studies.

Brenaut et al (16) investigated the contribution of mammary epithelial cells to the early stages of an immune response against Staphylococcus aureus in goats using the Agilent sheep gene expression microarray. Goyal et al (65) conducted microarray with advanced bioinformatic analysis on carotid arteries from the normoxic near-term ovine fetus at sea level and high altitude. These studies identified genes and networks involved in the investigated processes. Therefore, microarray technology may be used to examine various aspects of the sheep and goat species.

The present study used microarray methodologies to investigate hair and cashmere growth in the Laiwu black goat. Gene transcripts expressed in skin cells from hair-rich and hairless areas were compared, and hundreds of DE genes were identified. However, no keratin or keratin-associated protein genes were DE, which contradicts the results of a previous study (34). Whether this problem is the result of species differences in probe sets remains to be elucidated.

Gene annotation is incomplete and imperfect (certain DE genes remain to be identified), which may explain why certain factors, including keratin and keratin-associated proteins, were not identified in the present study. However, the repertoire of gene probes used (8×15K) confirmed the involvement of a number of growth factors and HOX genes. The application of an original algorithm to construct gene networks of temporal regulation of hair growth was not performed due to limited DE genes. The present study identified molecules in the cashmere-bearing skin area of the Laiwu black goat, which may contribute to hair and cashmere traits.

Acknowledgments

The present study was funded by the National Hair Sheep Industry Technology System (grant no. CARS-40), Projects of Qingdao People's Livelihood Science and Technology (grant nos. 13-1-3-88-nsh and 14-2-3-45-nsh) and the National Natural Science Foundation of China (grant no. 31301936).

References