产蛋初期和后期蛋鸡胫骨转录谱的构建及骨代谢相关基因分析

张寅梁,岳巧娴,黄晨轩,陈 辉,王德贺,周荣艳

(河北农业大学动物科技学院,保定 071001)

蛋鸡在开产后对钙有较高的需求,在产蛋初期,蛋鸡可利用性成熟前储存积累的钙和日粮中所摄入的钙用于产蛋,而到了产蛋后期,蛋鸡积累的钙消耗殆尽,衰老导致蛋鸡新陈代谢减慢使得机体处于缺钙状态。当机体钙缺乏时,则会分解皮质骨为蛋壳形成提供钙源进而导致骨质流失,并造成一定的骨骼损伤。而骨骼损伤(龙骨弯曲、骨质疏松、骨折等)是威胁蛋鸡动物健康生产的一个重要问题[1]。高强度的产蛋对蛋鸡骨质健康造成威胁[2],蛋鸡骨骼健康既影响蛋鸡福利,也影响产业发展[3]。随着产蛋周龄增加,蛋壳质量下降[4],是限制蛋鸡超长产蛋期(“蛋鸡100周龄生产500个蛋”计划)的主要问题之一[5],这对实现蛋鸡养殖业延长产蛋期的目标[6]是一重大挑战。蛋鸡骨骼健康与营养、环境和遗传因素有关[7],因此研究蛋鸡骨骼健康问题对改善蛋鸡的动物福利和蛋鸡产业发展有重要意义。

蛋鸡的钙沉积是骨代谢维持稳态后的结果,当参与骨代谢的破骨细胞骨吸收活动增强或成骨细胞骨生成活动减弱时,原有平衡会被打破从而导致骨量逐渐减少[8]。骨代谢异常会威胁到蛋鸡的骨骼健康[9-10],并对蛋壳质量造成一定的影响[11]。骨是蛋鸡生长和生产的重要器官,因为在蛋壳形成期间30%至40%的钙来自骨组织中的钙源沉积。骨组织中的髓质骨是蛋鸡“钙库”的主要部位[12]。在产蛋初期,骨吸收较多的矿物质,以备需要时更容易供给机体[13]。在产蛋后期,由于钙沉积减少,当髓质骨供给不足时更多依赖消耗皮质骨形成蛋壳[14]。

转录组测序技术已广泛应用于挖掘表型形成相关的重要基因上。de Oliveira Peixoto等[15]筛选出5个影响肉鸡骨骼代谢发育的重要候选基因。Rubin等[16]筛选到3个在蛋鸡骨代谢中发挥作用的基因。Yue等[17]对产蛋后期蛋鸡胫骨高断裂强度组和低断裂强度组进行转录组分析,筛选出与产蛋后期骨强度有关的5个基因。目前转录组测序技术已经在差异表达基因的研究中得到广泛应用,但鲜有关于蛋鸡不同产蛋期骨组织的转录组研究。本研究通过对海兰灰蛋鸡的胫骨转录组测序结果进行分析,探索产蛋初期和产蛋后期蛋鸡骨代谢的调控过程,为进一步研究产蛋期对蛋鸡胫骨骨代谢的分子调控机制提供理论依据。

1 材料与方法

1.1 试验样品

以海兰灰蛋鸡为试验动物,采样前连续观察7 d,选择每天正常产蛋的19和79周龄鸡各3只,在产蛋后3 h内采集胫骨组织,将胫骨切割成骨块,置于液氮中速冻,转移至-80 ℃超低温冰箱中冻存,用于转录组测序。以19周龄蛋鸡为对照组,以79周龄蛋鸡为处理组。

1.2 cDNA文库构建和转录组测序

