基于WGCNA方法识别影响猪肉质性状的关键模块和基因
Identification of Key Gene Co-expression Modules and Genes Associated with Pork Quality Based on WGCNA Method
-
Keywords:
- pig /
- pork quality trait /
- gene expression /
- weighted gene co-expression network /
- key gene
-
百合(Lilium spp.)是百合科(Liliaceae)百合属(Lilium)所有种类的总称,通常为具有根状茎、块茎或鳞茎的多年生草本。西伯利亚百合(Lilium ‘Siberia’)属于东方百合,由天香百合(L. auratum Lindl.)、药百合(L. speciosum var. gloriosoides Baker)等品种杂交而成,由于其花大而美丽,开放时气味芳香,是非常受消费者欢迎的鲜切花之一[1]。切花的生产需要优质的商品种球,鳞片繁殖成本低、操作简便、繁殖系数高,是目前种球商品化生产最常用的技术途径[2]。百合种球从鳞片扦插到成球的生产周期较长,小鳞茎培养成合格的种球一般需要2~3年,如何提高小鳞茎的繁殖系数、促进鳞茎膨大对百合种球的生产至关重要。因此,研究促进百合小鳞茎萌发的分子机制及挖掘调控萌发的潜在基因具有重要意义。
AP2/ERF转录因子是植物中的一大类因子,该家族参与多种生物学过程,包括植物形态建成、病菌防御、对各种胁迫的响应机制、激素信号转导等[3]。AP2/ERF家族转录因子主要分为脱水反应元件结合蛋白(dehydration-responsive element binding protein,DREB)、乙烯响应因子(ethylene response factor,ERF)、APETALA2 (AP2)和ABI3/VP1相关因子(related to ABI3/VP1,RAV)四大亚家族以及少数未分类因子(Soloists)[4]。ERF亚家族不仅参与植物发育和生理过程的多种功能,还参与种子萌发和芽休眠。SlPti4是番茄(Solanum lycopersicum L.)果实成熟、种子萌发和非生物胁迫响应的重要调节因子[5];油桃(Prunus persica var. nectarina)中的ERF家族成员PpeEBB1可调控桃芽萌发[6];在砀山酥梨中,PpyERF060、PpyABF3和PpyMADS71互作可整合乙烯与脱落酸的信号通路,进而调控梨芽休眠进程[7];樱桃(Prunus pseudocerasus)中的PpcERF5表达被脱落酸(abscisic acid,ABA)诱导上调,这可能是其参与樱桃花芽休眠调控的方式之一,过表达PpcERF5基因还可提高转基因拟南芥(Arabidopsis thaliana)种子的萌发率[8]。
褪黑素在调节植物生长、种子发芽、叶片衰老、光合效率、果实成熟等生理生长过程以及调控植物增强抗逆性、抵御细菌和真菌病害方面具有重要作用[9-11]。褪黑素是与植物激素相关的基因表达的重要调节剂,它与植物激素之间表现出复杂的交叉作用,如参与生长素、吲哚乙酸、赤霉素、细胞分裂素、脱落酸、乙烯、水杨酸、茉莉酮、油菜素类固醇的代谢[12-14]。拟南芥中的褪黑素可通过降低IAA17及衰老相关基因的表达抑制植物的衰老进程[15];褪黑素处理可以通过影响氨基酸代谢来提高冷藏果实的可溶性糖含量、香气水平和果实品质[16];番茄中低浓度的外源褪黑素处理可以促进不定根的形成[17];相反,高浓度的褪黑素处理会抑制生长素的生物合成以及转运途径中许多基因的表达[18];用褪黑素处理多年生黑麦草(Lolium perenne L.),增加了热应激下褪黑素的内源含量,从而降低ABA含量[19];外源褪黑素可促进Cr胁迫下百日草(Zinnia elegans)的种子萌发,提高其对Cr胁迫的抗性[20];在黄瓜(Cucumis sativus L.)幼苗中,褪黑素通过调节ABA和赤霉素(gibberellin,GA4)的生物合成和分解代谢以减轻NaCl的胁迫抑制作用,促进发芽[21]。
百合萌发的相关作用机制很复杂,高通量转录组测序技术为探究萌发分子机制提供了途径和方法,可用于系统研究百合萌发过程的相关代谢通路。AP2/ERF转录因子在植物生长发育、代谢、生物和非生物胁迫方面扮演着重要角色,而关于促进百合鳞茎萌发的报道较少。为深入研究百合萌发阶段的分子调控机制,通过对褪黑素处理下的百合小鳞茎进行转录组测序分析,并对差异表达的基因进行功能注释,筛选受褪黑素诱导的LoERF96基因;对该基因进行克隆及生物信息学分析,利用实时荧光定量PCR (RT-qPCR)分析其在小鳞茎4个不同生长发育时期的基因表达量,为探究百合AP2/ERF转录因子在萌发过程中的生物学功能提供理论基础和参考。
1. 材料与方法
1.1 植物材料
供试百合小鳞茎由西伯利亚百合大种球(云南沃德球根生物科技有限公司)鳞片包埋繁殖获得。将西伯利亚百合大种球中外层的鳞片剥下,包埋于椰糠基质中,在28 ℃恒温培养箱中遮光培养45 d,然后将获得的小鳞茎置于4 ℃环境中低温冷藏70 d,冷藏结束后将小鳞茎洗净、晾干待用。
1.2 样本处理
选取生长状态良好、大小一致、鲜质量为(0.20±0.10) g的百合小鳞茎,随机平均分为2组,每组72个(含3次重复),其中一组为处理组,用100 μmol/L褪黑素溶液浸泡2 h;另一组为对照组,用蒸馏水浸泡2 h。将浸泡处理后的小鳞茎放在滤纸上晾干,将小鳞茎栽培到椰糠为基质的营养钵中,每个营养钵栽培3个;将营养钵放入25 ℃的光照恒温培养箱中培养,小鳞茎破土长出芽即认为萌发,分别在第0、7、14和21天取样,液氮速冻后置于−80 ℃超低温冰箱保存,用于后续试验。
1.3 总RNA提取
按照RNAprep Pure多糖多酚植物总RNA提取试剂盒(天根生化科技有限公司,DP441)说明书,分别提取4个不同时期处理组和对照组的西伯利亚百合小鳞茎总RNA,用NanoDrop 2000分光光度计检测RNA的纯度和浓度,用1%琼脂糖凝胶电泳检测RNA的质量;检测合格后,取部分试验第21天的RNA进行转录组测序及目的基因克隆试验;4个时期的全部RNA材料均用于后续的基因表达分析。
1.4 转录组测序
RNA检测合格后进行转录组文库构建,用含有Oligo (dT)的磁珠富集真核生物mRNA将mRNA随机打断,以mRNA为模板合成cDNA第1链和第2链并进行cDNA纯化,然后进行末端修复、加A尾并连接测序接头,PCR富集得到cDNA文库。用N50数值对组装文库进行质量评估检测,Unigene N50越长,说明组装质量越好。 N50的定义为:将所有Unigene从长到短排序,并依次累加长度,当累加片段长度达到总片段长度的50%时,对应的片段长度即为N50。文库质检合格后,用Illumina NovaSeq6000测序平台进行测序(委托北京百迈客生物科技有限公司完成)。
1.5 Unigene功能注释
将Unigene序列与非冗余蛋白质序列(non-redundant protein sequence,NR)数据库、Swiss-Prot蛋白质序列 (Swiss-Prot protein sequence,Swiss-Prot)数据库、原核生物同源蛋白簇 (clusters of orthologous groups,COG)、真核生物同源蛋白簇 (eukaryotic orthologous groups,KOG)、 eggNOG4.5数据库、Pfam数据库(protein family,http://pfam.sanger.ac.uk/)、京都基因与基因组百科全书 (Kyoto encyclopedia of genes and genomes,KEGG)和(gene ontology,GO,http://www.geneontology.org/)数据库比对,获得Unigene的注释信息。
1.6 差异基因的筛选
将差异倍数|log2FC|≥2和P<0.01作为标准筛选差异表达基因,使用DESeq2软件进行差异分析。
1.7 LoERF96序列克隆
将褪黑素诱导21 d的百合小鳞茎总RNA用All-In-One 5× RT MasterMix (abm,G592)逆转录为cDNA,基于转录组数据获得差异基因LoERF96的序列,用Premier 6.0设计扩增全长序列的特异性引物(表1),委托擎科生物科技股份有限公司昆明分公司进行引物合成,然后进行LoERF96基因的cDNA扩增。PCR总反应体系为25.0 µL,包括2×T5 Super PCR Mix (Basic) 12.5 µL,正、反向引物各1.0 µL,cDNA模板1.0 µL,ddH2O 9.5 µL。PCR反应程序为:98 ℃预变性3 min,98 ℃变性10 s,63 ℃退火10 s,72 ℃延伸15 s,共35个循环,循环结束后72 ℃终延伸2 min,4 ℃保存。用1%琼脂糖凝胶电泳检测PCR产物,对扩增得到的目的片段进行切胶,再使用DNA凝胶回收试剂盒(擎科生物科技股份有限公司,TSP602-200)进行胶回收和纯化;将目的片段与pClone007-T载体连接,并转化大肠杆菌DH5α感受态细胞,涂到含Amp抗生素的LB培养基平板后,倒置于37 ℃培养箱过夜培养;次日挑取单菌落并对菌液进行PCR验证,将阳性克隆菌液送至擎科生物科技股份有限公司测序。
表 1 引物序列Table 1. Primer sequence基因名称
gene name正向引物序列 (5'→3')
forward primer sequence反向引物序列 (5'→3')
reverse primer sequence用途
usageLoERF96 ATGCACCAGCTCATGGAAGGAAGC TTATTTGTCCCTGTCGTCGCTCTCCT 基因克隆
gene cloningLoERF96 GGAGCCCAAACACCAACCCAA AAGTGCCGAGCCAGACACGTA 荧光定量PCR
real-time PCRERF105 CAGTTCTAGCACCTGCCTCCCT GTATCGCTCTGGCGTCCGTATCT CDC5 CCGTCTCAAAGGACACCCTACCA TGTGGAGTGGGAGTGTGAAGAGG MYB36 CTGCCGTTGCCGAAGATGAACT CACCTCCACCACCACTACCAGAA WRKY28 TGCGATGGCAACATCTGCATCT ATGCGGGCATTGCAGGAATTGA WOX8 TGGCTCCGAACACCTTGTGGA GCTCCTCGAACTTCCCTTGTGC LOB30 GCTCTCATCTCCGCCGAAACAC TCATCGTGCCTCGCTCTGTTCT LoFP TCGCCTACATCGCTAACC TTCCCAATAATCGCAAGACC 荧光定量内参
reference gene of real-time1.8 RT-qPCR验证转录组及不同时期LoERF96的表达分析
为了验证转录组测序结果的准确性,随机挑选CDC5、LOB30、LoERF105、MYB36、WOX8和WRKY28共6个差异基因进行RT-qPCR验证。用Premier 6.0设计RT-qPCR基因特异性引物(表1),按照BlasTaq™ 2× qPCR MasterMix (abm,G891)试剂盒说明书,在Applied Biosystems QuantStudio 5实时荧光定量PCR仪上进行实时荧光检测,选择百合F-box家族蛋白LoFP基因作为内参基因[22]。RT-qPCR反应体系为20.0 µL,包括 BlasTaqTM 2× qPCR MasterMix 10.0 µL,10 µmol/L正、反向引物各0.5 µL,cDNA 1.0 µL和ddH2O 8.0 µL。 RT-qPCR扩增程序为:95 ℃ 3 min, 95 ℃ 15 s,60 ℃ 1 min,共40个循环。每个样品设置3次生物学重复。采用2−ΔΔCt法计算基因相对表达量;采用Excel统计数据;采用Origin 2021绘制图表;采用SPSS 26进行显著性分析。
1.9 LoERF96的生物信息分析
使用Expasy (https://web.expasy.org/protparam/)分析LoERF96蛋白的序列长度、分子质量、等电点、不稳定指数、脂溶性指数、亲水性等理化性质;使用DeepTMHMM (https://dtu.biolib.com/DeepTMHMM)解析蛋白跨膜区域;使用NCBI中CD-Search (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi)工具预测氨基酸序列的保守结构域;使用SOPMA (https://npsa-prabi.ibcp.fr/cgi-bin/npsa_automat.pl?page=/NPSA/npsa_sopma.html)预测蛋白质的二级结构;使用SWISS-MODE预测蛋白质的三级结构 (https://swissmodel.expasy.org/);使用Plant-mPLoc (http://www.csbio.sjtu.edu.cn/bioinf/plant-multi/)在线分析工具预测蛋白质的亚细胞定位。使用NCBI数据库中的Protein Blast工具寻找LoERF96蛋白氨基酸的同源序列;使用MEGA 11.0进行多序列比对,并采用邻接法构建系统进化树,Bootstraps设置重复
1000 次。2. 结果与分析
2.1 序列的组装和质量控制
完成6个样品的转录组测序,获得的有效数据总量为38.95 Gb,各样品的有效数据量均达到6.03 Gb,GC碱基含量为50.51%,Q30碱基的百分比大于87.88%。组装共获得
132750 个转录本和78100 条Unigene,其中长度大于1 kb的Unigene有20431 条,Unigene的N50为1 355 bp,表明组装的完整性较高(表2)。表 2 转录本和Unigene序列的长度分布Table 2. Length distribution of transcripts and Unigenes序列长度/bp
sequence length转录本 transcripts Unigenes 数量 number 比例/% percentage 数量 number 比例/% percentage >200~300 18515 13.95 15442 19.77 >300~500 34749 26.18 23329 29.87 >500~ 1000 35451 26.70 18898 24.20 > 1000 ~2000 28747 21.65 13313 17.05 > 2000 15288 11.52 7118 9.11 总数量 total number 132750 78100 总长度/bp total length 131646756 66907113 N50长度/bp N50 length 1512 1355 平均长度/bp mean length 991.69 856.69 2.2 Unigene功能注释统计
在NR、Swiss-Prot、COG、KOG、eggNOG4.5、Pfam、KEGG和GO共8个数据库获得
42458 条Unigene的注释结果,每个数据库注释到的基因数量分别为40403 、21492 、10357 、21880 、29112 、26857 、24294 和31149 。2.3 差异基因的筛选结果
对褪黑素处理组和对照组的转录组数据进行差异表达分析,得到780个差异表达基因,其中上调表达基因392个,下调表达基因388个(图1)。
图 1 差异表达基因火山图(a)和聚类热图(b)注:图a)中,FC为差异倍数;图b)中,CK1、CK2和CK3为蒸馏水对照组的3个重复,MT1、MT2和MT3为褪黑素处理组的3个重复。Figure 1. Volcano map and cluster heatmap of differentially expressed genesNote: In Fig. a), FC indicates fold change; in Fig. b), CK1, CK2, and CK3 are three repeats of the distilled water control group; MT1, MT2, and MT3 are three repeats of the melatonin treated group.2.4 差异表达基因的GO分类富集分析和代谢通路KEGG富集分析
所有差异表达基因的GO富集分析结果(图2)显示:在生物过程中富集到差异表达基因主要包括代谢过程、细胞过程、生物调节、刺激响应、信号、生殖过程等;在细胞组分中富集到差异基因功能主要包括细胞、细胞器、膜、大分子复合物等;在分子功能中富集到差异基因功能主要包括结合活性、催化活性、转运活性、核酸结合转录因子活性、分子功能调节等。
图 2 差异表达基因的GO富集图注/Note:生物过程(biological process):代谢过程(metabolic process)、细胞过程(cellular process)、单体过程(single-organism process)、生物调节(biological regulation)、定位(localization)、刺激响应(response to stimulus)、细胞成分组织或生物发生(cellular component organization or biogenesis)、发育过程(developmental process)、信号(signaling)、多细胞生物过程(multicellular organismal process)、生殖(reproduction)、生殖过程(reproductive process)、多有机体过程(multi-organism process)、解毒(detoxification)、免疫系统过程(immune system process)、生长(growth)、节律过程(rhythmic process)、移动(locomotion)、生物黏附(biological adhesion)、行为(behavior)、细胞聚集(cell aggregation)。细胞组分(cellular component):细胞(cell)、细胞部分(cell part)、细胞器(organelle)、膜(membrane)、膜部件(membrane part)、大分子复合物(macromolecular complex)、细胞器部分(organelle part)、膜腔(membrane-enclosed lumen)、细胞外区域(extracellular region)、细胞连接(cell junction)、合胞体(symplast)、其他生物体(other organism)、其他生物体部分(other organism part)、超分子复合物(supramolecular complex)、胞外区部分(extracellular region part)、类核物质(nucleoid)、突触(synapse)、突触部分(synapse part)、病毒粒子(virion)、病毒粒子部分(virion part)。分子功能(molecular function):结合活性(binding)、催化活性(catalytic activity)、转运活性(transporter activity)、结构分子活性(structural molecule activity)、核酸结合转录因子活性(nucleic acid binding transcription factor activity)、分子功能调节(molecular function regulator)、电子载体活性(electron carrier activity)、抗氧化活性(antioxidant activity)、信号转运活性(signal transducer activity)、分子转运活性(molecular transducer activity)、转录因子活性和蛋白质结合(transcription factor activity and protein binding)、营养活性(nutrient reservoir activity)、蛋白质标签(protein tag)、金属伴侣活性(metallochaperone activity)、翻译调节因子活性(translation regulator activity)。Figure 2. GO enrichment map of differentially expressed genes对差异表达基因进行KEGG代谢通路的富集分析,筛选出显著富集的代谢途径50条,主要包括植物激素信号转导、植物MAPK信号途径、核糖体、苯丙烷生物合成、淀粉和蔗糖代谢、光合作用、碳代谢、内质网中的蛋白质加工、植物病原菌互作等(图3)。
图 3 差异表达基因的KEGG代谢通路富集注/Note:细胞过程(cellular processes):自噬其他(autophagy-other),内分泌(endocytosis)。环境信息加工(environmental information processing):植物激素信号转导(plant hormone signal transduction),植物MAPK信号途径(plant MAPK signaling pathway)。遗传信息加工(genetic information processing):核糖体(ribosome),内质网中的蛋白质加工(protein processing in endoplasmic reticulum),RNA转运(RNA transport),氨酰基tRNA生物合成(aminoacyl-tRNA biosynthesis),泛素介导的蛋白水解(ubiquitin mediated proteolysis),真核生物核糖体的生物发生(ribosome biogenesis in eukaryotes)。RNA降解(RNA degradation),RNA聚合酶(RNA polymerase),mRNA监测途径(mRNA surveillance pathway),剪接体(spliceosome)。代谢(metabolism):苯丙烷生物合成(phenylpropanoid biosynthesis),淀粉和蔗糖代谢(starch and sucrose metabolism),光合作用(photosynthesis),色氨酸代谢(tryptophan metabolism),半胱氨酸和蛋氨酸代谢(cysteine and methionine metabolism),糖酵解/糖异生(glycolysis/gluconeogenesis),二苯乙烯、二芳基庚烷和姜酚的生物合成(stilbenoid, diarylheptanoid and gingerol biosynthesis),类黄酮生物合成(flavonoid biosynthesis),碳代谢(carbon metabolism),氧化磷酸化(oxidative phosphorylation),酪氨酸代谢(tyrosine metabolism),氨基糖和核苷酸糖代谢(amino sugar and nucleotide sugar metabolism),缬氨酸、亮氨酸和异亮氨酸降解(valine, leucine and isoleucine degradation),光合生物的固碳作用(carbon fixation in photosynthetic organisms),氨基酸的生物合成(biosynthesis of amino acids),角质、木栓和蜡的生物合成(cutin, suberine and wax biosynthesis),甘油酯代谢(glycerolipid metabolism),异喹啉生物碱生物合成(isoquinoline alkaloid biosynthesis),氰基氨基酸代谢(cyanoamino acid metabolism),维生素B6代谢(vitamin B6 metabolism),脂肪酸降解(fatty acid degradation),肌醇磷酸代谢(inositol phosphate metabolism),玉米素生物合成(zeatin biosynthesis),谷胱甘肽代谢(glutathione metabolism),泛醌和其他萜类醌生物合成(ubiquinone and other terpenoid-quinone biosynthesis),甜菜碱的生物合成(betalain biosynthesis),果糖和甘露糖代谢(fructose and mannose metabolism),乙醛酸和二羧酸代谢(glyoxylate and dicarboxylate metabolism),丙氨酸、天冬氨酸和谷氨酸代谢(alanine, aspartate and glutamate metabolism),甘油磷脂代谢(glycerophospholipid metabolism),抗坏血酸和阿糖二酸代谢(ascorbate and aldarate metabolism),甘氨酸、丝氨酸和苏氨酸代谢(glycine, serine and threonine metabolism),糖基磷脂酰肌醇(GPl)-锚生物合成[glycosylphosphatidylinositol (GPl)-anchor biosynthesis],类胡萝卜素生物合成(carotenoid biosynthesis)。生物系统(organismal systems):植物病原菌互作(plant-pathogen interaction),植物昼夜节律(circadian rhythm of plant)。Figure 3. Enrichment of KEGG metabolic pathway of differentially expressed genes2.5 小鳞茎萌发相关基因LoERF96的克隆及测序结果
根据转录组数据分析结果,筛选出1个褪黑素诱导后显著上调表达的ERF转录因子,以西伯利亚百合小鳞茎cDNA为模板,通过PCR扩增得到ERF基因的cDNA全长,1%的琼脂糖凝胶电泳显示在250~500 bp之间有1条清晰的条带(图4a),命名为LoERF96,切胶回收后测序并与转录组数据的基因序列比对一致,全长为435 bp,编码144个氨基酸(图4b)。
2.6 LoERF96编码蛋白的生物信息学分析
LoERF96蛋白相对分子质量约为16.32 ku,理论等电点为4.94,分子式为C692H
1078 N218O232S5,共含有2225 个原子。该蛋白由20种不同的氨基酸组成,以精氨酸的含量最高,占总量的12.5%;半胱氨酸、组氨酸和赖氨酸的含量最低,分别占总量的1.4%。LoERF96蛋白有27个带负电荷的残基(Asp+Glu)和20个带正电荷的残基(Arg+Lys),不稳定系数为70.95,表明该蛋白属于不稳定蛋白;脂肪系数为53.61;亲水性系数为−1.026,表明该蛋白为亲水性蛋白。预测LoERF96蛋白不存在跨膜区,不属于跨膜蛋白(图4c)。LoERF96基因结构域分析结果(图4d)显示:其在氨基酸40和103之间有1个AP2保守结构域,是典型的AP2/ERF家族成员,属于ERF亚家族,其亚细胞定位于细胞核。由LoERF96蛋白的二级结构(图4e)可知:该蛋白编码的144个氨基酸中有4种二级结构,其中,α-螺旋占31.94%,β-折叠占4.17%,无规则卷曲占52.08%,延伸链占11.81%,蛋白质以无规则卷曲为主,说明其结构不稳定。以含有1个AP2结构域的A0A0D6RAM1_ARACU (A0A0D6RAM1.1.A)基因为模型预测LoERF96蛋白质的三级结构,其覆盖率为56.00%,序列相似性为70.37%,以无规则卷曲和α-螺旋为主(图4f),符合二级结构预测结果。
2.7 LoERF96蛋白的同源性与系统进化树分析
进化树分析结果(图5)表明:LoERF96氨基酸序列与番茄(Solanum lycopersicum A0A3Q7I5Y9.1)、拟南芥(Arabidopsis thaliana Q8VY90.1)、烟草(Nicotiana tabacum Q40478.1)、蒺藜苜蓿(Medicago truncatula A9Y0K9.1)、喜马拉雅凤仙花(Impatiens glandulifera XP_047312446.1)等物种的同源性较高。
2.8 LoERF96在小鳞茎不同生长发育时期的表达
由图6可知:褪黑素处理后,LoERF96的相对表达量总体呈先上升后下降的趋势。第7天时,该基因的表达量比对照组高2.08倍;第14天时,其表达量较对照组极显著高10.06倍,第21天时,其表达量比对照组高2.94倍。可见,LoERF96基因受褪黑素诱导表达,且可以快速响应并表达。
2.9 RT-qPCR验证转录组数据
6个差异表达基因的RT-qPCR试验结果如图7所示。RT-qPCR结果和第21天转录组数据(表3)的基因表达趋势完全一致,表明转录组测序获得的数据比较准确。
表 3 转录组数据的基因表达量Table 3. Gene expression levels of transcriptome data基因名称
gene name基因号
gene IDFPKM P值
P-value基因表达量差异倍数
log2FCCK1 CK2 CK3 MT1 MT2 MT3 LoERF105 TRINITY_DN1258_c0_g1 24.52 44.71 20.79 93.13 77.71 43.91 0.0026 1.4920 CDC5 TRINITY_DN10098_c0_g1 113.32 92.5 110.89 43.36 49.1 40.51 0.0006 − 1.0428 MYB36 TRINITY_DN11434_c0_g2 6.45 3.56 5.58 1.69 1.81 1.17 0.0029 − 1.4653 WRKY28 TRINITY_DN11579_c0_g1 33.2 8.76 35.61 4.39 10.35 6.76 0.0037 − 1.5668 WOX8 TRINITY_DN6545_c0_g1 33.92 16.71 29.25 7.78 6.44 14.02 0.0003 − 1.3179 LOB30 TRINITY_DN10132_c0_g1 7.19 5.81 13.36 2.37 3.67 3.67 0.0064 − 1.2463 注:FPKM. 每千个碱基的转录每百万映射读取的碎片数,用于表示基因表达量;FC. 差异倍数;CK1、CK2和CK3为蒸馏水对照组的3个重复,MT1、MT2和MT3为褪黑素处理组的3个重复。
Note: FPKM. fragments per kilobase of exon model per million mapped fragments, indicating gene expression levels; FC. fold change; CK1, CK2, and CK3 are three repeats of the distilled water control group; MT1, MT2, and MT3 are three repeats of the melatonin treated group.3. 讨论
本研究对外源褪黑素处理下西伯利亚百合小鳞茎生长发育时期的转录组数据进行分析,KEGG富集分析结果显示:差异表达基因显著富集在包括植物激素信号转导、植物MAPK信号途径、核糖体、苯丙烷生物合成、淀粉和蔗糖代谢、光合作用、碳代谢、内质网中的蛋白质加工、植物病原菌互作等途径。植物激素代谢和信号转导途径中的不同激素可以单独或相互调控发挥作用,调控植物生长发育、开花结果、衰老等过程[23]。褪黑素与大多数植物激素的变化有关,它在植物生命周期的每个发育阶段起关键作用,大量研究表明:植物内源激素水平与休眠解除和萌发密切相关,但其在百合萌发过程中的调控机制尚不明确。
通过筛选褪黑素处理下的西伯利亚百合转录组数据,从小鳞茎中成功克隆出可能与萌发有关的乙烯响应因子LoERF96基因, 其cDNA全长为435 bp,编码144个氨基酸。对其蛋白理化性质进行分析,推测该蛋白属于不稳定的亲水性蛋白,且无跨膜结构;其含有1个AP2保守结构域, 是典型的AP2/ERF家族成员, 这与杧果(Mangifera indica L.)的MiERF2基因也只存在单一AP2结构域[24]的研究结果一致。蛋白二级结构预测结果显示:LoERF96蛋白以无规则卷曲和α-螺旋为主,β-折叠和延伸链占比最少,无规则卷曲为主说明其蛋白质结构不稳定,这与苦荞[Fagopyrum tataricum (L.) Gaertn][25]FtERF基因以及木薯(Manihot esculenta Crantz)[26]MeERF5基因的编码蛋白结构类似,也符合LoERF96蛋白的三级结构预测结果。亚细胞定位预测LoERF96基因位于细胞核,氨基酸多序列同源比对发现LoERF96蛋白与番茄、拟南芥、烟草、蒺藜苜蓿、喜马拉雅凤仙花等物种的同源性较高。相关研究表明:小乙烯响应因子ERF96参与拟南芥ABA反应的调控,是ABA反应的正调节因子[27]。因此,推测LoERF96在百合中也具有类似的功能。
褪黑素主要参与植物发芽、生长、根茎生成、促进果实成熟、延缓衰老等过程,且褪黑素与大多数植物激素的变化有关[28]。外源高浓度的褪黑素可以抑制或延迟种子发芽,表明其在维持种子休眠和防止过早发芽方面发挥作用[29]。褪黑素可以调控植物信号转导通路中ABA和GA相关基因的表达,诱导种子萌发并解除休眠[30]。AP2/ERF转录因子参与植物生长发育、胁迫响应、防御反应、激素信号转导、代谢物调控等,还可以反馈调节乙烯、细胞分裂素、赤霉素、脱落酸等植物激素的生物合成[31]。杏仁(Prunus dulcis)中的PdERF24可能在休眠杏仁对超低温冷冻胁迫的响应中发挥重要作用;在蓝莓(Vaccinium spp.)花芽休眠解除过程中,多数VcERF基因均能表达[32]。本研究中,在褪黑素处理下,4个不同生长发育时期的西伯利亚百合小鳞茎的LoERF96基因表达量均高于对照组,且都高于第1天的基因表达量,说明LoERF96在褪黑素处理后可以快速响应并表达。该基因的相对表达量总体呈先上升后下降的趋势,第14天时达到峰值,且极显著高于对照组;虽然之后其基因表达量下降,但仍高于对照组,这与前人的研究结果[33]一致。本研究表明:乙烯响应因子LoERF96受褪黑素诱导表达,但该基因在萌发、生长发育等过程中的生物学功能还需进一步研究。
4. 结论
本研究对褪黑素处理下的西伯利亚百合小鳞茎进行转录组测序,通过转录组数据筛选得到780个差异表达基因,其中,上调表达基因392个,下调表达基因388个。这些差异表达基因的KEGG代谢途径主要包括植物激素信号转导、植物MAPK信号途径、苯丙烷生物合成、淀粉和蔗糖代谢、光合作用、碳代谢、内质网中的蛋白质加工等。LoERF96基因是典型的AP2/ERF家族成员,其编码蛋白为亲水性且无跨膜结构。本研究证明了LoERF96基因受褪黑素诱导表达,为探究百合AP2/ERF转录因子在萌发和生长发育过程中的生物学功能提供了理论基础和参考。
-
图 1 大蒲莲和长白猪差异表达基因分布及表达水平
注:a) 差异表达基因的分布,红色和绿色点分别代表高表达和低表达基因;b) 差异表达基因的表达水平,由红到蓝表示基因表达水平由高到低。
Figure 1. Distribution and expression level of differentially expressed genes in Dapulian and Landrace pigs
Note: a) distribution of differentially expressed genes, the red and the green dots indicate highly and lowly expressed genes, respectively; b) expression of differentially expressed genes, the color red to blue indicates high to low level of gene expression.
表 1 差异表达基因极显著富集的GO条目
Table 1 GO terms with extremely significant enrichment of differentially expressed genes
表达水平
expression levelGO登录号
GO accession条目名称
term description富集基因数
gene number条目类型
term typeP值
P-value高
highGO:0030199 胶原纤维组织 collagen fibril organization 5 BP 1.37E−4 GO:0007520 成肌细胞融合 myoblast fusion 3 BP 5.83E−3 GO:0014065 磷脂酰肌醇3-激酶信号 phosphatidylinositol 3-kinase signaling 3 BP 9.88E−3 GO:0005578 蛋白质细胞外基质 proteinaceous extracellular matrix 12 CC 4.14E−7 GO:0005604 基膜 basement membrane 4 CC 3.77E−3 GO:0005615 细胞外空间 extracellular space 17 CC 8.85E−3 GO:0050840 胞外基质结合 extracellular matrix binding 4 MF 1.25E−3 GO:0005509 钙离子结合 calcium ion binding 14 MF 4.77E−3 低
lowGO:0035914 骨骼肌细胞分化 skeletal muscle cell differentiation 7 BP 1.29E−6 GO:0006457 蛋白折叠 protein folding 8 BP 2.69E−4 GO:0043524 神经元凋亡过程负调控 negative regulation of neuron apoptotic process 7 BP 3.77E−4 GO:0043065 凋亡过程的正调控 positive regulation of apoptotic process 7 BP 2.57E−3 GO:0006351 转录,DNA模板化 transcription, DNA-templated 12 BP 2.65E−3 GO:0090084 包涵体组装负调节 negative regulation of inclusion body assembly 3 BP 2.85E−3 GO:0009408 热反应 response to heat 4 BP 4.46E−3 GO:0051091 序列特异性DNA结合转录因子活性正调控
positive regulation of sequence-specific DNA binding transcription
factor activity5 BP 5.65E−3 GO:0043154 参与凋亡过程的半胱氨酸型内肽酶活性负调控
negative regulation of cysteine-type endopeptidase activity involved
in apoptotic process4 BP 7.64E−3 GO:0042692 肌肉细胞分化 muscle cell differentiation 3 BP 9.96E−3 GO:0097193 内在凋亡信号通路 intrinsic apoptotic signaling pathway 3 BP 9.96E−3 GO:0001894 组织内稳态 tissue homeostasis 3 BP 9.96E−3 GO:0005634 细胞核 nucleus 50 CC 1.41E−5 GO:0005829 细胞溶质 cytosol 25 CC 1.71E−5 GO:0005737 细胞质 cytoplasm 43 CC 3.01E−3 GO:0031674 I带 I band 3 CC 5.64E−3 GO:0031012 细胞外基质 extracellular matrix 6 CC 6.19E−3 GO:0019888 蛋白磷酸酶调节活性 protein phosphatase regulator activity 3 MF 2.38E−3 GO:0043565 序列特异性DNA结合 sequence-specific DNA binding 10 MF 8.63E−3 注:BP、CC和MF分别表示生物学过程、细胞组分和分子功能。
Note: BP, CC and MF indicates biological process, cellular component and molecular function, respectively.表 2 差异表达基因极显著富集的KEGG通路
Table 2 KEGG pathways with extremely significant enrichment of differentially expressed genes
表达水平
expression level通路ID
pathway ID通路名称
pathway term富集基因数
gene numberP 值
P-value高
highssc04974 蛋白质消化吸收 protein digestion and absorption 6 1.02E−3 ssc04960 醛固酮调节的钠重吸收
aldosterone-regulated sodium reabsorption4 7.18E−3 ssc04261 心肌细胞的肾上腺素能信号传导
adrenergic signaling in cardiomyocytes6 8.37E−3 ssc00910 氮代谢 nitrogen metabolism 3 9.00E−3 ssc04068 FoxO信号通路 FoxO signaling pathway 6 9.20E−3 低
lowssc04010 MARK信号通路 MAPK signaling pathway 10 2.67E−3 ssc04141 内质网蛋白加工
protein processing in endoplasmic reticulum8 3.06E−3 ssc04915 雌激素信号通路 estrogen signaling pathway 6 6.03E−3 表 3 红色模块基因极显著富集的前5个GO条目和KEGG通路
Table 3 Top 5 GO terms and KEGG pathways extremely significantly enriched by the genes in the red module
富集方式
enriched methodGO登录号/通路ID
GO accession/pathway ID条目名称
term description富集基因数
enriched gene numberP 值
P-valueGO条目
GO termGO:0006457 蛋白质折叠 protein folding 19 2.72E−8 GO:0035914 骨骼肌细胞分化 skeletal muscle cell differentiation 8 7.19E−5 GO:0043086 催化活性负调控 negative regulation of catalytic activity 7 1.24E−4 GO:0005977 糖原代谢过程 glycogen metabolic process 6 2.70E−4 GO:0045444 脂肪细胞分化 fat cell differentiation 8 1.90E−3 KEGG通路
KEGG pathwayssc05164 甲型流感 influenza A 23 2.07E−6 ssc01130 抗生素生物合成 biosynthesis of antibiotics 24 5.76E−6 ssc04915 雌激素信号通路 estrogen signaling pathway 14 2.00E−4 ssc04141 内质网蛋白加工 protein processing in endoplasmic reticulum 17 1.18E−3 ssc04010 MARK信号通路 MAPK signaling pathway 22 1.65E−3 表 4 蓝色模块基因极显著富集的前5个GO条目和KEGG通路
Table 4 Top 5 GO terms and KEGG pathways extremely significantly enriched by the genes in the blue module
富集方式
enriched methodGO登录号/通路ID
GO accession/pathway ID条目名称
term description富集基因数
enriched gene numberP 值
P-valueGO条目
GO termGO:0055010 心室心肌组织形态发生 ventricular cardiac muscle tissue morphogenesis 6 1.18E−4 GO:0060048 心肌收缩 cardiac muscle contraction 6 2.20E−4 GO:0014883 快慢肌纤维间转换 transition between fast and slow fiber 4 3.31E−4 GO:0043588 皮肤发育 skin development 6 4.82E−4 GO:0007160 细胞基质黏附 cell-matrix adhesion 8 6.75E−4 KEGG通路
KEGG pathwayssc04974 蛋白质消化吸收 protein digestion and absorption 14 9.75E−8 ssc04512 ECM受体互作 ECM-receptor interaction 13 8.71E−7 ssc04510 黏着斑 focal adhesion 18 8.29E−6 ssc04151 PI3K-Akt信号通路 PI3K-Akt signaling pathway 22 8.59E−5 ssc04260 心肌收缩 cardiac muscle contraction 9 5.27E−4 表 5 红色模块基因PPI子网络基因中心性分析
Table 5 Centrality analysis of genes in the PPI sub-network of the red module
序号
rank基因
gene子图
subgraph特征向量
eigenvector信息
information本地平均连通性
local average connectivity介性
betweenness近性
closeness网络
network1 FBXL4 2002.74 0.3015 6.0499 9 0 1 10 2 RNF115 2002.74 0.3015 6.0499 9 0 1 10 3 FBXO40 2002.74 0.3015 6.0499 9 0 1 10 4 ASB4 2002.74 0.3015 6.0499 9 0 1 10 5 ASB15 2002.74 0.3015 6.0499 9 0 1 10 6 UBE2L6 2002.74 0.3015 6.0499 9 0 1 10 7 RNF41 2002.74 0.3015 6.0499 9 0 1 10 8 ASB5 2002.74 0.3015 6.0499 9 0 1 10 9 HERC5 2002.74 0.3015 6.0499 9 0 1 10 10 SOCS1 2002.74 0.3015 6.0499 9 0 1 10 11 ASB11 2002.74 0.3015 6.0499 9 0 1 10 表 6 蓝色模块基因PPI子网络基因中心性分析
Table 6 Centrality analysis of genes in the PPI sub-network of the blue module
序号
rank基因
gene子图
subgraph特征向量
eigenvector信息
information本地平均连通性
local average connectivity介性
betweenness近性
closeness网络
network1 COL3A1 5196.21 0.29290 6.2447 9.17 220 0.700 11 2 COL5A1 5015.27 0.28785 6.0254 10.00 0 0.525 11 3 COL1A2 5015.27 0.28785 6.0254 10.00 0 0.525 11 4 COL14A1 5015.27 0.28785 6.0254 10.00 0 0.525 11 5 CRTAP 5015.27 0.28785 6.0254 10.00 0 0.525 11 6 COL5A3 5015.27 0.28785 6.0254 10.00 0 0.525 11 7 LEPREL1 5015.27 0.28785 6.0254 10.00 0 0.525 11 8 DEPTOR 5015.27 0.28785 6.0254 10.00 0 0.525 11 9 COL4A2 5015.27 0.28785 6.0254 10.00 0 0.525 11 10 COL6A2 5015.27 0.28785 6.0254 10.00 0 0.525 11 11 COL6A3 5015.26 0.28785 6.0254 10.00 0 0.525 11 12 COL6A1 5015.26 0.28785 6.0254 10.00 0 0.525 11 13 FBN1 907.06 0.03820 5.7853 7.20 216 0.656 9 14 FN1 819.99 0.01181 5.5213 8.00 0 0.477 9 15 CST3 819.99 0.01181 5.5213 8.00 0 0.477 9 16 APOE 819.99 0.01181 5.5213 8.00 0 0.477 9 17 MXRA8 819.99 0.01181 5.5213 8.00 0 0.477 9 18 RCN1 819.99 0.01181 5.5213 8.00 0 0.477 9 19 FAM20C 819.99 0.01181 5.5213 8.00 0 0.477 9 20 LAMB1 819.99 0.01181 5.5213 8.00 0 0.477 9 21 BMP4 819.99 0.01181 5.5213 8.00 0 0.477 9 22 GPC3 819.99 0.01181 5.5213 8.00 0 0.477 9 -
[1] DAMON M, WYSZYNSKA-KOKO J, VINCENT A, et al. Comparison of muscle transcriptome between pigs with divergent meat quality phenotypes identifies genes related to muscle metabolism and structure[J]. PLoS One, 2012, 7(3): e33763. DOI: 10.1371/journal.pone.0033763.
[2] ZAPPATERRA M, GIOIOSA S, CHILLEMI G, et al. Muscle transcriptome analysis identifies genes involved in ciliogenesis and the molecular cascade associated with intramuscular fat content in Large White heavy pigs[J]. PLoS One, 2020, 15(5): e0233372. DOI: 10.1371/journal.pone.0233372.
[3] YU K F, SHU G, YUAN F F, et al. Fatty acid and transcriptome profiling of longissimus dorsi muscles between pig breeds differing in meat quality[J]. International Journal of Biological Sciences, 2013, 9(1): 108. DOI: 10.7150/ijbs.5306.
[4] YAO C G, PANG D X, LU C, et al. Data mining and validation of AMPK pathway as a novel candidate role affecting intramuscular fat content in pigs[J]. Animals, 2019, 9(4): 137. DOI: 10.3390/ani9040137.
[5] LANGFELDER P, HORVATH S. WGCNA: an R package for weighted correlation network analysis[J]. BMC Bioinformatics, 2008, 9(1): 559. DOI: 10.1186/1471-2105-9-559.
[6] XING K, LIU H T, ZHANG F X, et al. Identification of key genes affecting porcine fat deposition based on co-expression network analysis of weighted genes[J]. Journal of Animal Science and Biotechnology, 2021, 12: 100. DOI: 10.1186/s40104-021-00616-9.
[7] GONG Y, HE J, LI B, et al. Integrated analysis of lncRNA and mRNA in subcutaneous adipose tissue of Ningxiang pig[J]. Biology, 2021, 10(8): 726. DOI: 10.3390/biology10080726.
[8] KOGELMAN L, CIRERA S, ZHERNAKOVA D V, et al. Identification of co-expression gene networks, regulatory genes and pathways for obesity based on adipose tissue RNA Sequencing in a porcine model[J]. BMC Medical Genomics, 2014, 7: 57. DOI: 10.1186/1755-8794-7-57.
[9] HE Y Q, MALTECCA C, TIEZZI F, et al. Transcriptome analysis identifies genes and co-expression networks underlying heat tolerance in pigs[J]. BMC Genetics, 2020, 21(1): 44. DOI: 10.1186/s12863-020-00852-4.
[10] PONSUKSILI S, DU Y, HADLICH F, et al. Correlated mRNAs and miRNAs from co-expression and regulatory networks affect porcine muscle and finally meat properties[J]. BMC Genomics, 2013, 14: 533. DOI: 10.1186/1471-2164-14-533.
[11] PONSUKSILI S, SIENGDEE P, DU Y, et al. Identification of common regulators of genes in co-expression networks affecting muscle and meat properties[J]. PLoS One, 2015, 10(4): e0123678. DOI: 10.1371/journal.pone.0123678.
[12] CHOTIMA P, CHRISTELLE H A, KANJANA T, et al. Revealing pathways associated with feed efficiency and meat quality traits in slow-growing chickens[J]. Animals, 2021, 11(10): 2977. DOI: 10.3390/ani11102977.
[13] DINIZ W, MAZZONI G, COUTINHO L L, et al. Detection of co-expressed pathway modules associated with mineral concentration and meat quality in Nelore cattle[J]. Frontiers in Genetics, 2019, 10: 210. DOI: 10.3389/fgene.2019.00210.
[14] ZHAO X, WANG C, WANG Y, et al. Weighted gene co-expression network analysis reveals potential candidate genes affecting drip loss in pork[J]. Animal Genetics, 2020, 51(6): 855. DOI: 10.1111/age.13006.
[15] RITCHIE M E, BELINDA P, WU D, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies[J]. Nucleic Acids Research, 2015, 43(7): e47. DOI: 10.1093/nar/gkv007.
[16] WEI H D, SHERMAN B T, LEMPICKI R A. Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists[J]. Nucleic Acids Research, 2009, 37(1): 1. DOI: 10.1093/nar/gkn923.
[17] DAMIAN S, ANDREA F, STEFAN W, et al. STRING v10: protein-protein interaction networks, integrated over the tree of life[J]. Nucleic Acids Research, 2015, 43: 447. DOI: 10.1093/nar/gku1003.
[18] PAUL S, ANDREW M, OWEN O, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks[J]. Genome Research, 2003, 13(11): 2498. DOI: 10.1101/gr.1239303.
[19] BADER G D, HOGUE C W. An automated method for finding molecular complexes in large protein interaction networks[J]. BMC Bioinformatics, 2003, 4: 2. DOI: 10.1186/1471-2105-4-2.
[20] TANG Y, LI M, WANG J X, et al. CytoNCA: a cytoscape plugin for centrality analysis and evaluation of protein interaction networks[J]. Biosystems, 2015, 127: 67. DOI: 10.1016/j.biosystems.2014.11.005.
[21] HAUSMAN G J. Meat science and muscle biology symposium: the influence of extracellular matrix on intramuscular and extramuscular adipogenesis[J]. Journal of Animal Science, 2012, 90(3): 942. DOI: 10.2527/jas.2011-4616.
[22] LUCARINI M, DURAZZO A, SCIUBBA F, et al. Stability of the meat protein type I collagen: influence of pH, ionic strength, and phenolic antioxidant[J]. Foods, 2020, 9(4): 480. DOI: 10.3390/foods9040480.
[23] LEPETIT J. A theoretical approach of the relationships between collagen content, collagen cross-links and meat tenderness[J]. Meat Science, 2007, 76(1): 147. DOI: 10.1016/j.meatsci.2006.10.027.
[24] MCBRYAN J, HAMILL R M, DAVEY G, et al. Identification of suitable reference genes for gene expression analysis of pork meat quality and analysis of candidate genes associated with the trait drip loss[J]. Meat Science, 2010, 86(2): 436. DOI: 10.1016/j.meatsci.2010.05.030.
[25] NAKAJIMA I, MUROYA S, TANABE R, et al. Positive effect of collagen V and VI on triglyceride accumulation during differentiation in cultures of bovine intramuscular adipocytes[J]. Differentiation, 2002, 70(2/3): 84. DOI: 10.1046/j.1432-0436.2002.700203.x.
[26] 曾勇庆, 王根林, 魏述东, 等. 莱芜猪肌肉胶原蛋白的发育性变化及其与肉质的相关性[J]. 中国农业科学, 2008, 41(2): 619. DOI: 10.3864/j.issn.0578-1752.2008.02.043. [27] LIM K S, LEE K T, PARK J E, et al. Identification of differentially expressed genes in longissimus muscle of pigs with high and low intramuscular fat content using RNA sequencing[J]. Animal Genetics, 2017, 48(2): 166. DOI: 10.1111/age.12518.
[28] BAO X J, ZENG Y Q, WEI S D, et al. Developmental changes of Col3a1 mRNA expression in muscle and their association with intramuscular collagen in pigs[J]. Journal of Genetics and Genomics, 2007, 34(3): 223. DOI: 10.1016/S1673-8527(07)60023-X.
[29] FERNÁNDEZ-BARROSO M, CARABALLO C, SILIÓ L, et al. Differences in the loin tenderness of Iberian pigs explained through dissimilarities in their transcriptome expression profile[J]. Animals, 2020, 10(9): 1715. DOI: 10.3390/ani10091715.
[30] VALLI M, BARNES A M, GALLANTI A, et al. Deficiency of CRTAP in non-lethal recessive osteogenesis imperfecta reduces collagen deposition into matrix[J]. Clinical Genetics, 2011, 82(5): 453. DOI: 10.1111/j.1399-0004.2011.01794.x.
[31] ZOU Y L, LI Z Y, ZOU Y J, et al. An FBXO40 knockout generated by CRISPR/Cas9 causes muscle hypertrophy in pigs without detectable pathological effects[J]. Biochemical and Biophysical Research Communications, 2018, 498(4): 940. DOI: 10.1016/j.bbrc.2018.03.085.
[32] SHI J, LUO L Q, EASH J, et al. The SCF-Fbxo40 complex induces IRS1 ubiquitination in skeletal muscle, limiting IGF1 signaling[J]. Developmental Cell, 2011, 21(5): 835. DOI: 10.1016/j.devcel.2011.09.011.
[33] ALSINA D, LYTOVCHENKO O, SCHAB A, et al. FBXL4 deficiency increases mitochondrial removal by autophagy[J]. EMBO Molecular Medicine, 2020, 12(7): e11659. DOI: 10.15252/emmm.201911659.
[34] EHRLICH K C, LACEY M, EHRLICH M. Epigenetics of skeletal muscle-associated genes in the ASB, LRRC, TMEM, and OSBPL gene families[J]. Epigenomes, 2020, 4(1): 1. DOI: 10.3390/epigenomes4010001.
[35] VERBRUGGE S A J, SCHÖNFELDER M, BECKER L, et al. Genes whose gain or loss-of-function increases skeletal muscle mass in mice: a systematic literature review[J]. Frontiers in Physiology, 2018, 9: 553. DOI: 10.3389/fphys.2018.00553.
[36] WEI W P, LI Y Q, LI Y Y, et al. Adipose-specific knockout of ubiquitin-conjugating enzyme E2L6 (Ube2l6) reduces diet-induced obesity, insulin resistance, and hepatic steatosis[J]. Journal of Pharmacological Sciences, 2021, 145(4): 327. DOI: 10.1016/j.jphs.2020.12.008.
[37] MARCELIN G, LIU S M, SCHWARTZ G J, et al. Identification of a loss-of-function mutation in Ube2l6 associated with obesity resistance[J]. Diabetes, 2013, 62(8): 2784. DOI: 10.2337/db12-1054.
[38] 连林生, 鲁绍雄, 严达伟, 等. 撒坝猪配套系选育初报[J]. 猪业科学, 2006, 6: 58. DOI: 10.3969/j.issn.1673-5358.2006.06.018. [39] 郭建凤, 吕山东, 王继英, 等. 大蒲莲猪与长蒲二元猪生长育肥性能、胴体性能及肉品质比较[J]. 畜牧与兽医, 2010, 42(11): 24.