點源煉油廢水處理系統微生物多樣性及代謝功能研究

宋佳宇， 秦 榕， 杜顯元， 吳慧君， 梁立寶

1.石油石化污染物控制與處理國家重點實驗室， 北京 102206

2.中國石油集團安全環保技術研究院有限公司， 北京 102206

3.中國石油大學(北京)化學工程與環境學院， 北京 102206

4.中國石油新疆油田公司百口泉采油廠， 新疆 克拉瑪依 834011

煉油廢水是一種高濃度有機廢水，具有水量大、污染物種類多、生物降解難度大等特點[1-3]. 目前，石化企業普遍采用隔油—浮選—生化“老三套”工藝處理此類難降解點源廢水[4]，但保障其出水水質穩定達標GB 31570—2015《石油煉制工業污染物排放標準》尚有困難. 此外，原油重質化、劣質化情況加劇，導致點源煉油廢水處理系統下游生物處理單元沖擊風險加大，污水處理廠無法實現長周期穩定運行[5].

作為生物處理單元的保護屏障，上游一級處理單元對保障廢水處理系統穩定高效運行具有重要意義[6-7]. 目前，已報道的研究多集中于對煉油廢水處理系統生物處理單元活性污泥群落結構變化的評價[8-12]，如石油類物質對活性污泥沉降性能影響及其干擾處理效果評價[13]，以及硫化物對活性污泥中亞硝酸鹽氧化細菌(NOB)和氨氧化細菌(AOB)的抑制作用[14]等. 關于煉油廢水處理系統沿程微生物群落結構變化及代謝功能的研究卻鮮見報道.

目前，國內外研究者主要采用熒光原位雜交技術(FISH)、變性梯度凝膠電泳技術(DGGE)、末端限制性片段長度多態性技術(T-RFLP)、克隆文庫(Clone library)和高通量測序等技術[15-18]開展微生物群落結構研究，其中，高通量測序技術因具有測序通量高、準確、全面等優勢[19]，已成為微生物生態結構和功能解析的強大技術支持手段.

綜上，該文以某石化污水處理廠點源煉油廢水為研究對象，開展廢水微生物多樣性及群落結構研究，明確點源煉油廢水處理系統(簡稱“系統”)微生物群落動態行為，闡明環境因子與廢水微生物群落的響應關系，預測其代謝功能，以期為現行廢水處理工藝綜合診斷及優化奠定理論基礎.

1 材料與方法

1.1 廢水處理系統處理效果及采樣信息

某石化企業污水處理廠運行穩定，COD、TOC (總有機碳)、TN、NH3-N、石油類等水質排放指標達到GB 31570—2015《石油煉制工業污染物排放標準》排放標準. 該污水處理廠一級處理單元各出水口設置采樣點，采樣點布設如圖1所示，水樣依采樣順序用英文字母編號(分別為A、B、C、D、E). 每個樣品采集3個平行樣，共計45個樣品. 樣品采集后，0～4 ℃冷藏運輸至實驗室，48 h內完成水質指標檢測分析. 各單元主要污染物去除貢獻率如表1所示. 由表1可知，該系統沿程主要污染物的去除效果顯著，去除率均在90%以上.

表1 各單元主要污染物去除貢獻情況

圖1 點源煉油廢水處理系統采樣點分布

1.2 水質指標及測試方法

水質指標分為無機指標、有機指標和金屬指標3類(見表2)，其中，無機指標包括pH、DO、ORP (氧化還原電位)、電導率、TN、NH3-N、TP ；有機指標包括COD、TOC、石油類；金屬指標包括Cd、Cu、Ni、Ba、Cr、Mn、Fe、Mg、Se、As、Pb、Al、Ca、K、Zn、Na. 分析測試方法參照《水和廢水監測分析方法》(第四版)[20].

表2 水質指標及測試方法

1.3 高通量測序方法

1.3.1樣本前處理

采用無菌微孔膜(0.22 μm，美國密理博公司)及無油式真空泵(Rocker 300，臺灣洛科儀器有限公司)抽濾1 L水樣，所得濾膜即為微生物菌群樣品，立即放入-80 ℃冰箱保存.

1.3.2DNA抽提、PCR擴增及Pair End文庫構建

