侯冠群,來劍斌,李 靜,劉 振,公華銳,王 兵,孫志剛,歐陽竹,侯瑞星, *

1 中國科學(xué)院地理科學(xué)與資源研究所,生態(tài)系統(tǒng)網(wǎng)絡(luò)觀測與模擬重點(diǎn)實(shí)驗(yàn)室, 北京 100101

2 中國科學(xué)院大學(xué)資源與環(huán)境學(xué)院, 北京 100049

3 中國地質(zhì)調(diào)查局煙臺海岸帶地質(zhì)調(diào)查中心, 煙臺 264004

濱海新生土地演替的研究有助于濕地資源的保護(hù),為改善當(dāng)?shù)氐纳鷳B(tài)環(huán)境提供科學(xué)依據(jù)[1—2]。了解新生土地生境更新、演替規(guī)律,還有助于評估人類活動的影響[3],使開發(fā)速度低于生境更新、演替的速度,或者引導(dǎo)和控制新生土地生境更新、演替的速度和方向[4]。

演替的監(jiān)測需要頻繁連續(xù)的時間序列。但是,初生演替的研究往往受到時空因素干擾或限制,通常采用時空替代法擬合時間序列,如冰川退縮跡地研究中通常采用樣線法[5],而土壤年齡的驅(qū)動力研究多采用樣地法[6]。在濱海三角洲演替的研究中,采樣和分析方法通常為“空間換時間”[7—8]或針對植被演替階段的樣地法[9—11]。雖然前人很早關(guān)注到河口三角洲土壤年齡對生態(tài)系統(tǒng)演替階段可能帶來的影響[8, 11—12],但缺少成土?xí)r間的驅(qū)動力的定量研究,這可能造成演替時間和空間的不匹配,進(jìn)而形成對演替過程的描述產(chǎn)生偏差。

黃河三角洲自1855年黃河奪大清河流路入渤海以來,黃河尾閭多次改道,形成多個呈疊瓦狀復(fù)式結(jié)構(gòu)的亞三角洲堆積體[13]。其中有4個亞三角洲堆積體與其他洲體的重疊程度較低,在地表形成相對獨(dú)立的沖積扇區(qū),成土?xí)r期分別為1904—1929年,1929—1934年,1964—1976年和1976年至今[14—15]。亞三角洲堆積體重疊時,表層土壤成土?xí)r期由晚近流路的沖淤決定[16],故相同亞三角洲堆積體的土壤具有相對一致的發(fā)育年齡。各扇區(qū)均從內(nèi)陸向?yàn)I海延伸,均包含濱海典型植被群落和植被演替過程,并且均在內(nèi)陸深處形成了農(nóng)田土壤。

黃河三角洲土壤含鹽量較高且以裸地為起點(diǎn)[17]。隨著演替進(jìn)行,植物種類增加,物種多樣性上升,植被蓋度、物種豐富度增加[3, 8],土壤質(zhì)量呈現(xiàn)逐漸改善的趨勢[7—8, 18]。植被演替和土壤演替在濱海生態(tài)系統(tǒng)內(nèi)交替發(fā)生,且均受到成土?xí)r間和離海距離等因素影響。成土?xí)r間的不同帶來生態(tài)系統(tǒng)年齡的差異,且與離海距離共同影響土壤脫鹽程度和熟化特征[8, 18]。鹽度影響海岸和河口的碳源和匯功能[19],進(jìn)而影響演替過程。白茅是輕度耐鹽植物,可作為可墾殖的植被指標(biāo)[8],耕地主要由白茅草甸和蘆葦-荻草甸開墾而來[3]。黃河三角洲土壤和植被隨時間的演替過程也伴隨著微生物的演替,微生物群落的演替受土壤化學(xué)性質(zhì)和植物種類影響[10, 12, 17],而微生物群落作為生物地球化學(xué)循環(huán)的關(guān)鍵因素,又塑造了環(huán)境特征和功能[20—22]。

黃河三角洲多個亞三角洲堆積體可提供對濱海土壤向農(nóng)田土壤演替的時空要素開展研究的平臺。本文探究成土?xí)r間是否為黃河三角洲植物和微生物群落演替的生態(tài)驅(qū)動力,為提高評估濱海生態(tài)系統(tǒng)植物和微生物演替規(guī)律的準(zhǔn)確性提供科學(xué)依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

本研究采樣區(qū)域位于山東省東營市黃河三角洲地區(qū)(37°31.807′—38°07.664′N,118°15.858′—119°13.384′ E),地跨東營市河口區(qū)、墾利區(qū)、利津縣。黃河三角洲屬暖溫帶季風(fēng)氣候,年平均氣溫12.3℃,年平均降水量542.3mm。土壤以濱海鹽土和濱海潮土為主, 土壤鹽漬化嚴(yán)重, 主要是氯化物鹽土和氯化物潮化鹽土, 土壤組成以泥沙為主, 養(yǎng)分含量低[23]。自然植被以鹽生植物為主。植被帶垂直于河床分布,沿河流由海向陸發(fā)展[9]。

1.2 樣品采集

2020年5月,在各時期沖積扇區(qū)范圍內(nèi)進(jìn)行土樣和植被信息采集(圖1),按成土歷史分別命名為路線1扇區(qū)、路線2扇區(qū)、路線3扇區(qū)和路線4扇區(qū)(下文分別簡稱R1、R2、R3和R4)。分別由河口所在海岸沿河道或古河道方向設(shè)計采樣路線,每條采樣路線由沿海到內(nèi)陸對樣地編號為樣地1、樣地2……樣地6(下文分別簡稱S1、S2……R6)。采樣標(biāo)準(zhǔn)為自成土以來受人類活動干擾較小的、有代表性的原生景觀植被,包括典型的鹽灘、鹽生草地、中生草地和灌叢草地,并采集農(nóng)田樣地?；趯?shí)地調(diào)查,選取樣地21處,如表1所示。R1扇區(qū)受人類近海工業(yè)和養(yǎng)殖活動影響,未發(fā)現(xiàn)原生鹽灘；R4扇區(qū)成土?xí)r間較晚,未發(fā)育形成灌叢草地。R1S4、R2S4、R3S5、R4S4為各沖積扇區(qū)發(fā)現(xiàn)的內(nèi)陸頂級景觀植被,農(nóng)田樣地在以上樣地5—10km外的內(nèi)陸選取。

