基于不同人为干扰的土壤全量氮磷钾空间变异性研究

刘靖朝，熊黑钢，何旦旦，乔娟峰，郑曼迪

(1．新疆大学资源与环境科学学院，教育部重点实验室，新疆 乌鲁木齐 830046;2．北京联合大学应用文理学院城市系，北京 100083)

土壤的空间异质性是受气候、母质、地形、生物以及时间等成土因素的影响，而随空间位置发生变化[1]。作为一种区域化变量，土壤同时具有地质结构特性和统计学的随机特性[2]。地统计学是研究土壤属性空间变异定量化的有效方法[3],并在土壤学科领域得到广泛的应用[4-5]。利用地理信息系统(GIS)技术可将区域内的样点属性数据与地理数据结合起来，更好地分析研究对象的分布格局与变异规律[6]。

土壤全氮(TN)、全磷(TP)和全钾(TK)是土壤肥力的重要组成部分，其空间变异在一定程度上可以反映养分的状态，可为土壤养分管理、生态环境保护和治理提供依据[7]。近年来，国内外学者对土壤养分的空间变异研究逐渐增多。探讨了不同土地利用方式下的土壤空间异质性及其影响因素[8-9]；研究不同尺度土壤属性的变异性，如田间、县域尺度、区域尺度和单种复合类型[10-12]；也从生态化学角度，分析了碳、氮、磷特征空间变异的影响因子以及人类干扰对土壤的影响[13-14]。

干扰是自然界中最常见的一种现象。人类活动的方式、强度以及作用时间都会对土壤性质产生重要影响，也是土壤理化性质发生改变的重要原因[15]。目前，对比不同干扰下土壤养分空间变异特征研究的文章较少。因此，本文以无人为干扰区25个样本，人为干扰区30个样本为例，以传统统计学、地统计学和GIS相结合的方法，定量研究在不同人为干扰下，0～10 cm土壤养分的空间分布特征，探讨人类活动对土壤养分空间分布的影响，为不同干扰下土壤养分有效利用与管理提供依据。

1 研究区概况

研究区地处天山北麓与准噶尔盆地南缘间(87°44′-88°46′E，43°29′-45°45′N)，位于三工河流域山冲洪积扇下部细土平原区。侏罗系中下统三工河组阶段，湖面面积扩大，碎屑物颗粒表现出以粗粒碎屑流为特征的河湖相沉积模式，主要以粗砂、砾岩、含砾粗砂等沉积为主[16]。流域下游土壤类型以灰漠土为主，成土母质为第三纪红泥页岩，荒漠区土壤以风沙土为主，成土母质为风积沙[17]。该地区冬季时间长，春秋季节不明显，夏季炎热，昼夜温差大，属中温带大陆性干旱气候。年均气温6.7℃，最高气温42.6℃，最低气温-41.6℃。降水稀少且空间分布不均匀，年均降雨量164 mm，年均蒸发潜力2 000 mm左右，冬季积雪3～29 cm。植被以梭梭、红柳、杂草以及人工梭梭林、榆树林等为主。pH值为7.76～8.98。

在研究区内布设A、B 2个采样区。其中，A区远离人类活动区，且与人类活动区间以上宽24 m，下宽6 m，长为15.30 km的水渠相隔，较少受到人类活动的影响，生态环境基本保持原貌，植被以梭梭林、红柳和杂草为主，将其视为无人为干扰区。B区离新疆建设兵团102团较近，受人类活动影响强烈，近两年区内土地大部分被开发成人工育林地，主要包括梭梭林地、榆树林地、育苗地等。其中梭梭林地表被开挖成上宽60 cm，下宽30 cm，深约50 cm的沟渠，沟渠间以3 m为间隔；榆树育林地则被挖出宽1 m，深0.5 m，间隔2 m的沟，沟间种植榆树；育苗地则被修整为“田”字格的平地。将该区视为人为干扰区。由于是普通林地，故不施用肥料。采样前一个多月未进行人工灌溉，且无自然降水。

2 研究方法

2.1 样品的采集

2014年7月中旬，在2个采样区分别布设覆盖整个区域的采样线，其中无人为干扰区内(A区)布设5条采样线，采样线间距为600～800 m；人为干扰区内(B区)布设采样线6条，采样线间距为800～1 000 m。在每条采样线上选择5个能代表该区域土壤背景的采样点，样点间距为300～500 m，并对每个样点进行GPS定位，共获得样点55个(图1)。采集每个样点0～10 cm的土壤样品，将其编号入袋，带回实验室，经自然风干、剔除杂质，过1 mm孔径筛后，送至中科院新疆生态与地理研究所由专业人员测定土壤养分。在每个采样点，选择一块具有对周边环境有代表性的10 m×10 m样方，统计样方内植被类型、株数及其覆盖幅度。

