無乳鏈球菌感染卵形鯧鲹脾臟轉(zhuǎn)錄組學分析

卞澤昌，陳健安，王亞丹，黃壯，杜春民，蔡小輝，胥鵬*，關杰耀,

（1.廣西大學 動物科學技術學院 亞熱帶農(nóng)業(yè)生物資源保護與利用國家重點實驗室，廣西南寧 530004；2.北部灣大學 海洋學院 廣西北部灣海洋生物多樣性養(yǎng)護重點實驗室，廣西欽州 535011）

卵形鯧鲹（Trachinotus ovatus）因肉質(zhì)鮮美、繁殖周期短、經(jīng)濟價值高而在中國被廣泛養(yǎng)殖。但近年來，隨著卵形鯧鲹養(yǎng)殖規(guī)模的不斷擴大，病害頻繁暴發(fā)，鏈球菌病就是其中之一[1]。鏈球菌病的頻繁發(fā)生給世界各地的養(yǎng)魚業(yè)造成了嚴重的經(jīng)濟損失[2-4]，受鏈球菌影響的魚類包括牙鲆（Paralichthys olivaceus）[5]、尼羅羅非魚（Oreochromis niloticus）[6]、斑點叉尾鮰（Ictalurus punctatus）[7]、卵形鯧鲹[3]和觀賞魚[3,8]等，感染鏈球菌的魚通常表現(xiàn)出非常相似的癥狀，如眼球突出、游泳異常、腦膜炎和內(nèi)臟出血等[2-3]。無乳鏈球菌作為鏈球菌中的一種，不僅會影響哺乳動物[9]，還會對魚類造成感染，如齊口裂腹魚（Schizothorax prenanti）[10]、尼羅羅非魚（Oreochromis niloticus）[11]和尖吻鱸（Lates calcarifer）[12]等。研究顯示無乳鏈球菌引起的疾病可造成卵形鯧鲹的高死亡率，嚴重阻礙卵形鯧鲹養(yǎng)殖產(chǎn)業(yè)的快速發(fā)展[3]。利用分子生物學技術研究卵形鯧鲹感染無乳鏈球菌后的免疫應答機制是至關重要的。

轉(zhuǎn)錄組測序（RNA-seq）作為一種快速、方便的基因組研究方法，可以識別許多功能基因和分子標記[13]，近年來已發(fā)展成為研究魚類病毒性和細菌性疾病發(fā)病機制的有力工具。例如，利用RNA-seq 技術檢測了牙鲆淋巴囊腫病毒、病毒性出血敗血癥和細胞腫大病毒引起的免疫應答和感染[13-14]，以及通過RNA-seq 分析羅非魚和斑馬魚感染無乳鏈球菌后的免疫反應[11,15]等。本研究通過給健康的卵形鯧鲹腹腔注射無乳鏈球菌進行攻毒，然后對卵形鯧鲹脾臟進行RNA 測序。構建測序數(shù)據(jù)庫后進行分析，利用GO 和KEGG 富集分析對差異表達基因（DEGs）和信號通路進行富集分析，并通過STEM 軟件分析，確定DEGs 的表達趨勢。本研究為進一步篩選和研究卵形鯧鲹無乳鏈球菌感染后的免疫相關基因和通路提供了理論依據(jù)。

1 材料與方法

1.1 實驗魚和菌株

36 條健康的卵形鯧鲹購自廣西欽州市犀牛角鎮(zhèn)鴻鈺水產(chǎn)科技有限公司，體重為（50±0.7）g。并將其置于200 L 水箱中用處理過的海水在28 ～30 ℃下適應10 d。無乳鏈球菌由北部灣大學蔡小輝副教授惠贈。

1.2 試劑與儀器

TRNzol Universal 總RNA 提取試劑（100 mL），北京天根生化科技有限公司；腦心浸出液肉湯（BHI）（250 g），北京索萊寶科技有限公司；瓊脂粉（250 g），上海生工生物工程股份有限公司；無水乙醇、氯仿，分析純，成都科龍化工試劑有限公司；NaoDrop2000 核酸濃度檢測儀，美國Thermo Scientific 公司；DYY-6C凝膠電泳儀，北京六一生物科技有限公司；GDBL-1000 凝膠成像儀，美國Axygen 公司；THZ-92CS 氣浴恒溫振蕩器，杭州川恒實驗儀器有限公司；2-16KL臺式低溫離心機，德國Sigma 公司。

