高郵鴨卵巢中雙黃蛋性狀相關lncRNA表達差異分析

唐現文，張蕾，章敬旗，章瑋玥，張響英

(江蘇農牧科技職業學院，江蘇 泰州 225300)

長鏈非編碼RNA簡稱為lncRNA(long non-coding RNA)，為長度超過200 bp且不編碼蛋白質的RNA長鏈分子[1].隨著高通量測序技術的廣泛應用，對cDNA文庫的注釋不斷增加，越來越多的lncRNA被逐漸發現[2-4].研究表明，lncRNA在表觀遺傳調控、轉錄調控以及轉錄后調控等多層面調控基因的表達，它們參與細胞增殖、分化和凋亡等多種生物學過程[5]，lncRNA作為調控性非編碼RNA的重要成員，目前已成為生命科學研究的熱點.

近年來，lncRNA在畜禽動物中的研究越來越受到人們的重視.Arriaga等[6]對雞胚發育過程中的性腺進行研究，發現lncRNA-αGT在雞發育晚期上調表達，是雞胚發育成熟的重要參與者.Li等[7]通過新一代測序技術對雞肝臟脂質調控機制進行研究，發現lncRNA和miRNA這類非編碼RNA在肝臟脂質代謝生理過程中發揮著重要的作用.任晉東[8]綜合分析了青年期、初產期以及老齡期鴨卵巢組織lncRNA與編碼基因的差異表達情況，驗證了lncRNA作為高產蛋鴨分子標記物的潛在可能.以上研究表明，lncRNA參與家禽的發育及生理機能調控.

高郵鴨因擅產雙黃蛋而聞名海內外，針對其雙黃蛋性狀的研究在DNA和mRNA層面相對比較完善，而對lncRNA介導的基因表達調控研究還比較匱乏.與高郵鴨性狀相關的lncRNA研究是探索其卵巢和卵泡發育機理的重要途徑，對于完善高郵鴨雙黃蛋性狀的調控機制具有重要意義.本研究利用全轉錄組測序技術，分析雙黃蛋高、低產組高郵鴨卵巢組織中差異表達的lncRNA，為進一步探索lncRNA在高郵鴨雙黃蛋性狀中的調控機制提供參考依據.

1 材料與方法

1.1 試驗動物

本試驗所用的高郵鴨由江蘇省泰州市豐達水禽育種場提供，選擇同一來源(同批孵化、同批出雛、同一環境下個體籠養、予以全價日糧)的健康鴨群，籠養330日齡并記錄個體產蛋數據.

1.2 樣品采集

根據個體產蛋記錄，篩選出雙黃蛋高產組(樣品編號G1，G2，G3)、低產組(樣品編號D4，D5，D6)個體各3只.按照國家實驗動物處理行為準則進行屠宰，并立即采取卵巢組織，將卵黃釋放后整體勻漿，放入無酶凍存管，經液氮速凍后，于-80 ℃保存備用.

1.3 試驗方法

1.3.1 RNA提取及轉錄組測序 通過Trizol法進行高郵鴨卵巢組織總RNA提取，Nanodrop檢測提取的RNA純度，Qubit 2.0檢測提取的RNA濃度，Aglient 2100檢測RNA的完整性.所有樣品RNA提取并檢測合格后，送至上海歐易生物技術有限公司，應用Illumina HiSeq2500高通量測序平臺進行全轉錄組測序分析[9].

1.3.2 lncRNA表達水平分析 對測序慮得的clean data，通過TopHat2軟件進行序列比對，獲得高郵鴨卵巢組織樣本特異序列信息，并進行基因組定位分析；基于每個樣本的基因組比對結果，應用StringTie軟件[10]流神經網絡算法，重新組裝并篩除已知編碼轉錄本或已知基因座的新轉錄本，用Cuffcompare軟件[11]分析轉錄本的位置類型，篩選出“i”、“u”、“x”、“o”字樣、長度大于200 bp、不具有編碼潛能的lncRNA轉錄本.用FEELnc軟件[12]分析lncRNA與已知蛋白編碼轉錄本的位置關系，統計lncRNA類型.用Bowtie2[13]和eXpress[14]估算FPKM(fragments per kb per million fragments，FPKM)，獲得lncRNA在各樣本中的表達豐度.

