鵝等級前卵泡組織轉(zhuǎn)錄組測序分析

孫永峰，武惠巖，王子遒，路洪濤，李書哲，楊 威，隋玉健

（吉林農(nóng)業(yè)大學動物科學技術學院，吉林長春 130118）

鵝等級前卵泡組織轉(zhuǎn)錄組測序分析

孫永峰*，武惠巖，王子遒，路洪濤，李書哲，楊 威，隋玉健

（吉林農(nóng)業(yè)大學動物科學技術學院，吉林長春 130118）

為研究鵝等級前卵泡組織的轉(zhuǎn)錄組差異，豐富鵝卵泡轉(zhuǎn)錄組數(shù)據(jù)信息，實驗選用3只經(jīng)產(chǎn)處于高峰期的籽鵝等級前卵泡，通過Illumina HiSeq 2500 測序平臺進行測序，并進行序列分析和注釋。經(jīng)測序質(zhì)檢后共得到224 929 214條序列，31 722 729 162 bp堿基，且質(zhì)量大于20的堿基比例高于95 %，從頭拼裝后得到114 486 條Unigenes序列和145 020條Transcripts序列，通過比對分析共有42 453條序列注釋到NR數(shù)據(jù)庫，并進行物種同源性比對，發(fā)現(xiàn)綠頭鴨與其同源性最高，通過GO數(shù)據(jù)庫比對，發(fā)現(xiàn)共有1 045條Unigenes注釋到GO數(shù)據(jù)上，主要涉及生物過程、分子功能和細胞組分。結果表明，共有 11 949 條Unigenes 比對到KEGG數(shù)據(jù)庫，共涉及到341條通路，其中癌癥通路、PI3K-AKT信號通路和MAPK信號通路顯著富集。為驗證基因表達差異，以FPKM值作為基因的表達量，F(xiàn)DR＜0.05，︱logFC︱≥1為標準篩選差異基因，其中以SY 等級卵泡為參考，在 PE等級卵泡組織間中得到15 516條差異基因，其中上調(diào)基因8 293條，下調(diào)基因7 223條。并對篩選得到的所有差異基因進行GO和KEGG比對。本實驗條件下，測得各等級卵泡之間基因表達差異較大，且癌癥通路、PI3K-AKT信號通路和MAPK信號通路較為活躍，可間接證明這3條通路對鵝等級卵泡發(fā)育影響較大，這對進一步分析鵝等級前卵泡發(fā)育機制提供基礎。

鵝；等級前卵泡；轉(zhuǎn)錄組測序；差異表達基因；功能分類

轉(zhuǎn)錄組測序技術（Transcriptome Sequencing）是指利用第2代高能高通量測序技術記性cDNA測序技術，可全面快速地捕獲某一物種特定器官或組織在某一狀態(tài)下所有轉(zhuǎn)錄本?？蓮恼w水平研究基因功能以及基因結構，發(fā)現(xiàn)不同生理或病理狀態(tài)下細胞、組織或個體間差異表達基因。目前隨著高通量測序技術的發(fā)展，基因研究速度加快，對畜禽卵泡、卵巢組織研究逐漸多了起來。吳陽升等[1]對綿羊小卵泡與中卵泡轉(zhuǎn)錄組差異特征分析，篩選出1 073條差異基因，其中上調(diào)基因744個，下調(diào)基因959個。朱志明等[2]對山麻鴨開產(chǎn)期和產(chǎn)蛋高峰期卵巢組織轉(zhuǎn)錄分析，得到1 929條差基因，通過GO和KEGG富集分析GnHR信號通路和TGF-β信號通路在卵巢的生理活動中發(fā)揮重要作用。韓坤鵬等[3]對京海黃雞卵巢轉(zhuǎn)錄組研究后與所選參考基因組序列對比發(fā)現(xiàn)4 431個新基因，其中1 809個基因通過與各個數(shù)據(jù)庫比對得到功能注釋。Tariq等[4]對鵝全血進行Denovo 轉(zhuǎn)錄組分析共得到211 198條功能基因（Unigenes）。

