段曉燕，劉 宇

(1.河北北方學(xué)院 教務(wù)處，河北 張家口 075000 2.河北北方學(xué)院實驗動物中心，河北 張家口 075000)

卵巢是蛋雞重要的生殖器官和內(nèi)分泌腺，其發(fā)育程度直接關(guān)系著蛋雞產(chǎn)蛋性能的優(yōu)劣。隨年齡變化，蛋雞卵巢狀況與產(chǎn)蛋性能直接相關(guān)[1]。同時，每個卵泡由卵泡卵囊層、顆粒層和卵母細(xì)胞周圍的基質(zhì)組成，而卵泡細(xì)胞的結(jié)構(gòu)和數(shù)量在卵巢發(fā)育階段有顯著差異。研究發(fā)現(xiàn)，在不同產(chǎn)蛋周期，蛋雞在經(jīng)歷了密集的新陳代謝后，不同組織與某些功能分子表達量的改變有關(guān)[2-3]。因此，篩選不同生長發(fā)育階段蛋雞卵巢組織的DEGs能夠從分子層面幫助研究人員蛋雞卵巢組織特征，進而能夠分析影響蛋雞繁殖性能的關(guān)鍵調(diào)節(jié)基因。 轉(zhuǎn)錄組學(xué)因其高通量、高靈敏性的優(yōu)勢可用來獲取某物種細(xì)胞或組織的完整轉(zhuǎn)錄組數(shù)據(jù)，可以分析物種在不同條件下的基因表達變化[4]?；诖?，本研究利用RNA-seq分析不同產(chǎn)蛋周期蛋雞卵巢組織DEGs，以期為蛋雞分子育種提供新的理論資料。

1 材料與方法

1.1 樣品采集與文庫構(gòu)建

本研究以張家口地區(qū)飼養(yǎng)的京紅蛋雞為研究對象。100只蛋雞按照標(biāo)準(zhǔn)飼養(yǎng)管理程序飼養(yǎng)，分別于5月齡(產(chǎn)蛋初期，產(chǎn)蛋率46%，LH_5M)，6月齡(產(chǎn)蛋高峰期，產(chǎn)蛋率93%，LH_6M)采集蛋雞卵巢組織(各組n=3)。樣品采集至無RNA酶EP管后立即投入液氮，運回實驗室用于總RNA提取?？俁NA采用trizol一步法提取，提取后采用Nanodrop檢測濃度，以RIN值大于7作為提取成功的標(biāo)準(zhǔn)。選取10 μg檢測合格的RNA，采用Epicentre Ribo-Zero Gold Kit (Illumina, San Diego, USA)去除核糖體RNA(rRNA)后隨機打斷為小片段，利用逆轉(zhuǎn)錄法合成cDNA文庫，平均雙末端文庫長度為(300±50)bp。

1.2 轉(zhuǎn)錄組測序及序列比對

構(gòu)建好的cDNA文庫采用Illumina Hiseq 4000平臺進行雙末端測序(杭州聯(lián)川)。對測序儀下機的原始數(shù)據(jù)(Raw Data)進行預(yù)處理，去除低質(zhì)量片段。使用tophat[5]對Valid Data進行雞參考基因組(Gallus_gallus_5.0)比對，統(tǒng)計相關(guān)信息。

1.3 基因表達分析及DEGs篩選與鑒定

用轉(zhuǎn)錄本組裝軟件StringTie[6]對reads進行組裝和定量，基因表達量采用每百萬測序堿基中每千個轉(zhuǎn)錄子測序堿基中所包含的測序片斷數(shù)(Fragments Per Kilobase of exon model per Million mapped reads，F(xiàn)PKM)表示。DEGs的篩選采用R包Ballgown進行差異統(tǒng)計和可視化繪圖，由于本研究所使用的樣本含有生物學(xué)重復(fù)，因此DEGs的篩選標(biāo)準(zhǔn)為校正后的P值padj<0.05。根據(jù)樣品基因表達譜的相近程度，將基因進行聚類分析，繪制熱圖。以log2(foldchange)為橫坐標(biāo)，-log10(pvalue)為縱坐標(biāo)，對差異表達分析中所有的基因繪制火山圖。根據(jù)生物學(xué)意義、引物表達量等信息選取5個DEGs，利用實時熒光定量(qRT-PCR)驗證表達量差異性。引物采用Primer 5.0軟件設(shè)計，詳情見表1。以3-磷酸甘油醛脫氫酶基因(GAPDH)作為內(nèi)參。反應(yīng)體系包括95 ℃ 3 min、95 ℃15 s以及60 ℃和72 ℃各30 s，共計40循環(huán)，熔解階段由95 ℃ 30 s和60 ℃ 1 min組成，運用2-ΔΔCt法計算不同組別lncRNA的倍數(shù)變化，數(shù)據(jù)以差異倍數(shù)的log2 值(log2FoldChange)表示，使用SPSS26.0單因素方差分析進行分析，設(shè)定P<0.05為差異有計學(xué)意義。利用DAVID軟件對篩選到的DEGs進行功能富集分析。采用ggplot2對GO富集分析結(jié)果和KEGG通路富集分析結(jié)果以散點圖展示。

