李 青， 何 斌， 陳群利， 游 萍

(貴州工程應用技術學院生態(tài)工程學院，貴州畢節(jié) 551700)

鯽魚()隸屬于鯉形目鯉科鯽屬，由于其營養(yǎng)價值高、味鮮肉嫩、生長速度快、雜食性和適應性強等特點，是我國重要的淡水養(yǎng)殖經(jīng)濟魚類。鯽魚多生活于水域的中下層，冬天低溫條件下可以忍受缺氧數(shù)月，具有極強的低氧耐受性，被公認為是低氧耐受力極強的魚類之一。此外，鯽魚具有二倍體、三倍體和四倍體3種不同的染色體倍性類型，不同類型具有不同的繁殖方式。二倍體通過兩性繁殖方式產(chǎn)生雌雄同體后代；四倍體通過雌核發(fā)育產(chǎn)生全雌后代；而三倍體鯽魚同時存在這2種繁殖方式，即通過精子激活雌核發(fā)育成全雌后代和雌雄交配產(chǎn)生雌雄同體三倍體后代，鯽魚是研究魚類進化發(fā)育基因組學和低氧脅迫適應分子機制的適當模型。

近年來，隨著測序技術的不斷發(fā)展，在研究目標物種全基因組數(shù)據(jù)未知的情況下，轉錄組測序已成為一種短期、低成本得到海量基因數(shù)據(jù)的高效方法?；虻谋磉_具有組織特異性，研究生物過程相關組織中基因的表達模式，可以為進一步闡明其分子機制提供基礎依據(jù)。目前，關于魚類多倍體起源、遺傳多樣性和耐低氧機制尚未完全了解。血液可以反映低氧應激時機體代謝水平和器官機能狀態(tài)的變化，肝臟作為動物體最大的代謝器官，低氧環(huán)境會導致魚類肝臟中氧化氫酶、谷胱甘肽過氧化物酶和超氧化物歧化酶等活性上升。此外，研究表明，低氧可能導致魚類產(chǎn)生強烈的應激反應而影響其生殖、生長和發(fā)育等。

因此，本研究通過RNA-Seq測序技術對鯽魚肝臟、血液和卵巢組織進行轉錄組測序分析，旨在發(fā)掘鯽魚低氧脅迫相關調控功能基因，以期為水產(chǎn)養(yǎng)殖過程中魚類適應低氧脅迫策略提供參考；此外，開發(fā)一批分子標記，以期為鯽魚種質資源評價、分子育種和群體遺傳多樣性分析等方面提供基礎信息。

1 材料與方法

1.1 樣品采集

試驗用鯽魚為個體大且體表無損傷的成魚(體質量為120～150 g，體長為 20～25 cm)，于2019年購自貴州省畢節(jié)市水產(chǎn)品交易市場，實驗室充氣暫養(yǎng)1周后，尾部取血后，分別取肝臟和卵巢組織，經(jīng)液氮速凍后于-80 ℃保存、備用。

1.2 總RNA提取

分別取肝臟、血液和卵巢組織各150 mg，采用TRIzol法提取其總RNA。瓊脂糖凝膠電泳分析RNA降解程度及是否有污染，Nanodrop檢測RNA純度，Qubit2.0對RNA濃度進行精確定量，Agilent 2100精確檢測RNA的完整性。將5尾個體相同組織樣品mRNA等量混合，用于cDNA文庫的構建。

1.3 建庫測序和拼接組裝

使用NEBNext? UltraRNA Library Prep Kit方法構建鯽魚肝臟、血液和卵巢轉錄組文庫，庫檢合格后，進行Illumina HiSeq高通量測序。得到的原始數(shù)據(jù)經(jīng)CASAVA堿基識別分析轉化為原始測序序列(raw data)，對原始測序序列過濾，得到clean reads。采用Trinity對clean reads進行拼接，轉錄本序列信息以FASTA格式儲存。

1.4 基因功能注釋

