长毛兔粗、细毛毛囊转录组比较分析

Comparative analysis of transcriptomes between coarse and fine hair follicles in Angora rabbits.

摘要

为筛选出影响兔毛纤维直径的相关基因及调控通路，本研究采用RNA测序（RNA-seq）鉴定长毛兔粗、细毛毛囊中差异表达的基因，并进行差异基因（DEGs）功能富集分析。通过Cytoscape构建蛋白互作网络，筛选出核心基因，最后采用实时荧光定量PCR（qRT-PCR）对差异表达基因进行验证。结果表明：1）RNA-seq共鉴定了 412个DEGs，相比粗毛毛囊，160个DEGs在细毛毛囊上调表达，252个DEGs下调表达。2）GO功能富集分析显示，差异基因显著富集到330个条目中，包括蜕皮周期、发周期、毛发周期的调节、超分子纤维等与粗、细毛生长相关的条目；KEGG富集分析表明差异表达基因显著富集在雌激素信号通路、蛋白质消化和吸收、内质网中的蛋白质加工途径、PPAR信号通路、FoxO信号通路和基底细胞癌等与长毛兔粗、细毛生长分化有关的信号通路中。3）蛋白互作网络分析表明SHH、LGR5、CTSB、ACTA2和CRNN等是调控毛纤维结构的重要功能基因。综上，本研究鉴定了调控兔毛直径的差异表达基因，挖掘出5个与长毛兔粗、细毛分化相关的候选基因，为细毛型长毛兔的选育提供理论依据和数据支撑。

Abstract

To screen for genes and regulatory pathways that affect the diameter of Angora rabbit hair fibers， RNA sequencing （RNA-Seq） was used to identify expressed genes in the coarse and fine hair follicles of Angora rabbits， and differentially expressed genes （DEGs） functional enrichment was performed. Additionally， a protein-protein interaction （PPI） network was constructed using Cytoscape to identify hub genes involved in regulating the diameter of rabbit hair fibers. Finally， quantitative real-time PCR （qRT-PCR） was employed for validation of the DEGs. The results showed that： 1） A total of 412 DEGs were identified by RNA-seq.Compared to coarse hair follicles， 160 DEGs were up-regulated in fine hair follicles and 252 DEGs were down-regulated. 2） GO functional analysis showed that differentially expressed genes were significantly enriched in 330 biological processes， including molting cycle hair cycle， hair cycle regulation， supramolecular fiber formation， and other articles related to the growth of coarse and fine hair. Additionally， KEGG enrichment analysis demonstrated that differentially expressed genes were significantly enrichment in the estrogen signaling pathway， protein digestion and absorption pathway， endoplasmic reticulum protein processing pathway， PPAR signaling pathway， FoxO signaling pathway， and basal cell carcinoma pathway， which were related to the growth and differentiation of coarse and fine hair in Angora rabbits.3） Protein interaction network analysis highlighted SHH， LGR5， CTSB， ACTA2， and CRNN as important functional genes regulating the structure of hair fibers. In conclusion， this study elucidated the molecular mechanism underlying the regulation of rabbit hair fineness， identified DEGs involved in controlling rabbit hair diameter， and discovered 5 potential candidate genes associated with coarse and fine hair differentiation in Angora rabbits. These findings provide a solid theoretical foundation and empirical evidence to support the breeding program of fine-haired Angora rabbits.

Graphical abstract

关键词

长毛兔 / 毛囊 / 纤维直径 / RNA-seq / 蛋白互作

Key words

Angora rabbit / fiber diameter / RNA-seq / protein-protein interaction

近年来，在乡村振兴的政策背景下，各地政府积极引导兔毛加工企业在新产品开发上的加大投入、加强市场推广、开拓产品市场，已取得初步成效。据调查，2022年我国兔毛制品的出口量为12.53万件，相较于20021年，同比增长20.63%；与此同时，兔毛原料出口总量为138.71 t，同比下降5%^［1］，具有高附加值的兔毛制品受到越来越多的关注。这一转变，增加了高品质兔毛的需求量，也对毛兔品种改良提出更高的要求。兔毛以长、白、松、净、细、软成为纺织工业中最受欢迎的天然纤维之一，兔毛纤维细度在13~14 μm左右为细毛，适合作为精纺原料，更受市场欢迎^［1］。细度是影响兔毛纺织性能、决定兔毛市场价格等经济价值的重要指标^［2］。然而目前细毛产量及品质难以满足市场需要，为进一步改善兔毛品质，深入了解研究调控兔毛细度相关分子机制必不可少。

