高原雨点鸽与詹森鸽胸肌转录组差异表达基因筛选与功能分析

兰国湘，金思琪，李星润，刘喜雨，李国美，董新星，*

(1.云南农业大学 动物科学技术学院，云南 昆明 650201；2.大理农林职业技术学院 动物科学学院，云南 大理 671003)

随人民生活水平提高，对鸽子的利用逐渐从肉用向观赏、比赛方向发展，赛鸽运动逐渐兴起。赛鸽应具备坚固结实的骨骼、发达有力的肌肉、紧密坚挺的羽毛。胸肌在飞行中贡献95%以上的力量，是鸟类中最重要的飞行肌肉，与飞翔有关的胸大肌负责向下煽动翅膀，胸小肌负责扬起翅膀，一般胸大肌比胸小肌发达。著名的詹森鸽原产于比利时，其特点是爆发力强，飞行速度快，适宜中短距离赛段。高原雨点鸽是我国第一代军鸽品种，具有远翔的能力和抗高原强磁的特点，耐受力强，可负重飞行。本研究对高原雨点鸽和詹森鸽的胸肌转录组数据进行分析，挖掘调控鸽子与飞行能力有关的功能基因，为赛鸽培育提供参考。

1 材料与方法

1.1 样品采集

选取性别相同、10～12月龄的高原雨点鸽和詹森鸽各3只，屠宰后采集胸肌组织，-80 ℃保存备用。

1.2 纤维性状组织学测定

参考蒋明雅等的方法测定肌纤维横截面积和直径。

1.3 RNA提取与转录组测序

用Trizol法提取总RNA，Nano Photometer分光光度计检测RNA纯度，Qubit3.0Flurometer检测RNA浓度，安捷伦2100 RNA Nano 6000 Assay Kit检测RNA的完整性。RNA检测合格的样品各取3 μg构建转录组文库，用Illumina HiSeq 2500测序平台进行双端测序，得到150 bp的双端测序reads。

1.4 转录组数据分析

1.4.1 数据质控

原始数据去除接头和低质量reads后得到高质量测序数据(clean data)，将clean data通过HISAT2(v2.1.0)比对到参考基因组，获取序列在参考基因组上的位置与特征信息。

1.4.2 差异表达基因筛选

转录本的表达量计算用fragments per kilobase million reads (FPKM)法，以|logFC|>1且<0.05为阈值筛选差异表达基因。

1.4.3 差异表达基因GO与KEGG分析

利用基因本体(gene ontology，GO)富集分析计算每个GO Term的基因数目。用京都基因与基因组百科全书(kyoto encyclopedia of genes and genomes，KEGG)数据库对每个pathway应用超几何检验进行富集分析，确定差异基因参与的信号通路。

1.5 蛋白质互作网络分析

用STRING数据库进行蛋白质互作分析，采用R语言统计与其他蛋白质相互作用关系最多的节点，取前45个进行可视化作图。

1.6 实时荧光定量PCR(qRT-PCR)

提取总RNA，逆转录建cDNA文库，根据、3、1、1、1基因序列设计引物(表1)，引物由擎科生物有限公司合成。以反转录得到的胸肌cDNA为模板进行qRT-PCR，以β肌动蛋白基因(-)作为内参，验证基因表达量，每个样本进行3 次重复。PCR扩增体系为20 μL：cDNA模板1 μL，上、下游引物各0.8 μL，2×SYBRSelect MasterMix 10 μL，ddHO 7.4 μL。PCR反应条件为：95 ℃预变性1 min；95 ℃变性10 s，55 ℃退火5 s，72 ℃延伸10 s，40个循环；用2法计算基因的相对表达量。

表1 实时荧光定量PCR引物序列Table 1 Primer sequences of qRT-PCR

2 结果与分析

2.1 测序数据质量评估

6个样品的RNA测序文库共得297.16 M读长为150 bp的双末端原始读段(raw reads)，产生高质量碱基417.25 Gb，高质量序列所占比例均在90%以上，Q30均在93%以上，所有样品比对基因组上的序列比例均在83%以上(表2)，说明测序质量较高，可用于后续分析。

表2 测序输出结果Table 2 Output of sequencing data

2.2 高原雨点鸽与詹森鸽胸肌转录组差异基因筛选

以|logFC|>1、<0.05为筛选差异显著基因的阈值，共检测到75个显著差异表达基因，其中，高原雨点鸽相对于詹森鸽49个基因上调表达，26个基因下调表达(图1)。

