蛋雞繁殖性能相關長鏈非編碼RNA的篩選與鑒定

劉 宇,段曉燕

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

隨著測序技術的不斷發(fā)展，大量哺乳動物基因組已經(jīng)完成了測序并進行了組裝。在哺乳動物基因組中，mRNA被認為是最重要的元件，但是其僅占基因組或轉錄組的一小部分，其余為非蛋白編碼基因包括siRNAs、miRNAs、piRNAs以及l(fā)ncRNA等等，這些非編碼RNA也參與生命活動各個方面的調(diào)控過程[1-2]。

不同于小RNA分子，lncRNA是一種類似于mRNA長度大于200個堿基，在轉錄水平、轉錄后水平和染色體重塑等多種層面調(diào)控基因表達的的保守性較低的非蛋白編碼RNA分子[3-4]。非編碼RNA的作用包括且不僅限于X染色體失活、P53介導的細胞凋亡、癌細胞轉移、誘導干細胞的重編程以及其他的一些過程[5-6]。隨著第二代測序技術的發(fā)展以及l(fā)ncRNA數(shù)據(jù)庫的建立，為大規(guī)模篩選與性成熟調(diào)節(jié)過程相關的lncRNA提供了技術依據(jù)。

家禽的性成熟過程受到不同轉錄本構建的分子調(diào)控網(wǎng)絡共同作用，而卵巢是一個獨特又重要的生殖器官和內(nèi)分泌腺，蛋雞的卵巢發(fā)育程度直接關系著其產(chǎn)蛋性能的優(yōu)劣，卵巢組織有關基因或非編碼RNA表達的差異類型及其表達量的變化在分子水平上影響蛋雞繁殖性能[7]。但相關研究，特別是關于lncRNA仍報道甚少。

本研究基于二代測序技術對不同產(chǎn)蛋階段蛋雞卵巢組織差異表達lncRNA進行篩選與鑒定，進一步挖掘影響蛋雞繁殖性能的關鍵調(diào)節(jié)基因以及功能元件，研究成果有助于研究人員更好的理解蛋雞性成熟過程中的分子調(diào)控機制，從而為蛋雞育種工作服務。

1 材料與方法

1.1 樣品采集與測序

以張家口地區(qū)飼養(yǎng)的京紅蛋雞為研究對象。整個試驗周期內(nèi)，所有蛋雞均自由取食、飲水，執(zhí)行常規(guī)光照和免疫等程序。于20周齡、30周齡分別隨機選取產(chǎn)蛋雞，稱重后取最接近平均值的3羽產(chǎn)蛋雞，屠宰并取出卵巢組織，放入液氮中保存?zhèn)溆?20周齡：LH_5M，30周齡：LH_6M)。獲取的卵巢組織采用Trizol一步法提取總RNA，采用Bioanalyzer 2100及RNA 6000 Nano Lab Chip Kit RNA對RNA質(zhì)量與純度進行檢測。采用Ribo-Zero Gold Kit對10μg總RNA進行核糖體RNA去除，離子打斷法將總RNA裂解為小片段，逆轉錄法合成cDNA文庫，以PCR法檢測文庫大小及濃度，采用Illumina Hiseq 4000(杭州聯(lián)川)雙末端測序法上機測序。

1.2 序列比對、lncRNA鑒定與表達分析

對測序原始數(shù)據(jù)去除低質(zhì)量數(shù)據(jù)，得到的有效數(shù)據(jù)(Valid Data)采用Bowtie2[8]和物種的參考基因組(Gallus_gallus_5.0)進行比對，同時根據(jù)基因組注釋文件(gtf和gff)指定的基因位置信息分別進行統(tǒng)計。對獲得的轉錄本去除已知的mRNA和小于200 bp的轉錄本，再對剩下的轉錄本進行l(wèi)ncRNA預測。預測軟件為CPC(Coding Potential Calculator)和CNCI(Coding-Non-Coding Index)。

采用StringTie和Ballgown對測序獲得的轉錄本進行表達水平計算，基因的表達水平主要采用每百萬測序堿基中每千個轉錄子測序堿基中所包含的測序片斷數(shù) FPKM度量基因表達的豐度值。差異表達lncRNA使用R語言的Ballgown包進行差異分析，選擇閾值矯正P值Padj<0.05。

1.3 lncRNA靶基因預測及功能富集分析

對篩查到的差異表達lncRNA進行的功能富集分析依據(jù)其靶基因進行。根據(jù)lncRNA位置確定lncRNA靶基因。篩選獲得的靶基因用于功能富集分析，包括Gene ontology和KEGG通路分析。采用David軟件進行功能富集分析，ggplot2作散點圖展示。