转录组测序（RNA sequencing，RNA-seq）技术能够对基因结构与功能进行深层次地挖掘和分析，探索各种经济性状与基因表达之间的内在联系^［4］。研究表明，兔毛的生长受到营养、环境及遗传等诸多因素的影响，与角质形成细胞增殖、表皮发育及毛囊形态发生密切相关^［3］。目前，利用RNA-seq技术，部分与兔毛生长密切相关通路已被报道，如Wnt信号通路^［5］、TGF-β信号通路^［6］、MAPK信号通路等^［7］。同时，鉴定出了大量与兔毛生长过程相关的基因，包括TYR、TYRP1、ASIP、WNT2、CREB1、Oca2、Dct、Slc24a5、Slc7a11、Bcl2等与兔毛生长相关^［8-10］；SFRP2、FRZB、CACNG1、PTPLA、TTN、DUSP26等与家兔毛囊和皮肤发育相关^［11］；FGF5、WNT5A、BMP4和BMP7等参与调控兔毛长度^［12］。

兔皮肤有2种不同类型的毛囊，即产生粗毛的初级毛囊和产生细毛的次级毛囊。毛囊是毛发生长的基础^［14］，毛囊底部的真皮状细胞提供信号来调控毛发的平均纤维直径和特性^［15］。毛纤维由角质形成细胞形成，而角质形成细胞来源于毛囊上皮中角质形成细胞干细胞^［16］。因此，毛囊在毛发纤维生长和形态发生过程起决定性作用^［17］。兔毛是由毛囊产生的重要纺织原料，不同类型毛囊的发育过程决定了兔毛产量和品质。毛囊结构功能和形态发生是一个复杂的生物学过程，解析影响不同类型毛囊发育的分子调控机制对于长毛兔种质资源创新利用及兔毛品质的改善，提供有价值的参考。然而，目前大多数与毛发纤维直径的研究都集中在羊毛上，对兔毛的研究有限。因此，本研究利用RNA-seq技术对长毛兔粗、细毛毛囊的转录组进行分析比较，鉴定调控兔毛直径的差异表达基因，了解调控兔毛细度的分子机制，为细毛型长毛兔的选育提供理论依据。

1　材料与方法

1.1　试验动物

本试验在安徽省农业科学院长毛兔育种基地开展，所有试验兔饲养条件相同，单笼饲养，自然采光、通风，定时定量采食，自动饮水。

1.2　样品采集

试验选取5只皖系长毛兔作为试验动物，同时进行拔毛处理，经过73 d养毛期后，选取体侧部皮肤采样。采样前通过耳静脉注射0.7%戊巴比妥钠（6 mL/kg）对试验群体进行麻醉，于每只毛兔体侧分别采集粗毛及细毛样本并从毛根部拔出，剪取毛根保存于2.0 mL冻存管，立即置于液氮保存备用。

1.3　cDNA文库构建及测序

使用Trizol试剂（Invitrogen，Carlsbad，CA，USA）从毛根样品中提取总RNA。使用NanoPhotometer分光光度计（Implen，Westlake Village，CA，USA）和Agilent 2100生物分析仪（Agilent Technologies，Santa Clara，CA，USA）测定RNA的纯度、浓度和完整性。对RNA样品进行质量控制，获得构建cDNA文库所需的高质量RNA。使用NEBNext^®Ultra™RNA Library Prep Kit for Illumina^®（NEB，Ipswich，MA，USA）构建粗毛毛囊和细毛毛囊文库，并使用Qubit2.0荧光仪进行定量。使用Agilent 2100生物分析仪确定文库插入物的大小，并使用Illumina NovaSeq 6000平台（Illumina，USA）进行测序。

1.4　差异基因表达分析

使用fastp（版本0.19.7）从原始测序reads中去除低质量的reads，包含接头（adaptor）和包含N过多的reads。使用HISAT2^［18］将粗毛和细毛毛囊样本的clean reads与OryCun2.0 数据库^［19］进行比对。转录本使用stringtie（v1.3.3b）进行组装，基因表达使用RSEM^［20］进行量化。利用FPKM的值来估计组织样本中的基因表达水平。采用DESeq2进行差异表达分析^［21］。采用标准|log₂（FC）| > 1，FDR< 0.05筛选差异表达基因。并使用GO（Gene Ontology）（http：//geneontology.org）和KEGG （Kyoto Encyclopedia of Genes and Genomes）通路（http：//www.genome.jp/kegg/）对差异基因进行功能富集分析。

