虹鱒肝組織新轉(zhuǎn)錄本分析及基因結(jié)構(gòu)優(yōu)化

馬芳，劉哲，康玉軍，權(quán)金強(qiáng)

(甘肅農(nóng)業(yè)大學(xué)動(dòng)物科學(xué)技術(shù)學(xué)院，蘭州 730070)

轉(zhuǎn)錄組是指細(xì)胞在特定階段產(chǎn)生的全部轉(zhuǎn)錄本，包括mRNA、rRNA、sRNA和tRNA[1]。轉(zhuǎn)錄組學(xué)已經(jīng)被廣泛用來研究生物體對(duì)環(huán)境的各種復(fù)雜反應(yīng)，解釋基因組的功能元素。隨著RNA-seq變得越來越便宜，常常成為研究環(huán)境壓力的方法。豐富的RNA-seq數(shù)據(jù)可以構(gòu)建完整的轉(zhuǎn)錄組，提供豐富的差異基因表達(dá)信息，并可用于識(shí)別涉及熱應(yīng)激反應(yīng)的生物通路。在魚類中，通過利用RNA-seq識(shí)別不同魚類溫度適應(yīng)機(jī)制的研究正在迅速增加。但是，現(xiàn)有數(shù)據(jù)庫(kù)中對(duì)轉(zhuǎn)錄本的注釋還不全面，通過RNA-seq技術(shù)，還能檢測(cè)到新轉(zhuǎn)錄本。目前，越來越多的研究開始關(guān)注RNA-seq技術(shù)在新轉(zhuǎn)錄本預(yù)測(cè)和基因結(jié)構(gòu)優(yōu)化方面的應(yīng)用[2-3]。

虹鱒(Oncorhynchusmykiss)作為鮭科魚類的一員，正迅速的成為水產(chǎn)養(yǎng)殖中的重要魚類。作為典型的冷水魚，最適的生活溫度是12 ～ 18℃。對(duì)于高溫的耐受性低，當(dāng)溫度超過24℃時(shí)免疫功能嚴(yán)重下降，組織受損[4]。以前的研究利用微陣列技術(shù)驗(yàn)證了虹鱒對(duì)溫度變化的反應(yīng)[5]，并對(duì)虹鱒不同種類的熱應(yīng)激反應(yīng)進(jìn)行了比較[6]。

本課題組前期應(yīng)用RNA-seq技術(shù)對(duì)虹鱒熱應(yīng)激下肝組織中差異表達(dá)基因進(jìn)行了鑒定[7]，本研究中，在前期研究的基礎(chǔ)上運(yùn)用生物信息學(xué)方法鑒定新的轉(zhuǎn)錄本并對(duì)已注釋基因的結(jié)構(gòu)進(jìn)行優(yōu)化及，為深入理解虹鱒熱應(yīng)激的機(jī)制奠定基礎(chǔ)，同時(shí)為虹鱒基因組的進(jìn)一步完善提供數(shù)據(jù)基礎(chǔ)。

1 材料與方法

1.1 實(shí)驗(yàn)動(dòng)物及飼養(yǎng)

選擇身體健壯，平均體重為(400 ± 10.5)g的全同胞虹鱒200尾運(yùn)送至實(shí)驗(yàn)室，置于一個(gè)3000 L的水箱在18℃下訓(xùn)養(yǎng)7 d。試驗(yàn)前，隨機(jī)挑選120尾分為六組，每組20尾，分別置于6個(gè)300 L的室內(nèi)循環(huán)流水水箱中暫養(yǎng)一周。飼養(yǎng)期間嚴(yán)格按照虹鱒飼養(yǎng)標(biāo)準(zhǔn)飼喂，光照周期為12 h光照和12 h黑暗，嚴(yán)格按照虹鱒飼養(yǎng)標(biāo)準(zhǔn)飼喂。

1.2 熱處理及采樣

暫養(yǎng)結(jié)束后，選3組繼續(xù)18℃飼養(yǎng)，對(duì)其余3組進(jìn)行熱處理升溫，從18℃到24℃以恒定的速率每24 h升高1℃。然后隨機(jī)從各組取1尾魚采取肝組織，18℃作為對(duì)照組，24℃作為熱處理組。采樣時(shí)，用0.05 g/L的間氨基苯甲酸乙酯甲磺酸鹽(MS-222)對(duì)實(shí)驗(yàn)用魚進(jìn)行麻醉，采集肝組織，迅速貯存到液氮中，然后-80℃保存?zhèn)溆谩?/p>

