基于轉(zhuǎn)錄組測(cè)序挖掘廣西麻雞生長(zhǎng)發(fā)育相關(guān)基因

陳 婷，盧晟盛，楊麗麗，李添寶，張 敏，江明生

(1.廣西大學(xué)動(dòng)物科學(xué)技術(shù)學(xué)院，南寧 530004；2.廣西大學(xué)農(nóng)牧產(chǎn)業(yè)發(fā)展研究院，南寧 530004)

廣西麻雞具備優(yōu)質(zhì)肉雞的良好產(chǎn)肉性能、風(fēng)味佳、營(yíng)養(yǎng)高，近年來(lái)在市場(chǎng)占有份額越來(lái)越大，母雞的營(yíng)養(yǎng)優(yōu)于公雞[1]，胸肌和腿肌具有豐富的氨基酸，且腿肌比胸肌風(fēng)味更佳[2]。腿肌是雞主要的產(chǎn)肉部位之一，其在肌纖維生長(zhǎng)發(fā)育[3]、屠宰性能和肉品質(zhì)[4]、風(fēng)味[5]等方面有相關(guān)報(bào)道。廣西大學(xué)動(dòng)物科技學(xué)院動(dòng)物遺傳育種實(shí)驗(yàn)室在前期生產(chǎn)試驗(yàn)中發(fā)現(xiàn)，廣西麻雞母雞生長(zhǎng)速度慢、產(chǎn)蛋較少、母雞就巢性強(qiáng)、繁殖性能差等問(wèn)題，其生長(zhǎng)高峰期在60日齡左右。因此，加強(qiáng)對(duì)廣西麻雞的遺傳育種研究具有重要的意義。

雞的生長(zhǎng)發(fā)育是一個(gè)重要的經(jīng)濟(jì)性狀，在家禽生產(chǎn)中，隨著家禽選育程度的不斷提高，對(duì)雞生長(zhǎng)發(fā)育的調(diào)控機(jī)制應(yīng)該得到重視。生長(zhǎng)發(fā)育受遺傳和環(huán)境等多種因素的影響，在一系列細(xì)胞因子和轉(zhuǎn)錄因子的調(diào)控下按照嚴(yán)格的進(jìn)程完成，出生后肌肉的形成對(duì)生長(zhǎng)起決定作用。通過(guò)轉(zhuǎn)錄組學(xué)研究方法已有學(xué)者對(duì)雞的肉質(zhì)[6]、脂肪沉積[7]、肌內(nèi)脂肪代謝[8]、精子活力[9]、長(zhǎng)鏈脂肪酸變化[10]、繁殖性狀[11]、肝臟抗氧化酶活性[12]、多不飽和脂肪酸[13]等相關(guān)差異表達(dá)基因進(jìn)行研究以加速品種的選育和改良。在雞生長(zhǎng)發(fā)育的研究中，腿肌轉(zhuǎn)錄組測(cè)序的研究相對(duì)較少。本研究以廣西麻雞母雞出生時(shí)(1日齡)和生長(zhǎng)高峰期(60日齡)腿肌為研究對(duì)象，對(duì)2個(gè)不同生長(zhǎng)階段的腿肌組織進(jìn)行轉(zhuǎn)錄組測(cè)序，從mRNA水平上對(duì)差異表達(dá)基因進(jìn)行研究，對(duì)肌肉生長(zhǎng)發(fā)育過(guò)程中差異表達(dá)基因進(jìn)行篩選，并進(jìn)行功能注釋以挖掘影響生長(zhǎng)發(fā)育的關(guān)鍵基因，以期為闡明廣西麻雞生長(zhǎng)發(fā)育的分子機(jī)制提供參考。

1 材料與方法

1.1 樣品采集

1和60日齡健康廣西麻雞母雞各3只，由廣西南寧良慶區(qū)農(nóng)利來(lái)種禽有限公司提供。屠宰后，采集腿肌組織(1日齡腿肌組織分別記為1LM-1、1LM-2、1LM-3；60日齡腿肌組織分別記為60LM-1、60LM-2、60LM-3)，迅速放入―80 ℃保存?zhèn)溆谩?/p>

1.2 RNA提取

