仔猪生长发育关键基因筛选和功能分析

宋 阳 ，芦春雪，殷淑洁 ，冯 静 ，刘 欢 ，苏 辰 ，周百灵，沈丛丛，陈清帅，黄平平*，于 雪*

（1.德州学院生物物理研究院/山东省生物物理重点实验室，山东 德州 253023；2.德州学院生命科学学院，山东 德州）

断奶仔猪是指处在从依靠母猪乳汁到独立依靠固体饲料的仔猪[1]。断奶过渡期是仔猪生长过程中最重要的阶段之一，该阶段仔猪生长发育快，但消化系统尚未发育完善，因此容易出现采食量差、生长性能下降和胃肠功能受损等严重后果[2]。过去研究者们主要通过对营养、健康和管理等有关技术进行改进，尽量减少断奶应激对仔猪生长的不利影响[3]。另外，骨骼肌和脂肪组织不仅是重要的营养部位，更是生长的关键组成部分，且能够对各种胁迫进行生理适应，包括疾病和压力等[4]。因此，揭示在仔猪背最长肌和皮下脂肪组织生长发育过程起关键作用的基因，将有助于改善家畜肉质和产量，同时有利于治疗和预防相关疾病。

近年来随着高通量测序技术的发展，与仔猪生长发育相关基因的挖掘越来越受到国内外研究者们的关注[5]。如Pilcher等[5]和Chen等[6]分别基于猪Ensembl release 70和Ensembl release 93版本的基因组和注释文件分析背最长肌和皮下脂肪组织发育相关的关键基因。根据 Ensembl官方网站，截止 2021年底猪基因组和注释文件已经更新至Ensembl release 104版本。因此本研究在上述工作基础上从 GEO数据库下载背最长肌和皮下脂肪组织各18个转录组测序数据，采用Ensembl release 104版本的猪基因组文件和注释文件进行生物信息学分析，进一步挖掘与这两种组织发育相关的关键基因，对补充认识仔猪生长发育相关的分子机制具有重要意义。

1 材料与方法

1.1 材料

本实验原始数据来自美国国家生物技术中心（National Center for Biotechnology Information，NCBI）的基因表达数据库（Gene Expression Omnibus, GEO），共36个样本，ID为GSE 65983。其中，18个样本来自猪的背最长肌（LD M），18个样本来自皮下脂肪（BF）组织。猪的基因组和注释文件自 Ensembl数据库下载（http://ftp.ensembl.org/pub/release-104/fasta/sus_scrofa/；http://ftp.ensembl.org/pub/release-104/gtf/sus_scrofa/）。

1.2 方法

1.2.1 原始数据质控 为了获得高质量的 clean data以便增强后续数据分析的准确性，需要对获得的原始数据质控。采用 Trimmomatic软件将原始数据中包含的测序接头、低质量 reads以及 N（无法确定碱基信息的比例大于5 % 的reads）、较短序列进行过滤去除，最终获得高质量 clean data。

1.2.2 参考序列比对 目前研究表明，Hisat2（http://daehwankimlab.github.io/hisat2/）相较于TopHat2等对比软件具有运行速度快、精确度较好等优点。因此本研究使用Hisat2软件进行序列比对，为提高其对比速度，实验前需构建其基因的索引（index），具体操作如下：（1）下载猪（sus_scrofa）基因组文件。（2）索引构建（hisat 2-build–p 16 genome.fa genome）。（3）进行 Hisat2比对。比对之后利用 samtools软件进行数据格式转换，将sam文件转为排好序的bam文件。

1.2.3 转录本拼接、定量及差异表达分析 分别采用 Stringtie软件组装转录本、Samtools软件对所有转录本定量，转录本的 TPM 值至少在一个样本中大于 1被认为是表达的。将 36个样本根据组织部位（背最长肌和皮下脂肪组织）分为 2组，然后采用 DEseq2软件对所有转录本进行差异表达分析。满足以下条件的转录本认为是差异表达的：在两组样本中∣log2FC∣＞1且 padj.＜0.05。

1.2.4 功能富集分析 为筛选与背最长肌和皮下脂肪组织发育相关的关键基因，本研究利用DAVID对获得的差异表达的 mRNA进行基因本体论（Gene Oncology, GO）功能注释和京都基因与基因组百科全书（Kyoto Encyclopedia of Genes and Genomes, KEGG）通路富集分析。

2 结果与分析

2.1 转录本数据统计

转录本的 TPM 值至少在一个样本中大于 1认为是表达的，36个样本共获得 79 800个转录本，其中包括 29 118个 mRNA（11 858种基因），2 105个已知的lncRNA（2种基因，仅6个转录本含基因名），以及 47 142个新的转录本（表1）。

