高、低產(chǎn)蛋量新疆伊犁鵝性腺軸轉(zhuǎn)錄組分析：基因結(jié)構(gòu)分析與新基因發(fā)掘注釋

吳盈萍，古力娜·巴克，李海英，趙曉鈺，馬 晨，姚瑩瑩，曹 妍，梅志勇，潘 璐

（新疆農(nóng)業(yè)大學(xué)動物科學(xué)學(xué)院，新疆烏魯木齊 830000）

新疆伊犁鵝（Yili Goose）是新疆優(yōu)良、珍稀的地方禽種，是我國唯一由灰雁（Anser anser）馴化而來的中小型家鵝品種[1]。新疆伊犁鵝的馴化歷史較短、馴化程度較低，其具有一定的短距離飛行能力，耐粗飼、抗病力強、抗寒耐熱且肉質(zhì)優(yōu)良，于2000 年8 月被農(nóng)業(yè)部（現(xiàn)為“農(nóng)業(yè)農(nóng)村部”，下同）列入第一批國家級畜禽品種資源保護名錄[2]。但新疆伊犁鵝的產(chǎn)蛋量較低，成年母鵝平均年產(chǎn)蛋量僅為8～12 枚[3]，同時其就巢性很強，嚴(yán)重制約了地方鵝品種產(chǎn)蛋性能的提高。

轉(zhuǎn)錄組在研究畜禽生長發(fā)育及生理活動中發(fā)揮重要作用，也是表型關(guān)聯(lián)研究的重要方法。隨著轉(zhuǎn)錄組測序技術(shù)的快速發(fā)展，RNA-seq 技術(shù)已廣泛應(yīng)用于檢測新的轉(zhuǎn)錄本[4]、差異基因篩選、可變剪接[5-7]、基因結(jié)構(gòu)優(yōu)化[8-9]以及SNP 篩查[10-11]等相關(guān)研究。張蕾[12]等人利用轉(zhuǎn)錄組測序技術(shù)對高郵鴨雙黃蛋高、低產(chǎn)組卵巢組織進行測序，對測序數(shù)據(jù)進行質(zhì)控和注釋，篩選與高郵鴨雙黃蛋產(chǎn)蛋性能相關(guān)的SNP 位點，在所覆蓋的基因組數(shù)據(jù)中，檢測到位于基因外顯子區(qū)域的108 260～116 478 個SNP 位點，轉(zhuǎn)換類型總數(shù)極顯著高于顛換類型總數(shù)（P＜0.01）。韓坤鵬[13]等通過轉(zhuǎn)錄組測序技術(shù)對京海黃雞卵巢進行測序分析，基于所選參考基因組序列，共發(fā)掘4 431 個新基因，其中1 809 個新基因得到注釋。本研究應(yīng)用RNA-seq 技術(shù)對高、低產(chǎn)蛋量伊犁鵝性腺軸組織進行SNP 查找、可變剪切篩查、基因結(jié)構(gòu)優(yōu)化以及新轉(zhuǎn)錄本預(yù)測，為進一步解析產(chǎn)蛋量調(diào)控的分子機制提供理論依據(jù)。

1 材料與方法

1.1 樣品采集 本研究所用實驗動物由新疆額敏縣恒鑫實業(yè)有限公司提供。根據(jù)系譜資料的連續(xù)產(chǎn)蛋量記錄以及本次記錄的全期產(chǎn)蛋量，在體重相近、飼養(yǎng)條件相同的3 歲齡新疆伊犁鵝母鵝中篩選產(chǎn)蛋量具有極顯著差異的高、低產(chǎn)母鵝各4 只。高產(chǎn)母鵝的產(chǎn)蛋量分別為：14、16、18、19 個，低產(chǎn)母鵝的產(chǎn)蛋量分別為：5、6、7、5 個。屠宰后快速采集新疆伊犁鵝下丘腦、垂體、卵巢組織，用PBS 沖洗1～2 次之后，立即分裝放入含有RNA 保存液的凍存管中，做好標(biāo)記，樣品編號分別為X01-04（低產(chǎn)母鵝下丘腦）、X21-24（高產(chǎn)母鵝下丘腦）、C01-04（低產(chǎn)母鵝垂體）、C21-24（高產(chǎn)母鵝垂體）、L01-04（低產(chǎn)母鵝卵巢），L21-24（高產(chǎn)母鵝卵巢），4℃靜置過夜，隨后再放到-80℃冰箱保存，用于組織樣總RNA 提取。

