張曉鎮(zhèn)，陸志波，季樹勛，夏珊珊，王 穎

（1. 同濟大學環(huán)境科學與工程學院，上海 200092；2. 寧波市水利發(fā)展規(guī)劃研究中心，浙江 寧波 315000）

河流是淡水生態(tài)系統(tǒng)的重要組成部分，河道生態(tài)系統(tǒng)中的浮游微生物是河道生物多樣性的重要組成部分，且與河流中溶解氧、pH、硝態(tài)氮、總磷和重金屬等眾多環(huán)境因子息息相關[1-3]，做好河流微生物多樣性研究是河流健康評價、水系連通與生態(tài)修復、供水安全保障等工作的重要基礎。

近 30 年來，各國形成了側重點不同、各具特色的評價方法，如，美國環(huán)境保護局（US EPA）提出的快速生物監(jiān)測協(xié)議（Rapid Bio-assessment Protocols，RBP）[4]、歐盟的《水框架法令》(WFD)[5]等都已在河流水生態(tài)系統(tǒng)健康管理中得到應用。2020 年世界自然基金會聯(lián)合我國8 家科研院所共同發(fā)布《長江生命力報告2020》，提出長江生命力指數(shù)。這些評價方法均將水生態(tài)作為重要維度。

微生物作為有機化合物的主要分解和礦化者[6-7]，在食物鏈和食物網中所處的營養(yǎng)級比大型水生動植物更低[8]，水體微生物網絡的復雜性一定程度上能反應水體中大型生物的網絡復雜性與社區(qū)對擾動的抵抗力和恢復力[9]。微生物群落結構、豐度等也可以作為有機物分解、反硝化功能的指標[10]。

近一個世紀以來，地表水微生物安全主要是通過檢測糞便污染的細菌指標來確定[11]，但除致病菌外，在地表水微生物監(jiān)測中污染指示微生物更能反應水體受污染的現(xiàn)狀[12]。有研究者[13]提出可以創(chuàng)建可靠、標準化的參考數(shù)據(jù)庫，作為指示健康微生物群落的基線指標。通過相關性計算可以分析環(huán)境因子與物種分布間的相關關系[14]。已有研究發(fā)現(xiàn)微生物多樣性指數(shù)與總磷和磷酸鹽含量呈極顯著負相關關系[15]。SOUFFREAU et al[16]也發(fā)現(xiàn)總磷和氨氮對湖泊的微生物群落組成有顯著影響。

微生物多樣性研究的傳統(tǒng)方法是微生物分離純化培養(yǎng)，存在如培養(yǎng)周期長、分類學鑒定不準確等弊端[17]。分子生物學技術可以高效、準確地從遺傳物質水平上研究環(huán)境中微生物的群落結構[18]。高通量測序技術是近年來微生物群落分析的常用技術[19-21]，具備測序快速、通量高等優(yōu)勢，且可對多樣品的多個可變區(qū)同時測序[22]。通過高通量測序技術研究三江干流水體微生物群落結構及其多樣性，可以為構建三江流域生態(tài)系統(tǒng)管理方案提供科學依據(jù)。

寧波市的余姚江、奉化江和甬江合稱“三江”。三江口由于其獨特的自然地理條件，成為寧波市的重要標志[23-24]。三江河道的治理方案主要為水利工程建設[20]，需要對河道水生態(tài)進行長期綜合監(jiān)測以提供數(shù)據(jù)支撐，極有必要對其中微生物群落結構與多樣性展開全面研究。當前對于寧波三江流域微生物群落結構與多樣性的研究較為缺乏，本文開展的研究說明通過高通量測序的方法開展河流的生態(tài)指標監(jiān)測具有可行性。

1 材料與方法

1.1 樣品采集

此次三江干流巡航測量范圍為余姚江（牟山閘/永思橋～蜀山大閘～姚江大閘～三江口）、奉化江（三江口～澄浪堰）、甬江（三江口～入?？冢?。航測總長度約105 km。以下稱三江干流均為上述范圍。巡航測量時間為2021 年5 月11～14 日。測量方式為使用航測船只沿余姚江、奉化江和甬江深泓航道線單向巡航一次。途中在三江干流12 處點位停船進行水體取樣測量。測量斷面及取樣點位置分布，見表1 和圖1。