1.3 細菌培養(yǎng)及組織采集

將實驗菌株無菌接種到BHI 固體培養(yǎng)基中過夜，然后挑取單克隆至10 mL BHI 液體培養(yǎng)基中，在37 ℃下以200 r/min 振蕩孵育過夜，最后在50 mL BHI液體培養(yǎng)基中按1 ∶100（菌液∶液體培養(yǎng)基）體積比例繼續(xù)在37 ℃下以200 r/min 振蕩孵育10 h。4 ℃離心收集細菌菌體，用PBS 緩沖液（pH=7.4）重懸，通過梯度稀釋法將重懸菌液稀釋106倍，最后使用平板計數(shù)法進行計數(shù)。

將36 條魚分為0 h（對照組，C）、6 h（實驗組，W6）和24 h（實驗組， W24）3 組，每組包含12 條魚；對照組腹腔注射100 μL PBS 緩沖液（pH = 7.4），W6 組與W24 組腹腔注射100 μL 1×107CFU/mL 無乳鏈球菌懸液。注射后0 h、6 h 和24 h 收集每組9條魚的脾臟，每3個脾臟混合成一個樣品并標明名稱，置于-80 ℃保存至RNA 提取過程。

1.4 RNA 提取、文庫構建和測序

使用Trizol 提取總RNA，通過瓊脂糖凝膠電泳和核酸濃度檢測儀評估RNA 完整性和濃度，之后構建cDNA 文庫，送往廣州基迪奧生物科技有限公司測序。

1.5 組裝和注釋

為了保證數(shù)據(jù)質(zhì)量，對測序結果進行數(shù)據(jù)過濾，以減少無效數(shù)據(jù)所帶來的分析干擾，首先利用fastp對raw reads 過濾低質(zhì)量數(shù)據(jù)得到clean reads，再通過Trinity 軟件進行組裝獲得unigenes。將獲得的unigenes通過Nr數(shù)據(jù)庫、Swiss-Prot蛋白數(shù)據(jù)庫、KEGG數(shù)據(jù)庫、GO 數(shù)據(jù)庫、COG/KOG 數(shù)據(jù)庫進行注釋。

1.6 差異表達基因分析方法

使用DESeq2 軟件分析3 個不同組之間的差異表達（并在兩個樣本之間使用edgeR）。以錯誤發(fā)現(xiàn)率（FDR）＜0.05 和絕對折疊變化≥2 為參數(shù)，判斷存在差異表達的基因?；赪allenius 非中心超幾何分布，GOseq 包對差異表達基因進行GO 富集分析。通過KEGG 顯著富集分析找出基因顯著富集的通路，通過通路顯著富集確定基因參與的最主要生化代謝通路和信號轉(zhuǎn)導通路。

1.7 趨勢分析

為檢驗DEGs 的表達模式，先將每個樣本的表達數(shù)據(jù)歸一化為0、log2（v1/v0）、log2（v2/v0），再通過STEM 軟件進行聚類。隨后，對所有圖譜中的DEGs 進行GO 和KEGG 富集分析。通過顯著性分析和FDR 校正假設檢驗，將Q 值≤0.05 的GO 條目或通路分類為顯著富集的GO 條目或通路。

2 結果與分析

2.1 RNA 測序和功能注釋

測序結果顯示，raw reads 經(jīng)過過濾后，所有clean reads 被組裝成72 857 個unigenes，其中平均長度和N50 分別為1 121 bp 和2 553 bp。通過主要功能數(shù)據(jù)庫（NR、KOG/COG、Swiss-Prot 和KEGG）共對72 857 個unigenes 進行了注釋，獲得結果分別為25 089 個（占34.56%）、23 324 個（占32.01%）、18 500 個（占25.39%）和14 712 個（占20.19%）unigenes。通過DESeq2 分析共發(fā)現(xiàn)了10 014 個DEGs。相較于對照組，6 h 中發(fā)現(xiàn)了4 749 個DEGs，其中上調(diào)和下調(diào)DEGs 分別為1 812 和2 937 個（圖1a）；24 h 中共發(fā)現(xiàn)8 651 個DEGs，其中上調(diào)和下調(diào)DEGs 分別為4 418 和4 233 個（圖1b）。在6 h 和24 h 組中同時出現(xiàn)表達的DEGs 有3 386 個。

圖1 差異表達火山圖

2.2 GO 和KEGG 富集分析