用TRIzol法分別提取6只雞的腿肌組織總RNA，利用Bioanalyzer 2100和RNA 6000 Nano LabChip Kit分析總RNA質(zhì)量和純度，選取D260 nm/D280 nm>1.8、D260 nm/D230 nm>1.0、RNA濃度>50 ng/μL的樣品用于后續(xù)試驗(yàn)。

1.3 轉(zhuǎn)錄組測(cè)序分析

樣品檢驗(yàn)合格后，以mRNA為模板合成雙鏈cDNA，雙鏈cDNA經(jīng)過(guò)末端修復(fù)，加A尾修飾并連接測(cè)序接頭，再利用Oligo磁珠對(duì)其片段大小進(jìn)行篩選和純化，用UDG酶消化雙鏈。PCR反應(yīng)程序：95 ℃預(yù)變性3 min；98 ℃變性15 s，60 ℃退火15 s，72 ℃延伸30 s，共8個(gè)循環(huán)；72 ℃終延伸5 min，形成300 bp±50 bp的cDNA文庫(kù)。使用Illumina NovaseqTM6000進(jìn)行轉(zhuǎn)錄組測(cè)序，利用Cutadap軟件對(duì)原始數(shù)據(jù)進(jìn)行質(zhì)量控制，利用HISAT2軟件將 clean reads比對(duì)到基因組上，使用FPKM 值作為基因表達(dá)量。

1.4 差異表達(dá)基因分析

利用DESeq2軟件進(jìn)行基因的差異表達(dá)分析，設(shè)定閾值為log2|FoldChange|>0，校正P值(Padj)<0.05，當(dāng)FoldChange≥1.5或≤0.05且P值<0.05時(shí)被定義為差異表達(dá)基因(DEGs)。

1.5 差異表達(dá)基因GO功能和KEGG通路富集分析

利用GO數(shù)據(jù)庫(kù)對(duì)篩選出的差異表達(dá)基因進(jìn)行功能注釋，利用KEGG數(shù)據(jù)庫(kù)對(duì)差異表達(dá)基因的信號(hào)通路進(jìn)行富集分析，篩選出與生長(zhǎng)發(fā)育相關(guān)的功能基因。

1.6 實(shí)時(shí)熒光定量PCR驗(yàn)證

隨機(jī)選取視黃醇結(jié)合蛋白7(RBP7)、肌鈣蛋白T2(TNNT2)、谷氨酸脫羧酶1(GADL1)、細(xì)胞表明糖蛋白(CD44)、肌肉生長(zhǎng)抑制素(GDF8)、磷酸化酶激酶調(diào)節(jié)亞基β(PHKB)、熱休克蛋白H1(HSPH1)、RCSD結(jié)構(gòu)域包含1(RCSD1)共8個(gè)差異表達(dá)基因進(jìn)行實(shí)時(shí)熒光定量PCR驗(yàn)證。根據(jù)GenBank中各基因序列，以GAPDH作為內(nèi)參基因，利用Oligo 7.0軟件設(shè)計(jì)引物，引物信息見(jiàn)表1。引物均由南寧捷尼斯生物科技有限公司合成。使用實(shí)時(shí)熒光定量PCR驗(yàn)證轉(zhuǎn)錄組測(cè)序結(jié)果。PCR反應(yīng)體系10 μL：SYBR Green Ⅰ 5 μL，上、下游引物各0.4 μL，cDNA模板1 μL，ddH2O 3.2 μL。PCR反應(yīng)程序：95 ℃預(yù)變性30 s；95 ℃變性5 s，60 ℃退火30 s，65 ℃延伸30 s，共45個(gè)循環(huán)。采用2-ΔΔCt法計(jì)算目的基因相對(duì)表達(dá)量。用SPSS 25.0軟件進(jìn)行t檢驗(yàn)，用GraphPad Prism 8.0軟件進(jìn)行作圖。

表1 引物信息Table 1 Primer information

2 結(jié) 果

2.1 轉(zhuǎn)錄組測(cè)序數(shù)據(jù)分析

通過(guò)文庫(kù)構(gòu)建和Illumina NovaseqTM6000平臺(tái)測(cè)序，原始數(shù)據(jù)經(jīng)過(guò)質(zhì)量控制后過(guò)濾后，6個(gè)樣本分別獲得34 296 198～41 647 642條clean reads，有效比對(duì)率在88.33%以上，Q20條含量在99.90%以上，Q30含量在97.51%以上，GC含量在50.00%～51.50%，說(shuō)明測(cè)序數(shù)據(jù)可靠，可用于后續(xù)研究。

