煤矿微生物复垦区灌木林下土壤养分的空间异质性

毕银丽，胡晶晶，刘 京

(1.西安科技大学 地质与环境学院,陕西 西安 710054; 2.中国矿业大学(北京) 地球科学与测绘工程学院,北京 100083)

中国是世界上最大的煤炭生产国与消费国，煤炭资源开发引发了一系列的生态环境问题，如造成沉陷区地表塌陷、土壤结构受损、土壤肥力降低，植被和生态退化，土壤和生态系统可持续发展受到制约[1-3]，而传统工程复垦投资大、周期长、生态效果不稳定，微生物复垦技术以其在土壤改良方面独特的优势，成为了当前国内外矿区土地复垦的新热点[4]。丛枝菌根真菌(arbuscular mycorrhiza fungi，简称AMF)是土壤中普遍存在的一类微生物，能与80%以上的陆生植物形成互惠共生关系，在促进植物养分吸收、修复根系功能、改良土壤结构、修复生态等方面具有重要作用[5]。毕银丽和申慧慧[[6]采用丛枝菌根等微生物复垦技术能够促进采煤沉陷地植被恢复、改良土壤、加快生态复垦速度，微生物复垦多年后，生态系统种群实现了自我演替变化。微生物复垦区域的土壤理化性状经多年微生物累积作用而发生质变，土壤养分空间异质性变化更值得在空间尺度上实现由点到面的关注。

土壤养分的空间异质性是土壤的重要属性之一。土壤是一个自然连续体，土壤在自然和人为因素的共同影响下，形成了不均一和动态变化的时空连续体，导致其在不同范围内都存在空间变异性[7-8]。碳、氮和磷营养是表征土壤肥力的重要元素，掌握其空间特征是实现土壤与生态可持续利用与管理的前提[9-13]。从20世纪90年代开始，地统计学与GIS相结合的方法成为研究土壤性质时空变异并揭示变异本质和机制的重要方法和手段，地统计学是研究空间变异最有力的手段，是在经典统计学基础上发展起来的空间分析方法[14]。目前，空间插值法被广泛用于土壤养分元素的空间分布和地图制图[15-16]，该方法能有效地实现空间数据由点向面的尺度扩展。协同克里金是地统计学中较为重要的插值方法之一，其可以用空间中的一个或多个辅助变量对所感兴趣的主要变量进行插值估算，从而获取研究范围内未知点的主变量属性值[17]。实际上，很多变量之间存在一定的相关性，普通克里金插值法(Ordinary Kriging)很难达到预期的精度，而协同克里金插值(Cokriging)为利用可靠的辅助变量对目标变量进行最佳线性无偏估值，不仅考虑了变量的空间连续性，还考虑了变量之间的相关关系,可有效地提高精度估值[18-20]。高文武等[21]发现基于与土壤锰元素含量相关的多辅助变量的协同克里金插值法可以提高锰元素含量的精度。杨奇勇等[22]发现以土壤含水率为辅助数据的协同克里金插值法能有效提高土壤干容重的预测精度。协同克里金插值适用于耕作层土壤全氮空间估值研究，不仅能提供局部变异信息也可优化采样数量[15]。方学燕等[23]对茎柔鱼的空间分布进行分析，发现协同克里金插值对未知区域的预测具有一定的可靠性。

采用协同克里金插值方法来研究微生物复垦多年后人工灌木林下土壤碳、氮和磷含量的养分空间异质性，对于揭示微生物复垦技术的长效土壤改良作用提供了一种新的统计方法。因此，笔者利用野外采集土壤样点数据，结合地统计学、GIS空间分析以及传统统计方法，揭示出采煤沉陷地微生物复垦多年后，其灌木林下表层土壤有机质、全氮和速效磷含量空间分布规律，并分析该养分含量与丛枝菌根真菌关键因子间的相关性，对于揭示采煤沉陷地微生物复垦区以丛枝菌根真菌为辅变量的土壤碳、氮和磷含量预测研究与空间分布具有重要的理论和现实意义，阐明微生物复垦持续土壤改良与生态效应，可为采煤沉陷地土地复垦和生态修复提供理论与技术指导。

1 材料和方法

1.1 研究区概况

试验区位于陕西省神木县大柳塔镇东山煤矿开采沉陷区国家水土保持示范园区内微生物复垦基地(109.22°E～111.12°E，38.83°N～39.78°N)，该区域处于黄土高原沟壑区向毛乌素沙漠的过渡地带，典型的干旱半干旱高原大陆性气候，盛行偏西和西北大风，全年干旱少雨多风沙，年均大风日14 d，年均降水量441 mm，主要集中在7至9月，年均蒸发量1 337 mm，年均气温8.9 ℃，年均日照2 875.9 h。当地先锋植被主要以1年生的狗尾草(Setariaviridis)、猪毛蒿(Artemisiascoparia)等为主。该基地微生物复垦前土壤类型为风积沙土，土质疏松，蓄水保肥、抗蚀力差。土壤基本理化性质:pH值8.35，电导率164.36 μS/cm，有效磷0.008 g/kg，速效钾0.066 g/kg，碱解氮0.032 g/kg，全氮0.41 g/kg，最大持水量15.12%，按照全国土壤养分含量分级标准，养分属于极贫瘠的沙土。

1.2 样地设置

设置13块分别为200 m×200 m试验地(图1)。分别是:纯对照区，未人工干扰的区域，该区以油蒿草本为主，是无人工灌木的生态自修复区;人工修复区，分别是接种丛枝菌根真菌和未接种丛枝菌根真菌的沙棘、文冠果和紫穗槐灌木林区，3种灌木林又分别为复垦5 a(2012年复垦)和3 a(2014年复垦)两个复垦时段，各样地基本情况见表1。2012年微生物复垦基地建设时，人工种植小灌木并接菌，供试丛枝菌根真菌的菌种为摩西管柄囊霉(Funneliformismosseae，Fm)，由中国矿业大学(北京)微生物复垦实验室通过盆栽沙土扩繁培养而得的菌根菌剂，该菌剂每10 mL基质含孢子为300个。

1.3 土样采集与分析

样区采样时间为2017年7月，每个样地随机选取15个10 m×10 m的样方点并做标记，按照“S”布点采样法和多点混合取样法，挖开植物根系周围0～20 cm土层，轻轻抖落收集根际土壤，同时快速剪下部分细根，用于测定根系侵染率;收集的土样均匀混合后装入自封袋带回实验室内风干，仔细除去其中可见的植物残体及土壤动物，过2 mm土壤筛后风干备用测试土壤理化性质。采后回填表土，防止水分蒸发对植物造成损伤。利用手持GPS读取采样点的经纬度、海拔，观察记录样方周围的生态环境特征。

土壤养分含量的测定:全氮采用凯氏定氮法，有机质采用重铬酸钾外加热法，速效磷采用钼锑抗比色法。菌根侵染率测定用Phillips和Hayman的KOH脱色-曲利苯蓝染色法，玻片镜检测定侵染根段数，其计算公式为:菌根侵染率=菌根段数/被检根段数×100%。土壤根外菌丝密度采用真空泵微孔滤膜抽滤-网格交叉法测定[23]。

1.4 数据分析

利用统计分析软件SPSS17.0对样本数据进行一般性统计描述、K-S检验和Pearson相关性分析，对并不符合正态分布的数据进行对数变换，利用地理统计软件GS+7.0对处理后数据进行半变异分析，根据残差越小越好，决定系数越大越好的原则，确定最优拟合模型。根据确定的最优拟合模型，利用ArcGIS 10.3.1做协同克里金插值土壤养分的空间分布图。从采样数据中取出一点，用其他点预测该点，将预测值与实测值进行比较，再将该点放回原始采样数据，然后再取出一个未被取出过的点，用剩余点进行预测，如此往复，直到所有点都被取过为止[24]。对研究区域林下土壤全氮、有机质和速效磷以及其他辅助变量进行地统计学分析，并比较不同的半变异函数模拟模型和预测方法。偏度表征土壤养分含量数据分布偏斜方向和程度的度量，峰度表征土壤养分含量数据分布形态的陡缓程度，变异系数表征土壤养分特性的空间变异程度。分析土壤有机质、全氮和速效磷含量的极小值、极大值、均值、标准差、偏度、峰度和变异系数。变异系数可定义为

(1)

2 结果与分析

2.1 土壤养分经典统计描述

表2为土壤全氮、有机质和速效磷描述性统计特征及数据正态分布检验，经K-S检验，全氮和有机质服从正态分布，速效磷需经对数变换服从正态分布。土壤全氮、有机质和速效磷含量的均值为0.50 g/kg,11.93 g/kg,2.24 mg/kg，均属于中等偏低水平，而且3者范围变化幅度较大。变异系数反映数据间的离散程度，也反映了样本的空间变异性大小，土壤全氮、有机质和速效磷的变异系数分别为52.18%,29.96%,54.33%，呈现出中等强度变异。偏度和峰度是用来描述数据分布特征，分别用以衡量数据分配的集中程度和不对称程度，通过对偏度和峰度的观测，可知全氮、有机质和速效磷均呈现右偏态，同时呈现一定的峰态。经对土壤全氮、有机质和速效磷分布检验，土壤全氮和有机质均呈现正态分布，而速效磷服从经对数变换的正态分布。这可能与灌木林植被类型、是否接种丛枝菌根真菌等因素密切相关。全氮、有机质和速效磷的含量变化范围与杨之江等[8]研究稻田土壤养分空间变异的结果相接近，表明土壤养分空间分布受人为干扰和结构性因素共同作用的影响。

表2 土壤养分描述性统计特征以及正态检验Table 2 Descriptive statistical characteristics of soil nutrients and normal test

2.2 土壤养分空间自相关分析

利用GS+7.0做土壤全氮、有机质和速效磷空间相关图(图2)。Moran’sI的数值大小反应了空间自相关性大小[25]。Moran’sI>0时，表示空间正相关，其值越大空间相关性越明显;Moran’sI<0时，表示空间负相关，其值越小空间差异越大;Moran’sI=0时，空间呈随机性。

全氮的Moran’sI系数在-0.219～0.438，I系数达到0.438，表明全氮对空间依赖程度较大;有机质的Moran’sI系数在-0.226～0.271，具有一定空间相关性，其空间自相关性小于全氮。速效磷的Moran’sI系数在-0.157～0.313，I系数在滞后距较小时变化较大，随着滞后距的增加，逐渐具有规律性。全氮和有机质的Moran’sI系数均随滞后距的增大而减小，全氮在211.29 m左右出现拐点，Moran’sI下降为负相关并趋于平缓;有机质空间自相关整体上随着间隔距离的增加先呈正相关后转变为负相关，具有一定的空间自相关性;速效磷的Moran’sI系数在150.00 m开始出现先增加后减少，可进行地统计插值分析。综上所述，土壤全氮、有机质和速效磷含量存在空间自相关性，可以进行地统计插值，其空间相关性随着间隔距离的增大呈现波动性的变化，其变异可能是由采煤塌陷和丛枝菌根真菌共同作用所引起。

图2 研究区土壤全氮、有机质和速效磷空间相关Fig.2 Spatial correlation map of soil total nitrogen and organic matter and available phosphorus

2.3 土壤养分空间变异特征

为了研究微生物复垦3 a多后土壤养分的变化，对土壤全氮、有机质、速效磷与丛枝菌根真菌的特性(菌根侵染率和菌丝密度)指标进行Pearson相关性分析(表3)。结果表明，全氮、有机质和速效磷与菌根侵染率、根外菌丝密度呈极显著正相关(P<0.01)，说明丛枝菌根真菌菌丝量越大，土壤中的全氮和有机质含量越多，磷的有效含量高。因此，丛枝菌根真菌特性菌丝密度可作为协同变量进行土壤全氮、有机质和速效磷空间异质性分析。选择丛枝菌根真菌典型特征菌丝密度为协变量进行土壤全氮、有机质和速效磷的空间变异分析(图3～5)，各变异函数模型参数见表4。单一因素的土壤全氮、有机质和速效磷的最优模型分别为球状模型、指数模型和指数模型，以土壤根外菌丝密度为协变量的土壤全氮、有机质和速效磷的最优模型分别为高斯模型、球状模型和高斯模型。

表3 土壤养分与丛枝菌根真菌相关性分析Table 3 Correlation analysis among the nutrient of soil and AM fungi

单一因素的土壤全氮、有机质和速效磷的块金值为0.030 6,9.01,0.118 9，均大于0，表明全氮、有机质和速效磷均是由采样误差、随机和结构变异引起的空间异质性;单一因素的全氮、有机质和速效磷的变程范围为211～1 088 m，表明土壤全氮、有机质和速效磷在较大范围内呈现空间自相关，而且有机质的空间连续性较好，全氮和速效磷对空间依赖程度较高。块金效应揭示了土壤属性的空间依赖程度，单一因素的全氮、有机质和速效磷的块金系数分别为0.433,0.287,0.502，在0.25～0.75，均为中等强度的空间相关性，养分空间变异受结构因素和随机因素共同作用，与苏松锦等[26]研究自然林下土壤养分空间变异规律相一致。

图3 土壤全氮的变异函数及协同变异函数Fig.3 Variation function and synergistic variogram of total nitrogen in soil

图4 土壤有机质的变异函数及协同变异函数Fig.4 Variation function and synergistic variogram of organic matter in soil

图5 土壤速效磷的变异函数及协同变异函数Fig.5 Variation function and synergistic variogram of available phosphorous in soil

表4 土壤养分变异和协同变异函数模型类型及参数Table 4 Types and parameters of soil nutrient variation and synergistic variogram models

以土壤根外菌丝密度为协变量的全氮、有机质和速效磷的块金系数分别达到0.038,0.002,0.097，均为强空间相关性，表明辅助变量的土壤全氮、有机质和速效磷的空间变异性受结构性因素作用为主，空间变异受随机因素的影响低于单一因素的全氮、有机质和速效磷。以土壤根外菌丝密度为协变量的土壤全氮、有机质和速效磷的R2分别为0.709,0.757,0.467，高于土壤单一因素全氮和有机质的决定系数，可看出以菌丝密度作为辅助变量的模型比单一变量的模型拟合精度要高。可能是接种丛枝菌根真菌能显著提高采煤沉陷地土壤全氮、有机质和速效磷含量[27]，使复垦区土壤全氮、有机质和速效磷空间分布影响程度变大所致。

2.4 土壤养分空间分布格局

利用土壤根外菌丝密度为辅助变量的最优半变异函数模型在Arcgis 10.3.1中对土壤全氮、有机质和速效磷进行协同克里金插值，得到了土壤全氮、有机质和速效磷的空间分布图(图6)。土壤全氮呈块状分布，总体分布趋势为南北低，这是因为南和北两个方向是未接种丛枝菌根真菌的对照区，中间部分为接菌区，接菌微生物修复3 a后土壤养分含量有提高趋势。土壤全氮的含量表现为:复垦5 a的接菌区>复垦5 a的未接菌区≥复垦3 a的接菌区>自修复区>复垦3 a的未接菌区，随着复垦的年限增加，土壤全氮含量在增加。土壤有机质总体上接菌区明显优于对应年限的对照区，呈不规则斑块状，接菌区的西部和东南部含量较高，接种丛枝菌根真菌3 a后土壤有机质含量提高，顺序为:复垦5 a的接菌区>复垦3 a的接菌区≥复垦5 a的未接菌区>复垦3 a的未接菌区 >自修复区，复垦时间增加，土壤有机质的含量也随着增加。土壤速效磷分布呈现斑块状，破碎化程度高，含量高的斑块在接菌区分布多，复垦5 a(2012年复垦)紫穗槐接菌样地速效磷最高，复垦5 a区域高于复垦3 a。土壤全氮、有机质和速效磷空间分布特征与复垦年限、植物种类、本底土质及沉陷区塌陷微地形等因素有关。接菌促进了植物生长，人工修复后土壤养分的积累比自修复区高，与生物多样性高，加速枯落物的分解与有机质累积。

2.5 交叉验证

选择均方根和平均标准误差作为检验方程精度的评价指标，土壤全氮、有机质和速效磷的检验方程及评价指标结果见表5，根据均方根和平均标准误差越小，精度越好的原则，以土壤根外菌丝密度作为辅助变量的协同克里金插值比单一变量的全氮、有机质和速效磷克里金插值精度稍高。可知以菌丝密度为辅助变量的土壤全氮、有机质和速效磷协同克里金插值的预测值与真实值较为接近，精度较好，其中以菌丝密度为辅助变量的全氮的均方根最小，为0.217;平均标准误差最小，为0.254，表明预测精度最高。可见利用土壤根外菌丝密度作为土壤全氮、有机质和速效磷空间变异分析的协同变量更合理可行。

表5 土壤养分的检验方程及评价指标Table 5 Test equations and evaluation indexes of soil nutrient

3 讨 论

本研究区位于大柳塔采煤沉陷地微生物复垦区内，属于黄土高原沟壑区向毛乌素沙漠的过渡地带，主要是以风积沙土为主。煤炭开采引发局部地表塌陷，加剧了土质的不均一性和微地形起伏，损伤了生态环境。本研究采用协同克里金插值总体评价了微生物复垦3～5 a后土壤养分的空间异质性变化格局，能更全面、直观地反应土壤全氮、有机质和速效磷的空间分布。复垦区土壤中虽然存在着土著菌根真菌，但是与人工接菌效果相比，植物的菌根侵染率及根际土壤菌丝密度仍然很低。人工接菌的长期生态效应仍然持续产生影响[23]，丛枝菌根真菌菌丝增加了灌木林植被根系与土壤的接触面积，提高了豆科植物如紫穗槐根际土壤氮浓度，并吸收利用土壤中的无机氮和简单的氨基酸，加速矿化土壤中的有机氮，增加氮有效性[28]。丛枝菌根真菌是对土壤磷比较敏感的一类土壤微生物，在活化土壤磷有效性方面具有积极作用[5]，磷的有效性除了根的作用，还与土壤的性状、植物的种类、植物的生长量有关[29]。采煤沉陷地接种丛枝菌根真菌3 a后土壤开始表现出一定的养分积累，因植物种类及生长量不同，其土壤养分的空间变异性大，斑块破碎度大，其累积效应还需要一定积累期才能表现出来。碳、氮、磷等营养元素是表征土壤肥力的重要元素，土壤全氮、有机质和速效磷与菌根侵染率和土壤根外菌丝密度呈极显著正相关，存在协同反馈效应，丛枝菌根真菌在土壤养分转化过程中起着至关重要的作用[30]。

空间异质性导致空间格局的存在，协同克里金插值空间分布图能更全面、直观地反应土壤全氮、有机质和速效磷的空间分布。研究区土壤全氮和有机质呈现对角线状分布，速效磷呈现斑块状分布，这是由于土壤中的全氮、有机质和速效磷受丛枝菌根真菌的影响。复垦区的土壤的全氮、有机质和速效磷含量均表现为强空间自相关性，其变异可能是由结构性因素引起为主，比如土壤母质、气候、地形等因素;同时其变异也可能是由随机因素引起，而这一区域的随机因素主要是采煤塌陷导致地表土壤受到严重的扰动，土壤有机质原有的空间自相关性遭到破坏。在采煤沉陷地接种丛枝菌根真菌对土壤养分全氮、有机质和速效磷含量具有促进积累作用，采煤沉陷地微生物复垦区灌木林土壤全氮、有机质和速效磷存在不同的空间异质性和空间关联性。通过土壤养分空间分布格局可以看出，土壤全氮、有机质和速效磷的含量与复垦年限、是否接菌、植物种类、本底土质及沉陷区塌陷微地形等因素有关[31]。接菌区的土壤全氮、有机质和速效磷的含量的积累比自修复区高，可能是接菌促进了植物的生长及植被生物量，加速了复垦区枯落物的分解与有机质累积所致[32]。

采煤沉陷地微生物复垦区灌木林下土壤全氮、有机质和速效磷存在不同的空间异质性和空间关联性，这为采煤沉陷地微生物复垦后续生态可持续发展提供一定的理论依据，对维持矿区生态系统稳定性和可持续性提供了一种新的研究方法。本文主要是利用协同克里金插值法探究土壤全氮、有机质和速效磷在水平方向上的空间异质性，若进一步探究土壤养分在垂直方向上的异质性，并与水平方向上的异质性相结合，能更全面地揭示出采煤沉陷地全氮、有机质和有效磷空间异质性影响规律，其作用机理仍有待进一步深入研究。

4 结 论

(1)研究区域土壤全氮、有机质和速效磷存在空间自相关性，有机质的空间连续性较好，全氮对空间依赖程度较高。有机质空间格局呈现中等强度的空间变异，受随机因素和结构性因素的影响。接菌促进了土壤养分含量提高，复垦5 a后土壤养分产生明显的空间格局变异。

(2)林下土壤全氮、有机质和速效磷以菌丝密度为辅助变量的协同变异函数的最佳拟合模型分别为高斯模型、球状模型和高斯模型，决定系数分别为0.709,0.757和0.467。以菌丝密度为辅助变量的土壤全氮、有机质和速效磷均呈现强空间自相关，有机质的变程较长，空间连续性较好;土壤全氮空间分布呈块状分布，对空间依赖性较强;土壤速效磷呈不规则斑块状，破碎化程度高。

(3)协同克里格插值法能够有效地获取土壤全氮、有机质和速效磷的空间分布，以菌丝密度为辅助变量进行协同克里格插值能够显著提高土壤全氮、有机质和速效磷的预测精度。丛枝菌根真菌在土壤养分转化过程中起着至关重要的作用，丛枝菌根真菌的菌丝密度与土壤全氮、有机质和速效磷呈极显著正相关性。

(4)采煤沉陷区微生物复垦基地植物群落的演变与灌木林下土壤全氮、有机质和速效磷等养分累积的空间异质性和空间关联性，仍需进行继续监测和深入探讨。