鵝（Goose）屬于鳥綱雁形目鴨科動物，雜食性水禽，是一種重要的經(jīng)濟動物，品種約150多種，廣泛分布于世界各地。中國作為鵝品種資源最為豐富的國家，飼養(yǎng)著世界約90%的家鵝[5]。但鵝業(yè)發(fā)展一直落后于其他家禽業(yè)，這主要是由于鵝的卵泡發(fā)育與其他家禽不同，具有嚴格的等級性和選擇性[6]，產(chǎn)蛋性能低下。為研究鵝卵泡的發(fā)育機制，本實驗利用高通量RNA-seq，以鵝等級前卵泡組織為研究對象，對其轉(zhuǎn)錄組進行測序分析。通過分析鵝卵泡發(fā)育過程中的轉(zhuǎn)錄組基本信息，來探討研究鵝卵泡的發(fā)育機制，并進一步研究和發(fā)掘卵泡發(fā)育過程中的差異基因及通路，為鵝品種優(yōu)良性狀的改良提供理論依據(jù)。

1 材料與方法

1.1 實驗材料 實驗選擇經(jīng)產(chǎn)處于產(chǎn)蛋高峰期籽鵝3只（由吉林農(nóng)業(yè)大學鵝研中心種鵝繁殖基地提供）。實驗藥品及儀器有TRIzol?Reagent，RNA Purification Reagent，磁力架，TBS380 Picogreen，由Invitrogen公司生產(chǎn)；TruseqTMRNA sample prep Kit，cBot Truseq PE Cluster Kit v3-cBot-HS，Hisseq2000 Truseq SBS Kit v3-HS，Miseq Reagent Kit V2，由Illumina生產(chǎn)；Certifed Low Range UItra Agarose 由Bio-Rad生產(chǎn)。

1.2 實驗方法

1..2.1 實驗樣品的采集及等級劃分 將籽鵝用乙醚麻醉后放血處死，迅速將卵泡取出并用預冷生理鹽水清洗干凈，同時根據(jù)卵泡直徑大小對其進行等級劃分。同時去除卵泡外層血管膜，并釋放卵泡液。清洗干凈后迅速放入液氮中保存。

1.2.2 卵泡等級劃分 鵝卵泡可根據(jù)卵泡的發(fā)育將其分為等級卵泡和等級前卵泡。本實驗所需卵泡為等級卵泡，等級卵泡一般由大到小分為F1、F2、F3、F4、F5、F6共6 個等級。等級前卵泡根據(jù)其直徑大小可分為初級卵泡（Primary Follicles，PE）直徑小于2 mm、小白卵泡（Small White Follicles，S）直徑在2～4 mm、中白卵泡（Middle White Follicles，M）直徑在4～6 mm、大白卵泡（Large White Follicles，L）直徑在6～8 mm、小黃卵泡（Small Yellow Follicles，SY）直徑在8～10 mm[7]。

1.2.3 總RNA提取及cDNA文庫的建立 利用Trizol法將實驗所采集的卵泡組織進行總RNA提取，通過Nanodrop 2000 對所提RNA 的濃度和純度進行檢測，瓊脂糖凝膠電泳檢測RNA 完整性，Agilent2100 測定RIN 值。將所提取的RNA通過Oligo（dT）磁珠富集mRNA用于分析轉(zhuǎn)錄組信息，利用金屬離子將富集mRNA隨機斷裂成200 bp左右的小片段。在逆轉(zhuǎn)錄酶的作用下，以mRNA 為模板，利用隨機六堿基引物，反轉(zhuǎn)錄合成cDNA。加入End Repair Mix對雙鏈cDNA粘性末端進行修復，PCR擴增15個循環(huán)后，進行目的條帶回收，TBS380（Picogreen）定量，在cBot上進行橋式PCR擴增，生成Clusters后，利用Hiseq2500

1.3 數(shù)據(jù)分析