測序數(shù)據(jù)采用無參考基因組分析方法，為獲得全面的基因功能信息，進行七大數(shù)據(jù)庫的基因功能注釋，包括非冗余蛋白序列數(shù)據(jù)庫(Non-Redundant Protein Sequence Database，NR)，核酸序列數(shù)據(jù)庫(Nucleotide Sequence Database，NT)，蛋白質家族的集合數(shù)據(jù)庫(Protein Family，Pfam)，蛋白相鄰類的聚簇數(shù)據(jù)庫(Cluster of Orthologous Groups of Proteins，COG)，真核生物蛋白質同源簇數(shù)據(jù)庫(euKaryotic Orthologous Groups，KOG)，Swiss-Prot蛋白質序列數(shù)據(jù)庫(Manually Annotated and Reviewed Protein Sequence Database，Swiss-Prot)，京都基因和基因組百科數(shù)據(jù)庫(Kyoto Encyclopedia of Genes and Genomes，KEGG)和基因本體數(shù)據(jù)庫(Gene Ontology，GO)。

1.5 SNP、InDel和SSR分析

使用Samtools和PicardTools工具對比，將結果進行染色體坐標排序，去掉重復的reads，通過變異檢測軟件GATK2分別進行單個核苷酸變異多態(tài)性(single nucleotide polymorphisms，SNP)和核苷酸插入或缺失(insertion-deletion，InDel)，并對原始結果進行過濾，過濾掉質量值小于30，距離小于5的SNP。采用MISA 1.0對Unigene進行重復序列標記(simple sequence repeats，SSR)檢測，并對不同SSR類型在基因轉錄本的密度分布進行統(tǒng)計。

1.6 數(shù)據(jù)處理

采用RESM軟件，以Trinity拼接得到的轉錄組作為參考序列(ref)，將每個樣品clean reads往ref上做mapping，得到每個樣品比對到每個基因上的readcount數(shù)目，使用每百萬fragments中來自某一基因每千堿基長度的fragments數(shù)目(expected number of fragments per kilobase of transcript sequence per millions base pairs sequenced，F(xiàn)PKM)表示基因的表達水平。使用DEGSeq軟件分析不同組織樣品基因表達差異，從而篩選出差異基因，然后使用K-means和SOM進行差異基因聚類分析，獲得3種組織中差異基因表達量的聚類模式。為了更好地對差異基因的功能進行研究，分別使用GOSeq和KOBAS軟件對3個組織間差異基因進行GO富集和KEGG富集分析，同時將差異基因根據(jù)上調或下調分別進行富集分析。

2 結果與分析

2.1 數(shù)據(jù)質量評估與序列組裝分析

經(jīng)過數(shù)據(jù)質量控制，肝臟、卵巢和血液組織分別得到20 332 138、20 209 531、20 451 910個clean reads，平均GC含量為46.21%。將3個組織的reads片段拼接組裝，共得到219 192條轉錄本(transcript)，平均長度為1 088 bp，最長的片段長度為17 286 bp，最短的片段長度為181 bp，N50和N90(將拼接轉錄本按照長度從長到短排序，累加轉錄本的長度，到不小于總長50%～90%的拼接轉錄本的長度)分別為2 096 bp和397 bp。共得到127 801條unigene，平均長度為735 bp，最長的片段長度為17 286 bp，最短的片段長度為201 bp，N50和N90分別為1 288 bp和283 bp。轉錄組拼接組裝的數(shù)據(jù)已提交至BioProject數(shù)據(jù)庫(BioProject ID：PRJNA735422)。

2.2 功能基因注釋與分類

2.2.1 功能注釋與相似性 通過BLAST與相應數(shù)據(jù)庫比對，最終獲得有注釋信息的unigene數(shù)量為117 414個，約占總unigene數(shù)量的91.87%。比對到NT數(shù)據(jù)庫的unigene數(shù)量最多，為116 885條，占總unigene數(shù)量的91.45%，其次，為NR(39.44%)、KO(17.23%)、KOG(16.21%)、Swiss-Prot(6.52%)、GO(0.08%)，而沒有unigene比對到Pfam數(shù)據(jù)庫。在以上7個數(shù)據(jù)庫中至少1個數(shù)據(jù)庫注釋成功的unigene有117 414個，占總unigene數(shù)量的91.87%。

以NR數(shù)據(jù)庫為例，對鯽魚unigene序列的相似性進行分析，在斑馬魚()中相似序列匹配的比例(73.1%)最高，其次為墨西哥麗脂鯉()(52)、虹鱒()(3.4%)、鯉魚()(1.2%)和羅非魚()(1.2%)。匹配序列的相似度都大于40%，其中相似度在80%～95%占的比例(45.1%)最高，其次為相似度介于95%～100%(29.6%)、60%～80%(21.1%)、40%～60%(4.2%)。值在0區(qū)間內的unigene數(shù)量最多，占總體的19.2%，值介于10～10之間的unigene數(shù)量占總體的15.0%，值介于10～10之間的unigene數(shù)量最少，占總體的10.5%。