图1 采样点分布Fig.1 Sampling points distribution map

2.2 研究方法

地统计学是一种以区域化变量理论为基础，以半方差函数为主要工具，揭示属性变量在空间上的分布、变异和相关特征的空间分析方法[18]。半方差函数是用来描述区域化变量结构性和随机性并存的空间理论特征。作为半方差函数的重要参数，块金值、基台值和变程用来表达区域化变量在一定尺度上的空间变异和相关程度。变异函数理论模型的最优选择可综合考虑块金值和有效变程，用决定系数(R2)和残差平方和(RSS)决定。R2越大，模型拟合的曲线精度越高，RSS值越小，表明实测值与回归线越靠近，曲线拟合效果越好[9]。

Kriging插值方法，又称空间局部估计或空间局部插值法，是建立在变异函数理论及其结构分析基础上的，在有限区域内对区域化变量的取值进行无偏最优估计的一种方法。选用C-V交叉检验法[19]对插值结果进行验证，利用均值误差(ME)、均方根误差(RMSE)、平均标准差(ASE)和标准化均方根预测误差(SRMSE)对结果进行评价。均值误差越接近于0，说明预测值越是无偏的；平均标准差越接近于均方根误差越好；均方根误差越小，说明预测值与实测值偏差越小；标准化均方根预测误差越接近于1越好。

2.3 数据处理

特异值会使连续表面离散、半方差函数发生畸变[4]。采用阈值法处理特异值，并用平均值代替。之后对数据进行Kolmogorov-Smirinov 正态分布检验。采用SPSS 21.0软件计算数据描述性统计特征；GS+9.0拟合半方差函数；使用ArcGIS 10.2软件实现土壤养分Kriging空间插值分析。

3 结果与分析

3.1 土壤养分描述性统计分析