1.3.3 差異表達lncRNA篩選 利用DESeq[15]軟件對各個樣本lncRNA的counts數目進行標準化處理，計算差異倍數，并采用NB(負二項分布檢驗的方式)對reads數進行差異顯著性檢驗，以|log2Fold change|≥2(P<0.05)為閾值進行差異表達基因的篩選，結合GO(gene ontology，基因本體)和KEGG(Kyoto encyclopedia of genes and genomes，京都基因與基因組百科全書)數據庫對篩選獲得的差異表達基因進行功能注釋及通路富集分析.

1.3.4 qRT-PCR檢測 隨機選取1.2.2中獲得的4個差異lncRNA，采用Oligo 7.5軟件設計qRT-PCR引物(表1)，以反轉錄得到的開產前后高郵鴨卵巢組織樣本(樣本號1～6)的cDNA為模板進行PCR擴增，驗證上述基因的表達規律與全轉錄組測序是否一致.

表1 qRT-PCR引物信息

1.4 數據分析

2 結果與分析

2.1 lncRNA鑒別和表達水平分析

本試驗共構建雙黃蛋高、低產組高郵鴨6個卵巢組織轉錄組文庫，分別獲得了84 540 140～87 479 660個有效數據(clean data)[9].在此基礎上，對轉錄本重新組裝并篩選后，最終獲得新的lncRNA共計962條(表2)，lncRNA長度在201～10 125 bp之間.以100 bp為間隔區域，對上述獲得的新lncRNA進行序列長度分布統計(圖1)，本試驗獲得的lncRNA最多分布在200～500 bp以及≥2 000 bp的范圍內.

結合FEELnc軟件，將上述獲得的新lncRNA與已知蛋白編碼轉錄本進行位置對比分析，對lncRNA所在方向(正鏈、反鏈)，類型(基因內、基因間)和位置(上游、下游、內含子、外顯子)進行統計分析，結果顯示新lncRNA最多分布于反鏈基因內含子(126個)以及正鏈基因外顯子(193個)中，基因間下游區域分布最少(圖1).

通過轉錄本表達水平估算，對所有樣本lncRNA的FPKM作密度分布圖(圖3-A)，發現FPKM值在0～0.5之間的轉錄本最多，有7 796～9 643個；FPKM值≥10的轉錄本在992～2 278之間；FPKM在0.5～1之間轉錄本共計有6 199個；而FPKM值在1～10之間的轉錄本最少，僅占總數的 3.54 %.另外，還發現6個樣本的轉錄本表達豐度模式基本相同(圖3-B).

2.2 差異表達lncRNA篩選與分析

以|log2Fold change|≥2(P<0.05)為閾值對lncRNA進行顯著表達篩選，共獲得279條差異表達lncRNA，其中顯著性上調差異表達lncRNA 109個，顯著性下調差異表達lncRNA 170個(圖4).

表2 lncRNA序列信息統計表

圖1 lncRNA序列長度分布圖Figure 1 Length distribution of lncRNA

A：新lncRNA正反鏈分布圖； B：新lncRNA類型統計匯總圖.A：Novel lncRNA distribution between antisense and sense;B：Novel lncRNA classification summary chart.圖2 新lncRNA類型統計圖Figure 2 Novel lncRNA classification

A：FPKM表達量分布圖； B：豐度分布曲線圖.A：Transcript expression (FPKM) in each sample;B：Transcript expression density in each sample.圖3 FPKM表達量及豐度分布圖Figure 3 Transcript expression (FPKM) and density in each sample

A：差異表達lncRNA火山圖，灰色為非顯著性差異的lncRNA，紅色為顯著性上調lncRNA，綠色為顯著性下調lncRNA；B：差異表達lncRNA統計柱狀圖，紅色為顯著性上調lncRNA，藍色為顯著性下調lncRNA.A：Volcano plots of differently expressed lncRNA,grey colour represents non-differently expressed lncRNA,red colour represents up-regulated expressed lncRNA,green colour represents down-regulated lncRNA； B：Statistic of differently expressed lncRNA,red colour represents up-regulated expressed lncRNA,blue colour represents down-regulated lncRNA.圖4 差異表達lncRNA結果統計Figure 4 Statistic of differently expressed lncRNA