2.2.2 unigene的GO注釋分類 由圖1可知，鯽魚3個組織的轉錄組共有105條unigene在GO數(shù)據(jù)庫生物學過程(biological process)、細胞組分(cellular component)和分子功能(molecular function)三大類42個功能中找到對應。其中，生物學過程獲得的注釋信息最多，分子功能獲得的注釋信息最少。在生物學過程類中，獲得的注釋信息多集中于細胞學過程、單一有機體過程、生物調節(jié)、代謝過程、生物過程調節(jié)、應激反應和信號等功能。細胞組分類中，獲得的注釋信息多集中于細胞、細胞組件、大分子復合物和細胞器等方面。分子功能類中，在綁定和催化活性中分布較多，在酶調節(jié)器活性中分布最少。GO數(shù)據(jù)庫注釋分類有助于從整體上大致了解目的組織全部基因產(chǎn)物的功能，為進一步探究基因的功能提供參考。

2.2.3 unigene的KOG注釋分類 由圖2可知，KOG數(shù)據(jù)庫可以對基因產(chǎn)物進行直系同源分類，結果顯示共有20 725條unigene得到注釋，依據(jù)功能分為26個直系同源功能類型。其中，參與一般功能預測和信號傳導機制功能的最多，其次為翻譯后修飾、蛋白質折疊與分子伴侶，轉錄，細胞內轉運、分泌與囊泡運輸。

2.2.4 unigene的KEGG注釋 由圖3可知，通過與KEGG數(shù)據(jù)庫比對，鯽魚unigene歸屬于A(細胞過程)、B(環(huán)境信息處理)、C(遺傳信息處理通路)、D(代謝)和E(有機系統(tǒng)通路)五大類；有機系統(tǒng)通路注釋比例最高，而遺傳信息處理通路所占比例最低。其中注釋數(shù)量最多的代謝途徑是信號傳導途徑(3 856個)，注釋數(shù)量大于1 000的為內分泌系統(tǒng)(1 776個)，免疫系統(tǒng)(1 769個)，運輸和分解代謝途徑(1 426 個)，細胞群落(1 356個)，折疊、排序和降解途徑(1 191個)，神經(jīng)系統(tǒng)(1 190個)，信號分子與相互作用途徑(1 086個)。

2.3 不同組織差異表達基因分析

FPKM是每100萬個fragments中來自某一基因每千個堿基長度的fragments數(shù)目，其同時考慮了測序深度和基因長度對fragments計數(shù)的影響，是目前最為常用的基因表達水平估算方法。因此，筆者所在課題組將readcount數(shù)進行了FPKM轉換，以不同樣品作為橫坐標，不同樣品表達量(FPKM+1)的對數(shù)值作為縱坐標，繪制3個組織基因表達量的FPKM箱線圖，由圖4可知，3種組織表達量依次為血液>卵巢>肝臟，表明基因在不同組織中轉錄水平具有明顯差異。

進一步比較差異表達基因數(shù)目，由表1可知，卵巢vs血液差異表達基因數(shù)量最多(9 726個)，肝臟vs卵巢次之(8 390個)，肝臟vs血液最少(6 142個)，且卵巢vs血液上調基因數(shù)量最多，肝臟vs卵巢上調基因數(shù)量最少，與之相反，肝臟vs卵巢下調基因數(shù)量最多，卵巢vs血液下調基因數(shù)量最少。

表1 差異表達基因數(shù)目(DEGs) 個

由圖5可知，卵巢vs血液和卵巢vs肝臟差異表達基因模式相似，多集中于細胞組分中的胞內、細胞內組分、細胞器、細胞內細胞器和細胞質，和生物學過程中的初級代謝過程、細胞大分子代謝過程、細胞蛋白代謝過程、細胞組分組織或生物合成；在分子功能方面，僅在核糖體結構組成部分有少量差異表達基因，而肝臟vs血液差異表達基因數(shù)量最少，多集中于細胞組分中的蛋白質復合體、中間絲狀體和中間絲狀體細胞骨架，在生物學過程中的氣體運輸和氧氣運輸，分子功能中的氧氣轉運子活性、氧氣結合、鐵離子結合、亞鐵血紅素結合和四吡咯結合方面有少量差異表達基因(圖5-C)。差異基因的富集分析，為進一步篩選組織特異性表達基因及探究組織特異功能分子調控機制奠定了基礎。

