一株羅非魚無乳鏈球菌Streptococcus agalactiae CS2108的全基因組測序及功能分析

施金谷，馬振花，王志強，韓書煜，羅福廣，文衍紅，易 弋，蘇在權，佀再勇*

（1.廣西壯族自治區水產技術推廣站，廣西南寧 530199；2.廣西科技大學生物與化學工程學院，廣西柳州 545006；3.柳州市漁業技術推廣站，廣西柳州 545006）

羅非魚主要分布在熱帶、亞熱帶地區，對環境適應能力強，生長速度快，繁殖能力強，并且其肉味鮮美、肉質細嫩、食用方法簡單、營養豐富，深受廣大養殖戶的青睞[1-2]。發展羅非魚養殖產業對于增加農民就業機會和收入、促進國際貿易具有十分重要的意義[3-4]。

隨著羅非魚養殖模式向著高密度、集約化養殖發展，加上養殖水質、環境及氣候的惡化，漁民防控意識的欠缺，導致各類病害常大面積暴發[5]。其中鏈球菌病是羅非魚健康養殖的最大威脅，最主要的病原菌是無乳鏈球菌(Sreptococcus agalactiae)。近年來，羅非魚無乳鏈球菌病在我國廣東、廣西和海南等地時有發生，具有暴發范圍廣、死亡率高的特點[6]。無乳鏈球菌大多呈球形，鏈狀排列，革蘭氏染色陽性，具有莢膜和鞭毛，無芽孢，適宜在溫度28～37 ℃、pH 值7.4～7.6 條件下生長。在腦心浸液瓊脂固體培養基上呈卵圓形，菌落表面光滑。目前，利用Lance-field 法可以將其分成Ia、Ib、Ⅱ～Ⅸ共10 個血清型[7]。該病主要病理特征為眼球突出或渾濁發白[8]，眼眶充血，腹部膨大，肛門紅腫，少食或者絕食，游姿平衡失調，解剖病魚可見膽囊腫大，膽汁稀薄，色淺，腸內充滿淡黃色液體，肝臟增大[9-10]。組織病理學表現為鰓充血，腎臟受損、腫大、白細胞浸潤，肝臟顆粒變性，腸道粘膜上皮變性、壞死、脫落、固有膜上炎性白細胞浸潤，眼睛脈絡膜和眶骨膜組織炎性壞死等[11]。

目前細菌性疾病的診斷方法主要有生理生化鑒定、免疫學方法、分子生物學技術[12]。分子生物學診斷技術中常用的方法有PCR 技術、核酸分子雜交、多位點序列分型、環介導等溫擴增技術等。王均等基于羅非魚無乳鏈球菌的cfb和16S rRNA 基因序列建立了快速檢測無乳鏈球菌的雙重PCR 方法[13]；樊海平等根據羅非魚無乳鏈球菌16S rRNA 基因和sip基因序列用于羅非魚無乳鏈球菌病的流行病學調查[14]，這些分子生物學診斷方法都離不開無乳鏈球菌的基因組信息。

截至目前，NCBI 數據庫上有1 867 個無乳鏈球菌基因組信息，這些無乳鏈球菌基因組序列中人源和魚源的比較多，部分是牛源、蛙源和犬源的。前期實驗室從患病羅非魚中分離出1 株無乳鏈球菌，命名為S.agalactiaeCS2108，為深入了解該羅非魚無乳鏈球菌S.agalactiaeCS2108 的全基因組序列，分析該菌的基因組基本特征，為后續羅非魚無乳鏈球病的防治提供理論支撐。本研究采用Illumina Hiseq+PacBio 的測序方式對S.agalactiaeCS2108 菌株進行全基因組測序分析，并進行組裝和注釋、毒力基因和耐藥基因的預測，分析和預測結果將為羅非魚無乳鏈球菌的預防和防治提供理論支持。

1 材料與方法

1.1 菌株培養與基因組DNA提取

無乳鏈球菌S.agalactiaeCS2108 分離于患病的羅非魚肝臟組織，保存于廣西科技大學生物與化學工程學院。將S.agalactiaeCS2108 從甘油管接種到腦心浸液培養基（Brian Heart Infusion，BHI）固體培養基平板上，在37 ℃培養箱中過夜培養。挑取單菌落并接種至250 mL的BHI液體培養基中，在溫度37 ℃、轉速200 rpm的搖床中培養到對數期。然后在4 ℃、12 000 rpm 條件下離心10 min，收集菌體，送至上海美吉生物醫療科技有限公司進行全基因組測序。

1.2 系統進化樹的構建

基于S.agalactiaeCS2108 的16S rRNA 序列，通過與NCBI 數據庫對比，選擇在種屬水平上最接近的20株菌的16S rRNA 序列，通過MEGA 6.0 軟件選擇鄰近（Neighbor-Joining）法構建系統進化樹。

