Skip to main content

The fiber diameter traits of Tibetan cashmere goats are governed by the inherent differences in stress, hypoxic, and metabolic adaptations: an integrative study of proteome and transcriptome

Abstract

Background

Tibetan cashmere goats are served as a valuable model for high altitude adaptation and hypoxia complications related studies, while the cashmere produced by these goats is an important source of income for the herders. The aim of this study was to investigate the differences in protein abundance underlying the fine (average 12.20 ± 0.03 μm of mean fiber diameter) and coarse cashmere (average 14.67 ± 0.05 μm of mean fiber diameter) producing by Tibetan cashmere goats. We systematically investigated the genetic determinants of fiber diameter by integrated analysis with proteomic and transcriptomic datasets from skin tissues of Tibetan cashmere goats.

Results

We identified 1980 proteins using a label-free proteomics approach. They were annotated to three different databases, while 1730 proteins were mapped to the original protein coding genes (PCGs) of the transcriptomic study. Comparative analyses of cashmere with extremely fine vs. coarse phenotypes yielded 29 differentially expressed proteins (DEPs), for instance, APOH, GANAB, AEBP1, CP, CPB2, GPR142, VTN, IMPA1, CTSZ, GLB1, and HMCN1. Functional enrichment analysis of these DEPs revealed their involvement in oxidation-reduction process, cell redox homeostasis, metabolic, PI3K-Akt, MAPK, and Wnt signaling pathways. Transcription factors enrichment analysis revealed the proteins mainly belong to NF-YB family, HMG family, CSD family. We further validated the protein abundance of four DEPs (GC, VTN, AEBP1, and GPR142) through western blot, and considered they were the most potential candidate genes for cashmere traits in Tibetan cashmere goats.

Conclusions

These analyses indicated that the major biological variations underlying the difference of cashmere fiber diameter in Tibetan cashmere goats were attributed to the inherent adaptations related to metabolic, hypoxic, and stress response differences. This study provided novel insights into the breeding strategies for cashmere traits and enhance the understanding of the biological and genetic mechanisms of cashmere traits in Tibetan cashmere goats.

Peer Review reports

Background

The Tibetan cashmere goat (Capra hircus) is a critical source of cashmere and meat production and significantly contributes to the rural economy. The white cashmere goats of north-western Tibet trace their origins to cashmere goats in Kashmir valley. Cashmere derived from the secondary hair follicles, growing by hair development cycle, is a valuable commodity [1,2,3]. These goats are well-adapted to harsh climatic and hypoxic conditions of high altitudes owing to their dense coat. Due to these inherent characteristics, cashmere goats and other native species like Tibetan yaks and dogs are considered ideal models to study the genetic mechanisms underlying hypoxia related complications and high altitude adaptation [4,5,6,7]. The fiber diameter and production of cashmere determine the economic value. With an increase in demand, breeding cashmere goats for fineness and high quality has become a priority.

Several studies have been carried out on transcriptome and methylome to delineate the genetic basis of hair follicle development in cashmere goats and identified potential candidate genes [1, 8]. For instance, Wang et al. [8] revealed the regulatory mechanisms of hair follicle morphogenesis in Cashmere goats. Zhang et al. [9] investigated the hair follicle cycling in milk and cashmere goats. Meanwhile, our previous study reported potential mRNAs and lncRNAs associated with cashmere fineness [10]. These studies significantly contributed to our understanding of the genetic mechanisms underlying the complex biological landscape of cashmere.

From the biological perspective, transcriptome represents the intermediate state of gene expression and can reflect the mechanism such as transcriptional and post-transcriptional regulation, whereas methylome also plays a role [11]. However, the interpretation of this complex relationship of translation and post-translational modifications and its final outcome in the form of proteins, still remains unexplored. Due to protein being the direct functional executor of the organism, it has irreplaceable advantages through transcriptome and proteome data to obtain the characteristics of differences in gene expression level and protein level of samples.

Therefore, this study aimed to carry out integrative analysis of proteome with transcriptome (Fig. 1A) to provide novel insights about the final fate of differences observed earlier at the transcriptome level of the fine and coarse cashmere goats. It will be of great value to dissect the critical genes, signaling pathways, and the regulatory mechanism underlying the differences of cashmere fiber diameter in the cashmere goats.

Fig. 1
figure 1

Identification of proteins in the skins of Tibet. A The global study design. B Mean fibre diameter (MFD) of fineness type (F) cashmere and coarse type (C) cashmere. C The principal component analysis (PCA) plot shows separation patterns of skin samples from different mean fibre diameter (MFD) Tibet cashmere goats

Results

Protein identification and functional annotation at different mean fiber diameter (MFD) of cashmere using label-free proteomics

To comprehensively elucidate the proteome profiles that could be involved in the cashmere fiber differences, we identified proteins of skin tissues from Tibetan cashmere goats using label-free quantification proteomics. A total of 12,712 unique peptides were detected (Additional file 1: Tables S1) and 1965 proteins were identified in these skin tissues (Additional file 2: Tables S2). The molecular weight of most the proteins (> 50%) ranged from 0 to 60 kDa (Additional file 3: Fig. S1). In addition, most proteins (> 50%) had high peptide coverage (Additional file 4: Fig. S2). Among the identified proteins, 70% were represented by fewer than five peptides (Additional file 5: Fig. S3), indicating good sequence coverage. Principal component analysis (PCA) demonstrated the difference between the skin tissues of fine-type and coarse-type cashmere goats, indicating our newly generated protein datasets are reliable for further analyses (Fig. 1B, C).

