王龍斌， 泮淼軍， 王明陽， 王潤澤， 李 麗，2， 董雙林，2， 李衛(wèi)東， 田相利 ，2??

(1. 海水養(yǎng)殖教育部重點實驗室(中國海洋大學)， 山東 青島 266003； 2. 青島海洋科學與技術試點國家實驗室， 山東 青島 266237； 3. 唐山海都水產食品有限公司， 河北 唐山 063000)

新疆是中國鹽堿化土壤分布面積較廣、土壤積鹽較重的地區(qū)，其中南疆鹽堿水具有堿度高、離子比例不平衡、水質變動大等特點，合理有效地開發(fā)利用鹽堿地已成為該地區(qū)緊迫而艱巨的任務[1]。堿性水通常具有較高的堿度和pH，嚴重影響水生動物的正常生理過程，例如直接抑制魚類的氨排泄并增加CO2排泄[2]。盡管如此，仍有一些魚類可以在高堿的環(huán)境中生存[3]，如肯尼亞馬加迪湖的馬加迪羅非魚(Alcolapiagrahami)[4]、美國金字塔湖的美洲鮭(Oncorhynchusclarkihenshawi)[5]、我國達里湖的瓦氏雅羅魚(Leuciscuswaleckii)[6]以及青海湖的裸鯉(Gymnocyprisprzewa-lskii)[7]等。這些魚類已經(jīng)進化出特殊的生理機制以適應堿性環(huán)境[7-9]。

轉錄組分析是解釋功能基因組學元素和揭示細胞與組織分子機制的有力工具之一[10]。近年來，RNA測序技術(RNA-seq)已廣泛應用于基因表達譜的研究，包括水生動物響應鹽堿脅迫的通路研究。如馬加迪羅非魚[11]、瓦氏雅羅魚[6,12]、尼羅羅非魚(Oreochromisniloticus)[13]、奧尼羅非魚(Oreochromismossambicusfemale×O.urolepishornorummale)[14]等，這些研究提供了大量與堿度脅迫相關的候選通路和基因。

花鱸(Lateolabraxmaculatus)屬鱸形目(Perciformes)鮨科(Serranidae)花鱸屬(Latelabrax)，具有生長快、病害少、對鹽度適應廣等特點。近年來，關于花鱸對鹽度的適應能力[15-16]以及耐受機制[17-18]已得到了很好的研究?；|對堿度有較高的耐受性，在鹽度為10 時，堿度96 h半致死濃度為98.51 mmol/L[19]。但目前國內外對于花鱸耐堿機制的研究很少，僅見于堿度脅迫下Slc4基因家族[20]和熱休克蛋白 (Hsp)超家族[21]的相關研究。

在本研究中，我們利用RNA-seq技術探究了不同堿度下花鱸腎組織的基因表達情況，對其響應堿度脅迫相關的通路及基因進行了綜合分析，為了解花鱸適應堿性環(huán)境的分子生物學機制提供有價值的信息。

1 材料與方法

1.1 樣本來源及實驗設計

實驗所用花鱸是由唐山由海都水產食品有限公司(河北，唐山)提供的人工繁育苗種。暫養(yǎng)一周后隨機選取135尾健康、個體大小均勻的花鱸進行養(yǎng)殖實驗。

本研究參照南疆阿拉爾地區(qū)鹽堿水的鹽度水平(10左右)，通過添加NaHCO3，設計了不同的堿度梯度。將天然海水與淡水混合成鹽度為10的海水，再添加NaHCO3(分析純)，全部溶解并穩(wěn)定24 h后使用。根據(jù)前期急性脅迫的結果，設置了0 mmol/L(con組)、15 mmol/L(AW15組)和30 mmol/L(AW30組)的堿度梯度，具體實測堿度如表1所示。每組設置3個重復，每個重復15尾魚，實驗魚體長(15.21±1.56) cm，體質量(60.23±5.54) g，實驗時間為2020年9月15日至10月15日，實驗周期為 30 d。

表1 碳酸鹽堿度設置