表1 引物信息

2 結(jié)果與分析

2.1 轉(zhuǎn)錄組測序結(jié)果

各組平均測序數(shù)據(jù)產(chǎn)出如表2所示，共得到 57.95Gb Clean Data，樣品GC含量，Q30比值符合標(biāo)準(zhǔn)，表測序數(shù)據(jù)獲取完整，可用于后續(xù)分析。參考基因組比對結(jié)果顯示檢測樣品的Reads與參考基因組的平均比對效率在兩個組分別為78.4%和76.90%，大部分轉(zhuǎn)錄本在參考基因組中均是唯一比對位置(大于67%)。

表2 轉(zhuǎn)錄組測序及序列比對結(jié)果

依據(jù)比對結(jié)果，大部分reads比對到外顯子區(qū)域，其次為內(nèi)含子區(qū)域，而基因間區(qū)域最少(圖1)。此外，本研究對比對數(shù)據(jù)染色體分布進行了統(tǒng)計，結(jié)果顯示染色體內(nèi)部定位的reads總數(shù)與染色體大小呈正相關(guān)，1號染色體定位reads總數(shù)最多而W染色體最少。

圖1 比對區(qū)域統(tǒng)計

2.2 基因表達分析

通過與參考基因組比對結(jié)果，本研究在5月齡和6月齡京紅蛋雞卵巢組織共篩查到39985個轉(zhuǎn)錄本表達，其中包含19586個蛋白質(zhì)編碼基因表達，其余轉(zhuǎn)錄本可能為新轉(zhuǎn)錄本或蛋白質(zhì)編碼基因可變剪接產(chǎn)物。本研究主要對蛋白質(zhì)編碼基因表達進行分析。各樣本基因表達值分布統(tǒng)計如表3所示，基因FPKM大部分位于0.3～60之間，同組間樣本重復(fù)性良好。以不同樣本各自表達量處理后數(shù)值 log10(FPKM)做表達值密度圖，各生物學(xué)重復(fù)樣本的表達趨勢趨于一致(圖2)。

圖2 基因表達密度分布

表3 各樣本基因表達量分布統(tǒng)計

2.3 DEGs分析

以R語言的Ballgown包對StringTie組裝和定量完畢的基因進行差異分析，采用矯正P值≤0.05作為閾值，共篩選到63個DEGs，其中52個在6月齡齡蛋雞卵巢組織中高表達，11個在5月齡蛋雞卵巢組織中高表達(圖3A)。qRT-PCR結(jié)果顯示出與RNA-seq結(jié)果一致的趨勢(圖3B)。研究發(fā)現(xiàn)DEGs較少，但表達量聚類分析熱圖顯示DEGs仍可以反映兩個轉(zhuǎn)錄組的區(qū)別(圖3C)，但大量基因處于差異不顯著狀態(tài)。

圖3 DEGs分析

2.4 DEGs功能富集分析

對DEGs的功能富集顯示，63個DEGs顯著富集到273條GO條目以及2條KEGG通路中。GO和KEGG分析可視化結(jié)果以散點圖展示：Rich fac-tor表示位于該條目的差異基因個數(shù)/位于該條目的總基因數(shù)，Rich factor越大，富集程度越高(圖4)。GO富集分析結(jié)果顯示生物學(xué)功能方面富集在離子通道結(jié)合(ion channel binding, GO: 0044325)、糖酵解途徑正調(diào)控(positive regulation of glycolytic process, GO:0045821)以及鉀離子跨膜運輸負(fù)向調(diào)節(jié)(negative regulation of potassium ion transmembrane transport, GO:1901380)；在細(xì)胞組分方面富集在Z膜(Z disc, GO:0030018)、細(xì)胞骨架(cytoskeleton, GO:0005856)及細(xì)胞外基質(zhì)(extracellular matrix, GO:0031012)；在分子生物學(xué)途徑方面富集在固醇轉(zhuǎn)運ATP酶激活途徑(sterol-transporting ATPase activity, GO:0034041)以及金屬離子質(zhì)子泵逆向轉(zhuǎn)運途徑(metal ion:proton antiporter activity, GO:0051139)。在KEGG通路富集分析方面，本研究僅發(fā)現(xiàn)DEGs顯著富集在兩條通路中，分別為焦點粘連途徑(Focal adhesion)以及溶酶體途徑(Lysosome)。

3 討 論