1.2 組織樣總RNA 提取 采用Trizol 法分別提取24 個樣品的RNA，并分別用瓊脂糖凝膠電泳、Nanodrop、Qubit、Agilent 2100 4 種方法對其質(zhì)量進行檢測以保證使用合格的樣品進行轉(zhuǎn)錄組測序。

1.3 cDNA 文庫的構(gòu)建及Illumina 測序 樣品檢測合格后進行cDNA 文庫構(gòu)建，構(gòu)建完成后先使用Qubit 2.0進行初步定量，稀釋文庫至1 ng/μL，隨后使用Agilent 2100 對文庫的插入片段長度進行檢測，插入片段長度符合預(yù)期后，使用Q-PCR 方法對文庫的有效濃度進行準(zhǔn)確定量（文庫有效濃度＞2 nM），以保證文庫質(zhì)量。庫檢合格后，進行Illumina HiSeq 2500 高通量測序平臺測序，測序讀長為雙末端測序150 bp。

1.4 測序數(shù)據(jù)分析 高通量測序得到的原始圖像數(shù)據(jù)文件經(jīng)CASAVA 堿基識別（Base Calling）分析得到原始測序序列（Raw reads），對Raw reads 進行數(shù)據(jù)過濾，得到高質(zhì)量的過濾后測序數(shù)據(jù)（Clean reads），后續(xù)分析都基于Clean reads。使用序列高效比對軟件HISAT2.0.4 對各樣品測序數(shù)據(jù)與參考基因組（https://www.ncbi.nlm.nih.gov/assembly/GCF_000002315.6）進行序列比對，從而獲取Clean reads 在參考基因組上的位置及測序樣品特有的序列信息。

1.5 基因結(jié)構(gòu)分析 SNP 分析：首先通過samtools 和picard-tools 等工具對比對結(jié)果進行染色體坐標(biāo)排序、reads 去重復(fù)等處理，然后通過變異檢測軟件samtools分別查找NGS 數(shù)據(jù)與參考序列區(qū)別（SNP calling），最后進行過濾，得到SNP 信息。

可變剪切分析：以每個進行差異可變剪切分析的比較組為單位，本研究首先統(tǒng)計發(fā)生的可變剪切事件的種類及數(shù)量，然后分別計算每類可變剪切事件表達(dá)量，最后對每類可變剪切事件進行差異分析。

新轉(zhuǎn)錄本預(yù)測及基因結(jié)構(gòu)優(yōu)化：將所有測序reads數(shù)據(jù)的基因組定位結(jié)果放到一起，用Cufflinks 進行組裝，然后用Cuffcompare 和已知的基因模型進行比較，可以發(fā)現(xiàn)新基因（相對于原有基因注釋文件）和已知基因新的外顯子區(qū)域，以及對已知基因的起始和終止位置進行優(yōu)化。

2 結(jié)果

2.1 RNA 質(zhì)量檢測及測序數(shù)據(jù)質(zhì)量檢測

2.1.1 RNA 質(zhì)量檢測 由圖1 可知，RNA 分子28 s 和18 s 條帶清晰，且亮度比例約為2:1，表明RNA 質(zhì)量好、完整性好、RNA 降解率低。

圖1 總RNA 瓊脂糖凝膠電泳檢測結(jié)果

2.1.2 測序數(shù)據(jù)及其質(zhì)量控制 由表1 可知，通過測序，高、低產(chǎn)新疆伊犁鵝分別獲得了579 563 136 個reads和596 933 010 個reads；總堿基數(shù)為176.47 G，每個樣品均獲得了5.98 G 以上的堿基，Q20 在95.10%以上，Q30 在88.57%以上，GC 含量在48.28%～50.57%之間，說明轉(zhuǎn)錄組測序得到的結(jié)果質(zhì)量較高，滿足后續(xù)結(jié)果分析的需求。

表1 測序數(shù)據(jù)質(zhì)量評估

2.2 基因結(jié)構(gòu)分析