實驗所用水體的體積為200 L，鹽度為10，水溫控制在(20±2) ℃，溶氧(8.4±0.1) mg/L，pH=8.5±0.3。實驗期間，每天8:00和18:00各投喂花鱸專用顆粒飼料1次，投喂量約為魚體質量3%左右，每天換水50%，并檢測碳酸鹽堿度和pH。pH以DELTA320 型精密pH計測定，碳酸鹽堿度以酸堿滴定法(GB/T 9736—2008)測定[22]。

1.2 生長指標的測定

實驗開始和結束時測量花鱸的體質量，按以下公式計算體質量特定生長率(SGR)：

SGR=(lnW2-lnW1)/t×100%。

式中：W1為初始體質量(g)；W2為終末體質量；t為實驗時間(30 d)。

實驗數(shù)據(jù)以平均值±標準差表示,實驗結果采用SPSS 19.0軟件對實驗數(shù)據(jù)進行統(tǒng)計分析。先對數(shù)據(jù)做單因素方差分析，組間若有顯著性差異，再用Duncan’s多重比較分析，P<0.05時認為有顯著差異。

1.3 RNA提取、建庫和測序

實驗結束時，每個處理隨機取3尾魚(共9尾)使用MS-222麻醉，解剖以獲得腎組織，在液氮中速凍并儲存在-80 ℃用于RNA分離。使用Invitrogen公司的TRIzol?Reagent試劑盒提取組織RNA，然后使用Roche公司的DNaseI去除基因組DNA污染，采用Nanodrop 2000、Qubit 2.0和Aglient 2100方法檢測RNA樣品純度、濃度和完整性，檢驗合格的樣品送上海美吉生物醫(yī)藥科技有限公司，采用Illumina TruseqTM RNA sample prep Kit方法進行文庫構建，基于Illumina Novaseq 6000 測序平臺以Paired-end 150 bp雙末端測序模式分別對9個樣本進行轉錄組測序。

1.4 序列比對