1.3.1 原始數(shù)據(jù)質(zhì)控及轉(zhuǎn)錄組拼接組裝 通過Illumina 測序平臺對卵泡組織進行測序后，因得到的數(shù)據(jù)量十分龐大，為宏觀上直接反映樣本的測序質(zhì)量和文庫的構建質(zhì)量，對所有的測序reads 進行堿基分布和質(zhì)量波動統(tǒng)計，主要包括堿基含量分布統(tǒng)計、堿基質(zhì)量分布統(tǒng)計和堿基錯誤率分布統(tǒng)計。同時由于鵝屬于無參考基因組動物，針對無參考基因組的轉(zhuǎn)錄組研究，將測序數(shù)據(jù)進行質(zhì)量剪切，將測序數(shù)據(jù)中測序接頭序列末端質(zhì)量低于20 的堿基、含N比率超過10 %和質(zhì)量修剪后長度小于20 bp的序列去除，獲得RNA-seq高質(zhì)量測序數(shù)據(jù)。得到高質(zhì)量的測序數(shù)據(jù)后將所有測序讀段通過Trinity軟件進行從頭組裝，并對組裝后得到的數(shù)據(jù)進行統(tǒng)計分析[8]。

1.3.2 序列分析 序列分析可對組裝后得到的轉(zhuǎn)錄本進行功能注釋，來探究這些轉(zhuǎn)錄本的生物功能，在注釋前需利用Trinity軟件提供的ORF預測流程對組裝得到所有轉(zhuǎn)錄本序列進行核苷酸序列和蛋白質(zhì)序列預測，根據(jù)真核生物默認序列含有非編碼區(qū)原理，利用Markov模型預測該序列中的最佳ORF區(qū)域，并使用Pfam數(shù)據(jù)庫對預測結果進行校正，保留比對到數(shù)據(jù)可預測出ORF的核苷酸序列和蛋白序列，之后將拼接所得到的所有核苷酸序列別于NR、GO、COG、KEGG數(shù)據(jù)庫進行比對，獲得相應的注釋信息對其進行功能注釋[9-10]。

1.3.3 差異因的篩選 為實現(xiàn)所得到的基因表達量可視化，通過RSEM對不同樣品基因的表達量進行計算，并以FPKM值作為衡量基因表達水平的標準，使用EdgeR軟件進行差異基因計算[11]，本次實驗中，差異表達基因的篩選標準為FDR＜0.05，︱logFC︱≥1，同時將所得到的差異基因GO數(shù)據(jù)庫和KEGG數(shù)據(jù)庫富集進行注釋與統(tǒng)計分析。

2 結果與分析

2.1 原始數(shù)據(jù)分析 實驗采用Illumina Hiseq 測序平臺對鵝等級前卵泡組織完成轉(zhuǎn)錄組測序，并構建Illumina PE文庫進行2×101 bp測序，為了保證后續(xù)試驗質(zhì)量，對獲得的測序數(shù)據(jù)進行RNA質(zhì)量驗證并對原始測序數(shù)據(jù)質(zhì)控。

2.1.1 RNA質(zhì)量驗證 將從鵝等級前卵泡組織中提取的總RNA進行質(zhì)量驗證，其中由SY到PE各等級RNA 量為17.7、23.3、19.6、17.6、8.2 μg，且RNA濃度大于200 ng/μL，符合單次建庫RNA質(zhì)量要求。經(jīng)瓊脂糖電泳檢驗，5S、18S、28S條帶清晰（圖 1）。驗證結果表明RNA質(zhì)量良好，未發(fā)生降解，符合下一步建庫要求。

圖1 鵝5個等級前卵泡電泳結果圖

2.1.2 原始測序數(shù)據(jù)質(zhì)控 經(jīng)測序共得到2億多條序列，為宏觀反映出每條序列的質(zhì)量，通過統(tǒng)計學的方法對A/T/G/C堿基含量、堿基質(zhì)量、堿基錯誤率分布統(tǒng)計。經(jīng)堿基含量分布統(tǒng)計結果表明，在整個測序過程中，AT、GC含量相同，且除序列起始和測序接頭位置由于6 bp隨機引物引起A、T、C、G 出現(xiàn)正常波動外，其余呈現(xiàn)穩(wěn)定狀態(tài)，且模糊堿基N比率在合理范圍內(nèi)。堿基質(zhì)量檢測是對所有測序循環(huán)下的測序讀段的質(zhì)量波動進行統(tǒng)計分析，對于RNA-seq技術，根據(jù)下述公式進行計算質(zhì)量值Q。通過統(tǒng)計分析發(fā)現(xiàn)隨測序進行，Q值隨著測序試劑不斷消耗而下降，這屬正?，F(xiàn)象。且整體測序錯誤率較低，讀段堿基質(zhì)量較高，所獲得測序數(shù)據(jù)達到后續(xù)分析要求。