图1 高原雨点鸽与詹森鸽胸肌差异表达基因火山图Fig.1 Volcano plots of DEGs in breast muscle between Plateau raindrop pigeon and Janssen pigeon

2.3 差异显著基因GO富集分析

GO分析结果(表3)显示：差异基因主要富集在肌肉生长相关条目(JAK-STAT级联、骨骼肌细胞分化、肌管分化调控等)、能量利用相关条目(支链氨基酸代谢过程、脂质分解过程、支链氨基酸代谢过程等)、肌肉生长与能量利用相关条目(细胞代谢过程调节、初级代谢过程调节等)。

表3 肌肉生长、能量利用相关的差异表达基因主要显著富集GO条目Table 3 Main GO items of DEGs enrichment analysis related to muscle growth and energy utilization

2.4 差异基因KEGG富集分析

KEGG通路富集分析表明，差异表达基因主要富集在22条信号通路，其中，胰岛素信号通路、胰岛素抵抗等7条通路与肌肉生长、能量利用均有关(表4)，细胞因子信号传导抑制因子3基因(suppressor of cytokine signaling 3，3)、磷脂酶Cε1基因(phospholipase cepsilon-1，1)、信号转导和转录活化因子1基因(signal transducer and activator of transcription 1，1)、乙酰辅酶A羧化酶β基因(acetyl-CoA carboxylase beta，)、蛋白磷酸酶1调节亚基3C基因(protein phosphatase 1 regulatory subunit 3C，13)、低氧诱导因子1α基因(hypoxia inducible factor 1 alpha，1)等与肌肉生长、能量利用有关的基因显著富集于以上通路。

表4 肌肉生长与能量利用相关的差异表达基因主要显著富集信号通路Table 4 Pathways of DEGs enrichment analysis related to muscle growth and energy utilization

2.5 高原雨点鸽与詹森鸽胸肌差异表达基因编码的蛋白质互作网络分析

为进一步研究差异表达基因编码的蛋白质之间的调控关系，用STRING数据库进行蛋白互作分析。从图2可以看出，成纤维细胞生长因子6 (fibroblast growth factor 6，6)基因位于网络起始端，1、生肌决定因子1(myogenic differentiation 1，1)、血管内皮细胞生长因子A(vascular endothelial growth factor A，)位于网络核心位置；互作网络分为3大支，第1支由3、1、低氧诱导因子-1α(hypoxia inducible factor 1 subunit alpha，1)、1、SET和MYND结构域蛋白1(SET and MYND domain containing 1，1)等基因组成，通过自身编码蛋白激活转录调控肌肉细胞分化，或编码蛋白与受体激活作为转录激活因子参与调控肌细胞的分化；第2支由Mg/Mn依赖的蛋白磷酸酶1K(protein phosphatase，Mg/Mndependent 1K，1)、丙酮酸脱氢酶激酶4(pyruvate dehydrogenase kinase 4，4)、丙酮酸脱氢酶磷酸酶催化亚基1 (pyruvate dehydrogenase phosphatase catalytic subunit 1，1)等基因组成，通过磷酸化过程，参与调节葡萄糖代谢，增加能量生成；第3支由1、单酰甘油脂肪酶(monoglyceride lipase，)、脂肪酰胺水解酶(Fatty acid amide hydrolase，)等基因组成，编码蛋白与酶类通过催化反应转化为甘油，参与脂肪代谢。

圆圈代表蛋白质，其颜色深浅表示在网络中的权重，颜色越深表示权重越大。直线代表蛋白质之间的相互作用，线条越粗表示相互作用越大。Circles represented proteins，and their color shades represented their weights in the network，with darker colors representing greater weights.Straight lines represented interactions between proteins，with thicker lines representing grater interactions.图2 蛋白质互作网络Fig.2 Protein interaction network

2.6 高原雨点鸽与詹森鸽胸肌差异基因qRT-PCR

为检验RNA-seq数据可靠性，从与肌肉生长和能量利用有关的差异表达基因中随机选择、3、1、1、1基因进行qRT-PCR验证。转录组测序结果表明：高原雨点鸽胸肌1、1、1基因相对表达量高于詹森鸽，而、3基因相对表达量低于詹森鸽。qRT-PCR结果表明：高原雨点鸽胸肌1、1、1基因相对表达量高于詹森鸽，、3基因相对表达量低于詹森鸽，与转录组测序结果一致(图3)。

2.7 高原雨点鸽和詹森鸽胸肌组织形态比较