利用TRIzol法提取胫骨组织总RNA,Nanodrop检测RNA的纯度(OD260 nm/OD280 nm)、浓度、核酸吸收峰是否正常,Agilent 2100精确检测总RNA完整性。样品检测合格后,用带有Oligo(dT)的磁珠,通过A-T互补配对与mRNA的ployA尾结合的方式富集mRNA。随后加入fragmentation将mRNA打断成短片段,以mRNA为模板,用随机引物合成一链cDNA,后加入缓冲液、dNTPs和DNA polymerase I合成二链cDNA,随后利用AMPure XP beads纯化双链cDNA。最后用PCR法富集cDNA,构建海兰灰蛋鸡胫骨cDNA文库。文库构建完成后,使用Agilent 2100对文库的insert size进行检测。在Illumina HiSeq平台上对6个胫骨组织cDNA文库进行测序。

1.3 胫骨组织差异表达基因筛选与分析

对原始序列进行FastQC质控,用HISAT2软件将过滤后的Clean-Reads比对到鸡基因组(版本号:GRCg6a),运用DESeq2以|log2(Fold Change)|>1且FDR<0.05作为筛选条件筛选差异表达基因。利用clusterProfiler对差异表达基因进行GO功能富集分析与KEGG富集分析。

1.4 蛋白互作网络分析

利用STRING数据库预测蛋白质互作关系,用Cytoscape对蛋白互作网络进行可视化处理。

1.5 差异表达基因荧光定量PCR验证

随机选择5个差异表达基因,根据NCBI数据库中的基因序列,利用Primer 5.0软件进行引物设计(表1)。将提取的19周龄与79周龄组胫骨组织的总RNA样本反转录成cDNA,进行荧光定量PCR。PCR反应体系(20 μL):Forget-Me-Not EvaGreen qPCR Master Mix 10 μL,ROX 3 μL,上、下游引物(10 μmol·L-1)各0.4 μL,cDNA 2 μL,ddH2O 4.2 μL。反应条件:95 ℃预变性2 min;95 ℃变性5 s,60 ℃退火延伸30 s,共40个循环。每个样品3次重复。内参基因LBR与目的基因在同等条件下进行。采用2-ΔΔCt法计算其相对表达量。基因相对表达量结果用SPSS进行单因素方差分析,并用GraphPad Prism绘图。

表1 引物信息

2 结 果

2.1 RNA测序数据分析

转录组数据经质控后,分别从6个海兰灰蛋鸡胫骨cDNA文库中得到19 034 718～22 668 315条clean reads,所有样本Q30值至少为93.68%(表2),表明测序数据可用于后续分析。将clean reads与鸡的参考基因组进行比对,结果显示每个样品所产生的测序读数在鸡参考基因组上的比对率为83.78%～85.72%(表2)。

表2 转录组数据质量评估及Clean reads与参考基因组比对结果

2.2 蛋鸡胫骨差异表达基因筛选

蛋鸡产蛋初期和后期胫骨的转录组测序共鉴定到17 097个基因,利用DESeq2以|log2(Fold Change)|>1且FDR<0.05为条件共筛选到1 211个差异表达基因,其中645个为下调,566个为上调(图1)。

图1 蛋鸡产蛋初期和产蛋后期胫骨差异表达基因火山图Fig.1 Volcano map of differentially expressed genes in the tibia of laying hens during the early and late stages of laying

2.3 差异表达基因GO功能富集分析

为探究蛋鸡产蛋初期和产蛋后期胫骨差异表达基因的功能,利用clusterProfiler对1 211个差异表达基因进行GO功能富集分析,分别富集至生物过程(biological process,BP)、细胞组分(cellular component,CC)和分子功能(molecular function,MF)。为进一步研究差异表达基因在蛋鸡胫骨代谢中的调节作用,筛选排在前10的显着功能亚类。在BP方面,差异表达基因主要富集在细胞分化(cell differentiation)与细胞发育过程(cellular developmental process);在CC功能方面,差异表达基因主要富集在质膜(plasma membrane)、胞外区(extracellular region)和细胞外间隙(extracellular space);在MF中,差异表达基因主要富集在信号受体结合(signaling receptor binding)、信号受体活性(signaling receptor activity)和分子传感器活性(molecular transducer activity)等方面(图2)。

图2 19和79周龄蛋鸡胫骨差异表达基因TOP10 GO富集分析Fig.2 Analysis of TOP10 GO enrichment of differentially expressed genes in tibia of laying hens aged 19 and 79 weeks old