通過堿基錯誤率分布分析，發(fā)現(xiàn)堿基錯誤率會隨著測序序列長度增加和測序化學試劑消耗而增加，且建庫過程中反轉(zhuǎn)錄所使用隨機6堿基引物與RNA模板結合不完全導致測序前6個堿基位置錯誤率較高。但所有測序單個堿基位置平均錯誤率低于0.1%，符合實驗要求，可用于后續(xù)實驗。為保證后續(xù)的信息分析的準確性，對進行測序的數(shù)據(jù)進行質(zhì)量剪切，原始測序數(shù)據(jù)過濾，經(jīng)質(zhì)量剪修后統(tǒng)計共得到224 929 214條序列。具體結果見表 1。

表1 不同胎次繁殖母豬預定的淘汰率和產(chǎn)活仔數(shù)

2.2 測序數(shù)據(jù)分析 鵝屬于無參考基因組動物，通過Trinity軟件將所有測序讀段從頭組裝成重疊峰（Conting）和單一序列（Singleton），得到 277 077 條Unigenes序列和314 320條Transcripts序列，并對其長度分布進行統(tǒng)計，其中長度在1～600的Unigenes序列和Transcripts序列所占的總體比例較高，具體見表2。

表2 拼裝結果統(tǒng)計

2.3 序列分析

2.3.1 NR庫注釋結果 為對所組裝的轉(zhuǎn)錄本進行功能注釋，注釋前需先利用Trinity軟件對組裝得到的所有轉(zhuǎn)錄本進行基因預測，經(jīng)預測共得到74 405條核苷酸序列。將預測所得的核苷酸序列與NR數(shù)據(jù)庫進行比對，可得到本物種轉(zhuǎn)錄組序列與相近物種的相似情況及同源序列的功能信息。NR數(shù)據(jù)庫為NCBI非冗余蛋白庫，是一包含SwissProt、PIR、PRF、PDB蛋白質(zhì)數(shù)據(jù)庫中非冗余的數(shù)據(jù)以及從GenBank和RefSeq的CDS數(shù)據(jù)庫中翻譯蛋白質(zhì)數(shù)據(jù)的綜合數(shù)據(jù)庫。通過比對后統(tǒng)計共有42 453條序列注釋到NR數(shù)據(jù)庫，且進行物種同源性比對，具體結果見圖 2。圖2反映了測序樣本比對近源物種的情況，圖中每一扇形表示一個物種，扇形面積越大表示比對上序列越多，同源性越高。其中綠頭鴨比對序列的數(shù)目最多，因此可認定綠頭鴨與鵝的同源性最高。

圖2 對比物種分布

圖3 GO二級分類統(tǒng)計圖

2.3.2 GO 注釋結果 GO數(shù)據(jù)庫是基因本體論聯(lián)合會對不同數(shù)據(jù)庫中關于基因和基因產(chǎn)物的生物學術語進行標準化及基因和蛋白功能進行統(tǒng)一限定和描述的綜合數(shù)據(jù)庫。將所得Unigenes通過GO數(shù)據(jù)庫按照其參與的生物過程（Biological Process，BP），分子功能（Molecular Function，MF）和細胞組分（Cellular Component，CC）進行分類注釋，其三大分支下又有很多小的分級，級別數(shù)越多，功能注釋越詳細。經(jīng)比對共有1 045條Unigenes對比到GO數(shù)據(jù)上，其中在BP分類中細胞過程、代謝過程和單有機體過程類別所占比例最多；在CC分類中細胞、細胞部分分類和細胞器所占比例最多；在MF分類中，分子綁定功能和催化活性功能類別所占比例最多，具體結果見圖3。