表1 三江干流航測斷面及取樣點位置分布

圖1 寧波三江干流采樣點分布

取樣過濾裝置含真空泵或蠕動泵、六聯(lián)過濾器和直徑為47 mm、孔徑0.45 μm 的混合纖維素酯膜。

1.2 生物信息提取與分析

DNA 提取和測序在Personal Biotechnology Co. Ltd.（中國上海）進行。按照制造商的說明，通過DNeasy PowerSoil Kit（QIAGEN，Hilden，Germany）提取水樣中總微生物基因組 DNA。通過瓊脂糖凝膠電泳和 NanoDrop ND-1 000 分光光度計 (Thermo Fisher Scientific, Waltham, MA, USA) 評估提取的DNA 的質量和數(shù)量。

使用引物對515F_926R 和18SV4F_18SV4R來分別對 16S rRNA 基因、18S rRNA 基因的可變區(qū)進行PCR 擴增。

采用Illumina HiSeq X-10（Illumina，San Diego，CA，USA）平臺對群落DNA 片段進行雙端（Pairedend）測序。使用DADA2 方法進行去引物，質量過濾，去噪(denoise)，拼接和去嵌合體等步驟。用QIIME2 (2 019.4)進行物種分類學注釋。將有效數(shù)據(jù)在97%的相似水平上劃分為不同的運算分類單位(operational taxonomic unit, OTU), 將OTU 代表序列與相應的微生物數(shù)據(jù)庫比對, 利用Greengene數(shù)據(jù)等數(shù)據(jù)庫對物種進行注釋, 從而得到每個樣本所含的物種信息。通過對OTUs 進行Alpha 多樣性指數(shù)、豐度等分析, 并對物種注釋在各個分類水平上進行群落結構的統(tǒng)計分析, 得到細菌群落的組成信息, 再通過Beta 多樣性分析, 系統(tǒng)分析三江干流不同流域樣本間的微生物群落結構差異。

使用Alpha 多樣性來衡量局域均勻生境下的物種數(shù)目。α 多樣性使用Hill’s number family 多樣性計算[25]，主要包括Chao 指數(shù)；Faith’s pd 指數(shù)；Simpson 指數(shù)；Shannon 指數(shù)；Good’s coverage 指數(shù)；Pielou’s evenness 指數(shù)；Observed species 指數(shù)。

Chao1 是基于豐度的估計量，計算見式（1）[26]：

式中，Schao1為預測豐富度，Sobs為觀測到的物種數(shù)量，n1為單序列OTUs 數(shù)目，n2為雙序列OTUs 數(shù)目。

Faith’s pd 指數(shù)反映了譜系多樣性，等于生命發(fā)育樹（The Tree of Life Web）的分支數(shù)量[27]。

Simpson 指數(shù)是衡量個體在分類時的集中程度，見式（2）[25]：

式中，Sobs為觀測到的OTUs 數(shù)量，ni為第i個OUT 對應的個體數(shù)目，N為群落中所含的總個體數(shù)目（所有 OTU 的豐度之和）。

Shannon 指數(shù)量化了預測從數(shù)據(jù)集中隨機抽取的個體物種身份的不確定性，見式（3）[28]：

式中，Sobs為觀測到的OTUs 數(shù)量，ni為第i個OUT對應的個體數(shù)目，N為群落中所含的總個體數(shù)目。

Good’s coverage 指數(shù)用于評估測序結果代表的community 物種的總數(shù)[29]，見式（4）：

式中，F(xiàn)是單序列 OTU 的數(shù)量，N是群落中所含的總個體數(shù)目。

Pielou’s evenness 指數(shù)用于描述物種中的個體的相對豐富度或所占比例[30]，見式（5）：

式中，H’為Shannon 指數(shù)的值，H’max為H’的最大值。

Observed species 指數(shù)為群落中豐度＞0 的物種數(shù)之和，可以反映物種種類信息。