圖1 研究區(qū)域和采樣地點(diǎn)示意圖

表1 樣地信息

每處樣地采用五點(diǎn)取樣法,樣點(diǎn)間距5—10m。使用帶刻度的土鉆在每個樣點(diǎn)采集20cm深度土樣,將5個樣點(diǎn)各深度的土樣分別于自封袋中混合均勻。在每個樣地的中心樣點(diǎn)劃出1m×1m樣方,在自然群落采集所有草本植物地上部莖葉用于生物量測定,并記錄樣方內(nèi)物種信息；記錄該處自然群落的各植物類群的蓋度信息。

各樣地的5個樣點(diǎn)的土壤樣品混合,分成兩部分。其中一部分風(fēng)干并過20目標(biāo)準(zhǔn)篩,再從其中取約20g過100目標(biāo)準(zhǔn)篩,用于土壤理化性質(zhì)測定。另一部分立即過20目篩,填充5mL離心管并標(biāo)號,埋藏于干冰中冷凍,每個樣地設(shè)三個重復(fù),送至上海美吉生物醫(yī)藥科技有限公司測序。

1.3 土壤理化性質(zhì)測定

取過20目篩的干土,使用去離子水按水土比5∶1制備浸提液,以PHS—3E酸度計測定pH,以Mettler Toledo電導(dǎo)儀(FiveEasy Plus FE20)測定電導(dǎo)率(EC)。土壤含水量(SWC)用鋁盒烘干法測定；土壤有機(jī)質(zhì)(SOM)依據(jù)NY/T1121.6—2006標(biāo)準(zhǔn)滴定；全氮(TN)使用濃硫酸消解—凱氏定氮法測定；全磷(TP)和全鉀(TK)依據(jù)DZ/T0219.1—2006標(biāo)準(zhǔn)使用ARL 9900X射線熒光光譜儀測定；堿解氮(AN)依據(jù)LY/T 1228—2015標(biāo)準(zhǔn)滴定；速效磷(AP)依據(jù)LY/T1232—2015標(biāo)準(zhǔn)使用TU—1810紫外可見光分光光度計滴定；速效鉀(AP)依據(jù)NY/T890—2004標(biāo)準(zhǔn)使用ICE3500原子吸收分光光度計測定。以上10個指標(biāo)連同碳氮比(CN.ratio)作為本研究的土壤理化性質(zhì)因子集合。

1.4 植被多樣性和時空因子計算方法

植被樣品取回后,于烘箱中105℃處理24h,稱得生物量干重,記為plant.bioM。

植物物種豐富度為樣地中所能識別的物種種類數(shù),記為plant.S。統(tǒng)計植物物種數(shù)時,考慮到5月許多植物處于苗期或營養(yǎng)生長期而不易辨別,故將同一科的相同生活型的部分物種進(jìn)行了合并。樣地區(qū)域植被總蓋度記為plant.cov。

植物物種多樣性計算方法為Shannon-Wiener指數(shù)[24],計算公式為：

H=-∑PilnPi

式中,Pi指第i個物種的出現(xiàn)頻率。下文簡稱Shannon指數(shù)。

以樣方內(nèi)植被各物種的多度與植株個體總數(shù)的比值為Pi計算物種組成Shannon多樣性,記為Shannon.abun?？紤]到作為環(huán)境因子時,1m2的樣方對樣地植被信息代表性可能不足,故另以樣地植被各類群的蓋度占比為Pi計算類群蓋度Shannon多樣性,記為Shannon.cov。以上5個指標(biāo)作為本研究的植物多樣性因子集合。

各沖積扇區(qū)成土?xí)r間計算方式為亞三角洲體形成年份區(qū)間的中值與采樣時間(2020年)的年代差；河口海岸距離計算方法為樣地與黃河河口或古河口所在海岸方向的距離。

1.5 微生物群落組成測序方法

