基于RNA-Seq 高通量測(cè)序技術(shù)對(duì)肉牛下丘腦基因表達(dá)分析及新基因挖掘研究

雷 虹，孫 鏖，胡雄貴，李付強(qiáng)，張 爽，易康樂(lè)*

(1.湖南省畜牧獸醫(yī)研究所，湖南 長(zhǎng)沙 410130； 2.湖南天華實(shí)業(yè)有限公司，湖南 漣源 417100； 3.吉林省永吉縣口前鎮(zhèn)畜牧獸醫(yī)站，吉林 永吉 132100)

在家畜動(dòng)物發(fā)情周期和生殖活動(dòng)中，下丘腦-垂體-卵巢(the hypothalamic-pituitary-ovary，HPO)軸起著十分重要的作用[1-3]。經(jīng)典的繁殖免疫學(xué)理論中，雌性動(dòng)物的繁殖活動(dòng)過(guò)程要受到下丘腦、垂體前葉和性腺之間復(fù)雜的相互作用和調(diào)控。一般認(rèn)為，下丘腦部分泌促性腺激素釋放激素(GnRH)并形成脈沖，驅(qū)動(dòng)了整個(gè)性腺軸由大腦到卵巢的信號(hào)調(diào)節(jié)通路。 垂體和卵巢組織中的配子發(fā)生和類(lèi)固醇生成都要受到GnRH的調(diào)節(jié)。 垂體促性腺激素上的受體與GnRH 具有極高結(jié)合親和力，從而刺激促性腺激素脈沖式的釋放。在GnRH 的影響下，卵泡刺激素(FSH)和黃體生成激素(LH)對(duì)卵泡的發(fā)育和卵母細(xì)胞的成熟進(jìn)行調(diào)控。與此同時(shí)，脈動(dòng)式的GnRH 釋放方式對(duì)于動(dòng)物的正常生殖功能起到了至關(guān)重要的作用[4-5]。因此，鑒于GnRH 在雌性動(dòng)物繁殖活動(dòng)中至關(guān)重要的作用，體現(xiàn)出下丘腦部在動(dòng)物的發(fā)情調(diào)控及生殖活動(dòng)中起到了關(guān)鍵主導(dǎo)作用。沒(méi)有下丘腦分泌的GnRH，促性腺和性腺就不起作用[6-7]。下丘腦的激素分泌活動(dòng)受多種因素影響，例如品種、營(yíng)養(yǎng)和季節(jié)，而這些因素是如何影響下丘腦生理活動(dòng)的分子機(jī)制仍不清楚。

在動(dòng)物研究領(lǐng)域中，基于高通量測(cè)序 RNA-Seq 技術(shù)的轉(zhuǎn)錄組學(xué)的運(yùn)用越來(lái)越深入[8-9]。在肉牛遺傳標(biāo)記基因功能注解數(shù)據(jù)庫(kù)中，大部分的基因功能注釋是通過(guò)人類(lèi)、小鼠和其他哺乳動(dòng)物基因庫(kù)中比對(duì)后預(yù)測(cè)得到的,而肉牛特有特異基因的挖掘和注釋可能還不全面或存在遺漏。因此，本研究擬對(duì)不同年齡階段的安格斯母牛下丘腦組織進(jìn)行采集后，進(jìn)行mRNA水平轉(zhuǎn)錄組分析。通過(guò)應(yīng)用RNA-Seq 技術(shù)對(duì)其進(jìn)行高通量轉(zhuǎn)錄組測(cè)序分析，篩查相關(guān)表達(dá)基因，篩選出與性腺軸調(diào)控繁殖過(guò)程相關(guān)的基因與信號(hào)通路，為進(jìn)一步了解肉牛的環(huán)境適應(yīng)性、分子育種及探析繁殖性能等諸多方面提供依據(jù)。

1 材料與方法