根據DNeasy?96 PowerSoil?Pro QIAcube?HT Kit 說明書進行微生物群落總DNA抽提，使用超微量分光光度計(NanoDrop2000，美國賽默飛世爾科技公司)測定DNA濃度和純度(1%的瓊脂糖凝膠). 確定DNA濃度及純度后，微生物多樣性分析通過引物338F (5′-ACTCCTACGGGAGGCAGCAG-3′)和806R(5′-GG ACTACHVGGGTWTCTAAT-3′)對16S rRNA基因V3～V4可變區進行PCR擴增；宏基因組測序通過非接觸式超聲波破碎儀(Covaris M 220，美國Covaris公司)將DNA片段化，篩選約400 bp的片段，用于構建PE文庫.

1.3.3上機測序

微生物多樣性采用Miseq PE300(Illumina，美國Illumina公司)測序平臺進行測序；宏基因組采用NovaSeq 6000 (Illumina，美國Illumina公司)測序平臺測序. 上述測序平臺由上海美吉生物醫藥科技有限公司提供.

2 結果與討論

2.1 系統沿程微生物多樣性變化情況

Alpha多樣性指數[21]是用于評價微生物群落多樣性和豐富度的綜合指標，其中，Shannon-Wiener和Simpson指數反映群落多樣性[22]，Ace和Chao指數反映群落豐富度. Shannon-Wiener指數越大，說明群落多樣性越高；與之相反，Simpson指數越大，說明群落多樣性越低.

測序分析結果如表3所示. 由表3可見，所有樣品Coverage均在0.998以上，表明測序樣品微生物覆蓋率高，滿足菌群多樣性分析需求. 沿程微生物多樣性分析表明，Shannon-Wiener指數基本呈逐漸上升趨勢，表現為A

表3 水樣Alpha多樣性分析

2.2 系統沿程微生物群落結構變化情況

為明確系統沿程微生物群落結構變化情況，開展各單元出水樣品組間相似性分析及群落組成變化分析. 由相似性分析(analysis of similarities, ANOSIM)結果可知，不同處理單元出水群落結構具有顯著性差異(見圖2, Bray-Curtis ANOSIM=0.388 1,P<0.01)，表明處理工藝的不同對微生物群落結構影響較大.

注： Ave表示組間距離排名的平均值.

系統各單元出水門水平微生物菌群結構分布情況如圖3所示. 進水中變形菌門(Proteobacteria)相對豐度最高(90.55%)，與Shen等[24]研究結果一致. 變形菌門屬革蘭氏陰性異養細菌，是石油化工廢水中最常見的微生物之一[25-26]，具有高效降解石油烴化合物的能力. 該研究中，Proteobacteria在整個處理過程中亦占主導地位.

圖3 門水平上微生物群落結構組成

系統沿程各單元出水優勢菌組成均發生變化，其中，調節池、渦凹氣浮池、斜板氣浮池出水中Proteobacteria相對豐度分別降至54.90%、62.59%、52.89%，Epsilonbacteraeota相對豐度分別增至17.03%、29.60%、41.53%，厚壁菌門(Firmicutes)相對豐度分別增至26.13%、6.02%、2.33%. 已有研究[27-29]表明，Epsilonbacteraeota和厚壁菌門具有深度降解廢水中難降解有機物的能力. A/O池出水中Proteobacteria相對豐度降至39.69%，髕骨細菌門(Patescibacteria)和擬桿菌門(Bacteroidetes)相對豐度分別增至23.95%、15.65%. 其中，髕骨細菌門為難降解有機廢水處理工藝中常見優勢菌[30-31]；擬桿菌門具有較強的營養物質代謝能力[32]，為污水脫氮的重要參與者. 生物單元對系統TN和NH3-N去除貢獻率(分別為72.202%和28.750%)均較高，表明擬桿菌門在該單元脫氮過程中發揮重要作用.

2.3 微生物群落多樣性與環境因子響應分析

微生物群落多樣性與環境因子密切相關[33-35]. 進行環境因子分析前，使用VIF(variance inflation factor analysis，方差膨脹因子分析)對3類環境因子指標(14種)進行環境因子篩選，保留多重共線性較小的環境因子進行后續研究. 篩選后的環境因子(VIF<10)分別為COD(2.954)、TN(4.983)、TP(1.769)、pH(1.491)、NH3-N(6.555)、Na(2.878)和Ca(2.339)(見表4).

表4 篩選后環境因子的VIF值