In order to dissect the function of identified proteins, we further annotated all proteins functions and their classification using three databases. GO annotations revealed that most of the proteins were related to the metabolic and cellular processes (Additional file 6: Fig. S4). Clusters of orthologous groups for eukaryotic complete genomes (KOG) clustering of proteins functional categories included signal transduction mechanisms, posttranslational modification, protein turnover, chaperones, cytoskeleton, cell cycle control, cell division, and chromosome partitioning, indicating these functional categories might be closely related to the regulation of cashmere fiber differences (Fig. 2). Kyoto encyclopedia of genes and genomes (KEGG) pathway analysis (www.kegg.jp/feedback/copyright.html) revealed that most of the proteins participated in one of the carbohydrate metabolism (n = 219), signal transduction (n = 427), immune system (n = 355) (Additional file 7: Fig. S5). A total of 1962 proteins were annotated in these three databases (Additional file 8: Fig. S6).

Fig. 2
figure 2

The KOG functional categories of all identified proteins

To search the specific structure and independent functional region of the proteins, we performed protein domain prediction (Additional file 9: Table S3). Overall, the proteins were engaged in the EGF-like domain, laminin EGF domain, immunoglobulin domain, which showed that these candidate genes might play important roles in cashmere traits. To investigate whether the protein is a transcription factor (TFs) and which transcription factor family belongs to, we compared the predicted protein sequence with the corresponding animal TFdb (Additional file 10: Fig. S7). It revealed that the proteins mainly belong to the NF-YB family (n = 21, ranked top1), HMG family (n = 19, ranked top2), CSD family (n = 8, ranked top3). We also investigated the specific location of the proteins in the cell by subcellular localization (Additional file 11: Table S4).

Identification and functional enrichment of differentially expressed protein (DEPs)

We identified 29 DEPs between the fine-type cashmere and coarse-type cashmere skin tissues in Tibetan cashmere goats, including 5 up-regulated and 24 down-regulated DEPs (Fig. 3A). The heatmap was visualized to investigate the expression profiles patterns of all DEPs, and the results indicated that each group was clustered separately (Fig. 3B). In order to better understand the biological function of these DEPs, GO enrichment and KEGG pathway analysis were performed. We found potentially relevant GO terms such as vitamin transmembrane transport (such as GC), regulation of protein binding (VTN), metabolic process (GLB1), regulation of collagen fibril organization (AEBP1), G-protein coupled receptor signaling pathway (GNG12 and GPR142), these GO terms might directly or indirectly contribute to composition and function of cashmere fiber (Fig. 3C). KEGG pathway analysis (www.kegg.jp/feedback/copyright.html) revealed DEPs were significantly (P < 0.05) enriched for the PI3K-Akt signaling pathway, MAPK signaling pathway, and inositol phosphate metabolism (Fig. 3D).

Fig. 3
figure 3

Identification of differential expressed proteins (DEPs) between the fine-type (F) cashmere and coarse-type (C) cashmere skin tissues in Tibetan Cashmere goats. A Volcano plot for DEPs. The y-axis corresponds to the mean expression value of log10 (P-value), and the x-axis displays the log2 fold change value. The red dots represent the significantly up-regulated DEPs (fold change > 1.2 and FDR < 0.05), and the yellow dots represent the significantly down-regulated DEPs (fold change < 1.2 and FDR < 0.05); the red dots represent the non-DEPs. B Heatmap showed the expression patterns for DEPs. C and D The top 10 enriched gene ontology (GO) terms (biological processes, BP) and kyoto encyclopedia of genes and genomes (KEGG) pathways on the DEPs

Gene set enrichment analysis (GSEA) of all proteins revealed significantly (P < 0.05) enriched GO terms, such as regulation of hormone levels (Fig. 4A), and KEGG pathways, such as oxidative phosphorylation (Fig. 4B). Protein-protein interaction (PPI) analysis were performed for the 29 DEPs (7 edges identified; PPI enrichment P-value = 0.0002) (Fig. 4C). The PPI network revealed that GC, VTN, APOH, CPB2, and C11orf58 were the core nodes.

Fig. 4
figure 4

Gene set enrichment analysis (GSEA) of gene ontology (GO) terms (A) and kyoto encyclopedia of genes and genomes (KEGG) pathways (B) for proteins in skin tissues of Tibetan Cashmere goats. C The protein-protein interaction (PPI) network with differentially expressed proteins (DEPs)

Integrated analysis of proteome and transcriptome

As the transcriptomic data was obtained from the same samples of our previous study [10], we characterized the differences in gene expression and protein levels of these samples. The number of associations observed between transcriptome and proteome is shown in Fig. 5A. The distribution of the corresponding mRNA and proteins were shown through a scatter plot (Fig. 5B). The Pearson correlation coefficient was extremely low (0.004), indicating that the expression levels of mRNA and protein were not always positively correlated in skin tissue of Tibetan cashmere goats. We found a total of 69 concordant dots, representing a corresponding expressed trend between protein abundance and transcript abundance (red dots). In addition, 6 green dots (transcripts only) and 46 blue dots (proteins only) were identified, indicating differential expression was either found on the transcript or the protein levels in our study. We further investigated the functions of proteins in quadrants 4, 5, 6 regions; where most of them are engaged in the oxidation-reduction process, cell redox homeostasis, metabolic and cellular-related processes (Fig. 5C).

Fig. 5
figure 5

Association analysis of transcriptome and proteome data. A Venn plot of the number of all detected genes, all detected proteins, differentially expressed proteins (DEGs), differentially expressed proteins (DEPs). B Comparison of the expression between transcriptomic (y-axis) and proteomic (x-axis) profiling. Log2 expression ratios were calculated from the expression of fineness type (F) cashmere and coarse type (C) cashmere. Significant changes in expression are color-coded: green, transcripts only; blue, proteins only; red, both. C The top 5 enriched gene ontology (GO) terms (biological processes, BP) for genes in quadrants 4, 5, 6 regions