2.3.3 GOG注釋結果 GOG數(shù)據(jù)庫是蛋白直系同源簇數(shù)據(jù)庫，通過比對可以進行功能注釋，歸類以及蛋白進化分析，從而可對轉(zhuǎn)錄本進行功能分類，其中功能預測所占的轉(zhuǎn)錄本數(shù)目最多，其次為信號轉(zhuǎn)導機制和翻譯后修飾。具體見表3。

2.3.4 KEGG通路注釋結果 與基因組破譯公共數(shù)據(jù)庫（KEGG）進行比對分析，可通過KEGG數(shù)據(jù)庫中豐富的通路信息系統(tǒng)的了解與分析基因的生物學功能。通過分析結果發(fā)現(xiàn)，共有11 949 條Unigenes 比對到KEGG數(shù)據(jù)庫，共涉及到341條通路，其中富集最多的通路為癌癥通路，其次為PI3K-AKT信號通路和MAPK信號通路，說明這3條通路在鵝等級卵泡組織中最為活躍，具體見圖4。與KEGG數(shù)據(jù)對比可獲得轉(zhuǎn)錄本對應的KO編號，根據(jù)KO編號可知其可能參與的具體生物學通路，可根據(jù)其KO注釋對其所參與的KEGG通路分類為細胞過程、環(huán)境信息處理、遺傳信息處理、代謝和有機系統(tǒng)5個分支，具體見圖 5。

2.4 基因的表達 通過計算可得到鵝等級前卵泡組織的FPKM值，其中各等級前卵泡中FPKM值大于1的基因占20 % 以上，大于10的基因占3 % 以上（表4）。同時通過FPKM值對所有基因的表達量概率密度做圖，其中橫坐標為log10FPKM，此數(shù)值越高說明基因的表達量越高，縱坐標對應基因的密度，密度曲線的峰值為整個樣本中基因表達量最集中的位置，具體結果見圖 6。

2.5 差異表達基因的分析

2.5.1 差異基因的篩選與注釋分析 通過EdgeR軟件對鵝等級前卵泡中的差異基因計算，篩選標準為FDR＜0.05，︱logFC︱≥1，即基因在2個樣本間的差異顯著性結果校驗結果小與0.05，且樣本間差異倍數(shù)以2為底的對數(shù)絕對數(shù)值大于1可認為此基因在這兩樣本間差異表達。SY與PE等級卵泡均為鵝等級卵泡，SY等級卵泡為等級前卵泡發(fā)育的等級卵泡中最后1個等級，在眾多的SY等級卵泡中只有1個可以發(fā)育為等級前卵泡，而PE等級卵泡是等級卵泡中的初始等級卵泡，通過研究分析發(fā)現(xiàn)這兩等級卵泡中所存在的差異基因數(shù)量最多，其中以SY等級卵泡為參照，在PE 等級卵泡組織中共獲得15 516條差異基因，其中上調(diào)基因8 293條，下調(diào)基因7 223條，具體見圖 7 。

表3 GOG數(shù)據(jù)庫注釋結果

圖4 KEGG pathway中前20個信號通路

圖5 KEGG注釋統(tǒng)計圖

表4 鵝等級前卵泡組織FPKM值

圖6 表達量分布

圖7 PE-SY 卵泡組織存在的差異表達基因可視化圖

圖8 PE-SY 卵泡組織、上下調(diào)基因注釋柱形圖

通過與GO數(shù)據(jù)庫比對篩選得到的15 516條差異基因進行分子功能分類，共比對注釋到3 626條基因，其中上調(diào)基因1 450條，下調(diào)節(jié)基因2 176條，共注釋到參與的生物過程，細胞組分和實現(xiàn)分子功能3大類的22個小類上，其中生物過程24類，細胞組分17類，分子功能14類，具體結果見圖 8。