1.3 文庫構建與測序

基因組測序使用Nanopore PromethION 和Illumina NovaSeq 6000測序平臺進行全基因組測序。Illumina測序數據被用來評估基因組的雜合度、重復度，評估基因組大小，以及是否存在質粒及污染，輔助后續組裝策略的選擇。同時由于三代測序具有較高的測序錯誤率，用二代測序數據對三代數據進行校正，保證組裝結果具有較高的準確度。

取不少于1 μg 基因組DNA，利用Covaris 進行基因組DNA 片段化，將DNA 樣本剪切成約400 bp 的片段，瓊脂糖凝膠鑒定片段大小、分布，富集300～500 bp范圍的基因組片段，使用NEXTflex?Rapid DNA-Seq試劑盒進行文庫制備。具體步驟如下：連接A&B 接頭；去除接頭自連片段；瓊脂糖凝膠電泳進行片段篩選，保留一端是A 接頭、一端是B 接頭的片段；氫氧化鈉變性，產生單鏈DNA片段；橋式PCR擴增。

制備的文庫在Illumina Novaseq 6000 儀器上進行雙端測序（2×150 bp）。具體步驟如下：加入改造過的DNA聚合酶和帶有4種熒光標記的dNTP，每次循環只摻入單種堿基；用激光掃描反應板表面，讀取每條模板序列第一輪反應所聚合上去的核苷酸種類；將“熒光基團”和“終止基團”化學切割，恢復3′端粘性，繼續聚合第二個核苷酸；統計每輪收集到的熒光信號結果，獲知模板DNA 片段的序列。Nanopore 測序是一種納米孔測序技術。具體步驟如下：在測序前，雙鏈DNA 末端需添加2 個攜帶運動蛋白的接頭（引導接頭和發夾接頭）。當測序開始，引導接頭（leader adapter）將雙鏈DNA 片段引導至納米孔附近，引導運動蛋白將雙鏈DNA 解螺旋成單鏈，使第一鏈（模板鏈）穿過納米孔。在模板鏈的末端，發夾運動蛋白再以類似的方式介導互補鏈通過納米孔，從而保證DNA的雙鏈測序。當DNA 鏈穿過納米孔時，傳感器以恒定的頻率測量離子電流變化，根據電流變化的頻譜，應用模式識別算法得到堿基序列。

1.4 基因組組裝、基因預測與注釋

利用Nanopore PromethION 和Illumina NovaSeq6000平臺生成的數據進行生物信息學分析。所有分析均在上海美吉生物云（http://cloud.majorbio.com）上進行。具體程序如下：Illumina 測序下機的原始數據（raw data）以Fast 格式儲存，為了使后續的組裝更加準確，使用軟件Fast v0.23.0 對原始數據進行質量剪切。利用Unicycler v0.4.8[15]對質控之后的Illumina 數據和Nanopore 數據進行混合組裝，最后使用Pilon v1.22 將Illumina 短序列mapping 到組裝好的基因組上進行矯正。