Validation of the expression levels of DEPs

To verify the accuracy of the DEPs detected by proteome, we selected VTN, GLB1, AEBP1, and GPR142 for western blot analyses (Fig. 6A, Additional file 12: Fig. S8). The results showed that the protein abundance changes were consistent with those obtained by proteomics, and the protein expression levels of GPR142, VTN, GLB1, AEBP1 in the C group were significantly higher than F group (P < 0.05) (Fig. 6B-E).

Fig. 6
figure 6

Western blotting and quantitative analysis in the skin tissues of Tibetan cashmere goats. A Western blot bands of GPR142, VTN, GLB1, AEBP1, GAPDH was used as control; the experiments were performed using skin tissues for each group and repeated three times. The F1, F2 and F3 correspond to fine type cashmere (F) samples, C1, C2 and C3 correspond to coarse type cashmere (C) samples, respectively. The same as below. BE Relative protein expression levels of GPR142, VTN, GLB1 and AEBP1. All data are presented as means ± SE. * represents P < 0.05, ** represents P < 0.01 and *** represents P < 0.001, which regarded as statistically significant, highly significant and extremely significant, respectively

Discussion

Goat cashmere fiber has unique attributes for a valuable model that allow detailed study for hair follicles developmental biology, while the MFD greatly influences the quality of cashmere. In our proteomics study, these 29 DEPs were primarily involved in biological processes related to metabolism, cell process, and environmental information processing. Interestingly, six DEGs members of the solute-carrier family were observed in our previous RNA-seq study [10], and only SLC2A1 encoding glucose transporter protein type 1 (GLUT1) was upregulated and. Given the importance of GLUT1 in transporting glucose to the cell and vitality for its growth and proliferation [12]; it could be assumed that metabolic level changes might be central to the cashmere phenotype difference in the F and C group. This assumption of metabolic differences between both groups may be further augmented by the fact that hypoxia induces glucose uptake and metabolism through up-regulation of GLUT1 protein [13]. Meanwhile, Tibetan cashmere goats are widely studied for hypoxia and high altitude biology studies [4, 5]. As it is well-known that the upregulation trend of genes is a response to certain external stimulate. Our previous DEGs analysis has not depicted such big difference between both groups [10]. Therefore, it can be assumed from our proteomics approach that the inherent resilience of the F group is pronounced compared to the C group. These changes may be attributed to the role of post-transcriptional, translational, and post-translational modifications and regulatory processes [14], in addition to the regulation of protein degradation by hidden proteases regulatory activity [15]. These phenomena signify the magnitude of variation underlying the phenotype difference between the two groups.

The KEGG enriched proteins majority falls in metabolic, focal adhesion, Wnt signaling pathway, MAPK signaling pathway, inositol phosphate metabolism pathway, and PI3K-Akt signaling pathways. As cashmere goats are used as a model of hypoxia related studies, it is interesting to note a study showing metabolic changes attributed to intrauterine hypoxia impacting the hair metabolome [16]. The focal adhesion pathway is important in cell adherence and its interaction with the extracellular matrix, and thereby has a demonstrated role in hair follicle growth [17, 18]. Wnt pathway involved in hair follicle induction and regulates its downstream morphogenesis and differentiation [19, 20], is upregulated by the hypoxic conditions and hypoxia inducible factor 1-alpha (HIF-1α) gene [21,22,23,24]. These mechanisms have been shown to regulate metabolic processes contributing to cell survival advantages. PI3K-Akt pathway has also been shown to be crucially involved in hair follicle regeneration, and is regarded as a potential regulating pathway for hair regeneration therapy [25,26,27]. MAPK pathway essentially regulates the hair cycle and the self-renewal of hair follicle stem cells [28]. Besides, potentially GO terms were identified in this study such as vitamin D metabolic process (GC), regulation of protein binding (VTN), vitamin transmembrane transport (GC), regulation of collagen fibril organization (AEBP1), G-protein coupled receptor signaling pathway (GNG12 and GPR142); these GO terms directly or indirectly regulate and contribute to composition and function of the hair fiber. For instance, vitamin D regulates epidermal differentiation, hair follicle cycling, and hair growth [29]; and G-protein coupled receptor signaling pathway has an important role in hair follicle stem cell activation [30].

Proteins transcription factors prediction also revealed important insights about the unique characteristics of Tibetan cashmere goats. HMG box family of reductase inhibitors actively downregulate the activation of transcription factors NF-kappaB (NF-κB), activating protein-1 (AP-1), and HIF-1α [31, 32]. This phenomenon is noteworthy since NF-κB and HIF-1α are typically involved in metabolic process activation and survival oriented changes during hypoxia [33]. While studies have also shown the downstream regulation of the Wnt signaling pathway through the HMG box family of proteins [34]. NF-YB transcription factor family is widely studied in plants, where it has been related to root growth, drought, and heat resistance [35, 36]; however, related studies about hair growth or hypoxia lacks in mammals. The induction effects of the CSD family on the expression of HIF-1α and cold shock proteins are well-established [37,38,39], which highlights the importance of hypoxia and harsh climate resistance in Tibetan cashmere goats.