2.4 差异表达基因KEGG功能富集分析

用clusterProfiler对差异表达基因进行KEGG富集分析。如图3所示,有114个差异表达基因显着富集在11条KEGG通路(P<0.05),其中有5条与骨代谢相关,分别为MAPK信号通路(MAPK signaling pathway)、TGF-β信号通路(TGF-beta signaling pathway)、Toll样受体信号通路(toll-like receptor signaling pathway)、糖胺聚糖生物合成-硫酸软骨素/硫酸皮肤素(glycosaminoglycan biosynthesis-chondroitin sulfate/dermatan sulfate)和糖胺聚糖降解(glycosaminoglycan degradation)。在这5条通路中涉及到59个差异表达基因,其中29个上调,30个下调。其中TGFB1、TGFB2和TGFBR2L均与MAPK信号通路和TGF-β信号通路相关,MYD88、FOS和IL1B都在Toll样受体信号通路和MAPK信号通路上(图4),而糖胺聚糖生物合成-硫酸软骨素/硫酸皮肤素与糖胺聚糖降解通路和上述5条通路间均不存在共有的基因(图5)。

图3 19和79周龄蛋鸡胫骨差异表达基因KEGG富集分析Fig.3 KEGG enrichment analysis of differentially expressed gene in the tibia of laying hens aged 19 and 79 weeks old

图5 糖胺聚糖生物合成-硫酸软骨素/硫酸皮肤素与糖胺聚糖降解通路差异表达基因图谱Fig.5 The differentially expressed genes map of glycosaminoglycan biosynthesis-chondroitin sulfate/dermatan sulfate and glycosaminoglycan degradation pathways

2.5 蛋鸡胫骨骨代谢相关差异表达基因编码蛋白互作

为进一步研究上述骨代谢通路中所涉及的59个差异表达基因的编码蛋白之间的互作关系,利用STRING数据库预测蛋白质互作关系,用Cytoscape对蛋白互作网络进行可视化处理。在这59个差异表达基因中,共筛选出53种编码蛋白。在网络图中(图6),白细胞介素6(IL6)的Degree值(Degree值越高说明节点核心度越高)最高,共筛选出5个hub基因的编码蛋白(Degree≥14)IL6、IL1B、FOS、TGFB2和BMP4。

图6 蛋鸡产蛋初期和后期胫骨骨代谢相关差异表达基因编码蛋白互作网络图Fig.6 The protein interaction network diagram of differentially expressed genes related to bone metabolism in the tibia of laying hens during the early and late stages of egg laying

2.6 差异表达基因荧光定量PCR验证

随机选择5个差异表达基因(ROR1、SOST、SPP1、TGFB、WIF1),以LBR为内参基因,对19和79周龄胫骨组织检测该差异表达基因的相对表达水平。结果如图7所示,所有基因79周龄组表达水平极显着高于19周龄组(P<0.01),与转录组测序结果保持一致(图8)。

*. P<0.05;**. P<0.01图7 随机5个差异表达基因qPCR验证结果Fig.7 qPCR validation results of randomly selecting 5 differentially expressed genes

图8 随机5个差异表达基因RNA-Seq测序结果Fig.8 RNA-Seq results of randomly selecting 5 differentially expressed genes

3 讨 论