1.5　蛋白互作网络的构建

蛋白互作网络图（PPI， Protein-Protein Interaction Networks）采用string数据库（版本11.0），以中等置信度0.4为参数构建。通过Cytoscape软件中cytoHubba插件的Degree中心性算法筛选出蛋白互作网络中Degree数值前10的中心基因。

1.6　实时荧光定量qRT-PCR

引物设计与合成：试验选择CW（粗毛）、FW（细毛）各5个样本，挑选8个差异表达基因，利用Oligo 6.0设计引物，引物序列见表1，内参基因为GAPDH，采用SYBR Green I法实时荧光定量PCR（相对定量）进行mRNA表达量分析，退火温度为61 ℃，采用2 ^-△△Ct法进行相对表达量计算。

1.7　数据处理分析

使用 Excel 2010软件将试验中qRT-PCR的结果进行整理，试验数据使用SPSS 27.0软件进行统计分析，数据以“平均值±标准差”表示。采用t检验法和单因素方差分析（One-way ANOVA）进行差异显著性检验。其中，P<0.05表示差异显著，P<0.01表示差异极显著。

2　结果与分析

2.1　毛囊RNA-Seq结果分析

本研究对长毛兔粗、细毛毛囊组织进行转录组测序，表2为转录组数据质控信息。共获得了484 990 274 条Raw Data。对测序原始数据进行质量控制（去除含有adapter、poly-N及低质量reads）。获得 482 530 874条高质量的Clean Data，其中粗毛（CW）组和细毛（FW）组得到的Clean Data分别为241 922 552和240 608 322条。所有样本过滤后的Q20均在97%以上，Q30均在94%以上，GC含量在50%~54%。所有的样本reads比对到参考基因组上的总映射比率（Total mapped ratio）均在82%以上，唯一比对率（Unique mapped ratio）均在76%以上。

2.2　差异表达基因筛选

本研究共检测到412个差异表达基因，其中，160个基因在细毛毛囊组中上调表达，252个基因下调表达。对样本和差异表达基因进行聚类分析（图1），粗细毛样本间聚类情况良好，呈现不同的基因表达模式。

本研究也检测到14个角蛋白（KRT）基因及4个角蛋白相关蛋白（KAP）基因在两组间显著差异表达，如：KRT38、KRT77、KRT82、KRT36、KRT32、KRT72、KRT15、KRT79和KRT40在细毛毛囊组织中上调表达，KRT16、KRT25、KRT27、KRT21和KRT6A在细毛毛囊组织中下调表达， 4个角蛋白相关蛋白基因KRTAP26-1、KRTAP24-1、KRTAP16-1和KRTAP13-2均在细毛毛囊组织中高表达（表3）。

2.3　差异表达基因的富集分析

对筛选到的412个差异表达基因进行GO富集和KEGG通路分析。

2.3.1　差异表达基因GO功能分析

GO富集分析主要分为细胞组分（Cellular component，CC）、分子功能（Molecular function，MF）及生物过程（Biological process，BP）。如图2所示，GO富集分析发现，差异表达基因显著富集（P <0.05）在330个GO条目，其中生物过程共富集到264个GO 条目，细胞组分共富集到39个GO term，分子功能27个GO 条目。其中生物过程显著富集的前10个GO 条目为发育过程（GO：0032502）、解剖结构发育（GO：0048856）、蜕皮周期（GO：0042303）、细胞分化（GO：0030154）、中间丝为基础的工艺（GO：0045103）、中间丝细胞骨架组织（GO：0045104）、药物反应（GO：0042493）、细胞发育过程（GO：0048869）、解剖结构形态发生（GO：009653）、神经元投射发育（GO：0031175）；在细胞组分中以面上泡（GO：0031988）、囊泡（GO：0031982）、细胞外区域（GO：0005576）、细胞外区部分（GO：0044421）、中间丝（GO：0005882）5个类别占比最高；分子功能中最显著富集的条目是钠离子跨膜转运活性（GO：0015081）。

2.3.2　差异基因KEGG通路分析

KEGG通路分析发现，差异基因显著富集（P<0.05）在雌激素信号通路（Ko04915）、蛋白质消化和吸收（Ko04974）、内质网中的蛋白质加工途径（Ko04141）、PPAR信号通路（Ko03320）、FoxO信号通路（Ko04068）及基底细胞癌（Ko05217），富集结果如图3所示。此外，差异表达基因KRT38、KRT32、KRT16、KRT25、KRT27、KRT15和KRT36参与雌激素信号通路；COL15A1参与蛋白质消化和吸收通路；TGFBR2、CAT参与FoxO信号通路；SHH、Fzd2参与基底细胞癌通路。COL15A1、TGFBR2和Fzd2在细毛组的毛囊中上调表达，CAT在细毛组的毛囊中下调表达。