Among the important significant DEPs, VTN is a glycoprotein, and has multifunctional roles in cell adhesion and migration via peri-cellular proteolysis and growth factors signaling. Studies have reported that VTN is a common constituent of the extracellular deposits, therefore it is assumed as a candidate gene and selected for validation [40, 41]. GLB1 is shown to be directly induced by hypoxic conditions and has a particular role in the context of Tibetan cashmere goats, which are also exposed to hypoxic conditions due to their high altitude habitat [5, 42]. Another example is AEBP1, a potential regulator of cholesterol efflux and MAPK signaling, having higher expression in the telogen of hair growth cycle [43,44,45]. Similarly, GPR142 is a member of the G-protein coupled receptor signaling pathway, has an important role in hair follicle stem cell activation [30]. The results of protein abundance changes in VTN, GLB1, AEBP1 and GPR142 indicated that the expression profiles determined by our current labeling free proteomics method was reliable, and these candidate genes might be contributing to the difference of MFD in Tibetan cashmere goats.

Integrated transcriptome and proteome analyses are advised as an efficient strategy to explore the mechanisms underlying a biological system [46]. However, the pearson correlation coefficient between DEGs and DEPs from our study was low (0.004). There are instances of previous studies reporting a lower correlation between mRNA and protein data [47, 48]. The reason might be that the expression levels of protein coding genes are not an end in themselves but rather it fulfills the protein synthesis need of the biological systems in its near future. Mammalian cells are shown to produce fewer copies of relatively unstable mRNA at a much lower rate than dozens of stable protein copies per mRNA [49, 50]. Therefore, relative abundances of proteins might or might not occur in proportion to their relative mRNA levels [51, 52].

Conclusions

In this study, we characterized the protein abundance changes of different cashmere fibers in Tibetan cashmere goats. We highlighted the functional and regulatory networks by integrating proteome with transcriptome. Finally, we provided novel insights into the biological mechanisms underlying the adaptation of Tibetan cashmere goats to hypoxic high altitude elevations and the dissection of the molecular mechanisms underlying the different MFD of cashmere. Differences between inherent metabolic adaptation to harsh and hypoxic conditions may be correlated to the observed phenotype differences. Further studies are warranted in this direction to clarify these complex relationships.

Methods

Animals and sample collection

Tibetan cashmere goats were obtained from the goat breeding farm located in Tibet (longitude: 91°12′-93°02′ E, latitude: 30°31′-31°55′ N, altitude: about 5000 m), China. The local environmental challenge to animals is the low oxygen availability at high altitudes. The partial pressure of oxygen at 5000 m is only about 50% of the value at sea-level, and the resultant hypoxia imposes severe constraints on aerobic metabolism [5]. We collected cashmere samples of 50 female goats (1.5 years old) in anagen stage from the left mid-side scapular of each animal body, these female goats were artificially inseminated with fresh sperm from a single ram. The cashmere samples were measured MFD values following the standardized methods set by China Fiber Inspection Bureau. In total, 3 female goats with the lowest MFD values (average 12.20 ± 0.03 μm, 1.5 years old) as the fine type (F) group, and another 3 female goats with the highest MFD values (average 14.67 ± 0.05 μm, 1.5 years old) as the coarse type (C) group, were selected and sampled (Fig. 1B). Scapular skin tissues were collected from these 6 Tibetan cashmere goats in vivo and immediately frozen by liquid nitrogen for further analysis.

Skin tissue protein extraction and digestion

Total proteins were extracted using the cold acetone method [53]. Protein quality was examined with SDS-PAGE. The results of SDS-PAGE showed that the separated bands were clear, abundant and non-degraded, and the bands were consistent among all samples (Additional file 13: Fig. S9). BCA protein assay kit (Pierce, Rockford, IL) was used to determine the protein concentration of the supernatant The protein concentration of all samples is greater than 1 μg/μL and total protein content is greater than 0.2 mg, meet the label-free proteomics-sequencing requirements (Additional file 14: Table S5). Proteins (50 μg) extracted from cells were suspended in 50 μL solution, reduced by adding 1ul of 1 M dithiotreitol at 55 °C for 1 h, alkylated by adding 5 μl of 20 mM iodoacetamide in the dark at 37 °C for 1 h. The samples were precipitated and washed in acetone and re-suspended in 50 mM ammonium bicarbonate. Finally, the proteins were digested with trypsin (Promega, Madison, WI) at a substrate/enzyme ratio of 50:1 (w/w) at 37 °C for 16 h to obtain peptide mixtures.

Nano-HPLC-MS/MS analysis

The peptides were re-dissolved in 30 μL solvent A (A: 0.1% formic acid in water) and analyzed by online nano-spray LC-MS/MS on an Orbitrap FusionTM LumosTM coupled to EASY-nLC 1200 system (Thermo Fisher Scientific, MA, USA). The peptide sample was loaded onto the analytical column (Acclaim PepMap C18, 75 μm × 25 cm) and separated with a 120-min gradient, from 5 to 35% B (B: 0.1% formic acid in ACN). The column flow rate was maintained at 200 nL/min with a column temperature of 40 °C. The electrospray voltage of 2 kV versus the inlet of the mass spectrometer was used. The mass spectrometer was run under data dependent acquisition mode and automatically switch between MS and MS/MS modes. The parameters was: (1) MS: scan range (m/z) = 350–1200; resolution = 120,000; AGC target = 400,000; maximum injection time = 50 ms; Filter Dynamic Exclusion: exclusion duration = 30s; (2) HCD-MS/MS: resolution = 15,000; AGC target = 50,000; maximum injection time = 35 ms; collision energy = 32.

Protein identification and quantification

PEAKS Studio Version X (Bioinformatics Solutions Inc., Waterloo, Canada) mass spectrometry was used for analysis. PEAKS DB was used to search the protein database. PEAKS DB were searched with a fragment ion mass tolerance of 0.05 Da and a parent ion tolerance of 10 ppm. Qualitative analysis of proteins is to determine whether a protein is present in a sample and to identify the type of protein. To ensure the reliability of the results, the peptide false discovery rate (FDR) is required to be less than 1% to evaluate the error discovery rate. In addition, a unique peptide refers to a peptide that had been identified and only comes from a single protein sequence or a sequence from the same group, requiring a protein’s unique peptide ≥1. Peptides and proteins that met these requirements were used for subsequent analysis.