1.1 樣品采集在湖南漣源國(guó)家安格斯牛核心育種場(chǎng)選取9頭半同胞(同父異母)、體型發(fā)育基本一致、健康的9頭安格斯母牛，6，18和30月齡各3頭。于2018年4月9日進(jìn)行屠宰。肉牛屠宰后立即采集卵巢組織樣本,置于液氮中速凍,并于-70℃保存?zhèn)溆?用于提取組織總RNA。本研究不涉及動(dòng)物倫理及福利問(wèn)題,所有樣本均采集自政府審查批準(zhǔn)的相關(guān)屠宰場(chǎng)。

1.2 RNA定量和鑒定測(cè)試RNA樣品的純度、濃度和完整性，以確保將合格的樣品用于轉(zhuǎn)錄組測(cè)序。每個(gè)樣品總共1μg RNA用作RNA樣品制備的輸入材料。使用Illumina的NEBNext UltraTMRNA庫(kù)制備試劑盒(美國(guó)NEB)生成測(cè)序庫(kù)，并將索引代碼添加到每個(gè)樣品的屬性序列中。對(duì)文庫(kù)質(zhì)量進(jìn)行評(píng)估。

1.3 聚類(lèi)和排序使用BMKCloud(www.biocloud.net)進(jìn)行RNA-Seq分析。根據(jù)制造商的說(shuō)明，使用TruSeq PE Cluster Kit v4-cBot-HS(Illumia)對(duì)樣品數(shù)據(jù)在cBot簇生成系統(tǒng)上進(jìn)行聚類(lèi)分析，并將文庫(kù)的制備物在Illumina系統(tǒng)平臺(tái)上測(cè)序，從而生成配對(duì)末端的讀取數(shù)據(jù)。

1.4 質(zhì)量控制Fastq格式的原始數(shù)據(jù)(原始讀取)首先通過(guò)內(nèi)部perl腳本進(jìn)行處理。在此步驟中，通過(guò)從原始數(shù)據(jù)中刪除包含適配器的讀取，包含ploy-N的讀取和低質(zhì)量的讀取來(lái)獲得clean data。同時(shí)，計(jì)算了clean data的Q20、Q30、GC含量和序列重復(fù)水平。所有下游分析均基于高質(zhì)量的clean data。

1.5 比較分析低質(zhì)量序列讀數(shù)將從數(shù)據(jù)集中刪除。原始序列在數(shù)據(jù)處理后被轉(zhuǎn)換為clean reads。然后將這些clean reads進(jìn)行參考基因組定位，進(jìn)一步分析和注釋僅具有完全匹配或一個(gè)錯(cuò)配的讀數(shù)。使用Hisat2工具軟件繪制參考基因組圖[10-11]。

1.6 基因表達(dá)水平的定量基因表達(dá)水平的定量通過(guò)每百萬(wàn)個(gè)作圖片段轉(zhuǎn)錄本的每千堿基片段來(lái)估算基因表達(dá)水平。

1.7 基因功能注釋用于基因功能的注釋源于以下幾類(lèi)數(shù)據(jù)庫(kù)：蛋白質(zhì)家族數(shù)據(jù)庫(kù)(Pfam)和蛋白質(zhì)直系同源簇?cái)?shù)據(jù)庫(kù)(KOG/COG)KEGG；Ortholog數(shù)據(jù)庫(kù)(KO)和 Gene Ontology 數(shù)據(jù)庫(kù)(GO)；人工注釋和審查的蛋白質(zhì)序列數(shù)據(jù)庫(kù)(Swiss-Prot)； NCBI非冗余蛋白質(zhì)序列數(shù)據(jù)庫(kù)(Nr)和非冗余核苷酸序列數(shù)據(jù)庫(kù)(Nt)[12-15]。

2 結(jié)果

2.1 總RNA提取及文庫(kù)構(gòu)建對(duì)9頭牛的下丘腦組織進(jìn)行總RNA提取。提取后的mRNA用核酸蛋白測(cè)定儀進(jìn)行測(cè)定，在總RNA 純度(D260/D280)大于2.0以上、質(zhì)量濃度大于10 mg/L的情況下，進(jìn)行下一步試驗(yàn)。應(yīng)用qPCR方法對(duì)文庫(kù)質(zhì)量進(jìn)行定量驗(yàn)證。