表1 转录本总体情况

2.2 差异表达分析

共获得 34 873个差异表达转录本，其中9 671个转录本在背最长肌组织上调；25 202个转录本在背最长肌组织下调。另外包括12 536个差异表达的mRNA，3 352个mRNA在背最长肌组织上调，9 184个 mRNA在背最长肌组织下调（表2）。

表2 转录本差异表达情况

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

将获得的在两种组织中差异表达的mRNA进行GO富集分析，结果分别获得168个BP（biological process, BP）通路、139个 CC（cellular component, CC）通路和 122个 MF（molecular function，MF）通路（P-value < 0.05），图1a、b、c中分别展示位于前 20位的 BP、CC和 MF通路。在 BP通路中，与肌肉和脂肪组织相关的通路包括肌节组织（sarcomere organization）、肌动蛋白细胞骨架组织（actin cytoskeleton organization）和细胞迁移的积极调解（positive regulation of cell migration）等。61个mRNA转录本富集在肌节组织，如基因MYPN和TNNT1；124个mRNA转录本富集在肌动蛋白细胞骨架组织，如基因ARHGAP26和PDLIM3；150个mRNA转录本富集在细胞迁移的积极调解，如基因PECAM1。在 CC通路中，与肌肉和脂肪组织相关的通路包括丝状肌动蛋白（filamentous actin）、肌动蛋白丝（actin filament）和脂质颗粒（lipid particle）等。38个 mRNA转录本富集在丝状肌动蛋白，如基因PDLIM3；82个 mRNA转录本富集在肌动蛋白丝，如基因TPM4；83个mRNA转录本富集在脂质颗粒，如基因PLIN1。在 MF通路中，与肌肉和脂肪组织相关的通路包括磷脂酰肌醇结合（phosphatidylinositol binding）、脂质结合（lipid binding）和对旋肌球蛋白结合（tropomyosin binding）等。150个mRNA转录本富集在磷脂酰肌醇结合，94个 mRNA转录本富集在脂质结合，如基因APOA2；31个mRNA转录本富集在对旋肌球蛋白结合，如基因TNNT2和TMOD4。

图1 差异表达mRNA功能富集分析结果

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

通过KEGG信号通路分析共获得187个通路（P-value < 0.05），图1d仅展示位居前20位的通路。其中肌动蛋白细胞骨架调节（Regulation of actin cytoskeleton）与肌肉组织相关。203个mRNA转录本富集在肌动蛋白细胞骨架调节，如基因CFL和FGFR4。另外，脂肪细胞因子信号通路（Adipocytokine signaling pathway）、鞘脂信号通路（Sphingolipid signaling pathway）、脂肪酸代谢（Fatty acid metabolism）等通路与脂肪组织发育相关。73个 mRNA转录本富集在脂肪细胞因子信号通路，如基因AKT2和PCK1；113个mRNA转录本富集在鞘脂信号通路，另外，63个mRNA转录本富集在脂肪酸代谢结合，如基因SCD5和FASN。

3 讨论

在猪生长发育的一系列阶段中，断奶会对猪的生长发育产生不良影响。骨骼肌和脂肪组织是猪的营养部位，在猪的发育中发挥重要的作用，同时也是仔猪生长的关键部分，骨骼肌与皮下脂肪组织能够使仔猪对各种疾病和压力进行生理适应。皮下脂肪组织主要与能量代谢途径有关，在猪不同的生长阶段，皮下脂肪组织的生长状况亦不同[7]。近年来与猪肌肉和脂肪组织发育相关的基因及其功能逐渐被揭示，如Wang等（2015）[8]对不同品种猪的背最长肌组织进行高通量测序和生物信息学分析，发现CAV2、MYOZ2和FRZB与肌肉生长紧密相关，而FASN、SCD和ADORA1与脂质沉积紧密相关。Li等（2016）应用 RNA测序鉴定皖南花猪与约克夏猪背最长肌差异表达基因并进行功能富集分析，发现差异表达的基因主要富集在肌肉发育的生物学过程和脂肪酸代谢等通路[9]。但断奶仔猪背最长肌和皮下脂肪组织发育相关的基因及其功能的研究仍处于初级阶段。

为了深入研究在仔猪背最长肌和皮下脂肪组织生长发育过程起关键作用的基因，本研究采用数据库数据进行生物信息学分析。本研究详细介绍了可能影响仔猪生长性能的基因及通路，为仔猪生长发育相关分子生物学的研究提供理论支撑。在未来研究中，将进一步结合分子生物学和细胞生物学实验，揭示关键通路中基因的功能以及发挥功能的分子机制，为减少仔猪断奶应激、改善生长发育提供理论依据，对生猪产业具有重要意义。