Functional annotation analysis

The protein functions and classification were analyzed based on searches against the following databases: GO, KOG, and KEGG database (www.kegg.jp/feedback/copyright.html) [54]. Significant GO terms and pathways were examined with a P value≤0.05. The predicted transcription factor (TF) sequences were compared by hmmscan with the animalTFDB database [55]. The prediction of the protein domain was used the Pfam_scan [56]. The protein sequence was compared with the Pfam database to obtain the relevant annotation information of protein structure. The software WoLFPSort was used to predict the subcellular location of the protein and study the function of the protein [57].

Differentially expressed proteins and genes

DEPs analysis was performed by R package edgeR between two groups and FDR method was applied for correction [58]. We identified DEPs with a threshold of fold change>1.2 and a FDR < 0.05. Expressed genes were derived from our previous study [10]. We identified mRNA with a fold change> 2 and a FDR < 0.05 in a comparison as significant DEGs.

Transcriptome and proteome association analysis

Correlation analysis of genes and proteins was performed by R function cor.test. A nine-quadrant map was drawn based on changes in the expression of the gene in the transcriptome and proteome. Quantitative analysis and enrichment analysis was performed in genes of each region of the nine-quadrant map.

GSEA and PPI analysis

GSEA was performed by R package GSVA [59] and MSigDB [60] to identify whether a set of genes in specific GO terms/pathways show significant differences in the two groups. Enrichment scores and P-value were calculated in default parameters mode. PPI network was identified using String [61], which determined genes as nodes and interaction as lines in a network.

Western blot analysis of DEPs

Protein samples were subjected to SDS-PAGE, transferred onto PVDF membranes, and blocked with 5% skimmed milk (Boster, Wuhan, China) in TBST for 1 h at room temperature. We used different antibodies for different genes, the spacing was smaller, and the exposure intensity was different, so the blots were cut prior to hybridization with antibodies. Primary antibodies and their dilution ratio were followed as, goat anti-rabbit GPR142 antibody (1:500), goat anti-rabbit VTN antibody (1:1000), goat anti-rabbit GLB1 antibody (1:1000), goat anti-rabbit AEBP1 antibody (1:1000), and goat anti-mouse GAPDH antibody (1:5000) (Proteintech, USA). The protein bands were washed and visualized by using an enhanced chemi-luminescence ECL kit (ThermoFisher Scientific). Band intensities were visualized using Image Lab software with a Bio-Rad system (Bio-Rad, Hercules, CA, USA). The target protein contents were normalized to GAPDH level in each lane. The experiment for each gene was performed in triplicates.

Statistical analysis

The statistical analyses between two groups were analyzed by Student’s t-test. All values are shown as mean ± SE. For all tests, * represents P < 0.05, ** represents P < 0.01, and *** represents P < 0.001; regarded as statistically significant, highly significant, and extremely significant, respectively.

Availability of data and materials