2.2 轉(zhuǎn)錄組測(cè)序質(zhì)量分析對(duì)牛下丘腦的9個(gè)樣品的轉(zhuǎn)錄組測(cè)定數(shù)據(jù)分析，共獲得59.74 Gb clean data，各組的樣品clean data均達(dá)到6.64 Gb，Q30堿基百分比平均在94.60%及以上(表1)。分別將各樣品的clean reads與指定參考基因組數(shù)據(jù)庫(kù)進(jìn)行比對(duì)，比對(duì)率從94.48%～94.68% 不等(表2)。統(tǒng)計(jì)了各樣品mapped reads在指定的參考基因組外顯子、內(nèi)含子和基因間區(qū)的數(shù)目，對(duì)各樣品基因組不同區(qū)域上mapped reads的分布情況進(jìn)行繪圖(圖1)。

表1 測(cè)序數(shù)據(jù)統(tǒng)計(jì)表

表2 樣品測(cè)序數(shù)據(jù)與所選參考基因組的序列比對(duì)結(jié)果統(tǒng)計(jì)表

圖1 不同樣本測(cè)序序列與參考基因組序列各位置比對(duì)百分比

2.3 表達(dá)基因分析

2.3.1表達(dá)基因篩選 對(duì)于9個(gè)樣品數(shù)據(jù)，采用并集的方法，對(duì)表達(dá)閾值≥0.1 的基因進(jìn)行篩選，并將相同的樣品進(jìn)行混池后進(jìn)行比較后，維恩圖繪制(圖2)。結(jié)果顯示，6月齡組共有21 089個(gè)基因表達(dá)，18月齡組共有20 627個(gè)基因表達(dá)，30月齡共有20 782個(gè)基因表達(dá)。有18 979個(gè)基因在所有月齡段中都有表達(dá)，有1 706個(gè)基因只在2個(gè)月齡段中共同表達(dá)，有2 149個(gè)基因只在其中的1個(gè)月齡段中表達(dá)。

圖2 不同組織表達(dá)基因維恩圖

2.3.2表達(dá)基因細(xì)胞注釋和富集分析 對(duì)表達(dá)基因進(jìn)行數(shù)據(jù)庫(kù)的功能注釋結(jié)果顯示，下丘腦組織共有2 621個(gè)基因在GO數(shù)據(jù)庫(kù)中得到注釋(圖3)，13 453個(gè)基因在KEGG數(shù)據(jù)庫(kù)中得到注釋?zhuān)?2 930個(gè)基因在KOG數(shù)據(jù)庫(kù)中得到注釋(圖4)。

圖3 表達(dá)基因GO富集柱狀圖

圖4 表達(dá)基因KOG富集柱狀圖

2.3.3表達(dá)基因的細(xì)胞功能與代謝通路分析 對(duì)所有篩選出表達(dá)基因進(jìn)行GO數(shù)據(jù)庫(kù)的標(biāo)注分析，富集較多的為細(xì)胞過(guò)程、細(xì)胞組成和附著黏連。通過(guò)GO數(shù)據(jù)庫(kù)二級(jí)分類(lèi)注解分析表明，共發(fā)現(xiàn)12類(lèi)與器官發(fā)育、性激素反應(yīng)、受精活動(dòng)和胚胎發(fā)育相關(guān)。其中新基因主要集中在動(dòng)物器官發(fā)育、單個(gè)組織繁殖過(guò)程和授精等細(xì)胞和生物過(guò)程中(表3)。

表3 卵巢發(fā)育、受精活動(dòng)和胚胎發(fā)育相關(guān)表達(dá) 基因GO分類(lèi)表

將所有表達(dá)基因標(biāo)注到KOG數(shù)據(jù)庫(kù)中，富集較多的為普通功能預(yù)測(cè)(general function prediction only)、信號(hào)傳導(dǎo)機(jī)制(signal transduction mechanisms)和翻譯后修飾、蛋白轉(zhuǎn)換、分子伴侶(posttranslational modification,protein turnover,chaperones)(表4、圖4)。

表4 相關(guān)表達(dá)基因KOG分類(lèi)表