Beta 多樣性是指沿著環(huán)境梯度變化的不同群落之間，物種組成的相異性或物種沿環(huán)境梯度的更替速率，因此也被稱為生境間多樣性(betweenhabitat diversity)。β多樣性分析所使用的Braycurtis 距離基于OTUs 的計數(shù)統(tǒng)計，比較2 個群落微生物的組成差異，見式（6）[31]：

式（6）中，SA,i和SB,i表示第i個OTU 分別在A 群落和B 群落中的計數(shù)。min 表示取兩者最小值。

2 結果與討論

2.1 三江干流微生物群落組成

通過Silva 數(shù)據(jù)庫對所得數(shù)據(jù)進行比對，根據(jù)序列物種分類學注釋的結果對物種組成進行分類，展示相對豐度最高的10 個分類，其余合并顯示為others。統(tǒng)計樣本的物種注釋結果中門(Phylum) 、綱(Class) 、目(Order) 、科(Family) 和屬(Genus)水平上含有的分類單元的數(shù)量。

對于原核微生物，在門水平上，變形菌門(Proteobacteria)的平均豐度最高（圖2a），在C1、C3、C4 和PC2 點位處相對豐度均超過60%，但在靠近入?？诘腃13 點位處擬桿菌門(Firmicutes)的相對豐度最高，達到42.8%。需要注意在C4、C5 和C7 等采樣水域厚壁菌門(Firmicutes) 表現(xiàn)出較高豐度，說明已經存在污水排放導致的DO 含量下降，存在厭氧環(huán)境[32]。在綱水平上（圖2b），γ-變形菌綱(Gammaproteobacteria) 的平均豐度最高，其中在C1-C5、PC1、PC2 的7 個點位處相對豐度均超過50%，在臨近入?？诘腃13 處γ-變形菌綱的豐度明顯下降，擬桿菌綱(Bacteroidia)、梭菌綱(Clostridia)的相對豐度呈現(xiàn)上升趨勢。在門水平與綱水平上，優(yōu)勢菌種與其他研究中的結果一致[33-34]。在目水平上（圖2c ），在余姚江河段，假單胞菌目(Pseudomonadales)平均相對豐度較高，達到46.4%，在奉化江和甬江河段異常球菌目(Deinococcales)相對豐度可達到20%以上，在C9 處達到最高值42.7%。在科水平上（圖2d），莫拉菌科(Moraxellaceae)在余姚江河段各采樣點處都表現(xiàn)出較高的豐度，在C2-C5 點位處豐度均高于40%，尤其在C4 處達到54.5%，異常球菌科(Deinococcaceae)在甬江流域也表現(xiàn)出較高豐度，在C9 處達到最高值42.7%。在屬水平上，不動桿菌屬(Acinetobacter)的平均豐度最高，異常球菌屬(Deinococcus)在甬江流域也表現(xiàn)出了較高的豐度，最高在C9 支流交匯處，達到42.7%。最豐富的細菌種類包括奇異球菌(Deinococcus grandis)(9.44%)、魯氏不動桿菌(Acinetobacter lwoffii) (1.18%)、波西米亞不動桿菌(Acinetobacter bohemicus) (1.09%)等。

圖2 原核生物群落結構組成