1.4 lncRNA差異表達驗證

對高通量測序結果篩查到的差異表達lncRNA，利用實時熒光定量驗證表達量差異性。qRT-PCR引物由Primer 5.0軟件設計(表1)。以3-磷酸甘油醛脫氫酶基因(GAPDH)作為內(nèi)參。反應體系包括95 ℃ 3 min、95 ℃15 s以及60 ℃和72 ℃各30 s，共計40循環(huán)，熔解階段由95 ℃ 30 s和60 ℃ 1 min組成，運用2-ΔΔCt法計算不同組別lncRNA的倍數(shù)變化，使用SPSS 26.0單因素方差分析進行分析，設定P<0.05為差異有統(tǒng)計學意義。

表1 引物信息Table 1 Primer information

2 結 果

2.1 轉錄組測序及序列比對結果

各樣品測序數(shù)據(jù)產(chǎn)出如表2所示，共得到57.95 Gb Clean data，單個樣品Clean data大于8.29 Gb。GC 含量在48%～51%范圍內(nèi)，Q30堿基百分比大于94.71%；說明測序數(shù)據(jù)獲取完整，可用于后續(xù)分析。參考基因組比對結果顯示檢測樣品的Reads與參考基因組的比對效率介于76.90%～79.61%之間，大部分轉錄本在參考基因組中均是唯一比對位置。

表2 轉錄組測序及序列比對結果Table 2 Result of RNA-seq and sequence align

2.2 lncRNA篩選鑒定及特征分析

對測序獲得的數(shù)據(jù)，按照篩選依據(jù)，在20周齡和30周齡京紅蛋雞卵巢組織篩選到1 283個lncRNA。對lncRNA的分類(圖1A)表明，篩選到的lncRNA多為未知或基因間轉錄本，其次為反義鏈外顯子重疊區(qū)域。對篩選獲得的lncRNA進行基因組定位，結果見圖1B。對lncRNA與mRNA的結構特征和表達模式分析發(fā)現(xiàn)，lncRNA長度與mRNA類似，多數(shù)分布在>1 000 bp范圍(圖1C)。而表達量上，lncRNA要明顯低于蛋白質(zhì)編碼基因(圖1D)，log10(FPKM)分布在0～1之間，同時數(shù)量也少于蛋白質(zhì)編碼基因。

圖1 lncRNA轉錄組特征A.lncRNA分類統(tǒng)計；B.lncRNA染色體分布統(tǒng)計；C.lncRNA長度統(tǒng)計；D.lncRNA表達量與數(shù)量統(tǒng)計Fig.1 Transcriptom characteristics of lncRNAA. Classified statistic of lncRNA;B. Chormsome distribution of lncRNA;C. Statistic of lncRNA length;D. Expression level and amount of lncRNA

2.3 差異表達lncRNA篩選

在本研究中，采用矯正P值≤0.05作為閾值，共篩選到10個差異表達lncRNA，其中7個在20周齡蛋雞卵巢組織中高表達，3個在30周齡蛋雞卵巢組織中高表達(圖2A)?；鹕綀D展示了兩組間最具顯著性的差異基因(|log2FoldChange|≥1且Padj≤0.05)(圖2B)，可以看到篩選到的10個差異表達lncRNA表達量均較低，符合lncRNA表達模式，但大量基因處于差異不顯著狀態(tài)。對篩選到的10個差異表達lncRNA的qRT-PCR結果顯示與轉錄組測序結果相同(圖2C)。

圖2 差異表達lncRNA分析A.差異表達lncRNA數(shù)量；B.火山圖；C.qRT-PCR結果Fig. 2 Analysis of differential expressed lncRNAsA. amount of differential expressed lncRNA;B. volcano map;C. qRT-PCR result

2.4 lncRNA功能富集分析

基于位置關系，1 283個lncRNA潛在的靶基因有2 960個，結合表達模式分析，10個差異表達lncRNA中有2個lncRNA具有順式調(diào)控靶基因，分別為MSTRG.14606對應的肌動蛋白γ2基因(actin gamma 2，ACTG2)和MSTRG.15311對應的肌間線蛋白基因(desmin，DES)，同時lncRNA與mRNA表達量存在負相關。對lncRNA的功能富集分析結果如圖3所示。由圖3可見，顯著富集的前5條GO通路包括：心率調(diào)節(jié)(GO:0086091)、結構分子激活(GO:0005198)、蛋白激酶結合(GO:0019901)、離子通道結合(GO:0044325)以及基底外側質(zhì)膜(GO:0016323)。顯著富集的前5條KEGG通路包括：焦點粘連通路(4510)、mTOR 信號通路(4150)、孕酮介導的卵母細胞成熟(4914)、卵母細胞成熟(4114)以及p53 信號通路(4115)。