獲得差異表達lncRNA后，為進一步了解每個差異lncRNA的功能，結合差異表達lncRNA的來源及臨近基因進行GO富集分析，以對其功能進行描述.利用超幾何分布計算，篩選出對應差異lncRNA數目大于2的GO條目，再按照每個條目對應的-lgP由大到小排列，根據差異表達lncRNA上調及下調分類，列舉出上調和下調各組前30項GO條目，如圖5所示.

為顯示差異表達lncRNA富集的GO節點及其層次關系，按照上述生物過程、細胞成分以及分子功能三大類，對富集到的GO條目進行有向無環圖繪制，以顯示差異表達lncRNA之間的上下調控關系和所涉及的功能描述.圖6中分支代表包含關系，從上至下所定義的功能描述范圍越來越具體，從生物過程有向無環圖(6-A)可以看出，所富集到的lncRNA最終主要匯集于5個生物過程，分別是：細胞質翻譯(GO：0002181)，神經酰胺代謝過程(GO：0006672)，一氧化氮生物合成過程正向調控(GO：0045429)，心臟隔緣形成(GO：0060347)，以及內皮細胞分化正向調控(GO：0045603)；細胞成分有向無環圖(6-B)最終富集到：粗面內質網膜的細胞質側(GO：0098556)，胞漿小核糖體亞單位(GO：0022627)和剪接體snRNP組分(GO：0097525)3個過程；分子功能有向無環圖(6-C)最終歸為4個主要功能：mRNA鳥苷基轉移酶活性(GO：0004484)，多核苷酸5′-磷酸酶活性(GO：0004651)，三磷酸酶活性(GO：0050355)以及磷脂酰肌醇-3-磷酸結合(GO：0032266).

A：調表達lncRNA GO注釋分析TOP30；B：下調表達lncRNA GO注釋分析TOP30.A：Top 30 GO term of up-regulated lncRNA；B：Top 30 GO term of down-regulated lncRNA.上調：神經管閉合，肌動蛋白細胞骨架形成，蛋白質穩定，信號轉導，細胞遷移的正調控，細胞周期，細胞遷移，多細胞生物發育，凋亡過程，細胞分化，肌動蛋白絲，線粒體膜，褶邊，核仁周圍，基底膜，GABA能突觸，細胞外基質，細胞骨架，膜，細胞內，轉錄調控區序列特異性DNA結合，絲氨酸/蘇氨酸激酶活性，轉錄因子結合，信號受體活性，肝素結合，鋅離子結合，肌動蛋白結合，ATP結合，蛋白激酶結合，金屬離子結合；下調：初生胚層的形成，管子成形，細胞質翻譯，成肌細胞增殖，內皮細胞分化的正調控，神經酰胺代謝過程，翻譯，內皮細胞增殖的負調控，細胞極性的建立，腎水穩態，粗面內質網膜的細胞質側，多核糖體，胞漿小核糖體亞單位，回收內膜，核斑點，胞質囊泡膜，核質，核基質，細胞膜結合細胞器，質膜，磷脂酰肌醇-3-磷酸結合，醛固醇：NADP+1-氧化還原酶活性，乙醇脫氫酶(NADP+)活性，核糖體的結構組成，激活轉錄因子結合，E-box結合，氧化還原酶活性，磷脂酰肌醇-4，5-二磷酸結合，蛋白質二聚活性，β連環素結合.Up-regulated：neural tube closure,actin cytoskeleton organization,protein stabilization,signal transduction,positive regulation of cell migration,cell cycle,cell migration,multicellular organism development,apoptotic process,cell differentiation,actin filament,mitochondrial membrane,ruffle,perikaryon,basement membrane,GABA-ergic synapse,extracellular matrix,cytoskeleton,membrane,intracellular,transcription regulatory region sequence-specific DNA binding,protein serine/threonine kinase activity,transcription factor binding,signaling receptor activity,heparin binding,zinc ion binding,actin binding,ATP binding,protein kinase binding,metal ion bindingx;Down-regulated：formation of primary germ layer,tube formation,cytoplasmic translation,myoblast proliferation,positive regulation of endothelial cell differentiation,ceramide metabolic process,translation,negative regulation of endothelial cell proliferation,establishment of cell polarity,renal water homeostasis,cytoplasmic side of rough endoplasmic reticulum membrane,polysomal ribosome,cytosolic small ribosomal subunit,recycling endosome membrane,nuclear speck,cytoplasmic vesicle membrane,nucleoplasm,nuclear matrix,intracellular membrane-bounded organelle,plasma membrane,phosphatidylinositol-3-phosphate binding,alditol:NADP+1-oxidoreductase activity,alcohol dehydrogenase (NADP+) activity,structural constituent of ribosome,activating transcription factor binding,E-box binding,oxidoreductase activity,phosphatidylinositol-4,5-bisphosphate binding,protein dimerization activity,beta-catenin binding.圖5 差異表達lncRNA GO注釋Figure 5 GO annotation of diferentially expressed lncRNA