對于真核微生物，在門水平上，綠藻門(Chlorophyta)的平均豐度最高（圖3a），在C2、C6 和PC2 點位處相對豐度均超過40%，在C6、PC1、PC2 和C12 等采樣點處子囊菌門(Ascomycota) 也表現(xiàn)出較高的相對豐度，在PC1 點位處子囊菌門(Ascomycota)的相對豐度最高，達到65.5%。在綱水平上（圖3b），綠藻綱 (Chlorophyceae)的平均豐度最高，為25.6%，其次為引藻綱(Cryptomonadales)(12.7%)、甲藻綱(Dinophyceae)(7.8%)等分類。在臨近入?？诘腃12、C13 處綠藻綱的豐度明顯下降，擬桿菌綱(Bacteroidia)、梭菌綱(Clostridia)的相對豐度呈現(xiàn)上升趨勢。在目水平上（圖3c），在C1～C6， 環(huán) 藻 目 (Sphaeropleales)、 海 鏈 藻 目(Thalassiosirales)相對豐度較高，達到46.4%，在C7、C10 處，隱藻目(Cryptomonadales)相對豐度分別達到85.0%、69.1%，表現(xiàn)出顯著優(yōu)勢。在科水平上（圖3d），隱鞭藻科(Cryptomonadaceae) 平均豐度達到13.2%，在C7、C10 處表現(xiàn)出顯著的高豐度，尤其在C7 處達到85.0%，但在C1～C5 各點位處其豐度并不顯著。在屬水平上，隱藻屬(Cryptomonas)的平均豐度最高，主要在C7、C10 處表現(xiàn)出顯著的高豐度，小環(huán)藻屬(Cyclotella)在C1 處相對豐度達到46.6%。最豐富的真核微生物種類包括彎隱藻(Cryptomonas curvata) (6.85%)、寄生索囊藻(Choricystis parasitica)(0.79%)和馬索隱藻(Cryptomonas marssonii) (0.62%)等。

圖3 真核生物群落結構組成

2.2 微生物群落Alpha 多樣性分析

以Chao1 和Observed species 指數(shù)表征豐富度，以Shannon 和Simpson 指數(shù)表征多樣性，以Faith’s PD 指數(shù)表征基于進化的多樣性，以Pielou’s evenness 指數(shù)表征均勻度，以Good’s coverage 指數(shù)表征覆蓋度，微生物群落Alpha 多樣性分析結果，見圖4。

圖4a 可知，對于原核微生物，YYJ 樣品的Chao1 和Observed species 指數(shù)最高， FHJ 和YJ 樣品的Chao1 和Observed species 指數(shù)值差異不大，說明余姚江流域的群落豐富度最高，Chao1 值在3 000 以上。YYJ 樣品的Shannon 指數(shù)值相對較高，但3 組樣品整體差異不顯著，Simpson 指數(shù)基本穩(wěn)定在0.98 以上，說明三江干流流域各采樣點原核微生物群落多樣性較豐富，穩(wěn)定性較好，其中余姚江流域的原核微生物多樣性相對更高。YYJ 流域樣品的Faith’s PD 指數(shù)較小，而YJ 流域樣品Faith’s PD指數(shù)較大，F(xiàn)HJ 居中，說明三江干流自上游至下游原核微生物的親緣關系逐漸復雜化，基于進化的多樣性較高。YYJ 流域樣品的Pielou’s evenness 指數(shù)較高，F(xiàn)HJ 其次，YJ 最小，說明3 個流域的物種均勻度自上游向下游逐漸減低。Good’s coverage 指數(shù)值均在0.97 以上，說明測序對群落中物種的覆蓋度較高，未被檢測出的物種所占的比例很少。與瀾滄江[35]、科希河[36]和博斯騰湖流域[37]等流域相比，寧波三江流域原核微生物群落豐富度與多樣性更高。

圖4b 可知，對于真核微生物，F(xiàn)HJ 樣品的Chao1和Observed species 指數(shù)最高， YJ 其次，YYJ 樣品最低，說明奉化江采樣監(jiān)測流域的群落豐富度最高。FHJ 樣品的Shannon 指數(shù)值相對較高，67%的采樣點的Simpson 指數(shù)在0.9 以上，3 組樣品Simpson 指數(shù)值無顯著差異，說明三江干流流域各采樣點真核微生物群落多樣性較豐富，穩(wěn)定性較好，其中FHJ 流域的真核微生物多樣性相對更高。YYJ 流域樣品的Faith’s PD 指數(shù)較小，而YJ 流域樣品Faith’s PD 指數(shù)較大，F(xiàn)HJ 居中，說明三江干流自上游至下游真核微生物的親緣關系與原核微生物的發(fā)展一致，也呈現(xiàn)逐漸復雜化的趨勢，基于進化的多樣性較高。FHJ 流域樣品的Pielou’s evenness 指數(shù)較高，YJ 其次，YYJ 最小，說明3 個流域的物種均勻度自上游向下游逐漸升高。Good’s coverage 指數(shù)值均在0.999 以上，說明測序對群落中物種的覆蓋度較高，未被檢測出的物種所占的比例極少。與九龍江[38]、西安渭河、灞河和浐河[39]等流域相比，寧波三江流域真核微生物豐富度更低，多樣性相近。