2.4　蛋白互作网络的构建

为进一步探讨兔毛细度的分子调控机制，建立差异表达基因互作网络（图4）。构建162个差异基因点及173条互作关系线的蛋白互作网络图。其中SHH、LGR5、ACE2、IL18、GCG、CTSB、CD68和CRNN为主要核心基因，且SHH、LGR5基因Degree分数为9，连接网络最多，ACE2、IL18、GCG、CTSB和CD68基因次之，Degree分数为7。SHH、LGR5和ACE2在长毛兔细毛毛囊组织中上调表达，IL18、GCG、CTSB、CD68和CRNN下调表达。

2.5　qRT-PCR结果验证

为验证RNA-seq数据的可靠性，8个差异表达基因CRNN、TCHHL1、FABP9、KRT71、KRT32、KRT36、KRT82和KRT38进行qRT-PCR验证 （图5）。CRNN、TCHHL1、FABP9和KRT71基因在粗毛毛囊组中上调表达，KRT32、KRT36、KRT82和KRT38在粗毛毛囊组中下调表达，且其qRT-PCR结果与测序结果趋势基本一致。

3　讨 论

本研究通过对长毛兔粗细毛毛囊的转录组进行分析，比较粗毛和细毛毛囊mRNA表达谱，共鉴定412个差异表达基因，其中78个差异表达基因在参考基因组未被注释，另外334个差异基因为已知基因，相比粗毛毛囊组，160个基因在细毛毛囊组中上调表达，252个基因下调表达。GO功能富集分析结果显示差异表达基因主要显著富集在蜕皮周期（GO：0042303）、发周期（GO：0042633）、毛发周期的调节（GO：0042634）、超分子纤维（GO：0099512）等与毛发生长有关的条目中，同时本研究发现了甘油三酯代谢（GO：0006641）、中性脂质代谢过程（GO：0006638）、酰基甘油代谢过程（GO：0006639）等与脂类代谢相关的条目。已有研究指出，脂质代谢对毛发的生长有影响^［22］且脂质组成的变化，可导致毛囊发育的改变^［23］，此外，还有对雌激素的反应（GO：0043627）与雌激素代谢有关的条目，雌激素代谢对毛囊发育周期的影响已有报道，雌激素受体反应的变化，导致毛囊发育周期的变化^［24］。这些条目所富集到的基因作用于2组样本间，导致粗、细毛生长分化的差异。进一步分析差异基因的KEGG通路发现，差异基因显著富集在雌激素信号通路（Ko04915）、蛋白质消化吸收等通路（Ko04974）中。有研究指出，雌激素信号通路（Ko04915）与毛囊生长发育的周期相关^［25-26］，KRT38、KRT32、KRT16、KRT25、KRT27、KRT15 和 KRT36 等角蛋白基因也参与到雌激素信号通路的调控过程中。已有大量研究证明，KRTs和KRTAPs基因富集于毛发角蛋白及其相关蛋白质，是毛发纤维的主要结构蛋白^［17］，与毛囊的发育和分化有关，决定毛发纤维的生长过程^［27］，且KRATPs基因在决定绒毛纤维结构方面可能更为重要^［28］。本研究总共检测到14个KRTs基因和4个KRTAPs基因，其中KRT40、KRT79、KRT15、Krt72、KRT32、Krt36、Krt82、KRT77和KRT38；KRTAP24-1、Krtap26-1、Krtap16-1和KRTAP13-2基因在细毛毛囊中高表达，表明这些基因或与长毛兔细毛的生长和分化有关，该研究结果也佐证了上述结论。蛋白质消化和吸收（Ko04974）被报道与山羊粗、细毛的生长有关^［29］，该通路的调控因子COL15A1可抑制真皮细胞衰老并诱导球状体形成^［30］。同样的内质网中的蛋白质加工途径（Ko04141）也与山羊粗、细毛的生长有关^［31］。PPAR信号通路（Ko03320）与纤维卵泡、纤维毛囊发育相关^［32］。FoxO信号通路（Ko04068）的调控因子TGFBR2是毛囊即将发生的标志物^［33］，调控因子CAT其毛囊中的活性水平与羊毛生长速率相关，是影响羊毛生产或性能的因素^［34］，可用来衡量毛发的生长^［35］，在转基因小鼠中也与毛发发生过程中的调控相关^［36］。CAT在细毛组的毛囊中下调表达，这可能与CAT仅在毛发生长的活跃阶段才在其毛囊中合成有关。基底细胞癌（Ko05217）通路也被认为与毛发生长有关^［37］，该通路的调控因子SHH是初级和次级毛发向下生长和完全次生毛发长度所必需的基因^［38］，Fzd2是毛囊功能和维护过程中的关键基因^［39］。在本研究中，与毛囊发育相关的基因COL15A1、TGFBR2和Fzd2在细毛组的毛囊中上调表达，表明它们在细毛的纤维结构中具有重要作用。