通過KEGG數(shù)據(jù)庫比對分析，共篩選出4 950條差異基因注釋到299條KEGG通路上。為了能夠識別與生物現(xiàn)象最相關的生物過程，提高研究的可靠性，使用KOBAS（http://kobas.cbi.puk.edu.cn/ home.do）對所得到的KEGG PATHWAY富集分析，經(jīng)計算校正后P≥0.05的KEGG通路定義為在差異表達基因中顯著富集的KEGG通路。其中顯著富集通路有174條。其中顯著富集的通路有PI3K-AKT通路，癌癥通路和黏著斑通路，具體結果見圖 9。

圖9 PE-SY 卵泡組織KEGG通路富集圖

3 討 論

鵝卵泡發(fā)育過程是較為復雜的生物過程[12]。為探究這一過程，通過RNA-Seq技術對產(chǎn)蛋高峰期鵝各等級前卵泡組織進行轉(zhuǎn)錄組測序分析，RNASeq技術作為第2代高通量測序技術，對揭示研究動植物差異基因表達及調(diào)控發(fā)揮了重要作用[13]。本實驗通過Illumina HiSeq 2500測序平臺對鵝等級前卵泡組織進行轉(zhuǎn)錄組文庫的構建，并獲得大量轉(zhuǎn)錄組信息，同時通過NR、GO和KEGG數(shù)據(jù)庫進行比對，進行功能注釋及分類。共得到224 929 214條序列，Q20比例為95.5%以上，GC含量為47%以上，說明文庫構建成功且質(zhì)量較好。經(jīng)dnovo組裝拼接后共得到114 486 條Unigenes序列，并預測出74 405條核苷酸序列，說明測得到基因74 405個，經(jīng)NR數(shù)據(jù)庫對比發(fā)現(xiàn)，綠頭鴨、原雞和原鴿與其同源性較高。通過KEGG通路比對分析，其中最為活躍的通路為癌癥通路，PI3K-AKT通路和MAPK通路，而PI3K-AKT通路和MAPK通路在卵泡發(fā)育過程中的作用及功能以廣泛的被驗證及證實。PI3K-AKT通路為動物體內(nèi)一基本信號通路，可靈活的參與體內(nèi)卵母細胞增長和原始卵泡發(fā)育過程，同時還涉及顆粒細胞增殖、分化和凋亡，繼而調(diào)節(jié)卵泡發(fā)育。MAPK通路在動物的胚胎發(fā)育和細胞增殖、分化、凋亡中同樣也具有重要的調(diào)節(jié)作用[14-15]。由此可以說明，卵泡的發(fā)育是一個動態(tài)的空間調(diào)節(jié)過程，并不受單一通路影響。同時對SY與PE等級卵泡組織間對比共獲得15 516條差異基因，其中上調(diào)基因8 293條，下調(diào)基因7 223條。通過GO數(shù)據(jù)庫比對發(fā)現(xiàn)差異基因主要集中在生物過程類中。通過KEGG富集分析，富集最多的通路為PI3K-AKT通路、癌癥通路和黏著斑通路，這與蘭道亮等[16]在對牦牛卵巢組織做轉(zhuǎn)錄組測序數(shù)據(jù)分析時所得結果相同，可表明這3條通路在卵泡發(fā)育過程中的活躍狀態(tài)和對卵泡發(fā)育中所起到的作用。目前研究還尚未發(fā)現(xiàn)有對鵝等級前卵泡進行等級分類測序研究實例，本實驗極大豐富了鵝等級前卵泡發(fā)育時期組織的基因數(shù)據(jù)，同時通過各等級間卵泡組織所得差異表達基因數(shù)據(jù)分析，可進一步對鵝等級前卵泡發(fā)育關鍵基因進行發(fā)掘。

4 結 論

本實驗通過RNA-Seq技術首次對鵝等級前卵泡進行高通量轉(zhuǎn)錄組測序和分析，研究了鵝不同等級下的卵泡組織的差異表達基因，同時對所獲得差異表達基因進行功能分類和代謝通路劃分，不僅豐富了鵝卵泡組織的轉(zhuǎn)錄信息，還為以后開展鵝遺傳育種研究及分子調(diào)控研究提供理論基礎。