通過db-RDA分析(基于Bray-Curtis距離)揭示了樣品與篩選后環境因子的相關性(見圖4). 在db-RDA 二維圖中，帶箭頭的線段表示環境因子，樣點到數量型環境因子的箭頭投影長度反映該環境因子對樣品群落組成的影響. 結果表明，篩選后的環境因子均與A/O出水呈負相關，COD、NH3-N、TN濃度均與進水微生物群落結構呈正相關. 其中，COD(r2=0.599,P<0.05)、NH3-N(r2=0.574,P<0.05)濃度是影響廢水中微生物群落變化的最顯著因素.

圖4 環境因子與微生物群落基于db-RDA分析的相關性分析結果

為明確煉油廢水處理效果與微生物群落結構的響應關系，選取顯著相關的環境因子(COD、NH3-N)去除率與門水平優勢物種(相對豐度前10)進行相關性分析(見圖5). 結果表明，NH3-N去除率與Bacteroidetes (r=0.588,P<0.05)呈顯著正相關，COD去除率與Patescibacteria (r=0.530,P<0.05)和Bacteroidetes (r=0.706,P<0.01)均呈顯著正相關，這與陳重軍等[36-37]研究結果一致.

注： *表示0.01≤P≤0.05，** 表示0.001≤P≤0.01.

2.4 微生物代謝功能預測分析

基于宏基因組測序結果，采用PICRUST軟件與KEGG基因數據庫[38-39](Kyoto encyclopedia of genes and genomes, http://www.genome.jp/kegg)和通路數據庫進行對比分析，預測各樣品功能基因組成情況. KEGG通路數據庫主要包括6類功能基因，分別為新陳代謝、遺傳信息處理、環境信息處理、細胞過程、人類疾病和機體系統. 該研究中，各樣品中相對豐度最高的功能基因為新陳代謝功能基因，占比均高于60%(A、B、C、D、E樣品分別為68.37%、63.67%、63.20%、63.97%、62.59%). 其中，代表有機物去除的新陳代謝功能基因為能量代謝和外源化合物生物降解[40]，因此將其作為該研究的重點.

從能量代謝角度分析，沿程各單元能量代謝基因相對豐度大小表現為A>C>B>D>E(見圖6)，推測由于易降解物質逐級耗盡導致后續單元能量代謝基因相對豐度降低. 從外源化合物生物降解角度分析，一級處理工藝相對豐度波動較小，均在3.00%以下；進入生物處理單元后，其相對豐度迅速提升至15.38%，其中，注釋到石油烴降解代謝通路的功能基因相對豐度顯著增加，如苯甲酸降解基因(相對豐度為3.15%)、二甲苯降解基因(相對豐度為1.18%)、硝基甲苯降解(相對豐度為2.15%)等. 另外，生物處理單元對系統各主要污染物去除貢獻率較大(見表1)，進一步證實煉油廢水一級處理工藝無法完全去除廢水中難降解有毒有害組分，有毒物質的長期累積可能增加生物處理單元沖擊風險. 據此，亟待開展點源煉油廢水處理系統沿程廢水生物毒性削減負荷核算及關鍵致毒物質識別工作.

圖6 KEGG外源生物降解和能量代謝

3 結論

a) 系統沿程微生物多樣性Shannon-Wiener指數和豐富度Ace指數基本呈逐漸上升趨勢，表現為進水<渦凹氣浮池<斜板氣浮池<調節池

b) ANOSIM分析揭示了系統沿程微生物群落結構差異顯著，各單元出水優勢菌組成均發生變化. 其中，調節池、渦凹氣浮池、斜板氣浮池出水中Proteobacteria占比分別由90.55%降至54.90%、62.59%、52.89%，Epsilonbacteraeota占比分別增至17.03%、29.60%、41.53%，厚壁菌門(Firmicutes)占比分別增至26.13%、6.02%、2.33%；A/O池出水中Proteobacteria占比降至39.69%，髕骨細菌門(Patescibacteria)和擬桿菌門(Bacteroidetes)占比分別增至23.95%、15.65%.

c) 環境因子與微生物群落相關性分析表明，COD、NH3-N在系統微生物群落形成過程中起重要作用(P<0.05)，其中，NH3-N去除率與Bacteroidetes呈顯著正相關(P<0.05)，COD去除率與Patescibacteria(P<0.05)、Bacteroidetes(P<0.01)均呈顯著正相關.

d) KEGG功能注釋分析可知，生物處理后外源化合物生物降解相關基因占比迅速增加，其中，與石油烴降解代謝通路相關的功能基因豐度增加尤為明顯.