表2 轉(zhuǎn)錄組測(cè)序質(zhì)控結(jié)果Table 2 Quality control results of transcriptome sequencing

2.2 差異表達(dá)基因分析

通過(guò)DESeq2軟件對(duì)1和60日齡廣西麻雞腿肌轉(zhuǎn)錄組測(cè)序數(shù)據(jù)進(jìn)行分析，結(jié)果發(fā)現(xiàn)，與1日齡相比，60日齡共獲得差異表達(dá)基因2 304個(gè)，其中998個(gè)基因上調(diào)，1 306個(gè)基因下調(diào)(圖1)。

圖1 差異表達(dá)基因火山圖Fig.1 Volcano map of differentially expressed genes

2.3 差異表達(dá)基因GO功能和KEGG通路富集分析

2.3.1 GO功能富集分析 GO功能富集分析發(fā)現(xiàn)，差異表達(dá)基因GO富集的條目分為生物過(guò)程(biological process，BP)、細(xì)胞組分(cellular component，CC)和分子功能(molecular function，MF)3類，顯著富集的GO條目有572條，其中生物過(guò)程381條、細(xì)胞組分70條、分子功能121條，生物過(guò)程占大多數(shù)。由表3可知，生物過(guò)程中差異表達(dá)基因主要參與肌動(dòng)蛋白成核的正調(diào)控、成肌細(xì)胞分化的負(fù)調(diào)控、成纖維細(xì)胞生長(zhǎng)因子刺激的細(xì)胞反應(yīng)、橫紋肌組織發(fā)育，差異表達(dá)基因的數(shù)目分別為5、7、5、3和3個(gè)；分子功能中差異表達(dá)基因主要參與肌動(dòng)蛋白結(jié)合和膠原結(jié)合，差異表達(dá)基因數(shù)目分別為45和12個(gè)；細(xì)胞組分中差異表達(dá)基因主要參與細(xì)胞質(zhì)基質(zhì)、肌鈣蛋白復(fù)合物和肌球蛋白絲，表達(dá)差異基因數(shù)目分別為28、4和3個(gè)(表3)。

表3 部分生長(zhǎng)發(fā)育相關(guān)GO條目Table 3 Partial growth and development related GO terms

續(xù)表

2.3.2 KEGG通路富集分析 KEGG通路富集分析發(fā)現(xiàn)，共獲得34個(gè)顯著富集的通路。由圖2可知，主要顯著富集的通路包括：氧化磷酸化、心肌收縮、Apelin信號(hào)通路、心肌細(xì)胞中的腎上腺素能信號(hào)轉(zhuǎn)導(dǎo)、MAPK信號(hào)通路、卡波西肉瘤相關(guān)皰疹病毒感染、脂肪酸降解、煙酸和煙酰胺代謝、晝夜夾帶、GnRH信號(hào)通路、多巴胺能突觸、胃酸分泌、糖尿病并發(fā)癥中的AGE-RAGE信號(hào)通路、血管平滑肌收縮、鞘脂代謝、糖酵解/糖異生、亨廷頓病、cAMP信號(hào)通路等。

圖2 差異表達(dá)基因KEGG通路富集分析Fig.2 KEGG pathway enrichment analysis of differentially expressed genes

