Deciphering the Transcriptomic Complexity of Yak Skin across Different Ages and Body Sites

Skin plays an important role in the adaptability of yaks (Bos grunniens) to harsh environmental conditions. There were obvious differences in skin and hair phenotypes between scapular and ventral regions of yaks, and the differences became more prominent with age. However, the genetic mechanism that causes the difference of yak skin at different age has not been reported. In this study, we investigated the transcriptomic profile of yak skin across different ages (0.5 years, 2.5 years, and 4.5 years) and body sites (scapular and ventral regions). Differential gene expression analysis was initially conducted to explore the transcriptomic differences of skin in different ages and different body sites. Subsequently, Weighted Gene Co-expression Network Analysis (WGCNA) was employed to analyze the transcriptomic data from different ages and body sites comprehensively. A total of 1509 differentially expressed genes (DEGs) were identified in all the comparison groups. Among these, the largest number of DEGs was observed in Y2.5_S vs. Y2.5_V group, with a total count of 491 (179upregulated, 312downregulated). Then, the number of DEGs was higher in Y2.5_V vs. Y0.5_V group with 370 DEGs (103upregulated, 267downregulated). DEGs such as IGF2BP3, ADCY8, FOSL1 and S100A8 were found in all comparison groups of different ages, and multiple members of the HOX gene family including HOXC10, HOXA9, HOXA6, HOXB9 and HOXB6 were differentially expressed in comparison groups at different sites. Functional enrichment analysis showed that there were more obvious differences in immune function between different ages of skin, and more obvious differences in endocrine function between different parts of skin. WGCNA analysis revealed that genes related with immunity such as GLYATL2, ACSL5, SPDEF and FOXA1 were the core genes of the co-expression module associated with scapula region, and multiple genes related to hair follicle development and keratin such as FOXN1, OVOL1, KRT27, KRT71, DLX3, HOXC13 and TCHH were found to be the hub genes of the co-expression module associated with the ventral region. Overall, our study provides valuable insights into the transcriptomic complexity of yak skin across different ages and body sites. The differential gene expression patterns and co-expression network modules identified in this study lay the foundation for further research on skin biology and adaptation mechanisms in yaks.