为进一步探究DEGs的功能，进行了PPI分析，过滤了关键通路中高度连接的枢纽节点和基因，根据Degree分数的大小确定核心基因。其中SHH基因可改善毛发生长周期显著性抑制剂成纤维细胞生长因子的功能^［40］是基底细胞癌通路的重要调控因子。G蛋白偶联受体5（LGR5）在活跃增殖和多能干细胞群中存在，是毛囊干细胞表面标志物^［41-42］，能够产生新的毛囊^［43］。有研究表明，敲除小鼠毛囊细胞中的LGR5时，抑制了毛发再生过程^［44］。SHH和LGR5基因在细毛毛囊中上调表达，表明它们在决定长毛兔细纤维结构中具有重要作用。组织蛋白酶B（CTSB）活性在控制毛囊形态发生中必不可少^［45］，缺乏时对小鼠皮肤蛋白质组和脱皮都有影响^［46］。CRNN是调控山羊毛囊和羊毛弯曲的发育的关键候选基因^［29］。CTSB和CRNN基因在粗毛毛囊中上调表达，表明它们可能与粗毛的生长分化有关。肌动蛋白α2（ACTA2）主要在生产细毛的次级毛囊真皮细胞中表达^［47］，在本研究中ACTA2基因在粗毛毛囊中上调表达，与前人的研究结果相悖。差异基因的上下调表达是一种高度组织依赖性的瞬时表达结果，可能由于本研究与先前的研究在品种、采样方式及时间的选择上不同所导致。综上所述，在Degree分数靠前的核心基因中SHH、LGR5、CTSB、ACTA2和CRNN基因与粗、细毛生长有关，可能是与粗、细毛生长有关的关键候选基因。此外，本研究使用qRT-PCR选择了8个基因进行验证，结果与RNA-seq数据相匹配，证明了测序数据的可靠性。

4　结 论

本研究通过长毛兔粗、细毛毛囊的转录组谱，结合生物信息学方法分析筛选出了与粗、细毛生长分化相关的潜在功能基因，“雌激素信号通路，蛋白质消化吸收和代谢通路”、“内质网中的蛋白质加工”、“FoxO信号通路”、“基底细胞癌”等通路是可能影响长毛兔粗、细毛生长的关键通路；SHH、LGR5、CTSB、ACTA2和CRNN基因可能是与粗、细毛生长分化有关的关键候选基因。本研究鉴定并筛选出调控兔毛直径的潜在功能基因及通路，但仍需进一步的研究来验证此发现，以期为细毛型长毛兔的分子育种研究提供重要参考。

参考文献

[1] 武拉平， 王建勋， 秦应和. 2022年兔产业生产概况、2023年发展趋势及政策建议［J］. 中国畜牧杂志， 2023， 59（3）： 348-352

[2] Wu L P， Wang J X， Qin Y H. General situation of rabbit industry production in 2022， development trend in 2023 and policy suggestions ［J］. Chinese Journal of Animal Husbandry， 2023， 59（3）： 348-352 （in Chinese）

[3] Yue Y J， Guo T T， Liu J B， Guo J， Yuan C， Feng R L， Niu C E， Sun X P， Yang B H. Exploring differentially expressed genes and natural antisense transcripts in sheep （ Ovis aries ） skin with different wool fiber diameters by digital gene expression profiling ［J］. PLoS One， 2015， 10（6）： e0129249

[4] Danilenko D M， Ring B D， Yanagihara D， Benson W， Wiemann B， Starnes C O， Pierce G F. Keratinocyte growth factor is an important endogenous mediator of hair follicle growth， development， and differentiation. Normalization of the nu/nu follicular differentiation defect and amelioration of chemotherapy-induced alopecia ［J］. The American Journal of Pathology1995， 147（1）： 145-154