2.3.3 生長(zhǎng)發(fā)育相關(guān)基因篩選 通過(guò)GO功能和KEGG通路富集分析，獲得5個(gè)生長(zhǎng)發(fā)育相關(guān)基因，分別為球蛋白重鏈10(MYH10)、肌球蛋白鏈15(MYH15)、成纖維細(xì)胞生長(zhǎng)因子10(FGF10)、成纖維細(xì)胞生長(zhǎng)因子16(FGF16)、GDF8基因(表4)，其中MYH10、MYH15和FGF10基因?yàn)樯险{(diào)表達(dá)差異基因，F(xiàn)GF16和GDF8基因?yàn)橄抡{(diào)表達(dá)差異基因。

表4 生長(zhǎng)發(fā)育相關(guān)基因篩選結(jié)果Table 4 Screening results of growth and development related genes

2.4 實(shí)時(shí)熒光定量PCR驗(yàn)證

為了驗(yàn)證轉(zhuǎn)錄組測(cè)序結(jié)果的可靠性，隨機(jī)選取8個(gè)在轉(zhuǎn)錄組測(cè)序中的差異表達(dá)基因進(jìn)行實(shí)時(shí)熒光定量PCR驗(yàn)證。由圖3可知，TNNT2、RBP7和RCSD1基因表達(dá)水平上調(diào)，GDF8、CD44、GADL1、PHKB和HSPH1基因表達(dá)水平下調(diào)，與轉(zhuǎn)錄組測(cè)序結(jié)果一致，證明轉(zhuǎn)錄組測(cè)序結(jié)果可靠。

圖3 實(shí)時(shí)熒光定量PCR驗(yàn)證轉(zhuǎn)錄組測(cè)序結(jié)果Fig.3 Validation of RNA-Seq results by Real-time quantitative PCR

3 討 論

本研究對(duì)廣西麻雞腿肌組織進(jìn)行轉(zhuǎn)錄組測(cè)序，共獲得2 304個(gè)差異表達(dá)基因，其中998個(gè)基因上調(diào)，1 306個(gè)基因下調(diào)，結(jié)合GO功能和KEGG通路富集分析對(duì)篩選出的差異表達(dá)基因進(jìn)行挖掘，獲得MYH10、MYH15、FGF10、FGF16、GDF8共5個(gè)生長(zhǎng)發(fā)育相關(guān)基因。He等[14]在生長(zhǎng)快、慢的16周齡邊雞腿肌組織的差異轉(zhuǎn)錄本中篩選出108個(gè)差異表達(dá)基因，獲得了17個(gè)上調(diào)基因和91個(gè)下調(diào)基因。Yin等[15]在海陽(yáng)黃雞12、16 d胚胎及1日齡和10周齡4個(gè)生長(zhǎng)時(shí)期腿肌中篩選出6 150個(gè)差異表達(dá)基因，肌肉抑制素(MSTN)、肌源性分化1(MYOD1)、肌球蛋白重鏈6(MYF6)、肌球蛋白重鏈5(MYF5)和胰島素生長(zhǎng)因子1(IGF1)等基因是雞生長(zhǎng)發(fā)育的關(guān)鍵基因，這與本研究篩選的基因有所不同，說(shuō)明不同生長(zhǎng)時(shí)期所受到的基因調(diào)控會(huì)有所差異，但都可以為雞生長(zhǎng)發(fā)育研究提供一定的參考。

MYH10和MYH15是肌球蛋白超家族的成員。MYH10具有多種功能，參與調(diào)節(jié)胞質(zhì)分裂、細(xì)胞運(yùn)動(dòng)和細(xì)胞極性，在細(xì)胞骨架重組中起主要作用，其在小鼠神經(jīng)系統(tǒng)和心臟的形成過(guò)程中至關(guān)重要[16]，在心臟發(fā)育和胚胎存活中起重要作用[17]。MYH15參與肌動(dòng)蛋白結(jié)合和細(xì)胞骨架運(yùn)動(dòng)活性，Zhang等[18]研究3個(gè)品種雞肌肉生長(zhǎng)發(fā)育的機(jī)制發(fā)現(xiàn)，與生長(zhǎng)較慢的品種相比，生長(zhǎng)較快的品種MYH15基因?yàn)轱@著上調(diào)基因，還有研究發(fā)現(xiàn)MYH15在胚胎期肌肉中不存在，并且在出生后肌肉中可檢測(cè)到[19]。本研究發(fā)現(xiàn)，MYH15基因?yàn)轱@著上調(diào)基因，說(shuō)明MYH10基因?qū)Τ錾箅u的生長(zhǎng)發(fā)育起關(guān)鍵調(diào)控作用。FGF10基因參與胚胎發(fā)育、細(xì)胞增殖和細(xì)胞分化，Rivetti等[20]發(fā)現(xiàn)，F(xiàn)GF10在小鼠中對(duì)誘導(dǎo)乳腺基板和白色脂肪組織的形成至關(guān)重要，F(xiàn)GF10基因?qū)Χ涞恼T導(dǎo)和生長(zhǎng)也至關(guān)重要，F(xiàn)GF10基因突變小鼠在內(nèi)耳分化過(guò)程中表現(xiàn)出輕微的缺陷[21]，并且可促進(jìn)雞靜聲神經(jīng)節(jié)神經(jīng)突生長(zhǎng)和神經(jīng)元存活[22]。FGF16基因參與胚胎發(fā)育、細(xì)胞生長(zhǎng)、形態(tài)發(fā)生，Ahsan等[23]研究表明，F(xiàn)GF16基因是雞生長(zhǎng)發(fā)育的候選基因，推測(cè)其對(duì)雞的生長(zhǎng)調(diào)控可能有重要作用，但是也有研究報(bào)道FGF9和FGF16基因可協(xié)同作為胚胎心肌細(xì)胞的生長(zhǎng)因子[24]。GDF8基因具有骨骼肌生長(zhǎng)負(fù)調(diào)節(jié)功能，有研究表明，GDF8基因可抑制肌源性衛(wèi)星細(xì)胞的活化、增殖和分化，有研究發(fā)現(xiàn)GDF8對(duì)兩種不同骨骼肌的雞成肌衛(wèi)星細(xì)胞的影響有所差異，胸大肌肌肉中分離的衛(wèi)星細(xì)胞比股二頭肌的細(xì)胞中肌肉抑制素的抑制作用更敏感[25]，并且GDF8基因在雞胚胎的各種組織中廣泛表達(dá)[26]。本研究發(fā)現(xiàn)GDF8基因?yàn)轱@著下調(diào)基因，進(jìn)一步表明GDF8基因的表達(dá)有負(fù)調(diào)控作用，可影響雞肌肉的生長(zhǎng)發(fā)育。

