花生轉錄組denovo測序揭示果針發育相關的功能基因

　　近日，由廣東省農業科學院作物研究所梁炫強研究員、陳小平博士，山東省農業科學院花生研究所所長禹山林研究員及國際半干旱熱帶作物研究所(ICRISAT)基因組中心負責人Rajeev Kumar Varshney等完成的花生轉錄組de novo測序的結果發表于Plant Biotechnology Journal。該雜志創刊于2003年，影響因子5.44，是專注于應用植物科學的期刊。本研究是中國花生基因組計劃取得的又一階段性成果，測序和部分信息分析工作均在Macrogen千年基因完成。

　　花生(Arachis hypogaea L.)是我國重要的油料和經濟作物，因其“地上開花，地下結莢”的特性被稱之為“落花生”。花生受精后子房柄(又稱果針)開始伸長并啟動向地性生長，從而將子房帶入地下，到達一定深度后子房開始膨大并發育成莢果。研究表明，黑暗和機械刺激是果針發育的必要條件，而溫度、水分、空氣、營養等外界因素也影響果針發育的過程。因此，果針能否入土關系到花生能否正常結莢。但由于目前尚未完成花生全基因組圖譜的繪制，關于果針發育的分子機制目前尚不明確，本研究將從轉錄組水平闡明影響果針發育的關鍵基因。

　　本研究應用Roche 454(GS FLX Titanium)測序平臺，分別對花生的地上果針(AP)和2個不同發育時期的地下果針(SP1和SP2)進行1/2 run轉錄組de novo測序，得到274M(704,738條reads)、290M(711,496條reads)和238M(609,841條reads)原始數據(raw data)，測序reads的平均讀長為396bp，去除接頭等序列和低質量值reads后獲得83%的有效數據(clean data)。相關測序結果已提交至NCBI的SRA數據庫，登錄號為SRA053198。

　　AP、SP1和SP2的測序數據合并后以 Newbler(v2.6)和TGICL2.0分別進行組裝，并比較兩個軟件的組裝性能。組裝結果與現有的花生轉錄組數據合并去冗余后構建花生轉錄組參考序列，序列大小114.9M，約是花生基因組大小的4.1%。組裝結果分別進行功能基因預測、GO(gene ontology)分類、KEGG(Kyoto Encyclopaedia of Genes and Genomes)代謝通路注釋和基因差異表達分析，結果表明在2194個差異表達基因中，地上果針和地下果針分別有859和1068個轉錄本的表達顯著性上調。在地上果針中，與光合作用和衰老相關的基因表達顯著上調，這2個基因很可能阻止了果針的膨大從而導致敗育。本研究的結果為探索果針發育的分子機制提供了充分的理論依據，并為花生轉錄組和功能基因組的研究提供了更多的數據資源。

　　文獻解析：

　　1、 實驗材料

　　本研究選取的花生栽培種為“粵油 7 號”，該品種是由廣東省農業科學院作物研究所選育的高產花生品種，目前已在南方大面積推廣種植。在花生受精完成后，通過人工干預的方法使部分果針未能伸入土壤中，并分別取開花后第8、10、12、16、20天的地上果針樣品混合后成為一個AP樣本。另外，對于正常發育的果針，分別取果針伸入土壤后第2、 4、6、8、10、12天的地下果針樣品，并將第2、4、6天的樣本混合成為SP1，第8、10、12天的樣本混合成為SP2。所有樣本均以改良的 CTAB法抽提總RNA，并以安捷倫2100檢測RNA的質量。

　　2、 454測序

　　RNA樣品的文庫構建和測序均在Macrogen總部最高質量標準的基因組學實驗室完成，實驗試劑均采用進口原廠試劑，實驗流程和質量控制均嚴格受Macrogen董事會監管。AP、SP1和SP2樣本均進行1/2run的454測序，分別得到 274M(704,738條reads)、290M(711,496條reads)和238M(609,841條reads)原始數據，測序reads的平均讀長為396bp，且大部分reads分布于500bp左右。將測序得到的原始數據首先通過SeqClean軟件與UniVec數據庫進行比對后，去除接頭等序列，然后設置嚴格的過濾標準去除低質量reads。如果一條read中30%的堿基質量值小于Q20或者該read長度小于50bp，則會被過濾掉。經過以上數據處理，共得到1683,455條有效reads，有效數據(clean data)占原始數據(raw data)的比例為83%。相關測序結果已提交至NCBI的SRA數據庫，登錄號為SRA053198。

　　3、 轉錄組de novo組裝