圖6 GO條目有向無環圖Figure 6 Directed acyclic graph of GO annotation

將上述生物過程、細胞成分以及分子功能三大類最終匯集到的12個GO條目進行匯總，獲得25個差異表達lncRNA轉錄本(圖7-A)，其中包含6個雙黃蛋低產組高郵鴨卵巢特異性lncRNA，1個雙黃蛋高產組高郵鴨卵巢特異性lncRNA，兩組共表達差異lncRNA 18個.基于lncRNA作用模式，對上述25個差異表達lncRNA進行染色體位置和靶基因預測，合并歸屬于同一基因位置的lncRNA，共獲得18條特異lncRNA及其對應的靶基因(表3)，并且其中有兩對差異表達lncRNA和靶基因之間存在著潛在調控關系(圖7-B)：XLOC_000924和LOC106018336分別位于CREK基因正義鏈的上游區域和反義鏈的下游區域，LOC106018964和LOC113843289同時位于LOC101797459正義鏈的下游區域.以上這些潛在調控關系為研究lncRNA對高郵鴨雙黃蛋性能的綜合調控作用提供了新的素材.

A：GO富集lncRNA韋恩分析；B：lncRNA和基因之間存在的潛在調控關系.A：Venn analysis table of lncRNA by GO enrichment；B：Potential regulatory relationships between lncRNA and genes圖7 GO富集lncRNA聚類分析Figure 7 Cluster analysis by GO enrichment

表3 GO富集lncRNA信息統計表

為進一步探究lncRNA對靶基因的潛在調控作用，通過轉錄組測序方法對本試驗中雙黃蛋高、低產組卵巢樣本中的差異表達基因進行檢測，篩選出上述lncRNA的靶基因，并檢測其是否在兩組兩本中呈現差異性表達.從表4可以看出，本試驗通過GO富集分析篩選出的18條特異lncRNA對應的靶基因均呈現差異性表達，以P≤0.05為條件進行二次篩選，共獲得5個差異顯著表達基因：ADAMTS1，ADGRG6，CERK，HEY2和LOC101797459，對應lncRNA對上述差異基因的調控機制還有待進一步驗證研究.

表4 lncRNA靶基因差異性表達

2.3 差異表達lncRNA熒光定量驗證

為檢驗轉錄組測序技術篩選出的差異表達lncRNA結果是否準確，篩選上述5個差異表達靶基因對應的lncRNA：XLOC_000440，LOC106017612，LOC106018336，LOC110352 238和LOC1060189 64，并以β-actin為內參基因來進行qRT-PCR驗證.qRT-PCR結果顯示(圖8)：上述5個lncRNA的表達規律與全轉錄組測序獲得的表達規律基本一致.所以本試驗全轉錄組測序獲得的lncRNA結果可信，可用于后續的功能驗證分析.

圖8 差異表達lncRNA qRT-PCR驗證Figure 8 qRT-PCR results of diferentially expressed lncRNA

3 討論

家禽的產蛋性能主要依賴于卵巢上卵泡的發育，而卵泡的發育是禽類機體高度協調的復雜生理過程.當家禽性成熟時，卵巢上的原始卵泡開始生長，逐漸向等級前卵泡緩慢生長，最后優勢卵泡破裂，卵子排出，進入輸卵管并最終形成禽蛋[16].雙黃蛋是由于家禽卵巢上的兩個卵泡同時發育成熟并排卵，由輸卵管接納后形成，一般家禽品種的雙黃率極低.雙黃蛋雖在營養價值上高于單黃蛋，但在育種過程中被認為是畸形蛋的一種，不具備孵化能力，因而研究相對滯后.Salamon等[17]研究表明，當雙黃的大小和位置恰當時可以孵化出健康的家禽后代，打破了雙黃蛋無法孵化的傳統觀點，進一步說明了雙黃蛋性具能穩定遺傳的可能.對于雙黃蛋的產蛋機理，目前研究多歸于生理、遺傳以及病理原因，系統的調控機理尚未見報道.