利用Prodigal v2.6.3[16]對基因組中的編碼序列（CDS）進行預測，質粒基因采用GeneMarkS[17]軟件預測，tRNAscan-SE v2.0[18]進行tRNA 預測，Barrnap v0.9(https://github.com/tseemann/barrnap) 進 行rRNA 預 測。利用BLASTP、Diamond、HMMER 等序列比對工具，從NR、Swiss-Prot、Pfam、GO、COG、KEGG 數據庫中對預測到的CDS進行功能注釋。

1.5 毒力和耐藥基因預測

將無乳鏈球菌S.agalactiaeCS2108 基因組信息比對毒力因子數據庫VFDB (Version: 2016.03，http://www.mgc.ac.cn/VFs/main.htm)獲得毒力因子基因注釋概況并進行統計；通過耐藥基因數據庫（CARD Version 1.1.3，DOI: 10.1128/AAC.00419-13）獲得每個基因組中包含的耐藥基因的信息。

2 結果與分析

2.1 S.agalactiae CS2108的16S rRNA進化樹分析

基于菌株S.agalactiaeCS2108 16S rRNA 基因序列與NCBI 數據庫中相似菌株序列比對結果顯示無乳鏈球菌菌株CS2108 和GCF_000186445.1 相似度最高（見圖1）。

2.2 基因組組裝與注釋

原始測序數據結果表明，無乳鏈球菌S.agalactiaeCS2108 基因組包含1 條染色體和1 個質粒，染色體大小為2 189 951 bp，GC 含量為35.82%；質粒大小為4 440 bp，GC含量為37.82%。基因組共編碼2 110個基因，所有編碼基因總長度為1 922 967 bp，基因平均長度911.36 bp，占整個基因組總長度的96%，含有80個tRNA、21個rRNA(5S RNA、16S RNA和23S RNA 各7 個)；含有8 個基因島，基因島基因總長度為163 023 bp，平均每個基因島長20 377.88 bp，8個基因島共有基因179個，平均每個基因島有22.38個基因。包含2 個前噬菌體，一個位于622685～672418處，序列總長度為49 734 bp，包含47個基因，另一個位于2135280～2145186 處，序列總長度為9 907 bp，包含17 個基因。含有11 個CRISPR-Cas，共含有重復序列107 個，平均重復序列長度35.67 bp（見表1）。Circos 基因組圈圖可以全面展示基因組的特征，從外到內圓圈對應的信息依次為基因組大小標識、正鏈、負鏈上的基因信息、ncRNA、GC 含量、GC-Skew。菌株S.agalactiaeCS2108 的基因組圈圖見圖2A(見封三)，質粒圈見圖2B(見封三)，基因組序列登錄號分別為CP123969、CP123970。

表1 無乳鏈球菌S.agalactiae CS2108基因組特征

2.3 基因組功能分析

2.3.1 GO注釋結果

菌株S.agalactiaeCS2108共有1 640個基因注釋到GO 功能，占基因組的77.73%，結果如圖3(見封三)，其中與生物過程相關的基因有843個，與細胞組成相關的基因有807 個，與分子功能相關的基因有1 369 個。在生物過程中，與轉錄和翻譯相關的基因排前2 位(共109 個)，在細胞組分中，排前2 位的基因是與質膜和細胞質成分相關的基因（共642 個）；在分子功能方面，排前2 位的基因是與ATP 結合、DNA 結合相關的基因（共440個）。

2.3.2 COG注釋結果

菌株S.agalactiaeCS2108共有1 820個基因注釋到COG 功能，占基因組的86.26%，結果如圖4(見封三)，其中有462個未知功能的基因，功能已知的基因有173個與復制、重組和修復相關，有157 個基因與碳水化合物運輸和代謝相關，有148 個基因與翻譯、核糖體結構和生物發生相關，有138 個基因與氨基酸運輸和代謝相關，有115個基因與轉錄相關。

2.3.3 KEGG注釋結果

菌株S.agalactiaeCS2108共有1 213個基因注釋到功能信息，占基因組的57.49%，結果如圖5(見封三)，參與細胞過程的基因有72個，參與環境信息處理的基因有197個，參與遺傳信息處理的基因有157個，參與人類疾病的基因有68 個，參與代謝的基因有959 個，參與生物系統的基因有35 個。在KEGG 代謝通路中代謝相關的基因最多，在代謝途徑中排前三的分別是全局圖與概覽圖、碳水化合物代謝和氨基酸代謝。

2.3.4 毒力因子預測

將菌株S.agalactiaeCS2108 的基因組與毒力因子數據庫（Virulence factors database，VFDB）進行比對，得到毒力因子基因注釋的概況，結果表明有242 個毒力因子基因得到注釋，分別是68個防御性毒力因子基因（Defensive virulence factors）（圖6-A）、48 個非特異性毒力因子基因（Nonspecific virulence factor）（圖6-B）、88 個攻擊性毒力因子基因（Offensive virulence factors）（圖6-C）和20個毒力相關的調控基因（圖6-D）。在防御性毒力因子基因中，最多的是抗吞噬基因有30 個，其次是血清抗性相關的基因22 個，脅迫蛋白相關的基因11個，最少的是補充蛋白酶基因1 個；在非特異性毒力因子中，最多的是鐵離子吸收系統基因有42個，胞外酶和鎂離子吸收相關基因分別有3 個；在防御性毒力因子基因中，最多的是毒素相關的基因有36 個，粘附相關的基因有31 個，分泌系統和侵染相關的基因分別是14 個和7 個；毒力相關的調控基因有20個。由上可知，影響毒力的因素是多類基因相互協調共同作用的結果。

圖6 S.agalactiae CS2108的毒力因子統計

2.3.5 耐藥基因預測

將菌株S.agalactiaeCS2108 的基因組與綜合的抗生素抗性基因數據庫（Comprehensive Antibiotic Resistance Database，CARD）進行比對，獲得基因組中包含耐藥基因的信息。由表2 可見，S.agalactiaeCS2108 的基因組中含有耐藥基因141 個，排在前五的耐藥基因是46 個大環內酯類抗生素基因（Macrolide antibiotic）、29 個四環素抗生素基因（Tetracycline antibiotic）、25 個氟喹諾酮類抗生素基因（Fluoroquinolone antibiotic）、17 個鹽酸克林霉素基因（Lincosamide antibiotic）和16 個糖肽抗生素基因（Glycopeptide antibiotic）。

表2 無乳鏈球菌S.agalactiae CS2108耐藥基因統計表

3 討論與結論

分子生物學技術的飛速發展及大量羅非魚無乳鏈球菌基因組序列的公布，促進了對無乳鏈球菌基因組特征、流行病學、病原學、檢測方法和防控措施等方面的研究，為有效預防和治理提供了理論依據。

無乳鏈球菌基因組平均大小在2.05 Mb 左右，S.agalactiaeCS2108 的基因組大小為2.19 Mb，略大于平均值。染色體基因編碼2 104個蛋白質，小的質粒編碼6 個蛋白質。相較于其他病原微生物的基因組和編碼的蛋白質數量都較少，其他基因序列的功能還需要進一步研究。基因組島是有橫向起源跡象的一部分基因組，每個基因組島與多種生物功能相關、與共生或病原機理相關、與生物體的適應性相關，S.agalactiaeCS2108 的基因組中含有8 個基因組島，包含基因數量在10～63 個不等，每個基因組島含有功能相似的基因；前噬菌體序列上往往存在一些功能基因（抗生素基因、毒力基因）增強細菌對環境的適應性或使細菌成為致病菌，S.agalactiaeCS2108 中包含2 個前噬菌體，使得成為溶源菌；CRISPR/Cas 系統是一種原核生物的免疫系統，用來抵抗外源遺傳物質的入侵，比如噬菌體病毒和外源質粒。它可以識別出外源DNA、沉默外源基因的表達，S.agalactiaeCS2108 中包含6 個CRISPR/Cas，在抵抗外來病毒或者DNA入侵過程中起著重要作用。

有1 820 個基因注釋到COG 功能注釋中，其中約25%（462 個）的基因功能是未知的，這表明在無乳鏈球菌-羅非魚互作中仍有較多的新基因功能有待深入研究；其他基因中與復制、重組、修復相關的占比最高，在宿主中病原菌要面臨宿主細胞的攻擊和吞噬，所以復雜的重組和修復系統非常重要；其次是與碳水化合物轉運和代謝相關的基因，因為碳水化合物能為病原菌提供碳源，保證病原菌在宿主中能夠存活下來[19]。

通過比較基因組、蛋白質組學和組織病理學等方法來研究強弱毒株的毒力差異，揭示了強毒株在羅非魚體內的定植和組織損傷規律[19]。影響無乳鏈球菌毒力強度的基因包括莢膜多糖、表面蛋白、溶血素、CAMP 因子、絲氨酸蛋白酶和透明質酸等相關基因[20-21]，但對毒力基因缺乏深入的研究，我們基因組測序中發現了20個與毒力相關的基因。

抗生素是防治魚類疾病的主要手段，隨著抗生素使用量的增多產生了越來越多的耐藥性，抗生素的耐藥性成為十分嚴峻的問題，不同的無乳鏈球菌對抗生素的敏感性不同。張新艷等從廣東分離的無乳鏈球菌的耐藥性相對較弱，對氨芐西林等28 種抗生素敏感，對復方新諾明等15 種抗生素不敏感[22]；霍歡歡等發現海南分離到的無乳鏈球菌耐藥性很強，大部分菌株對氨基糖苷類和喹諾酮類耐藥[23]。了解無乳鏈球菌的耐藥規律，對于指導合理使用抗菌藥、防止多重耐藥菌株流行和擴散、保障水產品質量安全具有重要意義。

無乳鏈球菌為了適應不同環境下宿主體內的生活環境，會朝著不同的方向變異，同時基因組也會發生相應的變化。環境因素在無乳鏈球菌致病過程中有促進作用，例如高溫、高亞硝酸鹽、高氨氮、低溶解氧及pH 值等[24]。無乳鏈球菌通過躲避吞噬細胞的免疫防御，部分細菌被黏附吞噬后可在巨噬細胞內存活、生長和繁殖，隨巨噬細胞在體內游走、侵染其他組織，還可透過血-腦屏障侵入腦組織，導致羅非魚呈現典型的臨床癥狀和組織病理學變化[25]。

病原、宿主與環境之間的平衡關系決定疾病發生與否，了解鏈球菌感染與傳播途徑，可以更深入地認識羅非魚、鏈球菌和環境之間的關系[12]。羅非魚攝食了攜帶無乳鏈球菌的食物進入胃腸組織后，以單個個體、成組或成鏈的形式附著在腸細胞的頂端邊緣，穿越上皮細胞，在胃腸道黏膜上富集和增殖，然后隨體液循環進入各組織內，最后突破血腦屏障，進入腦部，導致腦膜炎，破壞羅非魚的腦神經，進而出現轉圈圈的現象[25]。在實際的養殖環境下，高濃度的菌液環境難以實現，鏈球菌通過口腔進入羅非魚的胃和腸道，可能是自然條件下感染羅非魚的主要傳播途徑。對S.agalactiaeCS2108 菌株進行全基因組測序分析，將為羅非魚無乳鏈球菌的預防和防治提供理論支持。