GO功能富集分析發(fā)現(xiàn)，差異表達(dá)基因顯著富集到了572個(gè)GO條目，其中包含了生長(zhǎng)發(fā)育過(guò)程的非常重要的條目，如肌動(dòng)蛋白成核的正調(diào)控、成肌細(xì)胞分化的負(fù)調(diào)控、成纖維細(xì)胞生長(zhǎng)因子刺激的細(xì)胞反應(yīng)、橫紋肌組織發(fā)育、肌動(dòng)蛋白結(jié)合、膠原結(jié)合、細(xì)胞質(zhì)基質(zhì)、肌鈣蛋白復(fù)合物和肌球蛋白絲等，其中MYH10和MYH15基因富集在細(xì)胞質(zhì)基質(zhì)，F(xiàn)GF16和FGF10基因富集在成纖維細(xì)胞增殖的正調(diào)控，GDF8基因富集在成肌細(xì)胞分化的負(fù)調(diào)控，說(shuō)明雞生長(zhǎng)發(fā)育需要復(fù)雜的基因調(diào)控。KEGG通路富集分析發(fā)現(xiàn)，MYH10、MYH15、FGF10、FGF1基因富集在心肌收縮、MAPK信號(hào)通路、心肌細(xì)胞中的腎上腺素能信號(hào)傳導(dǎo)、緊密連接、肌動(dòng)蛋白細(xì)胞骨架信號(hào)通路，其中心肌收縮、MAPK信號(hào)通路、心肌細(xì)胞中的腎上腺素能信號(hào)傳導(dǎo)富集靠前。MAPK信號(hào)通路對(duì)細(xì)胞增殖、分化、遷移、衰老和細(xì)胞凋亡起著至關(guān)重要作用[27]，另外緊密連接是細(xì)胞連接的重要組成部分，它構(gòu)成離子和分子通過(guò)細(xì)胞的屏障，調(diào)節(jié)質(zhì)膜頂端和基底蛋白和脂質(zhì)的運(yùn)動(dòng)[28]，緊密連接參與控制細(xì)胞增殖和基因表達(dá)的作用[29]。本研究也發(fā)現(xiàn)了MAPK信號(hào)通路和緊密連接被顯著富集，說(shuō)明MAPK信號(hào)通路和緊密連接對(duì)雞生長(zhǎng)發(fā)育起到關(guān)鍵作用。

4 結(jié) 論

本試驗(yàn)通過(guò)對(duì)廣西麻雞2個(gè)生長(zhǎng)時(shí)期腿肌進(jìn)行轉(zhuǎn)錄組分析，共獲得2 304個(gè)差異表達(dá)基因，其中998個(gè)基因上調(diào)，1 306個(gè)基因下調(diào)。通過(guò)GO功能和KEGG通路富集分析共獲得MYH10、MYH15、FGF10、FGF16和GDF8共5個(gè)生長(zhǎng)發(fā)育關(guān)鍵基因。