[5] 李祥翔， 夏茂华， 普天春， 卢岩， 张成林. 转录组测序： RNA-Seq技术在动物疾病检测中的应用研究进展［J］. 畜牧兽医科技信息， 2022（1）： 14-16

[6] Li X X， Xia M H， Pu T C， Lu Y， Zhang C L. Progress in the application of transcriptome sequencing-RNA-SEQ technology in animal disease detection ［J］. Chinese Journal of Animal Husbandry and Veterinary Medicine， 2022 （1）： 14-16 （in Chinese）

[7] Stoick-Cooper C L， Moon R T， Weidinger G. Advances in signaling in vertebrate regeneration as a prelude to regenerative medicine［J］. Genes & Development， 2007， 21（11）： 1292-1315

[8] Kahata K， Dadras M S， Moustakas A. TGF-β family signaling in epithelial differentiation and epithelial-mesenchymal transition ［J］. Cold Spring Harbor Perspectives in Biology Biol， 2018， 10（1）： a022194

[9] Eckert R L， Welter J F. Transcription factor regulation of epidermal keratinocyte gene expression［J］. Molecular Biology Reports， 1996， 23（1）： 59-70

[10] 牛晓艳， 任克良， 曹亮， 李燕平， 郑建婷， 冯国亮， 黄淑芳. 利用转录组测序方法研究獭兔毛色相关基因［J］. 中国草食动物科学， 2016， 36（2）： 1-6

[11] Niu X Y， Ren K L， Cao L， Li Y P， Zheng J T， Feng G L， Huang S F. Transcriptome analysis on coat color related genes in rex rabbits ［J］. China Herbivore Science， 2016， 36（2）： 1-6 （in Chinese）

[12] 施丽娟. 獭兔青紫蓝毛色候选基因筛选的初步研究［D］.扬州：扬州大学， 2015

[13] Shi L J. Preliminary study of genes related to fur color in chinchilla rex rabbit［D］.Yangzhou：Yangzhou University， 2015 （in Chinese）

[14] Zhao B， Chen Y， Hao Y， Yang N， Wang M， Mei M， Wang J， Qiu X， Wu X. Transcriptomic analysis reveals differentially expressed genes associated with wool length in rabbit［J］. Animal Genetics， 2018， 49（5）： 428-437

[15] Zhao B H， Chen Y， Yan X R， Hao Y， Zhu J， Weng Q Q， Wu X S. Gene expression profiling analysis reveals fur development in rex rabbits （Oryctolagus cuniculus）［J］. Genome， 2017， 60（12）： 1060-1067

[16] Zhao B， Chen Y， Hao Y， Yang N， Wang M， Mei M， Wang J， Qiu X， Wu X. Transcriptomic analysis reveals differentially expressed genes associated with wool length in rabbit［J］. Animal Genetics， 2018， 49（5）： 428-437

[17] Huang D W， Ding H S， Wang Y L， Cheng G L， Wang X F， Leng T S， Zhao H L. Hair follicle transcriptome analysis reveals differentially expressed genes that regulate wool fiber diameter in angora rabbits［J］. Biology， 2023， 12（3）：445

[18] Fuchs E. Scratching the surface of skin development［J］. Nature， 2007， 445（7130）： 834-842

[19] Driskell R R， Clavel C， Rendl M， Watt F M. Hair follicle dermal papilla cells at a glance［J］. Journal of Cell Science， 2011， 124（Pt 8）： 1179-1182

[20] Chi W， Wu E， Morgan B A. Dermal papilla cell number specifies hair size， shape and cycling and its reduction causes follicular decline［J］. Development， 2013， 140（8）： 1676-1683

[21] Dong Y， Xie M， Jiang Y， Xiao N Q， Du X Y， Zhang W G， Tosser-Klopp G， Wang J H， Yang S， Liang J， Chen W B， Chen J， Zeng P， Hou Y， Bian C， Pan S K， Li Y， Liu Y X， Wang W L， Servin B， Sayre B， Zhu B， Sweeney D， Moore R， Nie W H， Shen Y Y， Zhao R P， Zhang G J， Li J Q， Faraut T， Womack J， Zhang Y P， Kijas J， Cockett N， Xu X， Zhao S H， Wang J， Wang W. Sequencing and automated whole-genome optical mapping of the genome of a domestic goat （Capra hircus）［J］. Nature Biotechnology， 2013， 31： 135-141