由全量氮磷钾描述性统计特征(表1)及全国第二次土壤普查养分分级标准(表2)分析可知，两区土壤表层全氮含量介于0.197～0.747 g·kg-1之间，属低、很低水平。无人为干扰区全磷含量较人为干扰区高，属于很高、极高水平(>0.8 g·kg-1)，人为干扰区为高、很高、极高水平(>0.6 g·kg-1)。全钾在无人为干扰区属中、高、很高水平(10～25 g·kg-1)，而在人为干扰区属中、高水平(10～20 g·kg-1)。偏度系数、峰度系数分别测度数据分布的不对称性和集中程度。从偏度系数来看，除两区全磷和人为干扰区全钾数据分布为负偏态外，其余数据均为正偏态分布。峰度系数表明，无人为干扰区全氮和人为干扰区全钾数据分布的集中程度高于正态分布，而其余指标的分布则低于正态分布。从变异系数来看，无人为干扰区全磷和人为干扰区全钾，变异系数分别为0.079和0.097，属弱变异性(CV<0.1)，其余指标均为中等变异性(0.1

相较于无人为干扰区，由于人类活动的影响，人为干扰区土壤全量氮磷钾含量降低(表1)，全氮和全磷数据间的差异性增大，全氮和全钾数据的分布状态也发生了变化。造成这种状况的原因无疑是，在人为干扰区，土地多被开发成育林地和育苗地，原有植被的枯枝落叶及其根系等被移除，有机物来源减少，土壤生土层上翻，透气性变好，养分转化速度加快，加之人工梭梭林、人工榆树林等消耗了一定量的土壤养分，且在人工管理过程中未施加任何肥料，使得土壤养分含量降低。在人为干扰区，由于人类对土地的开发程度、利用方式不同，造成样本数据间差异较大；而无人干扰区基本保持原貌，样本数据空间变异较小，土壤养分含量的不同取决于周围环境、植被类型及其覆盖度等。

表1 无人为干扰区(A)和人为干扰区(B)土壤全量氮磷钾含量描述性统计特征

注 Note：TN-全氮 total nitrogen, TP-全磷 total phosphorus,TK-全钾 total potassium, 下同 the same below.

3.2 土壤养分空间变异及结构特征

3.2.1 土壤养分趋势面特征 对土壤养分进行半方差分析及Kriging插值时，先分离趋势面。图2中每根竖棒代表一个样点的属性值，将其投影到EW-SN正交平面，由投影点做出最佳拟合线，模拟特定方向上存在的趋势，如果拟合线平直，说明无全局趋势。从东西方向来看，两区全钾均呈现平滑下降的趋势，且无人为干扰区下降更快；人为干扰区全磷和全钾表现出凹形特征，具有二阶趋势；无人为干扰区全磷呈逐渐上升的趋势。从南北方向来看，人为干扰区全氮和全磷表现出南北高，中间低的凹形特点，具有明显的二阶趋势；无人为干扰区全量氮磷钾和人为干扰区全钾均表现出逐渐下降的趋势，其中，全钾变化最明显。

3.2.2 半方差函数 半方差函数有3个重要的参数，即块金值、基台值和变程，用来表达区域化变量在一定尺度上的空间变异和相关程度。对两区土壤全量氮磷钾进行半方差函数拟合，并确定最优理论模型，结果见表3。

表2 全国第二次土壤普查养分分级标准[20]/(g·kg-1)

注：旋转角度：水平方向120°，垂直方向15°；绿色：东西方向，蓝色：南北方向Note: rotation angle: horizontal 120 °, vertical 15 °； green: east-west, blue: north-south图2 无人为干扰区(A)和人为干扰区(B)土壤全量氮磷钾分布趋势Fig.2 Trend of soil total N, P, and K in two types of study areas

区域Area项目Item样本Sample块金值Nugget Co基台值(Co+C)有效变程/mEffectiverange A0块基比/%Co/(Co+C)理论模型Theoreticalmodel决定系数Coefficient ofdeterminationR2残差和Sum ofresidualRSSATN255.30×10-41.21×10-21812.08 4.39高斯模型 Gaussian0.692.61×10-5TP259.98×10-46.09×10-3498.0016.39球状模型 Spherical0.906.44×10-8TK259.99×10-47.75×10-31615.0012.89高斯模型 Gaussian0.892.42×10-6BTN309.94×10-41.38×10-2758.007.20高斯模型 Gaussian0.663.42×10-5TP305.31×10-41.06×10-2829.005.01高斯模型 Gaussian0.847.01×10-6TK300.183.371543.005.34高斯模型 Gaussian0.790.52

无人为干扰区，土壤全氮和全钾符合高斯模型，而全磷符合球状模型；人为干扰区土壤全量氮磷钾均符合高斯模型。两区理论模型精度R2为0.66～0.90，且无人为干扰区全量氮磷钾的理论模型的R2均高于人为干扰区。RSS也以无人为干扰区更小(表3)。

块金值表示由测量误差和空间变异引起的随机变异，较大的块金值表明较小尺度上的某些过程不容忽视。两区全量氮磷钾的块金值均接近于0，表明由测量、实验误差或小于最小取样尺度引起的随机变异较小。

基台值的大小可以反映变量变化幅度或系统的总变异程度，其值越高，表明总的空间变异程度越高。虽然不同区域土壤养分变量的块金值和基台值无可比性，但可通过块基比来反映区域变量的空间结构特征。块金值与基台值之比Co/(Co+C)表示由随机因素引起的空间变异占系统总变异的比例，可用其反映变量的空间自相关程度。一般认为Co/(Co+C)<25%，变量具有强烈的空间自相关性；25%75%，变量空间自相关性很弱[21]。依据空间自相关性程度的划分及两区土壤全量氮磷钾的块基比值可知，两区土壤全量氮磷钾由随机变异引起的空间异质性占系统总变异的比例均小于25%，表明在研究尺度上均具有强烈的空间自相关格局，其空间变异主要是由内在的结构性因素引起，如气候、母质、地形、土壤类型等[22]。

有效变程用以表示变量的最大相关距离，其大小反映区域化变量影响范围的大小，反映空间变异程度或空间自相关尺度。无人为干扰区有效变程波动幅度较大，其值介于498.00～1 812.08 m，其中全氮的变程最大，为1 812.08 m，表明其空间连续性较好。人为干扰区有效变程波动幅度较无人为干扰区小，有效变程在758.00 m到1 543.00 m间变化，其中全钾变程最大，值为1 543.00 m。

3.3 土壤养分空间分布特征

3.3.1 Kriging插值精度验证 对土壤全量氮磷钾的Kriging插值后的交叉验证结果见表4。无人为干扰区，全量氮磷钾的均值误差较小，趋近于0；均方根误差除全钾相对较大，值为1.271外，其余指标均小于1;平均标准差亦趋近于均方根误差;标准化均方根预测误差趋近于1。人为干扰区，均值误差和均方根误差均小于1，且平均标准差趋近于均方根误差;标准化均方根预测误差值也在1左右。总体来看，两区全量氮磷钾插值精度较好，能用来进行空间分布预测。

3.3.2 土壤养分空间分布规律 运用ArcGIS 10.2中的地统计模块，结合半方差函数参数以及趋势分布图，采用Kriging插值法，得到两区土壤全量氮磷钾的空间分布插值图(图3)。

自然状态下，土壤中的氮主要来源于动植物残体、生物固氮，少量来自大气降水，磷主要来源于成土母质和动植物残体的归还[23],钾含量与土壤类型、母质风化程度和土壤质地等因素有关，基本上源于土壤母质中的含钾矿物[24]。无人为干扰区，全氮有3个高值区，且在南北方向上有先升高后降低的趋势，东西方向呈现减少趋势；全磷分布相对均匀，无明显的高值区和低值区，含量分布由中心向四周逐渐减少；全钾含量空间分布有2个高值区，且以其为中心，含量分布向南北递减。这与趋势图分析结果相吻合。由表5可知，全量氮磷钾的高值区在L2-L4采样线范围内，该区域植被覆盖度在25%以上，多为梭梭、琵琶柴和红柳等自然植被，凋落的植物残体和根系分泌物为该区域提供大量的腐殖质，地表腐殖质层较其他区域厚。随着腐殖质的增加，土壤中的氮素、磷素逐渐积聚，致使该区域氮磷含量高。土壤有机质(腐殖质)与粘粒组成的有机-无机胶体是交换性阳离子的基本载体，具有突出的保钾作用[25]，因此该区全钾含量也相对较高。

表4 无人为干扰区(A)和人为干扰区(B)插值精度交叉验证参数

图3 无人为干扰区(A)和人为干扰区(B)全量氮磷钾含量空间分布Fig.3 Spatial distribution of total N, P, and K contents in two types of study areas

张希彪等[26]分析了人为干扰对人工林土壤物理性质的影响，表明随着干扰程度的增加，砂粒增加，粉粒和粘粒减少，土壤向粗骨化方向发展；微团聚体数量上升，土壤结构变差；土壤容重增加，土壤变紧实。人为干扰除使土壤物理性质发生改变外，土壤类型也发生了变化，与无人为干扰的风沙土相比，人为干扰区土壤以母质为第三纪红泥页岩的灰漠土为主。采样线L7-L10内，大面积土地被平整成人工育林地，人类活动对该区域土地的影响高于人为干扰区其他样点。在土地平整过程中，自然植被被移除，有机质来源减少，加之由于人工翻动，土壤透气性变好，加速了有机氮的矿化过程。另外，有机质矿化提供部分无机磷外，其分解产生的有机酸和螯合剂将部分固态磷释放为可溶态[24]，这也是全磷含量降低的原因。在该区域全氮和全磷形成低值区。全钾在L7和L8采样线范围内形成2个高值区的原因可能是，该区域干扰程度高于其他样点，除成土母质影响外，人类对土层的翻动，破坏原有生物和物理结皮，使地表失去对土壤水分的保护能力，加之研究区位于西北干旱区，水分蒸发强烈，自然降水少，主要通过人工灌溉补充育林地植物水分需求，干湿交替可加剧岩石风化，使钾释放[27]，增加了该区域土壤钾含量。

4 结 论

1)由于人类活动的影响，人为干扰区全量氮磷钾含量降低，全氮和全磷数据间的差异性增大，全氮和全钾数据的分布状态发生变化。造成这种状况的原因可能是该区土地多被开发成育林地和育苗地，对土壤养分有一定的消耗，加之土壤翻动，透气性变好，养分分解速度加快，但无人工肥料补充，造成土壤养分的含量、数据间差异及分布状态发生变化。