圖3 功能富集分析A.GO功能富集分析；rich factor表示位于該GO的差異基因個數(shù)/位于該GO的總基因數(shù)；B.KEGG通路分析; rich factor表示位于該KEGG的差異基因個數(shù)/位于該KEGG的總基因數(shù)Fig.3 Functional enrichment analysisA. GO enrichment analysis; rich factor represents the DE gene number/total gene number in this GO term;B. KEGG pathway analysis; rich factor represents the DE gene number/total gene number in this KEGG pathway

3 討 論

高通量測序和組學技術研究方式的發(fā)展使得研究人員能夠?qū)Σ煌纳镞^程中的分子調(diào)控網(wǎng)絡進行篩查，這已經(jīng)成為解析生物學過程的金標準[9]。本研究利用RNA-seq技術對不同生長發(fā)育階段的蛋雞卵巢組織的lncRNA進行篩選與鑒定，在20周齡和30周齡京紅蛋雞卵巢組織中，共鑒定到1 283個lncRNA表達，并對lncRNA進行了轉錄組學特征分析，lncRNA在不同染色體分布不同，且表達量與mRNA存在較大差異，這與其他研究結果一致[10]。同時qRT-PCR技術也驗證了RNA-seq的結果，提示RNA-seq在lncRNA研究中可以發(fā)揮重要作用。

20周齡為京紅蛋雞的產(chǎn)蛋初期，30周齡京紅蛋雞達到產(chǎn)蛋高峰期，這期間伴隨著蛋雞生長與繁殖過程，存在著復雜的分子調(diào)控網(wǎng)絡，特別是基因表達的變化，往往會受到非編碼RNA的影響[11]。本研究在20周齡和30周齡京紅蛋雞卵巢組織中篩查到10個差異表達lncRNA，提示這兩個生長發(fā)育階段lncRNA表達譜差距較小，但在本研究中發(fā)現(xiàn)差異表達lncRNA MSTRG.14606和MSTRG.15311對ACTG2和DES可能存在順式調(diào)節(jié)作用機制，且lncRNA表達模式與mRNA存在負相關。Felicioni等[12]報道了ACTG2基因與豬繁殖性能有關，MSTRG.14606是否能夠通過影響ACTG2表達量來影響蛋雞繁殖性能，MSTRG.14606和MSTRG.15311能否發(fā)揮分子海綿作用機制或是存在其他調(diào)控機制，仍有待于深入研究。

GO將采用lncRNA靶基因進而對lncRNA功能進行預測[13-14]。本研究發(fā)現(xiàn)，lncRNA在心率調(diào)節(jié)、結構分子激活、蛋白激酶結合、離子通道結合以及基底外側質(zhì)膜等通路中富集。蛋雞繁殖性能與多種因素有關。蛋白激酶激活通路中AMP激活的蛋白激酶與蛋雞生產(chǎn)性能、血液生化指標等相關[15]，提示被研究發(fā)現(xiàn)的差異表達lnRNA可能通過下丘腦軸調(diào)節(jié)性激素分泌從而影響蛋雞繁殖性能。細胞組分方面，基底外側質(zhì)膜位于前列，該成分與細胞代謝有關，參與代謝物的轉運，對維持細胞穩(wěn)態(tài)發(fā)揮重要作用[16-17]。蛋雞繁殖過程中代謝活動旺盛，后續(xù)研究中，可圍繞lncRNA與蛋雞代謝過程的關系開展。KEGG通路側重分析基因產(chǎn)物及其在代謝途徑作用中的關系[18-19]，本研究發(fā)現(xiàn)差異表達lncRNA的靶基因富集在與卵母細胞成熟有關的通路中。蛋雞性成熟過程中，卵黃前體物：極低密度脂蛋白VLDL和卵黃生成素VTG能夠經(jīng)由卵泡膜受體介導卵母細胞成熟[20]。這一過程受到多因素調(diào)節(jié)，本研究發(fā)現(xiàn)的差異表達lncRNA可能通過調(diào)節(jié)相關基因表達參與到上述生物學過程中。兩種功能富集結果表明，lncRNA可能參與到蛋雞不同生長發(fā)育階段卵巢組織各種生物學途徑中，但是具體的作用機制仍有待于功能研究、體外試驗等進一步驗證。

4 結 論

本研究表明，從20周齡和30周齡蛋雞卵巢組織篩選到了lncRNA，其轉錄組特征與mRNA存在差異，篩選獲得的差異表達lncRNA顯著富集在蛋雞繁殖性能相關的生物學過程，細胞組分，分子功能以及多種信號通路中。