根據(jù) E.Z.N.A. soil DNA kit (Omega Bio-tek, Norcross, GA, U.S.)說明書進(jìn)行微生物群落總 DNA 抽提。使用1%的瓊脂糖凝膠電泳檢測DNA的提取質(zhì)量,使用NanoDrop2000測定DNA 濃度和純度；使用338F (5′-ACTCCTACGGGAGGCAGCAG-3′) 和806R (5′-GGACTACHVGGGTWTCTAAT-3′) 對 16S rRNA基因V3-V4可變區(qū)進(jìn)行PCR擴(kuò)增,使用SSU0817F (5′-TTAGCATGGAATAATRRAATAGGA-3′)和SSU1196R (5′-TCTGGACCTGGTGAGTTTCC-3′)對18S rRNA基因V4—V5可變區(qū)進(jìn)行PCR擴(kuò)增。隨后基于Illumina MiSeq進(jìn)行16S rRNA、18S rRNA擴(kuò)增子測序。使用Fastp(https://github.com/OpenGene/fastp,version 0.20.0)軟件對原始測序序列進(jìn)行質(zhì)控,使用FLASH(http://www.cbcb.umd.edu/software/flash,version 1.2.7)軟件進(jìn)行拼接：(1)過濾reads尾部質(zhì)量值20以下的堿基,設(shè)置50bp的窗口,如果窗口內(nèi)的平均質(zhì)量值低于20,從窗口開始截去后端堿基,過濾質(zhì)控后50bp以下的reads,去除含N堿基的reads；(2)根據(jù)PE reads之間的重疊關(guān)系,將成對reads拼接成一條序列,最小重疊長度為10bp；(3)拼接序列的overlap區(qū)允許的最大錯配比率為0.2,篩選不符合序列；(4)根據(jù)序列首尾兩端的條形碼和引物區(qū)分樣品,并調(diào)整序列方向,條形碼允許的錯配數(shù)為0,最大引物錯配數(shù)為2。DNA 抽提、熒光定量PCR、MiSeq 測序、序列質(zhì)控、拼接均由商業(yè)測序公司完成(美吉生物,上海)。

1.6 數(shù)據(jù)分析

以97%的相似度[25]對序列進(jìn)行操作分類單元(OTU)聚類并剔除嵌合體,以RDP classifier比對Silva 16S rRNA、18S rRNA數(shù)據(jù)庫(v138.1)來對每條序列進(jìn)行物種分類注釋,設(shè)置比對閾值為70%。獲得OTU數(shù)據(jù)表后,分別篩選原核生物域和真菌界,分別篩得24380、1493個OTU。保留至少3個樣本中reads數(shù)≥5的OTU,并按最小樣本數(shù)抽平。所得結(jié)果中,每個重復(fù)組細(xì)菌的OTU頻數(shù)合計為23698,測序深度指數(shù)≥96.0%；真菌的OTU頻數(shù)合計為29965,測序深度指數(shù)≥99.8%；測序深度指數(shù)表明處理得到的數(shù)據(jù)集足以代表該區(qū)域真實(shí)的微生物群落。對細(xì)菌群落和真菌群落,分別將其Sobs指數(shù)、Shannon多樣性指數(shù)以及豐度占比前8的門類的各樣點(diǎn)reads數(shù)等10個因子組合作為微生物群落組成因子(bacteria, fungi)。

細(xì)菌物種的功能屬性的文獻(xiàn)查閱以原核生物標(biāo)準(zhǔn)命名列表(LPSN)數(shù)據(jù)庫的信息為主,若一個細(xì)菌屬下所有公開發(fā)表(validly published)的物種均有耐鹽3%及以上或生長所需鹽濃度的下限>0.5%,則該屬被識別為耐鹽菌；若一個細(xì)菌屬所有公開發(fā)表的物種均具備反硝化能力,則該屬被識別為反硝化菌；反硫化菌和自養(yǎng)菌等鑒定方法同上。

各扇區(qū)頂級景觀植被群落的差異可能體現(xiàn)成土?xí)r間對演替階段的限制。本文對各扇區(qū)頂級景觀植被群落的微生物進(jìn)行PICRUSt2(Phylogenetic Investigation of Communities by Reconstruction of Unobserved States 2)分析,再將取得的基于KEGG(Kyoto Encyclopedia of Genes and Genomes)數(shù)據(jù)庫的功能通路信息結(jié)果進(jìn)行線性判別分析(LDA),識別對各樣地功能通路組成的差異效果影響較大的通路,判別標(biāo)準(zhǔn)為LDA得分(log10轉(zhuǎn)化)大于2,且克魯斯卡爾-沃利斯秩和檢驗(yàn)(Kruskal-Wallis rank sum test)為統(tǒng)計學(xué)顯著(LDA score>2,P<0.05)。

對成土?xí)r間的驅(qū)動力進(jìn)行偏曼特爾檢驗(yàn)(Mantel)時,將本文中除成土?xí)r間外的其他環(huán)境因子(distance、vegetation、soil properties、bacteria、fungi)構(gòu)建控制環(huán)境矩陣(算法為歐氏距離),檢驗(yàn)植被、細(xì)菌、真菌群落的物種矩陣(算法為布雷柯蒂斯距離的可變性能夠被成土?xí)r間解釋的部分。環(huán)境變量和生物變量相關(guān)分析(Bioenv)選取本研究除bacteria、fungi之外的18個環(huán)境因子。VPA分析時,由于微生物OTU矩陣0值較多,在計算Bray-Curtis距離時進(jìn)行了海林格距離轉(zhuǎn)化,繪圖方法如Delgado-Baquerizo等的研究[6]。

多樣性統(tǒng)計、Bioenv分析、Mantel分析、方差分解分析使用R包vegan完成。相關(guān)性分析使用R軟件包ggcor完成。LDA分析和LEfSe分析于Galaxy/Hutlab網(wǎng)頁完成。物種Venn圖和PICRUSt2在美吉生信云網(wǎng)頁繪制。其他統(tǒng)計圖使用ArcGIS 10.2、Adobe Illustrator 2021和R軟件包ggplot2、reshape2繪制。

2 結(jié)果與分析

2.1 土壤理化性質(zhì)和微生物OTU組成異質(zhì)性檢驗(yàn)

在百年尺度上,黃河三角洲流域不同時期扇區(qū)的土壤理化性質(zhì)相近(圖2)。河口海岸距離與有機(jī)質(zhì)(r=0.75)、全氮(r=0.65)、電導(dǎo)率(r=-0.82)、pH(r=-0.42)極顯著相關(guān)(P<0.001),且以上4個土壤理化指標(biāo)相互間均有線性耦合(P<0.001)。成土?xí)r間僅與速效鉀顯著相關(guān)(r=-0.27,P<0.05),但非影響速效鉀的首要環(huán)境因素,速效鉀主要受到土壤有機(jī)質(zhì)(r=0.40,P<0.01)和全氮(r=0.48,P<0.001)的影響。

圖2 成土?xí)r間、海陸距離與土壤理化指標(biāo)的相關(guān)性特征

各扇區(qū)微生物群落的OTU組成相近(圖3)。細(xì)菌群落共篩選出6374個OTU,其中4381個(68.7%)為4個扇區(qū)共有OTU,272個(4.3%)為單個扇區(qū)特有OTU(圖3)。真菌群落共篩選出391個OTU,其中260個(66.5%)為4個扇區(qū)共有OTU,15個(3.8%)為單個扇區(qū)特有OTU。1512個(23.7%)細(xì)菌OTU的豐度隨成土?xí)r間的變化顯著(P<0.05),67個(17.1%)真菌OTU的豐度隨成土?xí)r間的變化顯著(P<0.05)。 4個沖積扇區(qū)自然稟賦相似,土壤微生物OTU組成相近,生態(tài)系統(tǒng)的空間演替可充分體現(xiàn)成土?xí)r間對生物演替的驅(qū)動作用。