[22] Kim D， Langmead B， Salzberg S L. HISAT： A fast spliced aligner with low memory requirements［J］. Nature Methods， 2015， 12： 357-360

[23] Pertea M， Pertea G M， Antonescu C M， Chang T C， Mendell J T， Salzberg S L. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads［J］. Nature Biotechnology， 2015， 33： 290-295

[24] Li B， Dewey C N. RSEM： Accurate transcript quantification from RNA-Seq data with or without a reference genome ［J］. BMC Bioinformatics， 2011， 12： 323

[25] Love M I， Huber W， Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2［J］. Genome Biology， 2014， 15（12）： 550

[26] Stenn K S， Karnik P. Lipids to the top of hair biology［J］.The Journal of Investigative Dermatology， 2010， 130（5）： 1205-1207

[27] Zhang D W， Tomisato W， Su L J， Sun L， Choi J H， Zhang Z， Wang K W， Zhan X M， Choi M， Li X H， Tang M， Castro-Perez J M， Hildebrand S， Murray A R， Moresco E M Y， Beutler B. Skin-specific regulation of SREBP processing and lipid biosynthesis by glycerol kinase 5［J］. Proceedings of the National Academy of Sciences of the United States of America， 2017， 114（26）： E5197-E5206

[28] Rinaldi F， Trink A， Mondadori G， Giuliani G， Pinto D. The menopausal transition： Is the hair follicle “going through menopause”［J］. Biomedicines， 2023， 11（11）： 3041

[29] Sun H R， Meng K， Wang Y F， Wang Y Y， Yuan X C， Li X H. LncRNAs regulate the cyclic growth and development of hair follicles in Dorper sheep［J］. Frontiers in Veterinary Science， 2023， 10： 1186294

[30] Wang J F， Sui J， Mao C， Li X R， Chen X Y， Liang C C， Wang X H， Wang S H， Jia C L. Identification of key pathways and genes related to the development of hair follicle cycle in Cashmere goats［J］. Gene， 2021， 12（2）：180

[31] Rogers M A， Langbein L， Praetzel-Wunder S， Winter H， Schweizer J. Human hair keratin-associated proteins （KAPs）［J］. International Review of Cytology， 2006， 251： 209-263

[32] Wu D D， Irwin D M. Evolution of trichocyte keratin associated proteins［J］. Advances in Experimental Medicine and Biology， 2018， 1054： 47-56

[33] Liu Y， Ding Y Y， Liu Z F， Chen Q， Li X， Xue X L B， Pu Y B， Ma Y H， Zhao Q J. Integration analysis of transcriptome and proteome reveal the mechanisms of goat wool bending［J］. Frontiers in Cell and Developmental Biology， 2022， 10： 836913

[34] Kim H， Choi N， Kim D Y， Kim S Y， Song S Y， Sung J H. TGF-β2 and collagen play pivotal roles in the spheroid formation and anti-aging of human dermal papilla cells ［J］. Aging， 2021， 13（16）： 19978-19995

[35] Ji D J， Yang B， Li Y J， Cai M Y， Zhang W， Cheng G H， Guo H Y. Transcriptomic inspection revealed a possible pathway regulating the formation of the high-quality brush hair in Chinese Haimen goat （Capra hircus）［J］. Royal Society Open Science， 2018， 5（1）： 170907

[36] Wu X M， Gu Y H， Li S Q， Guo S W， Wang J Q， Luo Y Z， Hu J， Liu X， Li S B， Hao Z Y， Li M N， Shi B G. RNA-seq reveals the roles of long non-coding rnaRNAs （lncRNAs） in Cashmere fiber production performance of Cashmere goats in China ［J］. Genes， 2023， 14（2）： 384

[37] Paus R， Foitzik K， Welker P， Bulfone-Paus S， Eichmüller S. Transforming growth factor-beta receptor type I and type II expression during murine hair follicle development and cycling［J］.The Journal of Investigative Dermatology， 1997， 109（4）： 518-526

[38] Damak S， Jay N P， Barrell G K， Bullock D W. Targeting gene expression to the wool follicle in transgenic sheep［J］. Bio/Technology， 1996， 14： 181-184

[39] Vogeli G， Wood L， McNab A R， Kaytes P， Wagner T E， Rea T J， Groppi V， Waldon D J， Kawabe T T， Buhl A E. High-sulfur protein gene expression in a transgenic mouse［J］. Annals of the New York Academy of Sciences， 1991， 642： 21-30