圖4 Alpha 多樣性分析

2.3 微生物群落Beta 多樣性分析

根據(jù)Bray-curtis 距離算法，對物種進行非度量多維尺度(non-metric multidi-mensional scaling，NMDS)分析(圖5)，發(fā)現(xiàn)YYJ 的樣品在維度NMDS1 和維度NMDS2 能夠很好地與FHJ、YJ 區(qū)分開，說明余姚江流域樣品的原核、真核微生物種類與FHJ、YJ 有較為明顯的區(qū)分性。其中Stress 值分別為0.057 1、0.107，均＜0.2，很好地證明了數(shù)據(jù)的可靠性。

圖5 基于相似性＞95%的OTU 水平上非度量多維度分析

采用層次聚類(hierarchical clustering)的分析方法，以等級樹的形式展示樣本間的相似度，通過聚類樹的分枝長度衡量聚類效果的好壞。對Braycurtis 距離矩陣采用UPGMA 算法（即聚類方法為average）在門水平上進行聚類分析，見圖6?？梢钥闯?，三江干流的YYJ、FHJ 和YJ 的3 部分采集的樣品中的原核微生物在不同河段呈現(xiàn)出鮮明的聚類效果，說明三江干流自上游向下游原核微生物在門水平上有比較明顯的變化，不同流域的原核微生物相似度較高，入?？贑13 處則展示出完全不同的群落組成。位置相近的采樣點的真核微生物的聚類也比較相似，如干流上游方位的C1～C5 點位處樣本之間的相似度很高，但與相對處于下游方向的點位之間的相似度較低。在支流匯入后的樣品的群落構成會出現(xiàn)比較明顯的變化，如C9 采樣點與其他點位的相似度就相對較低。

圖6 門水平上的層次聚類分析

2.4 物種組成熱圖

根據(jù)測序結果和注釋OTU 數(shù)據(jù)進行樣本間豐度和相似性聚類，構建寧波三江干流微生物屬層次物種分類Heatmap 圖譜，以顏色梯度變化表現(xiàn)豐度值的高低，見圖7。使用平均豐度前30 位的屬的豐度數(shù)據(jù)繪制熱圖。樣本按照物種組成數(shù)據(jù)的歐式距離(euclidean distance)進行UPGMA聚類。

圖7 三江干流微生物屬層次物種分類Heatmap

原核微生物在屬上的Heatmap 圖譜顯示，YYJ 流域樣品中微小桿菌屬(Exiguobacterium) 的豐度最高，次為馬賽菌屬(Massilia)、黃桿菌屬(Flavobacterium)、芽孢桿菌屬(Bacillus)和不動桿菌屬(Acinetobacter)等；FHJ 流域樣品中豐度最高的屬為新鞘脂菌屬(Novosphingobium)、砂單胞菌屬(Arenimonas)和短波單胞菌屬 (Brevundimonas)等；YJ 流域樣品中分支桿菌屬(Mycobacterium)、鞘脂菌屬(Sphingobium)和異常球菌屬(Deinococcus)等豐度更高。原核乳桿菌屬(Lactobacillus)、鼠桿菌屬(Muribaculaceae)等只在入?？贑13 處測得，該處組成與其余點位有明顯的區(qū)別。

真核微生物在屬水平上的Heatmap 圖譜顯示，YYJ 流域樣品中綠藻屬(Pseudomuriella)、麥可屬(Mychonastes) 和海鏈藻屬(Thalassiosira)的豐度最高， FHJ 樣品中囊藻屬(Tetracystis)、 扁藻屬(Tetraselmis)和直鏈藻屬(Aulacoseira)表現(xiàn)出較高的豐度，YJ 樣品中尖尾藻屬(Teleaulax)、隱藻屬(Cryptomonas)相對豐度較高。入?？贑13 采樣點處測得銅綠微囊藻屬(Hererocapsa)、小二仙草屬(Haloragis)、劍水蚤屬(Mesocyclops)和赤潮異彎藻屬(Heterosigma)等僅在此點位顯示出較高豐度的真核生物屬。