[1] 吳陽升, 林嘉鵬, 汪立芹, 等. 綿羊小卵泡與中卵泡轉(zhuǎn)錄組差異特征分析[J]. 江蘇農(nóng)業(yè)學報, 2016, (4):832-842.

[2] 朱志明, 陳紅萍, 林如龍, 等. 山麻鴨開產(chǎn)期和產(chǎn)蛋高峰期卵巢組織轉(zhuǎn)錄組分析[J]. 中國農(nóng)業(yè)科學, 2016, (5):998-1007.

[3] 韓昆鵬, 段煉, 李婷婷, 等. 京海黃雞卵巢轉(zhuǎn)錄組研究:基因結構分析與新基因發(fā)掘注釋[J]. 中國畜牧獸醫(yī), 2016, (4):854-861.

[4] Tariq M, Chen R, Yuan H Y, et al. De novo transcriptomic analysis of peripheral blood lymphocytes from the Chinese goose:gene discovery and immune system pathway description[J]. PLoS One, 2015, 10(3): e0121015.

[5] 王丹. 維生素A對朗德鵝肥肝抗氧化性能和血液生化指標影響的研究[D]. 石河子:石河子大學, 2011.

[6] 吳慧英, 賈汝敏, 劉鈾, 等. 鵝繁殖就巢機理研究進展[J].中國草食動物, 2009, (5):56-58.

[7] Johnson A. Ovarian follicle selection and granulosa cell differentiation [J]. Poult Sci, 2015, 94(4): 781-785.

[8] Grabherr M G, Haas B J, Yassour M, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome[J]. Nat Biotechnol, 2011, 29(7):644-652.

[9] Camacho C, Coulouris G, Avagyan V, et al. Blast+: architecture and applications[J]. Bmc Bioinformatics, 2009, 10(1):421.

[10] Conesa A, Cotz S, Terol J, et al. A universal tool for annotation, visualization and analysis in functional genomicsresearch.Bioinformatics[J]. 2005, 21:3674-3676.

[11] Robinson M D, Mccarthy D J, Smyth G K. A Bioconductor package for differential expression analysis of digital gene expression data[J]. Bioinformatics, 2010, 26(1):139-140.

[12] 張林媛. 鴻雁卵泡的組織學研究[D].哈爾濱:東北林業(yè)學, 2007.

[13] Zeng T, Zhang L P, Li J J. De novo assembly and characterization of muscovy duck liver transcriptome and analysis of differentilly regulated genes in response to beat stress[J].Cell Stress Chaperon, 2015,( 2):483-493.

[14] 甄艷紅. MAPK-CEBPβ對豬顆粒細胞凋亡、周期及類固醇激素分泌的調(diào)控作用和機制的研究[D]. 武漢:華中農(nóng)業(yè)大學, 2014.

[15] 吳艷青, 陳麗云, 張正紅, 等. 磷脂酰肌醇-3激酶/蛋白激酶B/骨形態(tài)發(fā)生蛋白-15通路與哺乳動物卵巢卵泡發(fā)育[J]. 中國醫(yī)學科學院學報, 2013, (2):224-228.

[16] 蘭道亮, 熊顯榮, 位艷麗, 等. 基于RNA-Seq高通量測序技術的牦牛卵巢轉(zhuǎn)錄組研究:進一步完善牦?；蚪Y構及挖掘與繁殖相關新基因[J]. 中國科學:生命科學, 2014, (3):307-317.

S835.2

：A

：10.19556/j.0258-7033.2017-06-063

2016-11-29；

2016-12-18

吉林農(nóng)業(yè)大學青年拔尖人才支持計劃項目（2016）；吉林省留學回國人員擇優(yōu)資助項目（2015）；吉林省科技廳重點科技攻關項目（20150204025NY）；國家科技部星火計劃項目（2015GA660001）；吉林省現(xiàn)代家禽產(chǎn)業(yè)技術體系項目（2016）作者簡介：孫永峰（1977-），男，吉林長春市人，博士，副教授，主要從事家禽遺傳育種工作研究，E-mail：sunyongfeng1977@

126.com

* 通訊作者