仔豬生長發育關鍵基因篩選和功能分析

宋 陽 ，蘆春雪，殷淑潔 ，馮 靜 ，劉 歡 ，蘇 辰 ，周百靈，沈叢叢，陳清帥，黃平平*，于 雪*

（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]。但斷奶仔豬背最長肌和皮下脂肪組織發育相關的基因及其功能的研究仍處于初級階段。

為了深入研究在仔豬背最長肌和皮下脂肪組織生長發育過程起關鍵作用的基因，本研究采用數據庫數據進行生物信息學分析。本研究詳細介紹了可能影響仔豬生長性能的基因及通路，為仔豬生長發育相關分子生物學的研究提供理論支撐。在未來研究中，將進一步結合分子生物學和細胞生物學實驗，揭示關鍵通路中基因的功能以及發揮功能的分子機制，為減少仔豬斷奶應激、改善生長發育提供理論依據，對生豬產業具有重要意義。