圖3 各沖積扇區(qū)細(xì)菌和真菌OTU數(shù)量分布特征

2.2 植物分布特征分析

伴隨成土?xí)r間的增加,植物物種數(shù)和科類豐富度整體呈增加趨勢(表2)。研究區(qū)域內(nèi),代表性的鹽生植物有蘆葦(phragmitescommunisTrin.)、檉柳 (TamarixchinensisLour.) 及翅堿蓬 (SuaedaheteropterKitag)等,優(yōu)勢植物類群為禾本科(Gramineae)、藜科(Chenopodiaceae)、菊科(Compositae)、檉柳科(Tamaricaceae)等。中生草地群落的指示種為白茅(Imperatacylindrica(L.) Beauv.)、羅布麻(ApocynumvenetumL.)、野大豆(GlycinesojaSieb. et Zucc.)等。成土年代最近的R4扇區(qū)尚未發(fā)育出灌叢草地的植被類型,而最早的R1扇區(qū)出現(xiàn)了自然生長的白蠟、洋槐等樹木。

R4S4、R3S5、R2S4、R1S4為各沖積扇區(qū)發(fā)現(xiàn)的內(nèi)陸頂級景觀植被(表1),建群種除檉柳外,從耐鹽的狗牙根、蘆葦向輕度耐鹽[8]的茵陳蒿和白茅過渡(表2)。R1S3、R2S4、R3S4和R4S3為各扇區(qū)中科類豐富度較高的樣地,且一致分布在距海20km附近。隨成土?xí)r間增加,距海20km附近依次呈現(xiàn)從鹽生草地、中生草地向灌叢草地或高禾草地景觀的過渡特征。

表2 各樣地植被物種多樣性

成土?xí)r間驅(qū)動植物生物量(r=0.53,P<0.001)和植物物種豐富度(r=0.29,P<0.05)增加,且植被蓋度、Shannon多樣性均隨成土?xí)r間增加而增大(圖4)。內(nèi)陸區(qū)域受人類活動干擾的影響,農(nóng)田植物物種多樣性較低,表現(xiàn)為植物Shannon多樣性與土壤有機(jī)質(zhì)和全氮呈負(fù)相關(guān)(P<0.01)。

2.3 細(xì)菌和真菌群落組成分布特征

研究區(qū)域內(nèi),各扇區(qū)微生物物種門類組成相似,但各門類相對優(yōu)勢度存在差異(表3)。細(xì)菌群落3個主要優(yōu)勢門類為變形菌門(Proteobacteria)、放線菌門(Actinobacteria)、綠彎菌門(Chloroflexi),總占比為59.7%。真菌群落3個主要優(yōu)勢門類為子囊菌門(Ascomycota)、擔(dān)子菌門(Basidiomycota)、毛霉門(Mucoromycota),總占比為94.9%。

表3 研究區(qū)域內(nèi)各沖積扇區(qū)細(xì)菌和真菌門類相對豐度/%

細(xì)菌群落的Shannon多樣性受到土壤pH限制(r=-0.26,P<0.05,圖5)。有機(jī)質(zhì)含量顯著促進(jìn)真菌群落物種數(shù)(Sobs指數(shù))增加(r=-0.62,P<0.001),電導(dǎo)率和pH則是限制真菌物種數(shù)的限制因子(r=-0.69,P<0.001,圖5)。植被蓋度和植物Shannon多樣性的提高促進(jìn)真菌物種多樣性增加(r>0.28,P<0.05)。

圖4 成土?xí)r間、海陸距離、土壤理化指標(biāo)與植物多樣性的相關(guān)性特征

圖5 細(xì)菌和真菌多樣性、門類相對豐度與環(huán)境因子的相關(guān)性熱圖

河口海岸距離、土壤有機(jī)質(zhì)、全氮、pH和電導(dǎo)率是影響微生物門類組成的主要因子,與多數(shù)主要門類的相對豐度有較強(qiáng)的相關(guān)性(圖5)。成土?xí)r間對微生物群落門類相對豐度影響則相對獨(dú)立。

隨著時間演替,部分門類的占比變化呈現(xiàn)相反趨勢。成土?xí)r間驅(qū)動硝化螺菌門(Nitrospirota)、酸桿菌門(Acidobacteriota)、Methylomirabilota相對豐度顯著增加(r>0.29,P<0.05,圖5),同時擬桿菌門(Bacteroidota)、芽單胞菌門(Gemmatimonadetes)、Patescibacteria相對豐度顯著降低(r<-0.27,P<0.05)。

成土?xí)r間還驅(qū)動Sumerlaeota、裝甲菌門(Armatimonadota)、異常球菌門(Deinococcota)、Acetothermia、Hydrogenedentes、鹽厭氧菌門(Halanaerobiaeota)、海微菌門(Marinimicrobia-SAR406-clade)、Calditrichota、Campilobacterota等次要門類相對豐度的降低(r<-0.26,P<0.05,圖4未列出)。

豐度與成土?xí)r間相關(guān)性顯著的真菌門類為捕蟲霉門(Zoopagomycota)(r=0.30,P<0.05)。研究區(qū)域內(nèi),子囊菌門(Ascomycota)和擔(dān)子菌門(Basidiomycota)在真菌群落的豐度占比達(dá)91.9%(表3)。對這2個門類的相對豐度影響最大的環(huán)境因子為pH(r=0.60,P<0.001,圖5)。

2.4 細(xì)菌和真菌群落差異物種分析

