影响羊肉品质的因素较多，羊肉脂肪酸的组成是影响肉品质的因素之一。营养调控和遗传改良是影响羊肉脂肪酸组成的两个途径。Gravador等[1]使用4种饲料饲喂羔羊，以大豆为主的饲料饲喂的羔羊屠宰后肌肉共轭亚油酸（C9,t11-C18∶2）、C18∶2n-6和n-6不饱和脂肪酸比例高于其余各组。Wachira等[2]研究发现，在饲喂条件相同的情况下，萨福克羔羊背最长肌中18∶3n-3的比例高于弗里斯兰羔羊。杂交改良可以进行优良性状的筛选。刘远等[3]进行努比亚黑山羊和福清山羊背最长肌差异表达基因转录组比较分析，挖掘707个新基因，筛选差异表达基因608个，上调61个，下调547个。路婉茹[4]对兰州大尾羊的不同部位脂肪组织的差异表达基因进行筛选，皮下脂肪组织与肌肉组织中共发现4 620个差异表达基因，3 642个基因上调，978个基因下调。赵珺等[5]进行内蒙古绒山羊骨骼肌肌肉差异研究，挖掘631个新基因，背最长肌和臂三头肌间发现790个差异基因，背最长肌中有342个基因上调，448个基因下调。目前，利用转录组测序进行同一品种羊的不同部位的同种组织、饲养模式对同一品种羊的某种组织的影响和不同生长期的组织差异表达基因的研究较多，但利用转录组测序对萨湖F1代与湖羊背最长肌关键差异表达基因及代谢通路进行发掘和分析的研究较少。本试验进行萨湖杂交F1代和湖羊纯繁羔羊背最长肌转录组测序分析，为深入探索影响萨湖杂交F1代肉质的重要基因和关键代谢通路提供参考。1材料与方法1.1试验材料选取10月龄萨福克与湖羊杂交F1代羔羊（SH）、湖羊纯繁羔羊（HH），两者的饲养条件、数量相同。进行屠宰，宰后采集背最长肌组织，剔除筋膜，装入冻存管（无酶），放入液氮，转存至-80 ℃冰箱保存。1.2测定指标及方法1.2.1创建文库和测序提取SH和HH背最长肌的RNA，建库测序，数据质控获得Clean Reads。1.2.2比对分析将Clean Reads比对到参考基因组得到定位信息。1.2.3差异基因筛选使用subread软件的feature Counts工具过滤掉比对不符合要求的reads，筛选SH和HH间的差异表达基因，获得SH和HH间的差异表达基因集（差异倍数设置为FoldChange≥2，padj值设置为≤0.05）。1.2.4聚类和富集分析根据基因的FPKM值对SH和HH间的差异表达基因集进行聚类分析，将表达模式相近的基因聚在一起。使用ClusterProfile软件分析GO功能和KEGG通路。1.2.5荧光定量PCR验证选出4个差异基因，进行实时荧光定量PCR引物验证。差异表达基因实时荧光定量PCR引物见表1。10.13557/j.cnki.issn1002-2813.2022.09.021.T001表1差异表达基因实时荧光定量PCR引物基因上引物下引物GAPDH（内参）5'ACCAGGGCTGCTTTTAATACT3'5'TTGACTGTGCCGTGGAACTT3'PLA2G55'TAAGCCAGTTGAGCGAGAAT3'5'GTGGATGGAAAGTGACTAAGAG3'AQP75'CCTTATTGTTATCATCGGAGTA3'5'GGTTCTAAAGACCTGTGAGCC3'HSPA65'GAGCAAGATGAAGGAGACGG3'5'GCTGCGAGTCGTTGAAGTAG3'ECI15'AGTCCAGTCTTTGCCCATTT3'5'GCCTAGCCAGTTTATTATCCC3'反转录后进行qPCR，22 μL体系：2×SYBR Green qPCR Master Mix 10 μL、上下引物各1 μL、cDNA 3 μL、ddH2O 7 μL。退火温度60 ℃、30 s。2结果与分析2.1RNA-seq产出数据质量评估（见表2）由表2可知，经Agilent 2100 bioanalyzer分析，提取的RNA符合建库测序要求。测序得到clean reads 41 478 440~47 549 168条，Clean bases均高于6.02 G，过滤后的clean reads数占raw reads比例均大于94.18%，Q30为94.18%~94.88%，碱基正确率达99.9%。10.13557/j.cnki.issn1002-2813.2022.09.021.T002表2RNA-seq产出数据质量评估样品原始reads/条过滤后reads/条过滤后碱基数/G错误率/%碱基占比（Phred值30）/%G和C占比/%SH141 511 73640 103 6866.020.0394.1853.86SH242 272 27840 236 8666.040.0294.5052.12SH341 478 44040 228 9706.030.0294.6153.62HH142 334 42240 981 2406.150.0294.3054.45HH244 062 61042 657 3626.400.0394.2253.58HH347 549 16846 215 7266.930.0294.8854.892.2clean reads与参考基因组比对情况统计（见表3）由表3可知，clean reads比对到参考基因组上的比例最低为88.43%，最高为90.05%，具有唯一位置的clean reads占比最低为79.69%，最高为80.94%。10.13557/j.cnki.issn1002-2813.2022.09.021.T003表3clean reads与参考基因组比对情况统计样品总reads可定位的reads唯一定位reads多点定位readsSH140 103 68635 742 357（89.12%）32 440 191（80.89%）3 302 166（8.23%）SH240 236 86636 234 136（90.05%）32 567 910（80.94%）3 666 226（9.11%）SH340 228 97036 043 882（89.60%）32 374 576（80.48%）3 669 306（9.12%）HH140 981 24036 601 245（89.31%）32 656 548（79.69%）3 944 697（9.63%）HH242 657 36237 722 505（88.43%）34 303 086（80.42%）3 419 41（8.02%）HH346 215 72641 386 716（89.55%）37 339 775（80.79%）4 046 941（8.76%）2.3差异基因统计（见表4）由表4可知，SH对比HH样品中，筛选出差异表达基因322个，114个上调，208个下调。10.13557/j.cnki.issn1002-2813.2022.09.021.T004表4差异基因统计样品总数/个上调/个下调/个阈值SH vs HH322114208DESeq2 P0.05，Fold Change≥22.4差异基因聚类分析（见图1）由图1可知，将SH和HH的差异基因取并集，均一化处理，根据基因表达方式聚类分析，横坐标为样品名，HH取3个重复，分别为H1、H2和H3，SH取3个重复，分别为S1、S2和S3，纵坐标为差异基因FPKM均一化后的数值，颜色越红，表达量越高，越绿，表达量越低。10.13557/j.cnki.issn1002-2813.2022.09.021.F001图1差异表达基因聚类分析2.5差异表达基因GO富集分析（见图2）由图2可知，GO注释278个差异基因，在B-p、M-f和C-c中各参与83个、140个和55个亚功能。生物学过程中蛋白质磷酸化、肽代谢和酰胺合成及代谢富集较多差异基因。分子功能中肽酶活性和作用于L-氨基酸肽的肽酶活性差异基因参与数各占11个，氧化还原酶活性有10个差异基因参与，辅因子结合和蛋白激酶活性各有7个差异基因参与。细胞组成中细胞器部分有20个差异基因参与，核过程、肌动蛋白细胞骨架、线粒体、细胞骨架部分和膜蛋白复合物富集差异基因较多。10.13557/j.cnki.issn1002-2813.2022.09.021.F002图2差异表达基因GO富集分析2.6KEGG通路分析肉质相关信号通路及其富集的差异基因见表5。10.13557/j.cnki.issn1002-2813.2022.09.021.T005表5肉质相关信号通路及其富集的差异基因序号信号通路上/下调差异基因1氨基酸的生物合成上调GPT2、IDH2、ENO3、PGK1下调2甘油酯代谢上调LIPG、LPIN1下调PNPLA3、DGAT23糖酵解/糖异生上调HK2、下调ENO3、PGK14肌动蛋白细胞骨架的调节上调ENAH、CFL2、ARHGAP35下调ITGAD、MYLK3、ITGB75甘油磷脂代谢上调LPIN1下调PLA2G56JAK-STAT信号通路上调CSF3R、IL2RG、JAK3下调7PPAR信号通路上调AQP7下调8Ras通路上调下调LAT、PLA2G59脂肪酸降解上调ECI1下调10脂肪细胞中脂解的调节上调AQP7下调11MAPK信号通路上调MAP2K3、MAP3K14下调HSPA6由表5可知，得到KEGG数据库注释的差异基因有147个，参与216条通路。KEGG通路的柱状图和散点图见图3、图4。10.13557/j.cnki.issn1002-2813.2022.09.021.F003图3差异基因显著富集的KEGG通路10.13557/j.cnki.issn1002-2813.2022.09.021.F004图4差异基因显著富集的KEGG通路散点图由图3、图4可知，人类疾病、免疫系统和氧化磷酸化富集较多差异基因。试验中筛选出可能与肉质相关的通路含有氨基酸的生物合成通路、甘油酯代谢、PPAR信号通路、MAPK信号通路和肌动蛋白细胞骨架的调节通路。2.7荧光定量PCR验证（见图5、表6）由图5、表6可知，Real-time PCR鉴定的表达量，与转录组测序数据表达趋势一致，说明转录组测序获得的差异表达基因信息准确、可靠。10.13557/j.cnki.issn1002-2813.2022.09.021.F005图54个差异基因的荧光定量PCR10.13557/j.cnki.issn1002-2813.2022.09.021.T006表6差异基因的RNA-Seq与Real-timePCR鉴定的表达量基因转录组测序表达倍数荧光定量PCR相对表达量PLA2G5-4.500.54AQP7-1.250.95HSPA6-1.850.65ECI1-1.300.923讨论评价羊肉质量的指标主要有肥度、肉脂硬度、肌肉发育程度和肉脂色泽等。粗蛋白含量、氨基酸含量、脂肪酸、肌肉剪切力等也是影响肉质的因素。影响肉质指标的途径主要有营养调控和遗传改良，杂交改良是从品种方面进行肉品质的提升途径。Jiao等[6]研究发现，n-3多不饱和脂肪酸在藏羊肉中含量高于小尾寒羊，而饱和脂肪酸、单不饱和脂肪酸和n-6多不饱和脂肪酸含量在藏羊与小尾寒羊中差异不显著，说明遗传品种对于羊肉品质的重要性。王杰[7]研究发现，适宜的酪氨酸含量，能够显著提高鸡胸黑色素的沉积量，对部分脂肪酸也具有显著影响。林金法等[8]研究发现，谷氨酰胺可以显著提高肌肉中胶原蛋白和肌苷酸的含量，达到改善肌肉品质的作用。姜锦鹏等[9]研究发现，烟酰胺可以较好地提高肉仔鸡屠宰性能，改善肉质性状，说明影响肉品质的部分因素，物质的含量是可以通过营养调控途径实现。本试验从萨福克羊和湖羊杂交F1代与湖羊转录组比较分析，得到差异基因322个，通过GO功能、KEGG通路等分析，结合肉质相关指标、影响因素，筛选与肉质相关的GO数据库的功能亚分类和KEGG通路以及富集在这些功能亚分类和KEGG通路上的肉质相关差异基因，为进一步研究肉质相关因素提供参考。GO注释278个差异基因，分析出参与生物过程（BP）、分子功能（MF）和细胞组成（CC）的83、140和55个功能亚分类。亚分类中蛋白质水解过程、蛋白质磷酸化过程、肽代谢过程、酰胺合成和代谢过程以及肌动蛋白细胞骨架可能与肉质相关。董玉影[10]研究发现，肌肉中蛋白质的结构与功能特性对肉品质影响显著，蛋白质水解与肉品质关系密切。蛋白质磷酸化主要发生在丝氨酸（包括苏氨酸）和酪氨酸上[11]。白稚子等[12]研究发现，苏氨酸可以改善肉的肌纤维特性及硬度品质。肌动蛋白在肌肉纤维中占细胞总蛋白质的20%。肌动蛋白基因与肌纤维发育、肌肉色泽、肉质嫩度等方面存在一定的关系[13]。通过上述研究基础，结合GO数据库，筛选出与肉质可能相关的亚分类功能，与肉质相关的差异表达基因，进行肉质相关的深入研究。本试验中，萨福克羊与湖羊杂交产生的F1代与湖羊之间的差异基因在KEGG通路分析时，氨基酸的生物合成中差异基因富集显著，富集在该通路中的ENO3基因，也在糖酵解通路中富集。赵玉强等[14]进行ENO3基因与大白猪的眼睑面积和眼肌厚之间的关系，结果二者显著相关。本试验中，差异基因LPIN1富集在甘油酯代谢通路，该基因与脂肪组织形成、分化和沉积密切相关[15]。段安琴等[16]研究发现，在水牛乳腺细胞中LPIN1蛋白表达量增加，会抑制甘油三酯的合成。程漫漫等[17]研究发现，LPIN1基因表达量与甘油三酯、高密度脂蛋白呈显著负相关，与肌内脂肪率呈显著正相关。段安琴等[18]研究发现，LPIN1基因过度表达会上调PPARα、PPARγ、SREBP1C基因，上调p-mTOR、PPARα蛋白，从而增加甘油三酯的表达量。在肌动蛋白细胞骨架的调节通路中富集的CFL2基因，主要编码肌动蛋白解聚因子，该蛋白可以调节肌动蛋白的装配，进而影响肌肉的发育和再生。CFL2基因可能与肉质相关。张继[19]发现，CFL2基因在背最长肌中的表达量和与肉色为显著负相关。张宇婷等[20]研究发现，CFL2基因可以通过调节MyHC类型转换达到调节红白肌肉比例的效果，进而影响肉品质。张亚楠等[21]对比过藏山羊的AQP7基因在各种组织中的表达量，结果发现，AQP7基因在脂肪组织中的表达量极显著高于其他组织。4结论本试验条件下，在转录组水平上筛选萨福克与湖羊杂交F1代、湖羊纯繁羔羊背最长肌差异表达基因322个，得到GO数据库注释的有278个，参与BP、MF和CC的分别为83个、140个和55个，得到KEGG数据库注释的有147个，参与216条通路，差异基因显著富集的通路24条。初步认为与肉质可能相关的通路有氨基酸的生物合成通路、甘油酯代谢通路、PPAR信号通路、MAPK信号通路和肌动蛋白细胞骨架的调节通路。筛选出与肉质相关的候选基因有LPIN1、CFL2、ENO3和AQP7。