3 結論與展望

3.1 結論

通過高通量測序，發(fā)現(xiàn)三江干流優(yōu)勢菌門主要有變形菌門、放線菌門、異常球菌-棲熱菌門和擬桿菌門，真核微生物優(yōu)勢門包括子囊菌門、甲藻門、綠藻門和硅藻門等，這些均為典型淡水微生物門類，與其他水域微生物群落組成檢測結構相似[40-42]。如劉峰等[43]發(fā)現(xiàn)汾河入黃口水體中變形菌門、放線菌門和厚壁菌門處于優(yōu)勢地位，在湖泊水體中變形菌門、藍細菌門、擬桿菌門、放線菌門和疣微菌門也是最主要的5 個原核生物門類[44]，這些典型的淡水菌門在淡水生態(tài)系統(tǒng)中的廣泛分布已經得到研究證實[45-46]。在三江干流下游異常球菌-棲熱菌門豐度的增加可能意味著該河段受到較多的熱污染影響[45]。在西安渭河、灃河、灞河和浐河的研究中，綠藻和真菌界子囊菌門也是在各類群中的優(yōu)勢物種[39]，與寧波三江流域真核微生物優(yōu)勢門一致，個別點位處（如C4、C9 等處）甲藻門豐度較高可能與溫度的升高或pH 的降低有關[46]，與原核微生物的優(yōu)勢菌門中異常球菌-棲熱菌門的出現(xiàn)共同說明寧波三江干流受到熱污染的可能性很高。

同時通過Beta 多樣性分析也注意到，不同流域微生物群落組成表現(xiàn)出明顯的區(qū)分性，如異常球菌-棲熱菌門、子囊菌門在三江干流的中下游水域豐度更高，而在上游的采樣點C1～C5 處的樣品中豐度較低。厭氧性厚壁菌的出現(xiàn)說明存在工業(yè)或生活廢水的排入，已經造成水生態(tài)環(huán)境的破壞[47]。入?？谔幩虻奈⑸锝M成也與其他采樣點處明顯不同。由物種組成熱圖發(fā)現(xiàn)，微生物組成在相鄰采樣點間相似性較高，YYJ、FHJ 和YJ 3 組樣品中也有不同的優(yōu)勢物種，可以作為流域水環(huán)境狀況的指示生物。

Alpha 多樣性分析反映了寧波三江干流流域微生物多樣性的綜合指標，自三江干流上游向下游進行采樣檢測，各采樣點水域處的原核微生物豐富度指數(shù)和多樣性指數(shù)均較高，YYJ 樣品與FHJ 樣品、YJ 樣品相比較，原核微生物豐富度差異不大，YYJ流域樣品多樣性最高。真核微生物方面，三江干流流域整體原核微生物豐富度與其他流域研究結果相比偏低，多樣性維持在較高的水平，其中FHJ 樣品的物種豐富度和多樣性相對更高。自上游至下游微生物的親緣關系呈現(xiàn)逐漸復雜化的趨勢。

3.2 展望

本實驗結果僅就寧波三江干流流域微生物多樣性進行表征，水體微生物多樣性受到多種環(huán)境因子的影響，一方面水體微生物多樣性與環(huán)境因子之間的關系，不同水域所處的地理位置對于水體微生物群落的影響都值得進一步探索。

同時，水體的微生物作為重要的環(huán)境指示生物也可以作為水生態(tài)環(huán)境現(xiàn)狀的指示，通過研究流域微生物群落結構多樣性與群落結構，可以為構建完善的監(jiān)測與管理體系提供理論支持。在常態(tài)化生態(tài)監(jiān)測工作中微生物多樣性的監(jiān)測與分析工作也具備可行性。通過對微生物群落結構多樣性與其驅動因子的研究，可以開發(fā)寧波三江干流流域特殊功能菌群，為區(qū)域水環(huán)境治理提供針對性支持。