在從門到屬的分類層級,部分微生物群落物種類群富集特征顯著(LDA score>2,P<0.05,圖6)。如,細(xì)菌厚壁菌門(Firmicutes)和真菌毛霉門(Mucoromycota)和捕蟲霉門(Zoopagomycota)在較早形成的R1扇區(qū)中可作為真菌群落的生物標(biāo)志物(biomarker),而擔(dān)子菌門(Basidiomycota)在R4扇區(qū)豐度最高(LDA score>2,P<0.05)。

圖6 各沖積扇區(qū)細(xì)菌和真菌LEfSe生物標(biāo)志物

富集于R3或R4扇區(qū)的細(xì)菌屬,包括Erythrobacter[26]、Fusibacter[27]、Gracilimonas[28]、Halanaerobium[29]、Limibaculum[30]、Marinobacter[31]、Palleronia-Pseudomaribius[32]、Pontibacter[33]、Salinimicrobium[34]、Zeaxanthinibacter[35]等,經(jīng)文獻(xiàn)分析發(fā)現(xiàn)除Erythrobacter[26]外均為耐鹽類群。富集于R4扇區(qū)的Campilobacterota為硫氧化化能自養(yǎng)細(xì)菌[36—38]。富集于R3或R4扇區(qū)的大型真菌為盤菌綱(Pezizomycetes)和煙管菌屬(Bjerkandera),而富集于R1的大型真菌為隸屬于蘑菇科(Agaricaceae)的側(cè)耳屬(Pleurotus)。以上物種均對各沖積扇區(qū)的微生物群落結(jié)構(gòu)組成的差異具有重要影響。

為控制海陸距離變量,對距海20km附近的R1S3、R2S4、R3S4和R4S3等4個樣地的細(xì)菌群落的屬進(jìn)行多級物種差異判別分析(Linear discriminant analysis Effect Size,LEfSe)。去除“norank”或“unclassified”開頭的沒有明確分類信息或分類名稱的類群,并去除平均豐度為0(Mean=0)的屬,最終得到554個細(xì)菌屬。215個(38.8%)屬對距海20km附近微生物群落結(jié)構(gòu)組成的差異具有重要影響(LDA score>2,P<0.05),其中73個屬的顯著富集區(qū)域?yàn)镽4S3樣地,其他68、21、53個屬分別為R3S4、R2S4和R1S3樣地的生物標(biāo)志物(biomarker)(表4)。在這215個屬中,至少100個屬為嗜鹽細(xì)群,其中有86個屬(86.0%)的顯著富集區(qū)域?yàn)槌赏凛^晚的R3和R4扇區(qū)；涉及反硝化、反硫化、光異養(yǎng)、嗜熱等功能性狀的細(xì)菌類群同樣集中于R3和R4扇區(qū)；R4扇區(qū)顯著富集的自養(yǎng)微生物多達(dá)13種,其中10種可營化能自養(yǎng)。

表4 距海20km 附近215個組間豐度差異顯著細(xì)菌屬的部分功能屬個數(shù)

2.5 各扇區(qū)頂級景觀群落細(xì)菌和真菌基于PICRUSt2通路豐度的LDA判別分析

針對細(xì)菌通路, PICRUSt2分析了各扇區(qū)頂級景觀群落功能通路的異質(zhì)性。各樣地細(xì)菌群落代謝通路相似, 但是代謝通路的豐度值不同,其中有12個level 2代謝通路存在組間顯著差異(圖7)。

R4S4樣地顯著富集的細(xì)菌群落功能通路較多涉及基礎(chǔ)生命活動。在該樣地的細(xì)菌群落,碳水化合物代謝(Carbohydrate metabolism)和氨基酸代謝(Amino acid metabolism)兩個level 2代謝通路的豐富度相對較高(LDA score>2,P<0.05,圖7)。存在組間顯著差異的level 3代謝通路中,所有從屬于遺傳信息處理(level 1: Genetic Information Processing)的9個通路均富集于R4S4扇區(qū),涉及氨酰tRNA生物合成(ko00970)、核糖體(ko03010)、RNA降解(ko03018)、DNA復(fù)制(ko03030)、蛋白質(zhì)運(yùn)出(ko03060)、核苷酸切除修復(fù)(ko03420)、錯配修復(fù)(ko03430)、同源重組(ko03440)和內(nèi)質(zhì)網(wǎng)蛋白質(zhì)加工(ko04141)等生命活動必需的復(fù)制和表達(dá)過程(LDA score>2,P<0.05,表5)。此外,涉及嘌呤代謝(ko00230)、嘧啶代謝(ko00240)、氨基酸生物合成(ko01230)等關(guān)鍵分子的通路的相對豐度均顯著高于其他扇區(qū)(LDA score>2,P<0.05)。光合作用(ko00195)和光合生物固碳(ko00710)通路在R4S4樣地細(xì)菌群落相對豐度較高(LDA score>2,P<0.05)。

R1S4樣地細(xì)菌群落富集的細(xì)菌群落功能通路則涉及多樣化的代謝過程。涉及細(xì)胞老化(Aging)、外源物質(zhì)的生物降解與代謝(Xenobiotics biodegradation and metabolism)、復(fù)制與修復(fù)(Replication and repair)等level 2通路豐度高于其他樣地(LDA score>2,p<0.05,圖7)。在營養(yǎng)代謝方面,R1扇區(qū)有關(guān)單糖的通路豐度較高,如糖酵解/糖異生(ko00010)、戊糖和葡萄糖醛酸轉(zhuǎn)化(ko00040)、果糖和甘露糖代謝(ko00051)、半乳糖代謝(ko00052)、丁酸(ko00650)等(LDA score>2,P<0.05,表5)。芳香族化合物降解(ko01220)的代謝功能也富集于R1S4樣地(LDA score>2,P<0.05)。

表5 細(xì)菌群落從屬基因信息處理和4類代謝通路的KEGG level 3差異通路