將所有表達(dá)基因標(biāo)注到KEGG數(shù)據(jù)庫(kù)中。下丘腦組織共富集到300個(gè)信號(hào)通路。對(duì)富集最顯著的前20條信號(hào)通路進(jìn)行繪圖(圖5)。在所有的信號(hào)通路中，篩選與繁殖相關(guān)的信號(hào)通路進(jìn)行歸類(lèi)整理。結(jié)果顯示，在cAMP信號(hào)通路、Jak-STAT信號(hào)通路、卵巢類(lèi)固醇激素信號(hào)通路等16個(gè)信號(hào)通路中有參與表達(dá)(表5)。

圖5 所有表達(dá)基因KEGG富集散點(diǎn)圖(前20)

表5 與繁殖相關(guān)KEGG信號(hào)通路統(tǒng)計(jì)表

3 討論

在雞、番鴨、綿羊等物種中均有利用轉(zhuǎn)錄組測(cè)序技術(shù)對(duì)下丘腦組織進(jìn)行測(cè)序的相關(guān)研究報(bào)道[16-19]。但以不同年齡階段母牛為研究對(duì)象，進(jìn)行系統(tǒng)得下丘腦轉(zhuǎn)錄組測(cè)序的相關(guān)報(bào)道不多。本試驗(yàn)以安格斯牛為研究對(duì)象，通過(guò)選取初情期左右的安格斯牛母牛下丘腦組織進(jìn)行RNA-Seq和生物信息學(xué)分析，希望深度挖掘安格斯牛與生長(zhǎng)發(fā)育、生殖激素調(diào)控和繁殖性能相關(guān)的分子標(biāo)記信息。

本研究對(duì)9頭半同胞母牛的下丘腦組織，按不同年齡段分組后進(jìn)行混池建庫(kù)測(cè)序。測(cè)序結(jié)果顯示，6月齡組共有21 089個(gè)基因表達(dá)，18月齡組共有20 627個(gè)基因表達(dá)，30月齡組共有20 782個(gè)基因表達(dá)。有18 979個(gè)基因在所有組織中都有表達(dá)，分別占到6，18和30月齡基因表達(dá)總數(shù)的89.99%，92.01%和91.32%。結(jié)果顯示，6月齡下丘腦的基因表達(dá)數(shù)量要多于其他月齡段，說(shuō)明肉牛在生長(zhǎng)發(fā)育的初期，下丘腦的基因調(diào)控活動(dòng)較其他時(shí)期更為活躍。

為了解各月齡段下丘腦組織表達(dá)基因相關(guān)功能，對(duì)測(cè)序篩選出的表達(dá)基因進(jìn)行了KOG、GO和KEGG數(shù)據(jù)庫(kù)比對(duì)注釋分析。

通過(guò)KOG數(shù)據(jù)庫(kù)注釋分析發(fā)現(xiàn)，有14 459個(gè)基因被注釋。其中肉牛各月齡段的下丘腦組織中富集最多的是普通功能預(yù)測(cè)，為3 029個(gè)，占到全部注解數(shù)的20.95%。

通過(guò)GO數(shù)據(jù)庫(kù)注釋分析發(fā)現(xiàn)，在6月齡組中富集到的GO條目為51個(gè)，18月齡組中富集到的GO條目為52個(gè)，30月齡組中富集到的GO條目為53個(gè)。在GO的分子功能分析中，3個(gè)年齡段富集最多的都是細(xì)胞過(guò)程、細(xì)胞組成和分子結(jié)合綁定。在對(duì)表達(dá)基因進(jìn)行GO數(shù)據(jù)庫(kù)二級(jí)分類(lèi)注解分析后，發(fā)現(xiàn)在所有涉及器官發(fā)育、生殖調(diào)控相關(guān)的調(diào)控分類(lèi)中，很大比例的基因(拷貝)涉及到動(dòng)物器官發(fā)育過(guò)程，而且這些基因大部分都是新拷貝基因。說(shuō)明下丘腦組織與肉牛生長(zhǎng)發(fā)育過(guò)程密切相關(guān)，且有大量新基因有待深入研究與探討。