詹森鸽和高原雨点鸽胸肌组织的石蜡切片(HE染色)见图4。詹森鸽胸肌肌纤维横切面能看到轮廓较为完整的肌纤维，肌纤维直径(22.77±2.05)μm，肌纤维密度为(1 443±256)根·mm，肌内脂肪沉积较充分；高原雨点鸽肌纤维直径(20.27±1.72)μm，肌纤维密度为(1 553±159)根·mm，肌内脂肪开始沉积。

A，转录组测序结果；B，qRT-PCR结果。A，Transcriptome sequencing results;B，qRT-PCR results.图3 高原雨点鸽和詹森鸽胸肌差异表达基因的相对表达水平Fig.3 Relative expression levels of DEGs in breast muscle between Plateau raindrop pigeon and Janssen pigeon

A，高原雨点鸽(40×40)；B，詹森鸽(40×40)。A，Plateau raindrop pigeon(40×40);B，Janssen pigeon(40×40).图4 高原雨点鸽和詹森鸽胸肌肌纤维横切面Fig.4 Breast muscle fiber transection of Plateau raindrop pigeon and Janssen pigeon

3 讨论

鸽子飞行依赖于发达的胸肌，飞行时胸肌收缩消耗大量能量，能量供给对赛鸽飞行至关重要，能量产生速率决定了骨骼肌的工作能力。本试验对高原雨点鸽、詹森鸽胸肌转录组进行测序分析，筛选到1、、1、3、1、1、1、等与肌肉生长和能量代谢有关的差异表达基因。

3.1 与肌肉生长相关的差异表达基因

1在心脏和肌肉特异表达，具有激活下游基因转录的功能，可促进肌肉细胞成熟、分化。王娟等研究表明，1可能是胰岛素样生长因子-1(insulin-like growth factor 1，-1)下游靶基因，1可能与存在协同激活效应，促进肌纤维分化。早期成肌细胞1缺失并不影响增殖，但抑制分化，导致肌纤维减少，骨骼肌特异性基因的表达减少。本试验中，高原雨点鸽1基因相对表达量比詹森鸽高，1基因富集于骨骼肌细胞分化、肌管分化、成肌细胞分化、横纹肌细胞分化等条目，可能导致高原雨点鸽胸肌肌纤维分化，形成更多的红肌纤维，耐力较好，更适宜中长距离飞行。

3.2 与能量利用有关的差异表达基因

基因可编码丝氨酸水解酶，催化单酰基甘油转化为游离脂肪酸和甘油，为能量消耗提供燃料。基因表达时，会降解2-花生四烯酰甘油，促进脂质储存。本试验中，基因在高原雨点鸽胸肌的相对表达量比詹森鸽低，主要富集于脂肪酸生物合成过程、二酰基甘油代谢过程、甘油三酯分解过程等条目，可能导致高原雨点鸽胸肌脂质减少，快速飞行时能量供应不足，爆发力不如詹森鸽。有研究表明，赛鸽飞翔骨骼肌运动能量消耗加大，主要靠水解游离脂肪酸、甘油提供能量。

1基因可编码支链α-酮酸脱氢酶，促进支链氨基酸分解。本试验中，高原雨点鸽1基因相对表达量比詹森鸽高，该基因富集于支链氨基酸分解过程、支链氨基酸代谢过程等条目，该基因可通过催化支链氨基酸分解，为高原雨点鸽长距离飞行提供足够的能量。

3.3 同时调控肌肉生长和能量利用的基因

本试验中，高原雨点鸽3基因表达量比詹森鸽低，3基因是细胞因子信号转导抑制因子家族的成员之一，影响肌肉组织中多种细胞因子的信号转导，可通过负反馈方式调节细胞因子与激素介导的酪氨酸蛋白激酶/信号传导子和转录激活子信号传导通路，促进成肌细胞分化。同时，3可正调控胰岛素样生长因子(IGF-1)信号通路，IGF-1可促进骨骼肌肥大、肌纤维体积增加。高原雨点鸽3基因表达量低，可能导致高原雨点鸽胸肌肌纤维直径变小，肌肉收缩时的力量不如詹森鸽，比赛时爆发力比詹森鸽差，组织学观测结果也表明，詹森鸽肌肉生长可能早于高原雨点鸽，詹森鸽胸肌肌纤维处于肥大期，肌纤维直径较大、发育较为充分，已经可以看到较为明晰的肌纤维轮廓，而高原雨点鸽则处于肌纤维生长期。另一方面，3基因是瘦素受体信号转导的负反馈调节因子，瘦素受体可激活JAK-STAT通路调控靶基因表达，调控脂肪沉积。高浓度瘦素可诱导3 mRNA和蛋白质在成熟脂肪细胞表达。当3基因表达时，抑制瘦素受体信号转导，导致长链脂肪酸合成与氧化速度加快，促黑激素分泌减少，食欲增加，脂肪沉积增加。高原雨点鸽3基因相对表达量下调，可能导致高原雨点鸽食欲降低、采食量减少、脂肪沉积减少，能量底物减少，从而导致爆发力不如詹森鸽。组织学观测表明，詹森鸽胸肌可以看到明显的脂滴，脂肪沉积较充分，而高原雨点鸽胸肌脂滴较少。