本研究对产蛋初期和后期蛋鸡胫骨组织进行转录组测序,共鉴定出1 211个差异表达基因,筛选出59个可能与骨代谢相关的差异表达基因,qPCR随机验证所选基因的mRNA相对表达量与转录组测序结果一致,证明转录组测序结果可靠。此外TGFB1、TGFB2和TGFBR2 L均在MAPK信号通路与TGF-β信号通路上,其中TGFB2和TGFBR2 L表现为上调,而TGFB1为下调。研究表明,TGFB1与下游信号的激活能阻止软骨细胞向脂肪细胞分化和组织退化[18-19],并且其基因表达的变化与该基因的甲基化修饰有关[20]。本研究发现,与产蛋后期相比TGFB1在产蛋初期表达下调,这可能说明该基因参与骨代谢相关细胞的分化调节过程。最近的研究表明,TGFB2对骨髓间充质干细胞成骨具有诱导分化的作用,且Elsafadi等[21-22]提出可将TGFB2作为成骨分化能力的鉴定指标。而TGFBR2 L基因功能未知,目前只发现miR-20a-5p/TGFBR2轴的调控与肝纤维化有关[23]。值得注意的是,糖胺聚糖生物合成-硫酸软骨素/硫酸皮肤素与糖胺聚糖降解相关通路与其他通路均没有共有基因相连。但是硫酸软骨素作为动物软骨基质[24],其合成与分解能够影响骨骼质量,还可作为药物用于治疗关节炎。在本研究中,糖胺聚糖生物合成-硫酸软骨素/硫酸皮肤素相关差异表达基因的差异倍数(log2FoldChange)均小于-1,糖胺聚糖降解相关差异表达基因的差异倍数均大于1,而HYAL3(-1.94)除外,这也说明产蛋后期蛋鸡糖胺聚糖的合成与分解产生了不平衡。因此,推测蛋鸡在产蛋后期的骨代谢变化很可能与糖胺聚糖的合成与分解有关。有研究表明,与HYAL1和HYAL2相比,HYAL3在骨关节炎和细胞外基质关节软骨降解过程中未起到关键作用[25],但本研究发现HYAL3作为调控透明质酸酶-3的差异表达基因表现为下调。Hadley等[26]发现,硫酸角质素(KS)可能作为一种髓质骨活动标记物,在本次试验中发现降解KS的β-氨基己糖苷酶亚单位前体的基因HEXA显着上调,表明产蛋后期髓质骨分解增多。此外,MYD88在Toll样受体信号通路和MAPK信号通路上,有研究发现细菌感染骨细胞的MYD88有调控骨吸收的作用[27],而在本试验中MYD88表现为下调,这说明正常产蛋情况下的骨吸收可能不受MYD88调节,而在机体缺钙情况下可能存在一种调控机制抑制了MYD88的表达进行自我保护。

本研究对差异表达基因进行了蛋白互作分析,发现IL6、IL1B、FOS、TGFB2、BMP4等蛋白节点,且其中的IL1B、FOS、TGFB2在KEGG通路分析中就被发现是多条功能通路的共有基因,这与蛋白结果保持了一致。IL6和IL1B对骨代谢的影响主要集中在人类医学和大鼠方面,研究发现IL6和IL1B的表达增加和血清中蛋白水平的增加易导致骨质疏松,认为是绝经后的人或大鼠雌激素分泌减少影响了IL6和IL1B的表达水平[28-29]。FOS与TGFB2都在骨髓间充质干细胞成骨方面起到重要作用[30-31]。BMP4的突变或缺失会导致小鼠或人类的骨骼发育异常[32]。然而本试验中上述5个基因与骨代谢相关差异表达基因的表达水平结果与他人研究结果均不一致,这可能是由于蛋鸡独特的骨代谢机制与哺乳动物的差异造成的,蛋鸡相较于哺乳动物,由于产蛋需要消耗更多的钙,蛋鸡可以在短时间内进行快速的骨代谢活动来动员钙形成蛋壳,因此蛋鸡与哺乳动物的骨代谢机制可能存在差异,具体原因有待进一步研究。此外,本试验还对这些编码蛋白进行了聚类分析,发现MAPK信号通路、TGF-β信号通路和Toll样受体信号通路的相关基因的编码蛋白互作紧密,而糖胺聚糖生物合成-硫酸软骨素/硫酸皮肤素与糖胺聚糖降解则各聚成两组。因此,本试验推测糖胺聚糖的合成与分解与上述3条通路之间可能不存在明显联系,而其合成与分解受何种机制的调控还有待进一步研究。

4 结 论

试验通过转录组技术构建了6个海兰灰蛋鸡胫骨cDNA文库,筛选出59个与骨代谢相关的差异表达基因,且显着富集到多个信号通路,筛选出5个(IL6、IL1B、FOS、TGFB2、BMP4)与蛋鸡骨代谢相关的重要调控基因。