The mass spectrometry proteomics data have been deposited to the ProteomeXchange Consortium (http://proteomecentral.proteomexchange.org) via the iProX partner repository with the dataset identifier PXD030146. ‍The RNA sequencing datasets used in the current study were obtained from our previous study [10] and has been deposited in National Center for Biotechnology Information Sequence Read Archive (NCBI SRA) database under BioProject No. PRJNA643003.

Abbreviations

DEPs:

Differentially expressed proteins

DEGs:

Differentially expressed proteins

KOG:

Clusters of orthologous groups for eukaryotic complete genomes

MFD:

Mean fiber diameter

FDR:

False discovery rate

GO:

Gene ontology

KEGG:

Kyoto encyclopedia of genes and genomes

TF:

Transcription factor

GSEA:

Gene set enrichment analysis

PPI:

Protein-protein interaction

PCG:

Protein-coding gene

PCA:

Principal component analysis

F:

Fine type

C:

Coarse type

References

  1. Wang S, Li F, Liu J, Zhang Y, Zheng Y, Ge W, et al. Integrative analysis of methylome and transcriptome reveals the regulatory mechanisms of hair follicle morphogenesis in cashmere goat. Cells. 2020;9:969.

    Article  PubMed Central  Google Scholar 

  2. Yang M, Song S, Dong KZ, Chen XF, Liu XX, Rouzi M, et al. Skin transcriptome reveals the intrinsic molecular mechanisms underlying hair follicle cycling in cashmere goats under natural and shortened photoperiod conditions. Sci Rep. 2017;7:13502.

    Article  PubMed  PubMed Central  Google Scholar 

  3. Liu GB, Liu RZ, Li QQ, Tang XH, Yu M, Li XY, et al. Identification of microRNAs in wool follicles during anagen, catagen, and telogen phases in Tibetan sheep. PLoS One. 2013;8:e77801.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  4. Kumar C, Song S, Jiang L, He XH, Zhao QJ, Pu YB, et al. Sequence characterization of DSG3 gene to know its role in high-altitude hypoxia adaptation in the Chinese cashmere goat. Front Genet. 2018;9:553.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Song S, Yao N, Yang M, Liu X, Dong K, Zhao Q, et al. Exome sequencing reveals genetic differentiation due to high-altitude adaptation in the Tibetan cashmere goat (Capra hircus). BMC Genomics. 2016;17:122.

    Article  PubMed  PubMed Central  Google Scholar 

  6. Qiu Q, Zhang GJ, Ma T, Qian WB, Wang JY, Ye ZQ, et al. The yak genome and adaptation to life at high altitude. Nat Genet. 2012;44:946–9.

    Article  CAS  PubMed  Google Scholar 

  7. Li Y, Wu DD, Boyko AR, Wang GD, Wu SF, Irwin DM, et al. Population variation revealed high-altitude adaptation of Tibetan mastiffs. Mol Biol Evol. 2014;31:1200–5.

    Article  CAS  PubMed  Google Scholar 

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

  9. Zhang YJ, Wu KJ, Wang LL, Wang ZY, Han WJ, Chen D, et al. Comparative study on seasonal hair follicle cycling by analysis of the transcriptomes from cashmere and milk goats. Genomics. 2020;112:332–45.

    Article  CAS  PubMed  Google Scholar 

  10. Fu X, Zhao B, Tian K, Wu Y, Suo L, Ba G, et al. Integrated analysis of lncRNA and mRNA reveals novel insights into cashmere fineness in Tibetan cashmere goats. PeerJ. 2020;8:e10217.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Zhao BR, Luo HP, He JM, Huang XX, Chen SQ, Fu XF, et al. Comprehensive transcriptome and methylome analysis delineates the biological basis of hair follicle development and wool-related traits in Merino sheep. BMC Biol. 2021;19:197.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Xiao HJ, Wang J, Yan WX, Cui YB, Chen Z, Gao X, et al. GLUT1 regulates cell glycolysis and proliferation in prostate cancer. Prostate. 2018;78:86–94.

    Article  CAS  PubMed  Google Scholar 

  13. Park HS, Kim JH, Sun BK, Song SU, Suh W, Sung JH. Hypoxia induces glucose uptake and metabolism of adipose-derived stem cells. Mol Med Rep. 2016;14:4706–14.

    Article  CAS  PubMed  Google Scholar 

  14. Buszczak M, Signer RAJ, Morrison SJ. Cellular differences in protein synthesis regulate tissue homeostasis. Cell. 2014;159:242–51.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Vilchez D, Boyer L, Morantte I, Lutz M, Merkwirth C, Joyce D, et al. Increased proteasome activity in human embryonic stem cells is regulated by PSMD11. Nature. 2012;489:304–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Yang J, Wei Y, Qi HB, Yin NL, Yang Y, Li ZL, et al. Neonatal hair profiling reveals a metabolic phenotype of monochorionic twins with selective intrauterine growth restriction and abnormal umbilical artery flow. Mol Med. 2020;26:37.

    Article  PubMed  PubMed Central  Google Scholar 

  17. Burridge K, ChrzanowskaWodnicka M. Focal adhesions, contractility, and signaling. Annu Rev Cell Dev Biol. 1996;12:463–518.

    Article  CAS  PubMed  Google Scholar 

  18. Yoon SY, Dieterich LC, Tacconi C, Sesartic M, He YL, Brunner L, et al. An important role of podoplanin in hair follicle growth. PLoS One. 2019;14:e0219938.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Rishikaysh P, Dev K, Diaz D, Qureshi W, Filip S, Mokry J. Signaling involved in hair follicle morphogenesis and development. Int J Mol Sci. 2014;15:1647–70.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Andl T, Reddy ST, Gaddapara T, Millar SE. WNT signals are required for the initiation of hair follicle development. Dev Cell. 2002;2:643–53.

    Article  CAS  PubMed  Google Scholar 

  21. Zhan L, Liu D, Wen H, Hu J, Pang T, Sun W, et al. Hypoxic postconditioning activates the Wnt/β-catenin pathway and protects against transient global cerebral ischemia through Dkk1 inhibition and GSK-3β inactivation. FASEB J. 2019;33:9291–307.

    Article  CAS  PubMed  Google Scholar 

  22. Varela-Nallar L, Rojas-Abalos M, Abbott AC, Moya EA, Iturriaga R, Inestrosa NC. Chronic hypoxia induces the activation of the Wnt/β-catenin signaling pathway and stimulates hippocampal neurogenesis in wild-type and APPswe-PS1ΔE9 transgenic mice in vivo. Front Cell Neurosci. 2014;10(8):17.

    Google Scholar 

  23. Zhang Q, Bai X, Chen W, Ma T, Hu Q, Liang C, et al. Wnt/β-catenin signaling enhances hypoxia-induced epithelial–mesenchymal transition in hepatocellular carcinoma via crosstalk with hif-1α signaling. Carcinogenesis. 2013;34:962–73.

    Article  PubMed  Google Scholar 

  24. Liu H, Liu D, Ding G, Liao P, Zhang J. Hypoxia-inducible factor-1α and Wnt/β-catenin signaling pathways promote the invasion of hypoxic gastric cancer cells. Mol Med Rep. 2015;12:3365–73.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Lu QM, Gao Y, Fan ZM, Xiao X, Chen Y, Si Y, et al. Amphiregulin promotes hair regeneration of skin-derived precursors via the PI3K and MAPK pathways. Cell Prolif. 2021;54:e13106.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Chen Y, Fan ZM, Wang XX, Mo MH, Zeng SB, Xu RH, et al. PI3K/Akt signaling pathway is essential for de novo hair follicle regeneration. Stem Cell Res Ther. 2020;11:144.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Huang HC, Lin H, Huang MC. Lactoferrin promotes hair growth in mice and increases dermal papilla cell proliferation through Erk/Akt and Wnt signaling pathways. Arch Dermatol Res. 2019;311:411–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Ozturk OA, Pakula H, Chmielowiec J, Qi JJ, Stein S, Lan LX, et al. Gab1 and Mapk signaling are essential in the hair cycle and hair follicle stem cell quiescence. Cell Rep. 2015;13:561–72.

    Article  Google Scholar 

  29. Bikle D, Christakos S. New aspects of vitamin D metabolism and action - addressing the skin as source and target. Nat Rev Endocrinol. 2020;16:234–52.

    Article  CAS  PubMed  Google Scholar 

  30. Miranda M, Avila I, Esparza J, Shwartz Y, Hsu Y-C, Berdeaux R, et al. Defining a role for G-protein coupled receptor/cAMP/CRE-binding protein signaling in hair follicle stem cell activation. J Invest Dermatol. 1800;142:53–64.e3.

    Article  Google Scholar 

  31. Dichtl W, Dulak J, Frick M, Alber HF, Schwarzacher SP, Ares MPS, et al. HMG-CoA reductase inhibitors regulate inflammatory transcription factors in human endothelial and vascular smooth muscle cells. Arterioscler Thromb Vasc Biol. 2003;23:58–63.

    Article  CAS  PubMed  Google Scholar 

  32. Pallottini V, Guantario B, Martini C, Totta P, Filippi I, Carraro F, et al. Regulation of HMG-CoA reductase expression by hypoxia. J Cell Biochem. 2008;104:701–9.

    Article  CAS  PubMed  Google Scholar 

  33. Bruning U, Fitzpatrick SF, Frank T, Birtwistle M, Taylor CT, Cheong AJC, et al. NFκB and HIF display synergistic behaviour during hypoxic inflammation. Cell Mol Life Sci. 2012;69:1319–29.

    Article  CAS  PubMed  Google Scholar 

  34. Bernard P, Harley VR. Acquisition of SOX transcription factor specificity through protein-protein interaction, modulation of Wnt signalling and post-translational modification. Int J Biochem Cell Biol. 2010;42:400–10.

    Article  CAS  PubMed  Google Scholar 

  35. Sato H, Suzuki T, Takahashi F, Shinozaki K, Yamaguchi-Shinozaki K. NF-YB2 and NF-YB3 have functionally diverged and differentially induce drought and heat stress-specific genes. Plant Physiol. 2019;180:1677–90.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Das S, Parida SK, Agarwal P, Tyagi AK. Transcription factor OsNF-YB9 regulates reproductive growth and development in rice. Planta. 2019;250:1849–65.

    Article  CAS  PubMed  Google Scholar 

  37. Yang CY, Wang LL, Siva VS, Shi XW, Jiang QF, Wang JJ, et al. A novel cold-regulated cold shock domain containing protein from scallop chlamys farreri with nucleic acid-binding activity. PLoS One. 2012;7:e32012.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Mihailovich M, Militti C, Gabaldon T, Gebauer F. Eukaryotic cold shock domain proteins: highly versatile regulators of gene expression. Bioessays. 2010;32:109–18.

    Article  CAS  PubMed  Google Scholar 

  39. Coles LS, Lambrusco L, Burrows J, Hunter J, Diamond P, Bert AG, et al. Phosphorylation of cold shock domain/Y-box proteins by ERK2 and GSK3beta and repression of the human VEGF promoter. FEBS Lett. 2005;579:5372–8.

    Article  CAS  PubMed  Google Scholar 

  40. Ferraris GMS, Schulte C, Buttiglione V, De Lorenzi V, Piontini A, Galluzzi M, et al. The interaction between uPAR and vitronectin triggers ligand-independent adhesion signalling by integrins. EMBO J. 2014;33:2458–72.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Ryoko N, Ohoshi M, Hironobu F, Gu J, Toru K, Saburo A, et al. Characterization of the ligand-binding specificities of integrin alpha3beta1 and alpha6beta1 using a panel of purified laminin isoforms containing distinct alpha chains. J Biochem. 2003;134:497–504.

    Article  Google Scholar 

  42. Hunt PW, Klok EJ, Trevaskis B, Watts RA, Ellis MH, Peacock WJ, et al. Increased level of hemoglobin 1 enhances survival of hypoxic stress and promotes early growth in Arabidopsis thaliana. Proc Natl Acad Sci U S A. 2002;99:17197–202.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Geyfman M, Gordon W, Paus R, Andersen B. Identification of telogen markers underscores that telogen is far from a quiescent hair cycle phase. J Invest Dermatol. 2012;132:721–4.

    Article  CAS  PubMed  Google Scholar 

  44. Morgan HJ, Benketah A, Olivero C, Rees E, Ziaj S, Mukhtar A, et al. Human basal cell carcinoma: the induction of anagen hair follicle differentiation. Clin Exp Dermatol. 2020;45:309–17.

    Article  CAS  PubMed  Google Scholar 

  45. Majdalawieh A, Ro H-S. PPARgamma1 and LXRalpha face a new regulator of macrophage cholesterol homeostasis and inflammatory responsiveness, AEBP1. Nucl Recept Signal. 2010;8:e004.

    Article  PubMed  PubMed Central  Google Scholar 

  46. Buccitelli C, Selbach M. mRNAs, proteins and the emerging principles of gene expression control. Nat Rev Genet. 2020;21:630–44.

    Article  CAS  PubMed  Google Scholar 

  47. Pascal LE, True LD, Campbell DS, Deutsch EW, Risk M, Coleman IM, et al. Correlation of mRNA and protein levels: cell type-specific gene expression of cluster designation antigens in the prostate. BMC Genomics. 2008;9:246.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Ghazalpour A, Bennett B, Petyuk VA, Orozco L, Hagopian R, Mungrue IN, et al. Comparative analysis of proteome and transcriptome variation in mouse. PLoS Genet. 2011;7:e1001393.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Sharova LV, Sharov AA, Nedorezov T, Piao Y, Shaik N, Ko MSH. Database for mRNA half-life of 19 977 genes obtained by DNA microarray analysis of pluripotent and differentiating mouse embryonic stem cells. DNA Res. 2009;16:45–58.

    Article  CAS  PubMed  Google Scholar 

  50. Eden E, Geva-Zatorsky N, Issaeva I, Cohen A, Dekel E, Danon T, et al. Proteome half-life dynamics in living human cells. Science. 2011;331:764–8.

    Article  CAS  PubMed  Google Scholar 

  51. Lundberg E, Fagerberg L, Klevebring D, Matic I, Geiger T, Cox J, et al. Defining the transcriptome and proteome in three functionally different human cell lines. Mol Syst Biol. 2010;6:450.

    Article  PubMed  PubMed Central  Google Scholar 

  52. de Godoy LMF, Olsen JV, Cox J, Nielsen ML, Hubner NC, Frohlich F, et al. Comprehensive mass-spectrometry-based proteome quantification of haploid versus diploid yeast. Nature. 2008;455:1251–U60.

    Article  PubMed  Google Scholar 

  53. Santa C, Anjo SI, Manadas B. Protein precipitation of diluted samples in SDS-containing buffer with acetone leads to higher protein recovery and reproducibility in comparison with TCA/acetone approach. Proteomics. 2016;16:1847–51.

    Article  CAS  PubMed  Google Scholar 

  54. Kanehisa M, Goto S. KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28:27–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Hu H, Miao YR, Jia LH, Yu QY, Zhang Q, Guo AY. AnimalTFDB 3.0: a comprehensive resource for annotation and prediction of animal transcription factors. Nucleic Acids Res. 2019;47:D33–8.

    Article  CAS  PubMed  Google Scholar 

  56. 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–D85.

    Article  CAS  PubMed  Google Scholar 

  57. Horton P, Park KJ, Obayashi T, Fujita N, Harada H, Adams-Collier CJ, et al. WoLF PSORT: protein localization predictor. Nucleic Acids Res. 2007;35:W585–W87.

    Article  PubMed  PubMed Central  Google Scholar 

  58. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.

    Article  CAS  PubMed  Google Scholar 

  59. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005;102:15545–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Liberzon A, Subramanian A, Pinchback R, Thorvaldsdottir H, Tamayo P, Mesirov JP. Molecular signatures database (MSigDB) 3.0. Bioinformatics. 2011;27:1739–40.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Szklarczyk D, Franceschini A, Wyder S, Forslund K, Heller D, Huerta-Cepas J, et al. STRING v10: protein-protein interaction networks, integrated over the tree of life. Nucleic Acids Res. 2015;43:D447–D52.

    Article  CAS  PubMed  Google Scholar 

Download references

Acknowledgements

We would like to thanks Tibet Cashmere Goat Breeding Farm (Tibet, China) for managing and collecting the experimental populations.

Funding

This work was supported by the Breeding and Healthy Farming of Sheep and Goats in Tibet (Grant No. XZ202101ZD0001N), the National Key Research and development Program (Grant No. 2021YFD1200902), and the China Agriculture Research System (Grant No. CARS-39). The funders played no roles in study design, in the collection, analysis, and interpretation of data, in the writing of the manuscript, and in the decision to submit the manuscript for publication.

Author information

Authors and Affiliations

Authors

Contributions

BZ, XF and YW conceived and designed the study. BZ conducted statistical and bioinformatics analyses. CW, ZM, LS, and YW contributed to sample collection and process. BZ wrote the manuscript and AS revised the manuscript. All authors have read and approved the final manuscript.

Corresponding authors

Correspondence to Yujiang Wu or Xuefeng Fu.

Ethics declarations

Ethics approval and consent to participate

All procedures pertaining to the handling of experimental animals were conducted in accordance with the ARRIVE guidelines. All sample collections were carried out with the approved guidelines of the Institutional Animal Science of Xinjiang Academy of Animal Science under Protocol No. 2019009 and all efforts were made to minimize discomfort and suffering. The skin wound was recovered in 2 weeks with carefully care.

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.

Supplementary Information

Additional file 1: Table S1.

Overview of peptides identification information.

Additional file 2: Table S2.

Overview of proteins identification information.

Additional file 3: Figure S1.

The number of proteins identified at various molecular weight ranges.

Additional file 4: Figure S2.

The distribution of protein’s sequence coverage.

Additional file 5: Figure S3.

The number of peptide fragments for protein identification.

Additional file 6: Figure S4.

The GO annotation of all expressed proteins.

Additional file 7: Figure S5.

The KEGG annotation of all expressed proteins.

Additional file 8: Figure S6.

Venn plot indicated the overlapped annotated proteins of GO, KEGG, and KOG database.

Additional file 9: Table S3.

Summary results of protein domain prediction.

Additional file 10: Figure S7.

The top 10 TF families detected in all expressed proteins.

Additional file 11: Table S4.

The specific location of the proteins in the cell by subcellular localization.

Additional file 12: Figure S8.

The target bands and western blots of GPR142, VTN, GLB1, AEBP1. GAPDH was used as control.

Additional file 13: Figure S9.

The results of SDS-PAGE showed that the separated bands were clear, abundant and non-degraded, and the bands were consistent among all samples. The F1, F2 and F3 correspond to fine type cashmere (F) samples, C1, C2 and C3 correspond to coarse type cashmere (C) samples, respectively.

Additional file 14: Table S5.

The quantitative results of protein samples.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. 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 in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Zhao, B., Wu, C., Sammad, A. et al. The fiber diameter traits of Tibetan cashmere goats are governed by the inherent differences in stress, hypoxic, and metabolic adaptations: an integrative study of proteome and transcriptome. BMC Genomics 23, 191 (2022). https://doi.org/10.1186/s12864-022-08422-x

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-022-08422-x

Keywords