在生物體內，不同基因相互協(xié)調行使其生物學功能，通過Pathway顯著性富集能確定差異表達基因參與的最主要生化代謝途徑和信號轉導途徑。由圖6可知，KEGG是有關Pathway的主要公共數(shù)據(jù)庫，依據(jù)KEGG代謝途徑數(shù)據(jù)庫，找出差異基因相對于所有有注釋的基因顯著富集的Pathway，繪制差異基因KEGG富集散點圖，KEGG富集程度通過Rich factor、Qvalue和富集到此通路上的基因數(shù)量來衡量。其中，Rich factor指差異表達的基因中位于該Pathway條目的基因數(shù)量與所有有注釋基因中位于該Pathway條目的基因總數(shù)的比值。Q-value是做過多重假設檢驗校正后的P-value，Q-value的取值范圍為[0，1]，越接近于0，表示富集越顯著，筆者所在課題組挑選富集最顯著的20條Pathway在圖6中展示，不足20條的則全部展示。由圖6可知，肝臟vs卵巢和卵巢vs血液中，差異表達基因富集數(shù)量在細胞周期和核糖體Pathway都較多；卵巢vs血液和肝臟vs血液中，差異表達基因富集數(shù)量在內質網(wǎng)蛋白質加工Pathway富集較多，肝臟vs卵巢在錯配修復Pathway富集數(shù)量最少，卵巢vs血液在脂肪酸延長Pathway富集數(shù)量最少，肝臟vs血液在藥物代謝-細胞色素P450 Pathway富集數(shù)量最少。

2.4 SNP和SSR統(tǒng)計結果分析

統(tǒng)計分析結果顯示，在鯽魚卵巢組織中檢測到的SNP位點最多(348 795個)，其次為血液(343 516個)，肝臟組織中最少(307 740個)。其中，肝臟和卵巢組織編碼區(qū)SNP數(shù)量多于非編碼區(qū)，而血液中非編碼區(qū)SNP數(shù)量多于編碼區(qū)數(shù)量。由表2可知，3個組織中同義突變數(shù)量都多于非同義突變數(shù)量。

表2 SNP統(tǒng)計

利用MISA從鯽魚127 801條unigene中共檢測到48 808個SSR位點，其中，具有SSR位點的unigene 32 769個，SSR位點大于1的unigene 10 260個。不同類型SSR出現(xiàn)頻率不同，由圖7可知，單堿基、二堿基和三堿基重復類型所占比例較高，其他3種重復類型所占比例較少。其中，單堿基重復類型中以9～12次重復數(shù)目的SSR最多，其次為 13～16次重復數(shù)目的SSR。二堿基重復類型中以 5～8次重復數(shù)目的SSR最多，9～12次重復數(shù)目的次之。三堿基、四堿基、五堿基和六堿基重復類型中均以5～8次重復數(shù)目的SSR最多。且發(fā)現(xiàn)二堿基重復基元中，AC/GT含量最高，其次為AT/AT，三堿基重復基元中，AAT/ATT含量最高，ATC/ATG次之。

3 討論與結論

鯽魚不僅是一種重要的水產(chǎn)養(yǎng)殖物種，而且還是研究進化基因組學和生理性適應機制的模式生物，近年來，國內外對其開展了大量研究，獲取了不同組織(嗅覺上皮、腦、肌肉、肝臟、腎臟和皮膚)的轉錄組信息。為進一步豐富其基因組數(shù)據(jù)，本研究對其肝臟、血液和卵巢組織進行了高通量轉錄組測序。

3.1 轉錄組質量

本研究利用Trinity軟件對所得的reads片段進行組裝，得到127 801條unigene，所得unigene的長度為201～17 286 bp，平均長度為735 bp。其N50為1 288 bp，N50值越大說明組織得到的長片段越多，組裝效果越好。Q30值在80%以上就認為測序質量非?？煽?，本研究中鯽魚肝臟、卵巢和血液轉錄組Q30值分別為90.42%、89.59%和90.21%，均大于80%，表明本研究構建的轉錄組數(shù)據(jù)庫準確可信，可以為后續(xù)鯽魚基因克隆及功能基因驗證提供基礎數(shù)據(jù)。

3.2 數(shù)據(jù)庫注釋與差異表達基因