1.3 RNA提取

利用TRIzol試劑盒 (Invitrogen, Carlsbad, CA, USA對(duì)肝組織的總RNA進(jìn)行提取，用NanoPhotometer? spectrophotometer (IMPLEN, CA, USA)和1%瓊脂糖凝膠電泳對(duì)總RNA的純度進(jìn)行檢測(cè)。用Qubit 2.0熒光光度計(jì)(LifeTechnologies, CA, USA)和Bioanalyzer 2100 (Agilent Technologies, CA, USA)檢測(cè)提取的總RNA的濃度和完整性。樣品檢測(cè)合格后進(jìn)行文庫(kù)構(gòu)建。

1.4 文庫(kù)構(gòu)建和轉(zhuǎn)錄組測(cè)序

應(yīng)用the NEBNext? UltraTMRNA Library Prep Kit構(gòu)建6個(gè)測(cè)序文庫(kù)。用帶有Oligo(d T)的磁珠純化mRNA，隨后加入NEBNext First Strand Synthesis Reaction Buffer(5X)使mRNA打斷成短片段，應(yīng)用六聚體引物和M-MuLV Reverse Transcriptase (RNase H-)合成一鏈cDNA，隨后用RNase H和DNA polymerase I合成二鏈cDNA。利用AMPure XP beads純化雙鏈cDNA選擇150～200 bp的cDNA片段。最后通過PCR擴(kuò)增得到cDNA文庫(kù)。構(gòu)建好的文庫(kù)用Agilent Bioanalyzer 2100檢測(cè)合格后，使用Illumina HiseqTM2500平臺(tái)進(jìn)行測(cè)序獲得150 bp的雙末端原始數(shù)據(jù)。

1.5 新轉(zhuǎn)錄本預(yù)測(cè)

原始數(shù)據(jù)(raw reads)去除帶接頭(adapter)的reads和含ploy-N和低質(zhì)量的reads后獲得clean data。同時(shí)計(jì)算clean data的Q20a、Q30和GC含量。后面的所有的分析都基于高質(zhì)量的clean data。應(yīng)用TopHat v2.0.12將clean data比對(duì)到參考基因組，隨后利用Cufflinks v2.1.1軟件對(duì)比對(duì)結(jié)果進(jìn)行組裝，用Cuffcompare和已知的基因注釋文件進(jìn)行比較，尋找潛在的新轉(zhuǎn)錄本。利用GOseq軟件對(duì)新轉(zhuǎn)錄本進(jìn)行GO功能注釋。KOBAS(2.0)軟件對(duì)KEGG注釋通路進(jìn)行分析。

1.6 新轉(zhuǎn)錄本熱應(yīng)激下基因表達(dá)分析

根據(jù)轉(zhuǎn)錄本的長(zhǎng)度和新轉(zhuǎn)錄本的比對(duì)結(jié)果計(jì)算每個(gè)轉(zhuǎn)錄本的FPKM(expected number of fragments per kilobase of transcript sequence per millions base pairs sequenced每百萬(wàn)片段中來自某一基因每千堿基長(zhǎng)度的數(shù)目)作為表達(dá)量的單位，歸一化處理后的數(shù)據(jù)用log2作為新轉(zhuǎn)錄本表達(dá)的數(shù)據(jù)。采用DESeqR package (v1.18.0)對(duì)熱處理組和對(duì)照組之間差異表達(dá)的新轉(zhuǎn)錄本進(jìn)行分析，P< 0.05的轉(zhuǎn)錄本為差異表達(dá)。

1.7 已知基因結(jié)構(gòu)優(yōu)化

組裝的轉(zhuǎn)錄本與虹鱒基因注釋信息進(jìn)行對(duì)比，如果在已注釋基因邊界外的區(qū)域有連續(xù)的匹配讀段，則將基因的5′和3′端進(jìn)行延伸，優(yōu)化已注釋基因的結(jié)構(gòu)。

2 結(jié)果

2.1 測(cè)序結(jié)果數(shù)據(jù)處理與分析

測(cè)序數(shù)據(jù)已提交NCBI數(shù)據(jù)庫(kù)(SRP092649)。6個(gè)文庫(kù)(CL1、CL2、CL3和HL1、HL2、HL3)總共產(chǎn)生了287 277 772條原始數(shù)據(jù)(raw reads)，去除帶接頭的reads，去除N的比例大于10%的reads和去除低質(zhì)量的reads后獲得277 680 702條clean reads。通過對(duì)堿基質(zhì)量進(jìn)行了評(píng)估和堿基組成的檢測(cè)，可以保證下游分析的準(zhǔn)確性。RNA-Seq測(cè)序的堿基質(zhì)量值是堿基識(shí)別出錯(cuò)概率的整數(shù)映射，使用Phred堿基質(zhì)量值公式計(jì)算。堿基質(zhì)量值越高表明堿基識(shí)別準(zhǔn)確度越高，例如堿基質(zhì)量值10 (Q10)、20 (Q20)、30 (Q30)和40 (Q40)分別表示堿基識(shí)別出錯(cuò)的概率為10%、1%、0.1%和0.01%。6個(gè)文庫(kù)中堿基質(zhì)量值≥Q30的堿基百分比分別為90.17%、89.77%、89.62%、91.25%、91.47%和91.30%，說明堿基質(zhì)量較高(表1)。6個(gè)樣品的堿基組成情況如圖1，各個(gè)堿基占的比例約為25%，G和C堿基及A和T堿基含量每個(gè)測(cè)序循環(huán)上分別相等，且整個(gè)測(cè)序過程穩(wěn)定不變，呈水平線，不存在堿基分離現(xiàn)象。利用Top Hat2軟件將clean reads與虹鱒參考基因組進(jìn)行比對(duì)，由表1可知，6個(gè)樣品中clean reads與虹鱒參考基因組進(jìn)行比對(duì)效率在66.17% ～ 68.61%之間，其中有單位點(diǎn)(uniquely mapped)比對(duì)率在64.83% ～ 67.31%之間，多位點(diǎn)(multiple mapped)比對(duì)率在1.23% ～ 2.28%之間，說明測(cè)序數(shù)據(jù)的比對(duì)率正常。

表1 clean data與參考基因組序列比對(duì)結(jié)果Table 1 Comparison of clean data with the reference genome sequences

圖1 原始數(shù)據(jù)堿基組成Figure 1 Base composition of raw data

2.2 新轉(zhuǎn)錄本的發(fā)掘

通過過濾掉少于50個(gè)氨基酸殘基的編碼肽鏈和只包含單個(gè)外顯子的序列，共獲得6555個(gè)新的轉(zhuǎn)錄本(表2)。其中表達(dá)量較低的(≤10)的基因?yàn)?991個(gè)，占30.4%；高表達(dá)的(>1000)的基因?yàn)?33個(gè)，占2.03%(圖2)。新轉(zhuǎn)錄本的長(zhǎng)度大都在500 bp以上，約占62.1%，說明新鑒定的轉(zhuǎn)錄本主要為蛋白質(zhì)編碼基因(圖3)。新轉(zhuǎn)錄本在染色體上的分布如圖4所示，在chrUn染色體上分布最多，有5411條；在染色體chrUn26上最少，有11條。