蛋雞的繁殖調(diào)控是一項復(fù)雜的生物學(xué)過程，涉及到不同分子的共同作用，通過基因互作，參與到不同的信號通路中[7-8]。本研究便是通過轉(zhuǎn)錄組測序手段檢測不同生長發(fā)育階段京紅蛋雞卵巢組織差異比到達基因，初步繪制了卵巢組織mRNA的表達譜，鑒定出部分與蛋雞繁殖性能相關(guān)的關(guān)鍵基因。本研究發(fā)揮了轉(zhuǎn)錄組測序高通量的優(yōu)點，在京紅蛋雞卵巢組織篩選到39985個轉(zhuǎn)錄本表達，其中包含19586個蛋白質(zhì)編碼基因表達，并對其染色體分布，基因表達量進行了分析。在本研究中，鑒定到了63個DEGs，盡管DEGs數(shù)目較少，但是某些高表達基因，可能在蛋雞生長發(fā)育過程中發(fā)揮基礎(chǔ)調(diào)節(jié)作用，其功能仍有待進一步研究。同時qRT-PCR技術(shù)也驗證了RNA-seq的結(jié)果，提示RNA-seq在不同生物學(xué)過程的研究中可以發(fā)揮重要作用。

細(xì)胞外基質(zhì)組分是維持細(xì)胞正常生態(tài)的重要成分，其對維持細(xì)胞環(huán)境穩(wěn)態(tài)具有重要作用[9]。膠原蛋白是細(xì)胞外基質(zhì)的主要成分之一，已有多個研究提示其不同類型膠原蛋白的含量變化與蛋雞繁殖性能具有明顯相關(guān)性[10]。在本研究篩選到的DEGs中，COL12A1，COL8A2在6月齡蛋雞卵巢組織中高表達，其余膠原蛋白基因在兩種卵巢組織中表達量差異不顯著，但也表現(xiàn)出類似趨勢。提示相關(guān)基因可能在蛋雞生長發(fā)育過程種發(fā)揮重要作用。原肌球蛋白1基因(tropomyosin, TPM1)并且可通過競爭性結(jié)合gga-miR-449調(diào)控VEGF表達進而影響血管再生過程，而血管再生與類固醇激素合成在卵泡形成、排卵周期等方面發(fā)揮重要作用[11-12]。本研究中，TPM1基因在6月齡蛋雞卵巢組織高表達，這與Li等[13]在海蘭褐蛋雞種的研究一致，提示TPM1可能通過調(diào)節(jié)相關(guān)基因表達影響蛋雞產(chǎn)蛋性能。此外，有研究報道，IGFs、TGF-β、SF-1以及FOXL2等[14]基因可通過多種信號途徑調(diào)節(jié)卵泡生成、成熟、細(xì)胞分化以及類固醇類激素合成。本研究發(fā)現(xiàn)IGF1基因在6月齡蛋雞卵巢組織中上調(diào)表達，log2FoldChange=1.19。值得關(guān)注的是，應(yīng)激相關(guān)基因HSPB7同樣在6月齡蛋雞卵巢組織高表達，而應(yīng)激反應(yīng)會誘導(dǎo)動物機體激素水平的變化，抑制蛋雞雌激素介導(dǎo)的相關(guān)信號通路，影響蛋雞蛋品質(zhì)[15]，應(yīng)激相關(guān)基因能否成為影響蛋雞產(chǎn)蛋性能的靶標(biāo)基因，還需要相關(guān)功能研究的驗證。本研究內(nèi)容，可為動物體抗應(yīng)激的分子遺傳機制的認(rèn)識提供新的思路。

功能富集分析能夠通過映射相關(guān)條目與通路，預(yù)測DEGs功能，為后續(xù)功能實驗提供思路。本研究發(fā)現(xiàn)的顯著富集GO通路有273條，涉及物質(zhì)合成、能量代謝、細(xì)胞周期調(diào)控等多個生物學(xué)過程，提示DEGs可在不同過程中發(fā)揮作用，后續(xù)的功能驗證可以圍繞相關(guān)內(nèi)容展開。而本研究僅發(fā)現(xiàn)DEGs富集在2條KEGG通路中，分別為焦點粘連途徑以及溶酶體途徑。盡管未達到顯著富集程度，仍發(fā)現(xiàn)一些DEGs參與到影響蛋雞繁殖性能的信號通路中，比如在6月齡蛋雞卵巢組織中下調(diào)表達的PLIN2基因位于PPAR信號通路中，PPAR信號通路是經(jīng)典脂代謝通路，調(diào)節(jié)脂肪酸代謝、脂質(zhì)合成等途徑[16]，研究結(jié)果提示蛋雞在不同生長發(fā)育階段可能存在脂質(zhì)代謝合成的差異性，這也可能影響蛋雞的產(chǎn)蛋性能。

4 結(jié) 論

本研究通過轉(zhuǎn)錄組學(xué)技術(shù)篩選了不同生長發(fā)育階段蛋雞卵巢組織DEGs，而這些DEGs參與了能量代謝、物質(zhì)合成、氧化應(yīng)激等多個生物學(xué)途徑，研究內(nèi)容反映了蛋雞卵巢組織在不同階段的差異性，同時為后續(xù)開展相關(guān)基因的功能研究以及蛋雞分子育種提供理論資料。