針對真菌通路, PICRUSt2分析得到了各扇區(qū)頂級景觀群落MetaCyc pathway信息及其豐度,共有25條MetaCyc代謝通路存在組間顯著差異(圖7)。R4S4樣地顯著富集的真菌代謝通路亦較多涉及基礎(chǔ)代謝,如TCA循環(huán)(PWY-5690)、丙二醇生物合成(PWY-7385)、有關(guān)細(xì)胞色素c的有氧呼吸I(PWY-3781、PWY-7279)和葡萄糖酵解(ANAGLYCOLYSIS-PWY)通路(LDA score>2,p<0.05,圖7)。R1S4樣地顯著富集的主要是有氧呼吸相關(guān)的泛醌類合成通路,如泛醌- 6(PWY-7235合成)、泛醌- 7(PWY-7235)和泛醌- 9合成(PWY-7235)(LDA score>2,P<0.05)。

圖7 細(xì)菌和真菌群落KEGG途徑豐度差異的線性判別分析(LDA)條形圖

2.6 成土?xí)r間驅(qū)動力的統(tǒng)計學(xué)檢驗(yàn)

成土?xí)r間對植物、細(xì)菌、真菌的群落組成變異的解釋能力分別為4.76%、2.98%、3.68%(圖8),均大于全球范圍土壤年齡對生態(tài)系統(tǒng)指標(biāo)的平均解釋能力2.1%[6],表明成土?xí)r間在黃河三角洲是重要的生態(tài)驅(qū)動力。成土?xí)r間對植物、細(xì)菌、真菌的群落組成變異的獨(dú)特方差解釋比例分別為3.41%、2.61%、1.12%,獨(dú)特解釋比例大于0,說明成土?xí)r間對生物群落組成變化的影響不能被其他狀態(tài)因子完全代替,因此成土?xí)r間是為黃河三角洲生物演替帶來額外解釋度的顯著驅(qū)動因子。

圖8 植被和微生物群落組成的方差分解分析

成土?xí)r間對細(xì)菌和真菌群落組成矩陣均顯著偏相關(guān)(Partial MantelP<0.05)(表6),表明成土?xí)r間單一因子對微生物演替過程中群落物種組成有極顯著驅(qū)動力。若不加入成土?xí)r間,細(xì)菌和真菌群落的Bioenv模型的相關(guān)系數(shù)相對最佳模型的擬合精度分別相對降低1.41%、5.78%(表7),表明在黃河三角洲微生物演替研究中將成土?xí)r間納入考慮可提高擬合精度。

3 討論

3.1 成土?xí)r間對黃河三角洲植物群落演替的影響

由于成土?xí)r間帶來演替時間的梯度,各扇區(qū)內(nèi)陸植物群落結(jié)構(gòu)呈由簡單向復(fù)雜的演變趨勢。研究區(qū)域內(nèi),植物分布和多樣性在隨成土?xí)r間的差異主要體現(xiàn)在生物量(r=0.53,P<0.001)和植物物種數(shù)(r=0.29,P<0.05)。距海約20km的自然植被樣地呈現(xiàn)了植被類型由低級向高級的過渡特征。生態(tài)系統(tǒng)的空間演替體現(xiàn)了成土?xí)r間對自然植被由低級向高級發(fā)展的驅(qū)動作用。劉志杰[39]的樣方調(diào)查同樣發(fā)現(xiàn),成土?xí)r間在1964年后的黃河三角洲葉瓣植物種類較少,植物種類不超過10科；灌草叢在成陸時期較長的區(qū)域成為主要植被類型之一。陳彭禎霓等[40]以冰川退縮區(qū)百年時間跨度的植被初生演替為研究對象,發(fā)現(xiàn)不同演替時間的樣地植物物種組成差異較大,木本植物生物量隨演替時間增加而增加；黃河三角洲植物群落結(jié)構(gòu)演替規(guī)律與之相符,豐富了百年尺度的植被初生演替規(guī)律的研究。

隨成土歷史的延長,白茅的分布向海岸推進(jìn),體現(xiàn)了成土?xí)r間對近海土壤生境條件改善的驅(qū)動作用。白茅的出現(xiàn)代表演替到達(dá)較高階段,種群結(jié)構(gòu)趨于穩(wěn)定[41]。本研究中,在白茅分布地點(diǎn),植被的共同點(diǎn)為植被蓋度在70%以上。在R1扇區(qū),白茅為常見種。在R4、R3和R2沖積扇區(qū),白茅分別在距入?？诩s20km、14km及9km進(jìn)入群落并成為伴生種。R4沖積扇區(qū)未演替出白茅占優(yōu)的草地群落,白茅僅在人類聚落附近偶見；R3S4樣地是在R3扇區(qū)發(fā)現(xiàn)的距海最近的白茅群落,白茅蓋度約25%；而在距海僅15.8km的R2S3樣地,白茅蓋度已達(dá)40%。前人研究證實(shí),白茅是黃河三角洲鹽土改良的指示植物,如侯本棟等[2]發(fā)現(xiàn),白茅群落的土壤有機(jī)質(zhì)、全氮含量顯著高于檉柳群落、鹽地堿蓬群落和裸地,且含鹽量相對較低。

3.2 成土?xí)r間對黃河三角洲微生物演替的影響

成土?xí)r間的增加可驅(qū)動細(xì)菌和真菌群落結(jié)構(gòu)(OTU的Bray-Curtis距離矩陣)的改變,這種改變具有一定的方向性(Partial Mantelr≥0.1377,P<0.05,表6)。除去土壤理化性質(zhì)、植物互作關(guān)系等的共線性干擾后,成土?xí)r間對微生物群落變化依然具有獨(dú)特影響；成土?xí)r間可為群落變化趨勢提供額外的解釋度(圖8,表7)。成土?xí)r間在黃河三角洲局地尺度是獨(dú)特而顯著的重要生態(tài)驅(qū)動力(圖8),這與Delgado-Baquerizo等[6]對土壤年齡驅(qū)動力的整合分析相符,即成土?xí)r間主要影響與生物活動相關(guān)的緩慢變化的生態(tài)指標(biāo)。