表2 利用RNA-seq技術(shù)鑒定的虹鱒新轉(zhuǎn)錄本Table 2 Novel transcripts in the rainbow trout identified by RNA-Seq technology

注：部分?jǐn)?shù)據(jù)未列出。

Note. Some data are not listed.

圖2 虹鱒肝新轉(zhuǎn)錄本表達(dá)量統(tǒng)計(jì)Figure 2 Statistical data of the expression of new transcripts in the liver of rainbow trout

圖3 虹鱒肝新轉(zhuǎn)錄本的長(zhǎng)度分布Figure 3 Length distribution of the new transcripts in the liver of rainbow trout

圖4 轉(zhuǎn)錄本在染色體上的分布Figure 4 Distribution of the transcripts on chromosomes

2.3 新轉(zhuǎn)錄本的注釋

利用Blast2Go軟件對(duì)篩選到的新轉(zhuǎn)錄本進(jìn)行GO富集分析，3097個(gè)新轉(zhuǎn)錄本注釋到細(xì)胞組成、生物學(xué)過程和分子功能。在分子功能中，綁定分子功能類別所占比例最多，其次是酶活性活性類別。在生物學(xué)過程中，代謝過程類別所占比例最多，其次是生物合成過程類別。在細(xì)胞組分中，胞外區(qū)類別所占比例最多(圖5)。