[40] McNab A R， Andrus P， Wagner T E， Buhl A E， Waldon D J， Kawabe T T， Rea T J， Groppi V， Vogeli G. Hair-specific expression of chloramphenicol acetyltransferase in transgenic mice under the control of an ultra-high-sulfur keratin promoter ［J］. Proceedings of the National Academy of Sciences of the United States of America， 1990， 87（17）： 6848-6852

[41] Ding H S， Zhao H L， Cheng G L， Yang Y X， Wang X F， Zhao X W， Qi Y X， Huang D W. Analyses of histological and transcriptome differences in the skin of short-hair and long-hair rabbits ［J］. BMC Genomics， 2019， 20（1）： 140

[42] Cui C Y， Kunisada M， Childress V， Michel M， Schlessinger D. Shh is required for Tabby hair follicle development［J］. Cell Cycle， 2011， 10（19）： 3379-3386

[43] Brennan-Crispi D M， Frankfurter M， Murphy C， Sheng E， Xu M G， Morrisey E E， Millar S E， Leung T H. FZD2 regulates murine hair follicle function and maintenance［J］. The Journal of Investigative Dermatology， 2022， 142（8）： 2260-2263.e2

[44] Zhang R， Li Y， Jia K， Xu X L， Li Y， Zhao Y， Zhang X S， Zhang J L， Liu G S， Deng S L， Lian Z X. Crosstalk between androgen and Wnt/β-catenin leads to changes of wool density in FGF5-knockout sheep ［J］. Cell Death & Disease， 2020， 11： 407

[45] Polkoff K M， Gupta N K， Green A J， Murphy Y， Chung J， Gleason K L， Simpson S G， Walker D M， Collins B， Piedrahita J A. LGR5 is a conserved marker of hair follicle stem cells in multiple species and is present early and throughout follicle morphogenesis［J］. Scientific Reports， 2022， 12： 9104

[46] Yan W， Hao F， Zhe X S， Wang Y M， Liu D J. Neural， adipocyte and hepatic differentiation potential of primary and secondary hair follicle stem cells isolated from Arbas Cashmere goats［J］. BMC Veterinary Research， 2022， 18（1）： 313

[47] Jaks V， Barker N， Kasper M， van Es J H， Snippert H J， Clevers H， Toftgård R. Lgr5 marks cycling， yet long-lived， hair follicle stem cells ［J］. Nature Genetics， 2008， 40： 1291-1299

[48] Hoeck J D， Biehs B， Kurtova A V， Kljavin N M， de Sousa e Melo F， Alicke B， Koeppen H， Modrusan Z， Piskol R， de Sauvage F J. Stem cell plasticity enables hair regeneration following Lgr5⁺ cell loss［J］. Nature Cell Biology， 2017， 19： 666-676

[49] Oortveld M A W， van Vlijmen-Willems I M J J， Kersten F F J， Cheng T， Verdoes M， van Erp P E J， Verbeek S， Reinheckel T， Hendriks W J A J， Schalkwijk J， Zeeuwen P L J M. Cathepsin B as a potential cystatin M/E target in the mouse hair follicle［J］. FASEB Journal： Official Publication of the Federation of American Societies for Experimental Biology， 2017， 31（10）： 4286-4294

[50] Tholen S， Biniossek M L， Gansz M， Gomez-Auli A， Bengsch F， Noel A， Kizhakkedathu J N， Boerries M， Busch H， Reinheckel T， Schilling O. Deletion of cysteine cathepsins B or L yields differential impacts on murine skin proteome and degradome［J］. Molecular & Cellular Proteomics： MCP， 2013， 12（3）： 611-625

[51] Wang Z Y， Wang Y R， Hui T Y， Chen R， Xu Y N， Zhang Y， Tian H， Wang W， Cong Y Y， Guo S P， Zhu Y X， Zhang X H， Guo D， Bai M， Fan Y X， Yue C， Bai Z X， Sun J M， Cai W D， Zhang X J， Gu M， Qin Y T， Sun Y G， Wu Y Z， Wu R N， Dou X T， Bai W L， Zheng Y Y. Single-cell sequencing reveals differential cell types in skin tissues of Liaoning Cashmere goats and key genes related potentially to the fineness of Cashmere fiber［J］. Frontiers in Genetics， 2021， 12： 726670

基金资助

国家兔产业技术体系(CARS-43-A-4)

安徽省重点研究与开发计划项目(202104a06020032)

安徽省农业科学院第一批学科领军人才计划(LJRC-202107)

安徽省特色农业产业技术体系