2.2.1 SNP 分析 由表2 可知，在24 個樣品上分別得到208 788～539 612 個可能的SNP 位點。在各樣品中，大多數(shù)堿基取代均為轉(zhuǎn)換大于顛換，其中轉(zhuǎn)換型SNP 所占百分比為72.22%～73.09%；顛換型SNP 所占百分比為26.9%～27.78%。同時，24 個樣品的雜合SNP 位點所占百分比為17.83%～35.75%，純合比例為64.25%～82.17%。

表2 SNP 位點統(tǒng)計

2.2.2 可變剪切事件分類和數(shù)量統(tǒng)計 性腺軸組織樣品的各種可變剪切方式的統(tǒng)計結(jié)果如圖2 所示。在定量過程中，rMATs 采取了2 種定量方式，分別為：Junction Count only、Reads on target and junction counts。本文只對其中5 類可變剪切事件包括外顯子跳躍（SE:Skipped exon）、外顯子選擇性跳躍（MXE:Mutually exclusive exon）、第一個外顯子可變剪切（A5SS:Alternative 5'splice site）、最后一個外顯子可變剪切（A3SS:Alternative 3'splice site）和內(nèi)含子滯留（RI:Retained intron）進行了鑒定。

圖2 可變剪切分類和數(shù)量統(tǒng)計

在樣品X2 vs X0（高產(chǎn)組下丘腦vs 低產(chǎn)組下丘腦）、C2 vs C0（高產(chǎn)組垂體vs 低產(chǎn)組垂體）和L2 vs L0（高產(chǎn)組卵巢vs 低產(chǎn)組卵巢）中同時使用Junction Counts 和reads on target 檢測到的可變剪切事件分別為20 531、21 634、22 044 個，表明許多基因都包含著多種可變剪切方式，從而造成基因在表達(dá)水平上的變化，由于基因剪切方式的不同，使得生物能夠表達(dá)多種蛋白，進而引起表型的多樣性。在5 類可變剪切事件中，X2 vs X0、C2 vs C0 和L2 vs L0 均為SE 事件比例最大，MXE 次之，RI 事件比例最小。

2.2.3 已注釋基因結(jié)構(gòu)優(yōu)化 由表3 可知，本研究共有11 383 個基因的起始位置和終止位置得到了優(yōu)化，其中正鏈5 484 條，負(fù)鏈5 899 條。

表3 部分已知基因結(jié)構(gòu)優(yōu)化

2.3 新轉(zhuǎn)錄本的發(fā)掘及注釋

2.3.1 新轉(zhuǎn)錄本的發(fā)掘 本研究共發(fā)掘新基因3 119 個，獲得注釋的基因902 個。本實驗垂體和卵巢組織差異表達(dá)基因篩選過程中將Q value＜0.05 作為篩選標(biāo)準(zhǔn)，篩選差異基因。當(dāng)Q value＜0.05 基因作為下丘腦差異表達(dá)閾值時，篩選到的差異新基因為0。因此，將下丘腦差異表達(dá)基因的篩選閾值調(diào)整為P value值＜0.005。差異表達(dá)新基因共64 個，其中49 個新基因表達(dá)上調(diào)，15個新基因表達(dá)下調(diào)（表4）。

表4 差異表達(dá)新基因數(shù)目統(tǒng)計表

新轉(zhuǎn)錄本在染色體上的分布如圖3 所示，新轉(zhuǎn)錄本分布前三的染色體分別為NW_013185656.1 染色體343條、NW_013185751.1 染色體321 條、NW_013185661.1染色體314 條。

圖3 轉(zhuǎn)錄本在染色體上的分布

2.3.2 新轉(zhuǎn)錄本的注釋 利用Blast2Go 軟件對篩選到的新轉(zhuǎn)錄本進行GO 富集分析，3 119 個新轉(zhuǎn)錄本注釋到細(xì)胞組成、生物學(xué)過程和分子功能。新轉(zhuǎn)錄本共富集到2 500 條GO Term，細(xì)胞組成占60.96%，生物學(xué)過程占26.16%，分子功能占12.88.%。其中，顯著富集的有多生物過程、DNA 拓?fù)洚悩?gòu)酶活性和DNA 拓?fù)渥兓▓D4）。在細(xì)胞組成方面主要富集到其他生物細(xì)胞、細(xì)胞外區(qū)域、線粒體膜部分、線粒體內(nèi)膜等；在生物學(xué)過程方面主要富集到調(diào)節(jié)免疫系統(tǒng)過程、ATP 生物合成過程、交配行為、繁殖行為、繁殖過程等；在分子功能方面主要富集到通道調(diào)節(jié)器活性、激素活性、結(jié)構(gòu)分子活性、雌激素受體活性等。