表6 成土?xí)r間對生物群落組成變化影響的偏Mantel檢驗(yàn)

表7 成土?xí)r間對微生物群落Bioenv最佳環(huán)境因子組合的擬合系數(shù)的影響檢驗(yàn)

微生物群落的功能特征也隨成土?xí)r間增加而發(fā)生顯著改變(圖7,表5)。在成土?xí)r期較近的區(qū)域,微生物群落涉及基因信息處理等基礎(chǔ)生命活動的功能通路相對豐度較高,與Yang等[42]關(guān)于黃河三角洲氨基酸代謝通路在高鹽土壤豐度更高的發(fā)現(xiàn)一致,推測在演替初期環(huán)境脅迫下,微生物群落依靠功能冗余提高抗逆性[43]及維持群落功能的穩(wěn)定性[44—45]；也可能是響應(yīng)于頻繁而幅度較大的環(huán)境變化,微生物需要能夠迅速利用周邊資源生存、生長和繁殖[21],故采取r-選擇競爭策略[46]。在生態(tài)系統(tǒng)年齡較長的R1扇區(qū),涉及細(xì)胞老化調(diào)控等功能通路的豐度顯著較高(LDA score>2,P<0.05)符合環(huán)境趨于穩(wěn)定帶來選擇過程優(yōu)勢的增加[47]的理論；單糖代謝、芳香族化合物降解等碳代謝相關(guān)基因通路的豐度發(fā)生顯著變化(表5),為前人關(guān)于微環(huán)境因素影響微生物碳源利用[42, 48]、土壤發(fā)育對生態(tài)位分化有促進(jìn)作用[21]等論點(diǎn)提供了支撐。

功能微生物的富集特征是微生物群落組成隨成土?xí)r間變化的原因之一。LEfSe分析表明,在各扇區(qū)距海20km附近的各樣地,對細(xì)菌群落組成變異貢獻(xiàn)顯著的215個可識別屬中,嗜鹽細(xì)菌、反硝化細(xì)菌、反硫化細(xì)菌等功能微生物的顯著富集區(qū)域均為成土較晚的R3和R4扇區(qū)(表4),其優(yōu)勢度隨演替時間有下降趨勢。對不同扇區(qū)細(xì)菌群落組成的整體差異起關(guān)鍵作用的細(xì)菌屬也以耐鹽菌為主,集中區(qū)域同樣為成土較晚的R3和R4扇區(qū)(圖6)。范運(yùn)梁[49]依據(jù)微生物培養(yǎng)實(shí)驗(yàn)發(fā)現(xiàn),伴隨鹽生植物的生長,耐鹽細(xì)菌的比例顯著下降。苗菁等[11]還發(fā)現(xiàn),隨著植被的覆蓋以及植物源有機(jī)物的輸入,鹽度對土壤細(xì)菌屬的作用逐漸減弱。耐鹽細(xì)菌類群的減少趨勢與鹽度隨演替降低的趨勢[8, 21]相適應(yīng)。隨鹽土向農(nóng)田土壤演替,根系生長可改善土壤結(jié)構(gòu)、降低土壤容重、增加土壤孔隙[50],通氣性的增加促進(jìn)好氧型生物活性的增加[51],而反硝化、反硫化均為厭氧過程[52—53],導(dǎo)致反硝化菌與反硫化菌競爭力隨成土?xí)r間增加的下降。

成土?xí)r間引起部分適應(yīng)高營養(yǎng)環(huán)境的類群相對豐度的增加,以及部分先鋒類群相對豐度的降低。酸桿菌門(Acidobacteriota)在環(huán)境土壤中一般占細(xì)菌總量20%左右[54],在土壤營養(yǎng)較高的黑土地或森林土壤,酸桿菌門占比更高[55—57]。酸桿菌相對豐度平均占比僅為8.26%(表3),苗永君[12]在黃河三角洲細(xì)菌多樣性的研究同樣發(fā)現(xiàn),酸桿菌在農(nóng)用地和自然植被土地類型的相對豐度占比均<10%。而且,酸桿菌隨土壤有機(jī)質(zhì)增加(r=-0.47,P<0.001)和電導(dǎo)率下降(r=-0.39,p<0.01)而增加(圖5),表明其在黃河三角洲屬于隨生境改良而增加的類群。擬桿菌門(Bacteroidota)、Sumerlaeota、Hydrogenedentes、海微菌門SAR406分支(MarinimicrobiaSAR406-clade)、異常球菌門(Deinococcota)和浮霉菌門(Gemmatimonadota)則隨土壤有機(jī)質(zhì)增加(r<-0.33,P<0.01)和電導(dǎo)率降低(r>0.33,P<0.01)而出現(xiàn)相對豐度的顯著下降,是優(yōu)勢度隨土壤改良而降低的先鋒類群。在成土早期相對貧瘠的極端生境,適應(yīng)寡營養(yǎng)生境的先鋒微生物容易取得較高優(yōu)勢度[10, 21],但是土壤營養(yǎng)限制的解除不利于這些種群的壯大[58],如李輝等[59]發(fā)現(xiàn),在赤泥自然成土的過程中,富營養(yǎng)細(xì)菌與貧營養(yǎng)細(xì)菌豐度比值隨演替明顯增加。

成土?xí)r間的增加驅(qū)動微生物生態(tài)位的分化。細(xì)菌群落功能通路的分化趨勢(圖6)體現(xiàn)了微生物生態(tài)位的分化。營養(yǎng)方式的多樣化同樣是微生物生態(tài)位分化的體現(xiàn),如寄生土壤線蟲、輪蟲和變形蟲的捕蟲霉門[60](Zoopagomycota)真菌(r=0.30,P<0.05),寄生昆蟲的黑僵菌[61](Metarhiziumanisopliae)真菌(r=0.32,P<0.05),寄生變形蟲的Neochlamydia[62]細(xì)菌(r=0.26,P<0.05),與共生固氮的Ensifer[63](r=0.31,P<0.05)和Azoarcus[64](r=0.33,P<0.01)細(xì)菌均為優(yōu)勢度隨成土?xí)r間增加而增大的類群,演替初期此類物種受制于宿主資源的貧乏。在對距海20km附近樣地細(xì)菌群落組成的對比中,還發(fā)現(xiàn)與反芻動物共生的瘤胃球菌屬[65](Ruminococcus)在R1S3樣地的顯著富集(LDA score=2.28,P<0.05),說明隨演替時間而增強(qiáng)[3]的人類活動能為微生物提供更多可選生態(tài)位。