圖5 虹鱒肝新轉(zhuǎn)錄本GO注釋結(jié)果Figure 5 GO annotation results for the new transcripts in the liver of rainbow trout

利用KOBAS(2.0)軟件對(duì)KEGG注釋通路進(jìn)行分析，3617個(gè)新轉(zhuǎn)錄本注釋到284條代謝通路。主要的10條代謝途徑見圖6，分別是代謝途徑(metabolic pathways)、粘著斑(focal adhesion)、內(nèi)吞作用(endocytosis)、PI3K-Akt信號(hào)通路(PI3K-Akt signaling pathway)、MAPK信號(hào)通路(MAPK signaling pathway)、胰島素信號(hào)通路(insulin signaling pathway)、碳代謝(carbon metabolism)、Rap1信號(hào)通路(Rap1 signaling pathway)、AMPK信號(hào)通路(AMPK signaling pathway)、細(xì)菌侵入上皮細(xì)胞(bacterial invasion of epithelial cells)。

2.4 新轉(zhuǎn)錄本表達(dá)譜分析

總共有30個(gè)新轉(zhuǎn)錄本在熱應(yīng)激下差異表達(dá)，參與了虹鱒熱應(yīng)激。其中15個(gè)顯著上調(diào)表達(dá)(Novel00236、Novel00736、Novel01309、Novel01495、Novel02292、Novel02550、Novel02698、Novel03125、Novel03334、Novel03377、Novel03766、Novel04249、Novel05645、Novel06326、Novel06367)(P< 0.05)，15個(gè)顯著下調(diào)(Novel00295、Novel00475、Novel00942、Novel01074、Novel01430、Novel03158、Novel03185、Novel03283、Novel03815、Novel04339、Novel05149、Novel05519、Novel05701、Novel06040、Novel06166)(P< 0.05)(圖7)。

圖6 虹鱒肝新轉(zhuǎn)錄本KEGG分析Figure 6 KEGG analysis for new transcripts in the liver of rainbow trout

2.5 已注釋基因結(jié)構(gòu)優(yōu)化

利用RNA-seq測(cè)序結(jié)果對(duì)已注釋基因的結(jié)構(gòu)進(jìn)一步優(yōu)化。如果在已注釋基因邊界之外的區(qū)域有連續(xù)的匹配讀段支持，則將基因的UTR區(qū)域向上游或向下游延伸，優(yōu)化基因邊界。基因結(jié)構(gòu)優(yōu)化結(jié)果顯示，本研究中共有19 424個(gè)已注釋基因5′或3′端在原有基礎(chǔ)上發(fā)生了延伸(表3)。其中5′端為14 719個(gè)延伸，3′端為14 796個(gè)延伸(表4)。

圖7 熱應(yīng)激后虹鱒肝新轉(zhuǎn)錄本的差異表達(dá)Figure 7 Differential expression of new transcripts in the liver of rainbow trout after heat stress

3 討論

虹鱒作為典型的冷水性魚，對(duì)高溫的耐受能力差，隨著全球氣候的變暖，對(duì)虹鱒的養(yǎng)殖造成了越來越嚴(yán)重的影響，因此了解虹鱒熱應(yīng)激下的生存機(jī)制，提高虹鱒的抗逆性至關(guān)重要。目前，越來越多的研究深入的探索魚類熱應(yīng)激的機(jī)制[8-10]。對(duì)于虹鱒，轉(zhuǎn)錄水平上進(jìn)行的熱應(yīng)激研究相對(duì)較少，一些研究采用活體[11-13]，或采用體外細(xì)胞培養(yǎng)的方法[14-15]，初步研究了熱應(yīng)激下虹鱒個(gè)別已知mRNA的表達(dá)水平變化，沒有系統(tǒng)研究熱應(yīng)激調(diào)控機(jī)理。目前發(fā)展的高通量轉(zhuǎn)錄組測(cè)序技術(shù)RNA-seq，在研究基因結(jié)構(gòu)和功能方面具有突出的優(yōu)勢(shì)，通過RNA-seq可以全面快速地獲得某一物種特定組織或器官在某一狀態(tài)下的轉(zhuǎn)錄本信息。