　　目前用于組裝的軟件種類較多，研究表明對于454測序數據的組裝，Newbler和TGICL的組裝性能高于其他軟件。為了獲得最佳的組裝效果，本研究先應用RepeatMasker軟件屏蔽簡單重復等低復雜性序列，然后運用Newbler(v2.6)和TGICL2.0分別進行組裝，并比較這兩個軟件對于花生轉錄組de novo組裝的性能。其中，Newbler是Roche 454平臺自帶的組裝軟件，本研究開展時已更新至v2.6版本，現已更新至v2.8版本，TGICL是專門用于轉錄組組裝的軟件。Newbler組裝得到近3萬個contigs，其中63%的contigs大于1kb，contig N50為1730bp，平均contig長度為1447bp。TGICL組裝得到近5萬個contigs，其中28%的contigs大于 1kb，contig N50為974bp，平均contig長度為859bp。另外，為了進一步驗證組裝效果，本研究以大豆基因組序列預測得到的蛋白質序列作為參考序列進行同源性比對，結果表明TGICL的組裝結果與大豆蛋白質序列的同源性更高。當單個contig與大豆蛋白序列比對的同源性設置為R50%時，Newbler 的組裝結果可比對到13,000個大豆蛋白，而TGICL可比對到15,000個大豆蛋白。因此，Newbler和TGICL在花生轉錄組數據的組裝過程中各有優勢，TGICL能夠得到更多的contigs數量，因而能夠反應轉錄本的多態性，而Newbler能夠得到更大的contigs長度。

4、 轉錄組參考序列的構建

將本研究的組裝結果與現有的花生轉錄組數據合并去冗余后，得到151,533個轉錄本，共計114.9M，約是花生基因組大小的4.1%。其平均contig長度為758bp，contig N50為866bp，平均GC含量為40.5%。現有版本的花生轉錄組組裝結果NCBI Peanut Unigene Build #3和PeanutDB分別是34.41和33.97M，將本研究的轉錄本分別與這兩個版本比對得知，共有5,286個轉錄本無法比對至參考序列，這些很可能是本研究發現的新轉錄本。

表3. 現有版本的花生轉錄組組裝結果比較

　　應用SSAHA2 v2.5.3將AP、SP1和SP2樣本測序得到的reads分別與構建的轉錄組參考序列比對，約99%的reads能夠匹配至參考序列，這表明該參考序列能夠較全面的體現三個樣本的轉錄情況。其中，三個樣本中平均比對至每個轉錄本的reads數為10.5、10.9和10.1，分別有73.79%、 76.30%和75.17%的reads與參考序列有唯一匹配，約23%-26%的reads有多個匹配，這部分reads在后續的分析中均被過濾掉了。根據比對結果得知，共鑒定出74,974個轉錄本，其中AP、SP1和SP2樣本中分別鑒定出49,632、49,952和50,494個轉錄本。雖然三個測序樣本分別取自不同的發育時期，但27,732個轉錄本是三個樣本共有的，且AP中約84%的轉錄本能夠在SP1或者SP2中檢測到。

　　通過RPKM值衡量三個樣本的基因表達水平，結果表明AP、SP1和SP2的平均表達量分別為RPKM15.2、14.3和15.6，大部分基因是低表達和中等表達水平(RPKM值介于3和50之間)，只有1.2%-1.4%的基因是高表達水平的(RPKM值大于100)。

　　5、 基因功能的注釋

　　為了鑒定相關的功能基因，將三個樣本中共74,974個轉錄本利用BLASTX與UniProtKB進行序列同源性比較，E值設定為1E-5。其中，59,342個轉錄本表現出較高的同源性，共匹配了28,633個功能蛋白和大量未知功能的蛋白。另外，約20%的轉錄本未能與數據庫中的任何蛋白匹配，很可能由于這部分轉錄本是花生特有的轉錄本。

　　對于已通過UniProtKB注釋的轉錄本，利用GO對已注釋基因的功能進行分類，共有42,995個轉錄本可比對至GO數據庫，且38,280個轉錄本與分子功能相關，31,152個轉錄本與生物過程相關，17,881個轉錄本與細胞組成相關。另外，分子功能分類中與代謝活性和結合功能相關的轉錄本所占比例最大，約20%以上。在生物過程分類中，細胞過程和代謝過程分別占13.39%和6.33%。在細胞組成分類中，11.42%和4.53%的定位在細胞內和細胞膜，只有很少一部分位于細胞外基質、細胞連接區。

　　為了進一步了解轉錄本的生物學功能，進行KEGG代謝通路注釋。利用KOBAS v2.0將蛋白數據庫注釋的轉錄本比對到KO數據庫，結果表明共有2,968個轉錄本比對到174個KEGG代謝通路中，其中168個KEGG通路是地上果針(AP)和地下果針(SP)共有的，這表明兩種樣本表達的大部分基因是相同的。

　　6、 差異表達分析

　　利用DEGseq和 RPKM分別對差異表達基因進行分析，在地上果針和地下果針樣本中分別鑒定出859和1068個表達上調的基因。其中，地上果針中與光合作用和衰老相關的基因表達顯著上調，這2個基因很可能阻止了果針的膨大從而導致敗育。另外，差異表達分析的結果也得到了RT-PCR的驗證。