利用KOBAS（2.0）軟件對KEGG 注釋通路進行分析，3 119 個新轉(zhuǎn)錄本注釋到49 條代謝通路。代謝途徑見圖5，顯著富集的有泛素介導(dǎo)的蛋白質(zhì)水解和核苷酸切除修復(fù)途徑。其中，核苷酸切除修復(fù)途徑富集因子最大。除顯著富集途徑外，富集排名前十的途徑有：RNA 轉(zhuǎn)運、精氨酸和脯氨酸代謝、組氨酸代謝、氨基酸的生物合成、氨?；?tRNA 的生物合成、β-丙氨酸代謝、鞘脂代謝、糖酵解/糖異生、抗壞血酸鹽和藻酸鹽代謝、核糖體。

圖5 性腺軸組織差異基因KEGG 富集散點圖

3 討 論

3.1 基因結(jié)構(gòu)分析 SNP 是指在基因組上由單個核苷酸變異形成的遺傳標(biāo)記，其數(shù)量很多，多態(tài)性豐富。從理論上看，每個SNP 位點都可以有4 種不同的變異形式，但實際上發(fā)生的只有2 種，即轉(zhuǎn)換和顛換。鄭杰[14]等研究顯示，犏牛新鮮囊胚和凍融囊胚兩樣本中轉(zhuǎn)換型SNP 比例都明顯高于顛換型，轉(zhuǎn)換/ 顛換比值分別為2.57 和2.45。黃育雯[15]對安格斯牛下丘腦和垂體的每個樣品的SNP 轉(zhuǎn)換和顛換進行計數(shù)，結(jié)果顯示，各樣品堿基轉(zhuǎn)換的可能性均超過70%，遠(yuǎn)遠(yuǎn)超過了顛換，同時，雜合SNP 位點超過37.68%。本研究在性腺軸組織24 個樣品中共獲得1 048 575 個SNP 位點，其中轉(zhuǎn)換型約占72.22%～73.09%，顛換型約占26.9%～27.78%，轉(zhuǎn)換與顛換的比值在2.60～2.72 之間。在高產(chǎn)伊犁鵝的垂體和卵巢組織中SNP 數(shù)目明顯高于低產(chǎn)伊犁鵝，提示SNP 在伊犁鵝產(chǎn)蛋量方面起重要的調(diào)控作用。

圖4 性腺軸組織差異基因GO 富集柱狀圖

可變剪切是調(diào)節(jié)基因表達(dá)和蛋白質(zhì)組多樣性的主要因素之一，在轉(zhuǎn)錄組中大量存在[16]。真核生物中存在著大量的可變剪切，動物細(xì)胞中也存在著組織特異性[17]。李慧鋒[18]等對白來航母雞首次產(chǎn)蛋前后的肝臟組織進行分析，得到可變剪切事件在產(chǎn)蛋前后分別為14 518個和14 683 個，其中外顯子跳躍事件最多，5' 可變剪切位點事件最少。Li 等[19]報道，20 周齡青年雞和30周齡產(chǎn)蛋雞肝臟的可變剪切類型主要為內(nèi)含子滯留（RI 32%）和外顯子跳躍（SE 24%）。本研究發(fā)現(xiàn)，在5類可變剪切事件中，X2 vs X0、C2 vs C0 和L2 vs L0均為SE 事件比例最大，MXE 次之，RI 事件比例最小，表明高、低產(chǎn)新疆伊犁鵝性腺軸組織的可變剪切類型中，SE 事件為主要的剪切方式，是導(dǎo)致蛋白質(zhì)形成差異的主要機制，也可能是表現(xiàn)產(chǎn)蛋量差異的主要機制。