表3 基因3′和5′端延伸情況Table 3 Extension of the 3′ and 5′ ends of genes

表4 部分3′或5′端延伸的基因Table 4 Partial 3′ or 5′ extension genes

隨著轉(zhuǎn)錄組測(cè)序技術(shù)的快速發(fā)展，越來越多的新轉(zhuǎn)錄本被發(fā)現(xiàn)，但是，在現(xiàn)有的數(shù)據(jù)庫(kù)中對(duì)新轉(zhuǎn)錄本的注釋還不全面。豬基因組自基因圖譜公布后，還有不少新的轉(zhuǎn)錄本被發(fā)現(xiàn)[16]。利用RNA-seq技術(shù)對(duì)綿陽(yáng)正常組合骨延遲愈合組進(jìn)行了測(cè)序，發(fā)現(xiàn)了12 431個(gè)新轉(zhuǎn)錄本[3]。

本研究應(yīng)用構(gòu)建虹鱒熱應(yīng)激下的6個(gè)轉(zhuǎn)錄本序列，將虹鱒熱應(yīng)激下肝RNA-seq結(jié)果中的原始數(shù)據(jù)，去除帶接頭的reads，去除N的比例大于10%的reads和去除低質(zhì)量的reads后獲得277 680 702條clean reads。然后對(duì)將clean data比對(duì)到參考基因組，隨后對(duì)比對(duì)結(jié)果進(jìn)行組裝，并與已知的基因注釋文件進(jìn)行比較，尋找潛在的新轉(zhuǎn)錄本。共發(fā)現(xiàn)6555個(gè)新轉(zhuǎn)錄本，表達(dá)量較低的(reads ≤ 10)的基因?yàn)?991個(gè)，說明虹鱒肝中至少表達(dá)了4564個(gè)新轉(zhuǎn)錄本，并且大多數(shù)是高表達(dá)[7]。熱應(yīng)激下新轉(zhuǎn)錄本的表達(dá)譜分析揭示了這些新轉(zhuǎn)錄本也參與了虹鱒抗熱應(yīng)激過程。盡管關(guān)于虹鱒對(duì)熱應(yīng)激脅迫的分子機(jī)制已經(jīng)有很多研究[14, 17-18]，但是還沒有對(duì)這些未知因子在熱應(yīng)激下的作用機(jī)制進(jìn)行深入分析，因此，本研究首次系統(tǒng)的分析了新轉(zhuǎn)錄本在熱應(yīng)激下的調(diào)控規(guī)律和作用。我們發(fā)現(xiàn)在肝中總共有30個(gè)新轉(zhuǎn)錄本受到熱應(yīng)激的調(diào)控，其中15個(gè)表達(dá)顯著上調(diào)，15個(gè)表達(dá)顯著下調(diào)，這給了解虹鱒熱應(yīng)激的機(jī)制進(jìn)一步奠定了基礎(chǔ)。

RNA-seq還在進(jìn)一步完善基因結(jié)構(gòu)信息方面發(fā)揮著重要的作用，將clean data比對(duì)到參考基因組后，發(fā)現(xiàn)共有19 424個(gè)已知基因的5′或3′UTR區(qū)在原有基礎(chǔ)上發(fā)生了不同的延伸。該結(jié)果表明，已知基因的5′或3′UTR區(qū)預(yù)測(cè)不完全，而這些延伸優(yōu)化了已知基因的結(jié)構(gòu)。

本實(shí)驗(yàn)通過對(duì)RNA-seq結(jié)果進(jìn)行生物信息學(xué)分析，發(fā)現(xiàn)了6555個(gè)新轉(zhuǎn)錄本，30個(gè)新轉(zhuǎn)錄本參與了虹鱒抗熱應(yīng)激過程。對(duì)已注釋基因的結(jié)構(gòu)進(jìn)行了優(yōu)化。這些結(jié)果使得虹鱒的全基因組更加全面，也為進(jìn)一步了解虹鱒熱應(yīng)激的機(jī)制提供更有力的理論基礎(chǔ)。