基于GO 注釋，進行GO 分析以得知DEGs 的功能，GO 分析的3 個主要類別是分子功能（Molecular Function）、細胞成分（Cellular Component）和生物過程（Biological Process）（圖2）。在生物過程類別中，6 h 和24 h 的DEGs 主要在細胞過程（Cellular process）、單有機體過程（Single organism process）和代謝過程（Metabolic process）中富集，其中免疫系統(tǒng)過程（Immune system process）在6 h 和24 h 組中同時出現(xiàn)。對于細胞成分類別，兩個時間點的DEGs 富集在細胞（Cell）、細胞部分（Cell part）和膜（Membrane）等條目中。在分子功能類別中，兩個時間點的DEGs 主要富集在結合（Binding）、催化活性（Catalytic activity）和分子傳感器活動（Molecular transducer activity）等條目中。

圖2 GO 富集分類柱狀圖

KEGG 富集分析顯示，在6 h 和24 h 分別確定了30 個和8 個顯著富集的通路（圖3）。在6 h 顯著富集的前3 條通路主要有細胞因子-細胞因子受體相互作用（Cytokine-cytokine receptor interaction）、瘧疾（Malaria）和炎癥性腸?。↖nflammatory bowel disease）等。在24 h 顯著富集的前3 條通路包括蛋白酶體（Proteasome）、溶酶體（Lysosome）和細胞因子-細胞因子受體相互作用（Cytokine-cytokine receptor interaction）。在6 h 和24 h 組中同時發(fā)現(xiàn)了細胞因子-細胞因子受體相互作用和類風濕性關節(jié)炎（Rheumatoid arthritis）兩個通路。

2.3 趨勢分析

利用STEM 分析，將20 691 個DEGs 聚類為8 個剖面圖，其中0、1、3 和7 這4 個剖面顯著聚類（P＜0.05），分別包含3 860 個、4 116 個、3 463 個和2 444 個DEGs（圖4）。通過KEGG 富集分析從這4 個剖面中發(fā)現(xiàn)了7 條主要免疫相關通路，其中產(chǎn)生IgA 的腸道免疫網(wǎng)絡、Th1 和Th2 細胞分化、白細胞經(jīng)內(nèi)皮遷移、造血細胞譜系、Th17 細胞分化和T 細胞受體信號通路6 條通路顯著下調(diào)，僅胞質(zhì)DNA 傳感通路顯著上調(diào)。

圖4 基因趨勢圖

3 討論與結論

本研究構建了卵形鯧鲹感染無乳鏈球菌后的脾臟轉(zhuǎn)錄組圖譜。研究結果顯示，無乳鏈球菌感染改變了10 014 個基因的表達，在感染6 h 處理組中發(fā)現(xiàn)4 749個DEGs，其中上調(diào)和下調(diào)基因數(shù)量分別為1 812個和2 937 個；在感染24 h 處理組中發(fā)現(xiàn)8 651 個DEGs，其中上調(diào)和下調(diào)基因數(shù)量分別為4 418 個和4 233 個。通過STEM 分析發(fā)現(xiàn)，DEGs 被聚類為8 個趨勢圖譜，其中0、1、3 和7 四個圖譜顯著聚集，在這4 張顯著富集的圖譜中發(fā)現(xiàn)了7 條主要與免疫相關的信號通路，其中產(chǎn)生IgA 的腸道免疫網(wǎng)絡、Th1 和Th2 細胞分化、白細胞經(jīng)內(nèi)皮遷移、造血細胞譜系、Th17 細胞分化、T 細胞受體信號通路6 條通路顯著下調(diào)，僅胞質(zhì)DNA 傳感通路顯著上調(diào)。在羅非魚無乳鏈球菌感染的研究中，早期時間點的DEGs 被快速誘導，晚期時間點的DEGs 表達水平逐漸下降，表明羅非魚可能具有快速識別入侵細菌并通過激活幾種免疫相關途徑來清除細菌的強大能力[11]。與羅非魚研究結果不同，本研究中大多數(shù)基因在早期迅速下調(diào)，在后期逐漸穩(wěn)定，這可能意味著卵形鯧鲹在早期缺乏快速識別和抵抗細菌入侵的能力，導致細菌在早期抑制了大多數(shù)免疫相關基因和通路，使其無法發(fā)揮免疫作用。本研究結果為進一步篩選卵形鯧鲹無乳鏈球菌病的免疫關鍵基因和通路提供了基礎信息，為卵形鯧鲹疾病預防和治療提供了新見解。