3.3 成土?xí)r間驅(qū)動力的微生物因素

植被是濱海生態(tài)系統(tǒng)發(fā)展其生態(tài)系統(tǒng)功能最關(guān)鍵的過程[10]。鹽生植物蓋度提高可減弱土壤的自然蒸騰,抑制鹽分在土壤表層的持續(xù)累積；翅堿蓬等真鹽生植物能大量吸收鹽分,使土壤脫鹽[66]。鹽生植被改善土壤質(zhì)量的作用為微生物提供良好的生存條件,植物根系的分泌物及植物殘體是土壤微生物的主要可利用的碳源[58, 67]。

但是,成土初期鹽灘光板地的存在表明先鋒植物的定居需要等待土壤的形成和改良過程。在極端生境下,高等植物無法存活,但細(xì)菌群落物種數(shù)和多樣性較少受環(huán)境因子影響(圖5,表明某些適應(yīng)力較強(qiáng)的菌種可定居。在微生物定居于原生裸地的同時,即開始了成土進(jìn)程[68],成土初期的功能微生物促進(jìn)了生境的改良。

自養(yǎng)微生物促進(jìn)土壤碳氮積累。光合作用相關(guān)通路富集在成土相對較早的扇區(qū),這與Wang等[10]光合作用相關(guān)基因顯著富集于新生濕地的發(fā)現(xiàn)一致。藍(lán)藻細(xì)菌門能夠進(jìn)行產(chǎn)氧的光合作用, 部分物種能進(jìn)行生物固氮,還可以通過光合作用釋放特異的蛋白質(zhì)信號分子來調(diào)控植物的生長和發(fā)育[67]。但是,本研究中藍(lán)藻門(Cyanobacteria)相對豐度與成土?xí)r間無顯著相關(guān)性(P>0.05,圖4),暗示藍(lán)藻以外的光合生物貢獻(xiàn)于細(xì)菌群落光合作用表型的富集。本研究在16S rDNA克隆文庫中發(fā)現(xiàn)了CandidatusChloroploca[69]、Rhodobacter[70]等光合硫氧化細(xì)菌以及光異養(yǎng)的Rubribacterium[71]、Rubellimicrobium[71]等紫色非硫細(xì)菌,還發(fā)現(xiàn)了距海20km附近發(fā)現(xiàn)的光異養(yǎng)細(xì)菌在R3、R4扇區(qū)的富集(表4),此類光合生物很可能對成土早期的細(xì)菌光合作用起到了關(guān)鍵作用。在R4S3扇區(qū)中還發(fā)現(xiàn)了化能自養(yǎng)細(xì)菌的富集現(xiàn)象(表4),包括3種耐鹽的自養(yǎng)硝化菌和6個可硫氧化自養(yǎng)的細(xì)菌屬。前人強(qiáng)調(diào)成土早期以藍(lán)藻為主的光合自養(yǎng)微生物對土壤碳氮的積累起到關(guān)鍵作用[59, 72],本研究推測濱海鹽土化能自養(yǎng)菌和藍(lán)藻以外的光合生物在土壤碳氮積累方面的作用同樣具有生態(tài)學(xué)意義,這與近海河流沉積物硫元素的富集和活躍的元素化學(xué)循環(huán)[73]相適應(yīng),有待實(shí)驗(yàn)驗(yàn)證。除藍(lán)藻和自養(yǎng)硝化細(xì)菌外,耐鹽固氮菌也具有促進(jìn)土壤氮積累的功能潛力。16S rDNA克隆文庫中發(fā)現(xiàn)的耐鹽的共生固氮菌Ensifer[63]和Azoarcus[64]分別可以與豆科植物和禾本科植物形成根瘤；耐堿的中慢生根瘤菌屬[74—75](Mesorhizobium)也有為先鋒植物提供了早期土壤所缺乏的氮源的潛力。

李輝等[59]還發(fā)現(xiàn),微生物可通過分泌有機(jī)酸、菌絲改良土壤結(jié)構(gòu)等方式參與土壤的改良。微生物能夠釋放土壤酶, 催化氧化還原、有機(jī)質(zhì)礦化、腐殖物質(zhì)在土壤中的合成、以及生長活性物的釋放等多種生物化學(xué)反應(yīng)[22], 改善植物的生長環(huán)境。鹽生植被與微生物循環(huán)反復(fù)作用[18], 最終實(shí)現(xiàn)濱海鹽土向農(nóng)田的轉(zhuǎn)化。

4 結(jié)論

成土?xí)r間對細(xì)菌和真菌群落變化具有顯著的偏相關(guān)性(MantelP<0.05),對植物、細(xì)菌和真菌群落組成變化的解釋能力高于全球尺度土壤年齡對生態(tài)系統(tǒng)指標(biāo)的平均解釋能力。成土?xí)r間驅(qū)動黃河三角洲植物、細(xì)菌和真菌群落組成變化,也造成黃河三角洲微生物群落的功能特征的顯著異質(zhì)性。在百年尺度上,成土?xí)r間是黃河三角洲植物、細(xì)菌和真菌群落演替的重要而顯著的驅(qū)動力,將成土?xí)r間納入研究變量能夠?yàn)樯锶郝溲萏娴难芯刻峁╊~外的解釋度。該結(jié)論為提高評估濱海生態(tài)系統(tǒng)植物和微生物演替規(guī)律的準(zhǔn)確性提供了科學(xué)依據(jù)。