此外，通過RNA-seq 越來越多的新轉(zhuǎn)錄本被發(fā)現(xiàn)，顧麗紅[20]等通過轉(zhuǎn)錄組測序技術(shù)，對180 日齡文昌雞和隱性白羽克洛雞胸肌進行測序，在4 個樣品中發(fā)現(xiàn)了289～339 個新轉(zhuǎn)錄本。葉保國[21]等對北京鴨6 個組織的RNA-seq 數(shù)據(jù)進行新轉(zhuǎn)錄本分析，預(yù)測得到了8 141個新轉(zhuǎn)錄本，單個樣品可以預(yù)測得到3 636～4 918 個新轉(zhuǎn)錄本。裴星朝[22]對番鴨下丘腦、垂體和卵巢3 種組織進行轉(zhuǎn)錄組深度測序，對基因結(jié)構(gòu)優(yōu)化分析，發(fā)掘了1 427 個新基因。本研究顯示，轉(zhuǎn)錄組測序結(jié)果中新轉(zhuǎn)錄本的表達(dá)譜分析揭示了這些新轉(zhuǎn)錄本也參與了伊犁鵝的繁殖過程。本研究發(fā)現(xiàn)在性腺軸組織中總共有3 119個新轉(zhuǎn)錄本，獲得注釋的基因902 個。差異表達(dá)新基因64 個，其中49 個差異新基因表達(dá)上調(diào)，15 個差異新基因表達(dá)下調(diào)。

3.2 功能分類分析 為揭示新轉(zhuǎn)錄本的生理功能，本實驗對篩選到的新轉(zhuǎn)錄本進行了GO 富集分析，其中細(xì)胞組成所占比例最大。此結(jié)果表明，參與伊犁鵝繁殖生理活動的分子組分廣泛分布于細(xì)胞組分，這提示細(xì)胞組分在繁殖生理活動中扮演著重要角色。此外，參與伊犁鵝性腺軸生理活動的新轉(zhuǎn)錄本顯著富集在多生物過程、DNA 拓?fù)洚悩?gòu)酶活性和DNA 拓?fù)渥兓?，說明這些生物學(xué)過程對于伊犁鵝繁殖性能具有重要調(diào)控作用。

KEGG 信號通路分析表明，這些新轉(zhuǎn)錄本共涉及到49 條代謝通路，其中泛素介導(dǎo)的蛋白水解和核苷酸切除修復(fù)途徑通路顯著富集。泛素-蛋白酶體途徑是一種依賴ATP 進行的、具有高度特異性和選擇性的蛋白質(zhì)降解途徑，主要負(fù)責(zé)降解細(xì)胞內(nèi)超過80% 的正?；虍惓５鞍踪|(zhì)，泛素介導(dǎo)的蛋白質(zhì)降解系統(tǒng)涉及細(xì)胞的多種重要生理功能，參與了基因轉(zhuǎn)錄、蛋白質(zhì)翻譯、信號傳導(dǎo)、細(xì)胞周期控制以及生長發(fā)育等幾乎所有的生命活動過程[23-26]。核苷酸切除修復(fù)為DNA 修復(fù)途徑之一，且核苷酸切除修復(fù)基因具有高度多態(tài)性。修復(fù)基因突變后可能引起DNA 修復(fù)能力的改變，從而導(dǎo)致機體免疫系統(tǒng)的改變[27]。生理上，DNA 修復(fù)過程需要多種酶作用，與修復(fù)蛋白表達(dá)水平密切相關(guān)，而后者受轉(zhuǎn)錄水平調(diào)控。因此，DNA 修復(fù)基因轉(zhuǎn)錄水平可以反映機體受損時細(xì)胞DNA 修復(fù)受損的能力[28]。新轉(zhuǎn)錄本富集到的信號通路分析只能為性腺軸生物學(xué)功能的信號調(diào)控途徑提供一個方向，但是新疆伊犁鵝繁殖性能的具體調(diào)控過程還需要對涉及到相關(guān)通路的對應(yīng)基因進行更深入的研究。

4 結(jié) 論

在高、低產(chǎn)蛋量新疆伊犁鵝性腺軸中共篩選出1 048 575個SNP位點，64 209個可變剪切，優(yōu)化了11 384個基因，篩選到3 119 個新轉(zhuǎn)錄本，獲得注釋的轉(zhuǎn)錄本902 個。這些新轉(zhuǎn)錄本的發(fā)現(xiàn)和功能得以注釋為進一步探索挖掘新基因奠定基礎(chǔ)，也為伊犁鵝產(chǎn)蛋量相關(guān)基因的篩選提供科學(xué)依據(jù)。