數(shù)據(jù)庫注釋結果顯示，鯽魚3個組織轉錄組拼接組裝的127 801條unigene在NR、NT、KO、Swiss-Prot、GO、KOG和KEGG公共數(shù)據(jù)庫中均得到注釋，有助于進一步深入了解基因的功能。其中GO、KOG和KEGG注釋結果表明，肝臟、卵巢和血液組織中除了參與一般功能預測外，參與信號傳導機制功能的基因數(shù)量最多，其次是內分泌系統(tǒng)和免疫系統(tǒng)，這可能與肝臟肩負免疫防御功能、血液是體液調節(jié)的聯(lián)系媒介，機體分泌的激素、酶和維生素等物質也是依靠血液傳遞才能發(fā)揮對代謝的調節(jié)作用，及與類淋巴細胞的免疫防御功能、卵巢主要承擔雌性生殖與內分泌兩大重要功能等密切相關。

不同組織的轉錄組表達基因量往往存在差異，如紅鰭東方鲀鰾的轉錄組表達基因量高于鰓的轉錄組表達基因量，日本七鰓鰻肝臟組織轉錄組表達基因量則低于血液組織轉錄組表達基因量，本研究中鯽魚肝臟、血液和卵巢組織的基因表達量分布顯示，血液的表達量最高，這可能與血液是體液調節(jié)的聯(lián)系媒介，參與機體許多生命活動的調節(jié)過程有關。差異表達基因數(shù)量顯示，卵巢和血液差異表達基因數(shù)量最多，肝臟和血液差異表達基因數(shù)量最少，這可能是由于卵巢主要參與雌性生殖的調控，與血液在功能上重合的較少，二者表達的基因存在較大差異。肝臟是機體的代謝中心，而調節(jié)機體代謝過程的一些物質主要依靠血液傳遞才能發(fā)揮對代謝的調節(jié)作用，所以肝臟和血液在功能上有較多的重疊，二者表達的基因也相似。一些差異表達基因則與組織的特異功能息息相關，如血液除了傳遞調節(jié)代謝的物質外，還主要承擔著氣體運輸、氧氣結合和運輸、鐵離子結合和亞鐵血紅素結合的特異功能。

KEGG注釋結果顯示，肝臟vs卵巢和卵巢vs血液差異主要體現(xiàn)在細胞周期和核糖體Pathway，推測這些活動相關Pathway主要與卵巢組織中雌性生殖細胞的產(chǎn)生密切相關。卵巢vs血液和肝臟vs血液差異主要體現(xiàn)在內質網(wǎng)蛋白質加工Pathway，推測其主要與血液中運輸?shù)牡鞍踪|物質加工過程相關。此外，肝臟和卵巢可能共同參與錯配修復Pathway，卵巢和血液共同參與脂肪酸延長Pathway，肝臟和血液共同參與藥物代謝-細胞色素P450 Pathway。

3.3 分子標記

SNP和SSR是利用轉錄組數(shù)據(jù)開發(fā)最多的2類標記，本研究依據(jù)轉錄組數(shù)據(jù)發(fā)現(xiàn)肝臟組織中307 740個、卵巢組織中348 795個和血液組織中343 516個SNP位點，此外，利用MISA軟件在鯽魚3個組織中共檢測到48 808個SSR位點。與傳統(tǒng)方法相比較，轉錄組數(shù)據(jù)可以直接反映基因的表達情況，利用轉錄組測序數(shù)據(jù)挖掘和開發(fā)與性狀相關的優(yōu)勢基因型與分子標記更省時、高效，為選育提供遺傳依據(jù)。因此，本研究將為鯽魚今后的多態(tài)性檢測、群體遺傳多樣性分析及分子鑒定等方面打下基礎。

本研究采用高通量測序技術，對鯽魚的肝臟、血液和卵巢轉錄組進行測序和分析。經(jīng)組裝后最終獲得有注釋信息的unigene共117 414條，不同組織間比較結果顯示，卵巢與血液差異表達基因數(shù)量最多，肝臟與血液的差異表達基因數(shù)量最少。此外，在鯽魚肝臟、血液和卵巢組織中分別檢測到 307 740、343 516、348 795個SNP位點，共檢測到 48 808 個SSR位點。研究結果為進一步克隆和挖掘鯽魚功能基因、多態(tài)性檢測及群體遺傳多樣性分析以及探究鯽魚耐低氧分子機制等方面研究奠定了基礎。