本试验中，高原雨点鸽1基因相对表达量比詹森鸽高，1基因编码蛋白为磷脂酶，催化磷脂酰肌醇-4，5-二磷酸生成三磷酸肌醇，其受体是肌肉特异性表达转录因子，可介导Ca释放，对成肌细胞分化和肌肉收缩至关重要。原代培养的大鼠骨骼肌成肌细胞肌管中含有大量的磷酸肌醇。成肌细胞分化为肌管时，磷脂酶活性增加。高原雨点鸽1基因相对表达量比詹森鸽高，可能导致高原雨点鸽胸肌肌管数量增加、肌纤维增多、肌纤维密度增加，耐力更好，更适合长距离飞行。另外，1基因编码的磷脂酶Cε1能水解磷脂生成二酰甘油和三磷酸肌醇，二酰甘油作为重要的第二信使分子，可激活蛋白激酶C，进一步激活一系列信号传递级联途径，上调胆固醇酰基转移酶-1，促进脂质沉积。高原雨点鸽1基因相对表达量比詹森鸽高，可能通过磷脂酰肌醇特异性磷脂酶C水解磷脂生成二酰甘油和三磷酸肌醇，激活蛋白激酶C，进一步激活一系列信号传递级联途径，上调胆固醇酰基转移酶-1，促进脂质沉积，为高原雨点鸽长距离负重飞行提供足够能量底物。

蛋白质互作网络分析显示，1、1、处于蛋白质互作网络的核心位置。1基因是信号转导子和转录激活子STATs蛋白家族成员之一，当1基因缺失时，骨骼肌生长加速，骨骼肌中甘油三酯含量降低。1基因在高原雨点鸽胸肌中上调，可能导致高原雨点鸽胸肌生长慢于詹森鸽，胸肌中甘油三酯含量高于詹森鸽。1基因是主要的生肌因子，是骨骼肌生长中基因表达的正调节剂。1基因转染到鼠源成纤维细胞，可促进60%～90%的成纤维细胞启动肌分化。1基因缺失可将卫星细胞来源的原代成肌细胞转化为棕色脂肪细胞。1基因在高原雨点鸽胸肌中下调，可能导致高原雨点鸽胸肌肌纤维细胞的肌分化慢于詹森鸽、肌肉生长变慢，脂肪细胞分化加快。基因可通过调节生长因子、转录因子、黏连蛋白基因的表达来调节脂肪组织的结构及功能，过表达会增加血管密度以抵制肥胖，使白色脂肪棕色化，基因在高原雨点鸽胸肌中上调，可能导致高原雨点鸽体内白色脂肪转化为棕色脂肪，分解白色脂肪为二氧化碳、水和能量，有利于远距离飞行。

4 结论

高原雨点鸽胸肌1、1、、1、1基因相对表达量上调。1可能导致高原雨点鸽胸肌肌纤维分化形成更多的红肌纤维；1可能导致高原雨点鸽胸肌肌管数量增加、肌纤维增多、肌纤维密度增加，同时促进脂质沉积；1可能导致高原雨点鸽胸肌生长慢于詹森鸽，胸肌中甘油三酯含量高于詹森鸽；可能导致高原雨点鸽体内白色脂肪转化为棕色脂肪；1可能催化支链氨基酸分解，最终导致高原雨点鸽飞行时的肌肉爆发力不如詹森鸽，同时为高原雨点鸽长距离负重飞行提供足够能量，耐力更好。

高原雨点鸽胸肌1、3、基因相对表达量下调。1可能导致高原雨点鸽胸肌肌纤维细胞的肌分化慢于詹森鸽，肌肉生长变慢；3可能导致高原雨点鸽胸肌肌纤维直径变小，爆发力比詹森鸽差；可能导致高原雨点鸽胸肌脂质减少，能量供应不足，爆发力不如詹森鸽，同时催化支链氨基酸分解，为高原雨点鸽长距离飞行提供足够的能量，更适宜远距离飞行。