2)两区除无人为干扰区全磷符合球状模型分布外，其余养分指标均以高斯模型拟合效果好。理论模型的R2介于0.66和0.90，拟合精度较好。且均以无人为干扰区全量氮磷钾的理论模型的R2较高，RSS较小。两区土壤养分指标的空间变异主要由结构性因素引起。

表5 无人为干扰区(A)和人为干扰区(B)采样线植被情况

3)无人为干扰区全量氮磷钾的高值区在植被覆盖度25%以上的区域内，较高的植被覆盖度为该区提供大量的腐殖质。伴随着腐殖质的增加，土壤中的氮素、磷素含量升高；全钾含量高则与有机-无机胶体的固钾作用有关。人为干扰区，尤其是受人为干扰程度较高的L7-L10采样线范围内，土地大面积平整为人工育林、育苗地，消耗一部分土壤养分，且在人工翻动过程中，土壤透气性变好，加速了土壤养分的分解，但由于没有人工肥料的补充，全氮和全磷含量在该区空间分布相对减少。全钾的2个高值区则是由于人类对土层的翻动，加快水分蒸发，频繁的干湿交替促进了含钾矿物的风化，使钾含量增加。因此，人为干扰区，尤其是人类活动影响较大的区域内，要适时补充一定量的氮肥和磷肥，以满足育林、育苗养分需求。