GnRHR(gonadorelin receptor，促性腺激素釋放激素受體)參與調控FSH(follicle-stimulating hormone，卵泡刺激素)的合成與釋放，從而調控卵巢的繁殖性能.有研究結果表明GnRHR基因對白來航雙黃蛋產蛋性能有一定的加性效應，但其最高雙黃蛋率在開產初期也僅為1.06%[18].Vieaud等[19]利用引自法國的白來航雙黃蛋品系WLDJ(white Leghon double-yolked Jade)，進行了雙黃蛋產蛋率與BMP15基因(bone morphogenetic protein 15，骨形態發生蛋白15)的相關性研究，未發現兩者之間有任何關聯.以上研究均從已知的單個或多個基因多態性入手，通過基因與雙黃蛋性狀的遺傳效應進行分析，尚未取得有效進展.

高郵鴨屬于蛋肉兼用型麻鴨品種，是我國的三大名鴨之一，原產于江蘇高郵地區，在我國江淮地區的飼養歷史逾百年，因善產雙黃蛋而聞名海內外，是研究雙黃蛋產蛋機理的優秀品種素材[20].本研究結合個體籠養方法，篩選出高郵鴨雙黃蛋高、低產組，并對其卵巢進行RNA-Seq全轉錄組測序分析，以此篩選與高郵鴨雙黃蛋性狀相關的lncRNA.結合lncRNA注釋和差異基因篩選，初步獲得了279條差異表達lncRNA，這些lncRNA廣泛參與了生物過程、細胞成分以及分子功能這三類調控功能，涉及范圍廣泛.為了進一步縮小與雙黃蛋性狀相關的候選lncRNA范圍，結合生物學功能分析和向無環圖繪制，最終獲得了18條候選lncRNA，它們參與細胞質翻譯、內皮細胞分化正向調控和磷脂酰肌醇-3-磷酸結合等多個生物過程.其中，磷脂酰肌醇-3-磷酸結合功能中主要涉及PI3K基因(phosphoinositide-3-kinase，磷脂酰肌醇 3激酶)，PI3K作用于多種細胞外因子的信號轉導，在細胞的生長、增殖、分化和凋亡等多種生命活動中起著重要作用[21]，顆粒細胞通過PI3K/AKT/mTOR的途徑介導自噬，進一步的調整控制卵泡發育和卵泡閉鎖[22]，因此，磷脂酰肌醇-3-磷酸結合生物過程中的lncRNA可作為高郵鴨雙黃蛋性狀的候選lncRNA進行深入驗證和研究.

本研究還通過靶基因預測，對上述18條候選lncRNA進行分析，其中有兩對差異表達lncRNA和靶基因之間存在著潛在調控關系：XLOC_000924和LOC106018336分別位于CREK基因正義鏈的上游區域和反義鏈的下游區域，LOC106018964和LOC113843289同時位于LOC101797459正義鏈的下游區域.另外，我們還關注到TCONS_00001086的靶基因ADAMTS1，其是I型血小板結合蛋白基序的解聚蛋白樣金屬蛋白酶[23](A disintegrin and metalloproteinase with thrombospondin motifs，ADAMTS)家族成員之一.有研究表明，ADAMTS1基因與卵巢功能關系密切：在卵泡的發育階段，ADAMTS1基因參與卵巢顆粒細胞功能調控，并可進一步促進卵母細胞的成熟[24]；在卵泡的排卵階段，ADAMTS1基因通過影響卵丘細胞的擴展，從而調控排卵的順利進行.本研究表明lncRNA可能參與了高郵鴨雙黃蛋性狀的調控，為研究高郵鴨雙黃蛋性能的綜合調控機制奠定了基礎.

4 結論

本試驗通過高郵鴨雙黃蛋高、低產組卵巢組織全轉錄組測序，初步篩選出了18個差異表達lncRNA和對應的靶基因，可能參與高郵鴨雙黃蛋性狀的調控.試驗結果為進一步了解高郵鴨雙黃蛋性能的綜合調控機制提供了新的參考依據.