利用fastp對測序獲得的原始數(shù)據(jù)進行質量評估，然后用SeqPrep(https://github.com/jstjohn/SeqPrep)和Sickle(https://github.com/najoshi/sickle)去除測序接頭序列、低質量讀段、N率較高序列及長度過短序列。通過TopHat2[23]將質控后的數(shù)據(jù)與花鱸參考基因組(http://gigadb.org/dataset/100458)進行比對，獲得用于后續(xù)分析的mapped data。

1.5 表達量差異分析

使用RSEM V1.2.15[24]對mapped data進行計數(shù)，然后通過FPKM進行數(shù)據(jù)標準化。利用DESeq2 R包[25]篩選差異表達基因，將P-adjust<0.05且|log2FC|≥1作為篩選標準。將花鱸腎組織在三個處理下兩兩對比得到的基因集(con_vs_AW15，con_vs_AW30，AW15_vs_AW30)中進行比較分析。

1.6 GO富集和KEGG通路富集分析

使用軟件Goatools(https: //github.com/tanghaibao/GOatools)對差異基因進行GO富集分析，采用Bonferroni方法對P值進行校正，當經(jīng)過校正的P值(FDR)小于0.05 時，認為此GO功能存在顯著富集。利用KOBAS軟件(v2.0.12)對差異基因進行KEGG通路富集分析，以判定差異基因主要影響的生物學功能和通路。

1.7 qRT-PCR驗證

選擇9個差異表達基因，使用Primer Premier 6.0軟件[26]設計引物(見表2)，并由生工生物工程(上海)股份有限公司合成，進行實時熒光定量PCR驗證。使用諾唯贊公司HiScript?III All-in-one RT SuperMix Perfect for qPCR試劑盒將總RNA逆轉錄合成cDNA，再以逆轉錄產物為模板、以18S ribosomal RNA為內參基因，使用諾唯贊公司的ChamQTMSYBR?Color qPCR Master Mix試劑盒，在ABI PRISM 7500型熒光定量PCR儀上進行表達量分析。總體系20 μL的反應體系如下：稀釋10倍的cDNA模板2.0 μL，2×ChamQ SYBR Color qPCR Master Mix 10 μL，正反向引物各0.4 μL，RNase-Free ddH2O 7.2 μL。反應程序：95 ℃ 30 s，95 ℃ 10 s，60 ℃ 30 s，40個循環(huán)。最后在60～95 ℃制作熔解曲線。每組樣品的每個基因均重復3次，采用2-ΔΔCt法[27]分析基因的相對表達量。使用SPSS 19.0軟件計算RNA-seq結果和qRT-PCR結果之間的相關性系數(shù)。

表2 RT-qPCR引物信息

2 結果

2.1 生長指標

養(yǎng)殖實驗結果見表3。由表3可以看出，堿度嚴重影響花鱸生長，隨著堿度升高，花鱸的終末體質量、特定生長率顯著降低。在15 mmol/L的中堿度下，花鱸的生長速率顯著降低。而在30 mmol/L的高堿度下，花鱸的終末體質量顯著低于初始體質量。

表3 不同堿度對花鱸生長性能的影響

2.2 測序數(shù)據(jù)統(tǒng)計與參考序列比對分析

所有樣品的RNA濃度均高于200 ng/μL，RNA總量8 μg，均符合建庫要求。對9個轉錄組測序結果進行分析可知，所有樣品獲得的Clean reads數(shù)量在46 266 138～55 723 922之間，利用SeqPrep軟件過濾掉不合格的序列后，所有樣品得到的有效讀數(shù)均占原始讀數(shù)總數(shù)量的97%以上(見表4)。表4中，Q20、Q30分別指測序質量在99%和99.9%以上的堿基占總堿基的百分比，在得到的9個轉錄組中，Q20均在96%以上，Q30均在91%以上。利用HISAT2軟件將有效讀數(shù)比對到參考基因組，結果發(fā)現(xiàn)有79.71%～82.12%讀數(shù)成功映射到參考基因組。這些結果表明，花鱸腎組織各處理組測序量大，獲得的有效讀數(shù)數(shù)量多，轉錄組質量高，符合進一步生物信息學分析的要求。

表4 花鱸腎組織轉錄組高通量測序基本數(shù)據(jù)統(tǒng)計

2.3 表達量差異分析

使用基于負二項分布的DESeq2軟件對raw counts進行分析，基于一定的標準化處理和篩選條件獲得比較組間表達差異的基因，參數(shù)如下：P-adjust<0.05 & |log2FC| ≥1。分析得到的差異表達基因結果表明，在不同的堿度環(huán)境下，花鱸發(fā)生顯著差異表達基因的數(shù)量也不同(見表5)。與con組相比，AW15組有142個基因上調，104個基因下調；AW30組有2 037個基因上調，1 882個基因下調；與AW15相比，AW30組有1 318個基因上調，964個基因下調。對不同堿度下差異表達基因進行綜合分析，結果如圖1所示，與con組相比，AW15組與AW30組分別有68和2 189個特有的差異表達基因；與AW15組相比，AW30組有577個特有的差異表達基因；3個基因集共有的差異表達基因為49個。

表5 花鱸腎組織差異表達基因統(tǒng)計結果

圖1 花鱸腎組織con_vs_AW15組、con_vs_AW30組和AW15_vs_AW30組顯著差異表達基因的維恩圖

2.4 GO和KEGG富集分析

為確定可能參與堿度脅迫的生物過程或通路，分別對3個基因集(con_vs_AW15，con_vs_AW30，AW15_vs_AW30)中的顯著差異基因進行GO和KEGG富集分析。

GO富集結果顯示，差異基因被富集到生物過程、細胞組分和分子功能3個類別體系中。其中最豐富的類別是生物過程，其次是分子功能和細胞組分。圖2顯示了腎中3個基因集前20個顯著GO術語?？梢钥闯?，con_vs_AW15中氣體運輸、氧結合、血紅蛋白復合物富集程度最高，con_vs_AW30中水解酶活性的正調控、GTPase活性的正調控、GTPase活性的調節(jié)富集程度最高，而AW15_vs_AW30中DNA復制啟動、蛋白質磷酸化、絲氨酸家族氨基酸代謝過程富集程度最高。其中，GTPase活性的正調控、水解酶活性的正調控、水解酶活性的調節(jié)、GTPase活性的調節(jié)、絲氨酸家族氨基酸代謝過程、α-氨基酸代謝過程、細胞氨基酸代謝過程等在con_vs_AW30和AW15_vs_AW30中均顯著富集。

圖2 花鱸腎組織con_cs_AW15(A)、con_vs_AW30(B)和AW15_vs_AW30(C)組差異表達基因的GO富集分析

KEGG分析表明，在3個基因集中，分別有254、344、343條通路與堿度脅迫有關。圖3顯示了腎con_vs_AW30和AW15_vs_AW30兩個基因集前20個顯著KEGG通路。可以看出，con_vs_AW30中細胞周期、葉酸一碳庫、過氧化物酶體顯著富集，AW15_vs_AW30中DNA復制、癌癥中的蛋白多糖、軸突導向顯著富集。其中，細胞周期、DNA復制、葉酸一碳庫、上皮細胞的細菌侵襲、趨化因子信號通路、軸突導向、肌動蛋白細胞骨架的調節(jié)、錯配修復、PPAR信號通路和癌癥中的蛋白多糖在con_vs_AW30和AW15_vs_AW30均顯著富集。

圖3 花鱸腎組織con_vs_AW30(A)組和AW15_vs_AW30(B)組差異基因KEGG散點圖

2.5 花鱸腎組織堿度脅迫響應的關鍵基因分析

圖1表明con_vs_AW15、con_vs_AW30、AW15_vs_AW30這3個基因集中有49個共有差異基因，這些基因是花鱸腎組織堿度脅迫響應的關鍵基因，其中大多隨著堿度升高而降低(見表6)。

表6 花鱸腎組織堿度脅迫響應的關鍵基因

2.6 qRT-PCR驗證

為了確認Illumina RNA-seq結果的準確性和可重復性，選擇9個代表性基因，通過實時定量PCR(qRT-PCR)驗證不同處理組基因的表達水平。進一步將所選基因的相對表達量與RNA-seq測序結果進行了相關性分析(見圖4)。可以看出，雖然所選基因在表達變化幅度上存在一定差異，但qRT-PCR檢測結果與RNA-seq測序結果中基因的表達趨勢基本一致(R2=0.959 4)，這表明了Illumina 測序獲得的轉錄組數(shù)據(jù)以及生物信息學分析的可靠性。

圖4 花鱸腎組織RNA-seq產生的9個差異表達基因進行qRT-PCR驗證

3 討論

作為水生動物的重要器官，魚類腎臟主要具有維持滲透壓和離子平衡的功能，在魚類的滲透調節(jié)過程中起著重要作用[28]。研究發(fā)現(xiàn)，高堿度會使魚類腎組織發(fā)生不同程度的退化，而堿度越高腎萎縮程度越大[29]，腎功能被明顯抑制[30]。在本研究中，我們首次研究了不同堿度脅迫下花鱸腎組織的有參轉錄組數(shù)據(jù)，探討了堿度對花鱸腎組織的影響，獲得了9 602個差異基因，富集分析后得到大量堿度脅迫相關通路以及對堿度敏感的基因，為花鱸耐堿機制研究提供了新視角。

3.1 GO富集分析

GO數(shù)據(jù)庫是一個國際化的基因功能標準分類體系，通過GO功能顯著性富集分析能確定差異表達基因行使的主要生物學功能[31]。本研究對花鱸腎臟轉錄組數(shù)據(jù)GO富集分析表明，相較于對照組，堿度15 mmol/L處理組差異基因GO術語顯著富集于氧氣和血紅蛋白相關術語，說明中堿度可能主要影響腎組織氧氣的運輸。相類似，何強等[32]發(fā)現(xiàn)，隨著堿度上升，瓦氏雅羅魚的耗氧率呈先下降后上升的趨勢，Lykkeboe等[33]發(fā)現(xiàn)，羅非魚(Tilapiagrahami)血紅蛋白的氧親和力隨著pH升高而降低，這些研究結果與本研究基本相一致。但高堿度對魚類的影響機制與中堿度可能有所差異。例如，Xu等[6]發(fā)現(xiàn)高堿度下瓦氏羅雅魚腎組織GO術語顯著富集于轉運蛋白活性和酶調節(jié)活性，高堿度導致腎組織跨膜轉運體活性明顯增加。在本研究中，高堿度(30 mmol/L)對花鱸腎的影響主要集中在蛋白質水解(水解酶活性的正調控、催化活性的正調控)和氨基酸代謝過程(絲氨酸家族氨基酸代謝過程、α-氨基酸代謝過程、細胞氨基酸代謝過程等)中，這可能表明腎組織在高堿度下能量消耗巨大，能量需求的增加導致腎臟蛋白質水解和氨基酸代謝過程活躍。

3.2 KEGG富集分析

KEGG是一個集成了基因組、化學和系統(tǒng)功能信息的生物系統(tǒng)數(shù)據(jù)庫，KEGG提供了一個參考知識庫，通過通路映射過程將基因組與生命聯(lián)系起來[34]。在本研究中，不同堿度下花鱸腎臟KEGG富集分析表明，相較于對照組，中堿度的15 mmol/L處理組對花鱸腎組織的影響主要集中在氨基酸合成與代謝相關通路(苯丙氨酸、酪氨酸和色氨酸生物合成、氨基酸生物合成、半胱氨酸和蛋氨酸代謝等)，而高堿度(30 mmol/L)對花鱸的影響則更加廣泛，主要包括細胞有絲分裂(細胞周期、DNA復制)、免疫與疾病相關通路(例如上皮細胞的細菌侵襲、癌癥中的通路)以及代謝(葉酸一碳庫)。從現(xiàn)有的研究看，Cheng等[35]發(fā)現(xiàn)，高堿度會使尼羅羅非魚肌肉中不飽和脂肪酸含量升高，Zhao等[36]發(fā)現(xiàn)不飽和脂肪酸的生物合成通路在尼羅羅非魚肝組織AW40_vs_AW60中顯著富集。除了脂肪酸，氨基酸代謝在滲透調節(jié)中也起著關鍵作用[37]，例如丙氨酸、天冬氨酸和谷氨酸代謝[38-40]以及纈氨酸、亮氨酸和異亮氨酸代謝等[41]。這些研究結果表明，脂肪酸與氨基酸的合成與代謝可能是堿度脅迫影響的關鍵通路之一。而在堿度脅迫下，奧尼羅非魚鰓組織富集通路與免疫、錯配修復、核苷酸切除修復、DNA復制以及細胞周期顯著相關[14]，這與本研究結果基本一致。因此，推測高堿度可能損傷花鱸腎細胞，抑制細胞有絲分裂。值得注意的是，本研究發(fā)現(xiàn)花鱸腎中軸突導向和肌動蛋白細胞骨架調節(jié)相關基因在高堿度下顯著富集，這在其他魚類相關研究中從未發(fā)現(xiàn)。已有研究表明，在腎病中腎小球軸突導向通路顯著富集，肌動蛋白細胞骨架失調隨著足突消失和粘連喪失而從基底膜脫離而動態(tài)重組[42]，說明這兩條通路的富集可能與腎組織疾病相關。相類似，Chen等[41]發(fā)現(xiàn)癌癥中的通路在鹽度和堿度脅迫下葉爾羌高原鰍(Triplophysayarkandensis)鰓中均富集，這也與本研究結果相一致。本研究發(fā)現(xiàn)花鱸腎臟癌癥中的蛋白多糖在con_vs_AW30和AW15_vs_AW30兩個基因集中均顯著富集，而該通路的富集被認為是腎細胞癌的標志[43]。上述結果表明，中、高堿度均會不同程度地影響花鱸免疫能力，增加了花鱸患病風險。

3.3 花鱸腎組織堿度脅迫響應關鍵基因分析

環(huán)境脅迫會導致組織細胞氧化應激、蛋白質損傷、DNA雙鏈斷裂和氧化堿修飾[44-45]。本研究發(fā)現(xiàn)，隨堿度升高花鱸腎組織生長相關基因表達顯著下調，例如fyn、ccna2、cks1、arhgap19、melk、cenpo、pbk、map3k15、rrm1等。fyn是Src家族的成員，fyn的抑制與細胞生長減少有關[46]。ccna2基因編碼細胞周期蛋白，有細胞周期調節(jié)劑的作用[47]。cks1在調節(jié)細胞周期中發(fā)揮作用，Cks1蛋白的缺失可導致細胞周期停滯在G2期，有絲分裂受阻[48]。arhgap19的低表達水平則會影響早期有絲分裂期間的細胞形狀，導致嚴重的染色體分離缺陷[49]。而melk基因在細胞周期控制方面發(fā)揮著重要作用，melk的表達水平與影響細胞周期進程[50]。cenpo基因編碼間期著絲粒復合體，其編碼的蛋白質在整個細胞周期中定位于著絲粒，是有絲分裂期間雙極紡錘體組裝、染色體分離和檢查點信號傳導所必需的[51]。pbk基因編碼的絲氨酸-蘇氨酸激酶在有絲分裂中活躍[52]。map3k15基因編碼MAP3K蛋白，在細胞凋亡、應激反應和各種疾病中起重要作用，從魚類相關研究看，三刺魚(Gasterosteusaculeatus)腎臟中MAP3K15的表達與鹽度呈顯著負相關[53]。rrm1基因編碼核糖核苷酸還原酶的大型催化亞基，對細胞周期S期的DNA復制以及多個DNA修復過程非常重要[54]。本研究發(fā)現(xiàn)上述相關基因在中高堿度下顯著下調，表明堿度嚴重抑制了腎細胞的有絲分裂。因此，堿度作為一種應激源，可能抑制了花鱸腎細胞的有絲分裂，影響了細胞的生長與分裂，甚至可導致細胞凋亡。

一般認為，在鹽度脅迫下水生生物為了維持體內平衡和滲透平衡，需要合成轉運蛋白和酶，這一過程會消耗大量能量[55]。本研究發(fā)現(xiàn)隨著堿度升高，花鱸腎組織能量代謝通路相關基因(acsl、gaphd)表達顯著上調。?；o酶A合成酶(ACSL)是負責脂質代謝初始步驟的關鍵酶[56]，gapdh基因編碼的3-磷酸甘油醛脫氫酶是能量代謝中的一種重要酶，在脊椎動物中具有多種細胞調節(jié)作用，增強糖酵解對于提供升高水平的 ROS 清除劑丙酮酸可能很重要[57]。我們推測隨著堿度升高，花鱸腎組織離子轉運與免疫調節(jié)會消耗大量能量，因此導致能量代謝相關基因表達上調。

不同類型的應激會導致魚類免疫和抗病能力發(fā)生變化，但具體應激反應取決于應激源的強度及其持續(xù)時間。如果應激源是慢性和長期的，免疫反應顯示出抑制作用，魚類患病的幾率就可能增加[58]。本研究發(fā)現(xiàn)，隨著堿度升高，花鱸腎組織免疫相關通路基因(ceacam1、clec4m和b4galt5等)表達顯著下調。ceacam1基因編碼免疫球蛋白相關糖蛋白，在先天性和適應性免疫反應的調節(jié)中發(fā)揮作用[59]。clec4m基因編碼的C型凝集素是固有免疫系統(tǒng)中的重要模式識別受體[60]。b4galt5基因編碼β-1,4半乳糖基轉移酶，參與多種免疫細胞表型修飾、募集和遷移[61]。白細胞介素 4 誘導1(il4i1)編碼的分泌型 L-苯丙氨酸氧化酶是一種免疫調節(jié)劑，主要與機體免疫調節(jié)功能有關[62]。crp基因編碼C-反應蛋白，參與膽固醇結合過程，功能與宿主防御相關[63]。已有研究報道顯示，堿度會導致尼羅羅非魚免疫基因表達下調[36]和大鱗鲃(Luciobarbuscapito)免疫抑制[64]。本研究結果表明，中高堿度導致花鱸腎組織免疫能力明顯下降，可能導致患病風險增加。當然，堿度下調花鱸免疫反應的詳細機制尚不清晰，需要進一步研究確認。

4 結語

本研究通過高通量測序探討了花鱸腎組織轉錄組對堿度脅迫的應答，獲得了大量堿度脅迫相關通路，并在花鱸腎組織中鑒定出49個持續(xù)變化基因。研究發(fā)現(xiàn)，中、高堿度對花鱸腎組織影響有所不同，其中，中堿度主要影響花鱸脂肪酸與氨基酸的合成與代謝，高堿度主要影響其有絲分裂過程及免疫能力。隨堿度升高，花鱸生長受到抑制，能量代謝顯著增強，機體免疫能力下降。本研究獲得了大量花鱸腎組織與堿度脅迫相關的通路和基因，為理解水生動物堿度適應的分子機制提供了新見解，同時也為花鱸在南疆鹽堿水域漁業(yè)中的潛在應用提供了科學依據(jù)。