對(duì)所有篩查出的表達(dá)基因在進(jìn)行了GO數(shù)據(jù)庫(kù)細(xì)胞功能注釋分析后又進(jìn)行了KEGG信號(hào)通路分析。其中富集最多的是癌癥通路。出現(xiàn)這一結(jié)果的原因可能有在所有基因組功能注釋的功能基因中，與癌癥和疾病相關(guān)的基因數(shù)量和信號(hào)通路相對(duì)較多。其次參與到癌癥通路的大部分基因可能同時(shí)涉及和參與了細(xì)胞的分化、增殖和器官組織的發(fā)育。但癌癥通路在性腺軸生理活動(dòng)中確切功能及作用尚不明確，有待進(jìn)一步的研究和驗(yàn)證[20]。

與此同時(shí)發(fā)現(xiàn)，在3個(gè)月齡段的下丘腦組織都有表達(dá)的并與動(dòng)物繁殖性狀相關(guān)信號(hào)通路有17條。其中GnRH 信號(hào)通路是性腺軸主導(dǎo)的調(diào)控動(dòng)物生殖活動(dòng)的經(jīng)典傳導(dǎo)通路[21-22]。同時(shí)本研究注意到，Notch信號(hào)通路、神經(jīng)營(yíng)養(yǎng)因子信號(hào)通路和血管內(nèi)皮生長(zhǎng)因子信號(hào)通路的參與基因在性腺軸中檢出頻率最高，說(shuō)明這些信號(hào)通路可能也通過(guò)性腺軸對(duì)雌性動(dòng)物生殖活動(dòng)產(chǎn)生影響。

研究表明，血管內(nèi)皮生長(zhǎng)因子(VEGF)家族及其受體在血管生成過(guò)程中起到重要的調(diào)控作用。VEGF家族共有7個(gè)因子，分別是VEGF-A、VE-GF-B、VEGF-C、VEGF-D、VEGF- E和胎盤(pán)生長(zhǎng)因子(PIGF)-l 、PIGF-2，同時(shí)還有2個(gè)受體分別為KDR和FLT1[23]。研究已證實(shí)，VEGF家族成員在雌性哺乳動(dòng)物的卵泡發(fā)育、卵母細(xì)胞的成熟、黃體的形成與功能的維持、子宮的周期性變化等生理過(guò)程中都起著重要的作用[24-26]。Notch 信號(hào)通路參與了生物的很多生理活動(dòng)，并通過(guò)與細(xì)胞間的互作方式來(lái)調(diào)控機(jī)體各組織器官的分化和發(fā)育[27]?，F(xiàn)已發(fā)現(xiàn)Notch 信號(hào)通路的主要參與因子由4種受體因子和5 種配體共同組成[28]。在成體干細(xì)胞的不對(duì)稱(chēng)分裂過(guò)程，與細(xì)胞凋亡和衰老過(guò)程中，Notch 信號(hào)通路起到了關(guān)鍵作用[29]。其中在卵巢生殖干細(xì)胞增殖分化及卵巢衰老過(guò)程中Notch 信號(hào)通路所起到的調(diào)控作用已被證實(shí)[30]。近年來(lái)的研究結(jié)果證實(shí)，Hippo 信號(hào)通路在機(jī)體組織和器官發(fā)育過(guò)程中發(fā)揮重要的調(diào)控作用[31-32]。現(xiàn)有的研究表明Hippo 信號(hào)通路參與到了動(dòng)物生殖生理活動(dòng)過(guò)程中，在卵泡的調(diào)控生長(zhǎng)、原始卵泡的啟動(dòng)、卵巢生殖干細(xì)胞功能調(diào)節(jié)和卵巢組織延緩衰老等過(guò)程中起到了重要的調(diào)控作用[33-35]。

本研究通過(guò)對(duì)不同月齡的肉牛下丘腦組織進(jìn)行轉(zhuǎn)錄組測(cè)序，對(duì)數(shù)據(jù)進(jìn)行比對(duì)、聚類(lèi)等分析后得出以下結(jié)論：血管內(nèi)皮生長(zhǎng)因子信號(hào)通路、Notch信號(hào)通路和Hippo 信號(hào)通路可能通過(guò)下丘腦和HPO軸參與了肉牛繁殖相關(guān)調(diào)控活動(dòng)。