土壤有机质是土壤养分的重要组成部分，有机质不但能够增加土壤的抗风蚀能力^[1]，而且能够改善土壤的物理和化学性质^[2]，能够促进土壤良好发育，并提高土壤肥力^[3]，同时土壤有机质也能促进植物生长^[4]，是各种植物生长的营养来源^[5]。土壤有机质空间异质性受到很多因素的影响^[6]，这些因素是生态系统碳氮循环的重要组成部分^[7]，对全球碳循环也有重要意义^[8]。由于土壤特性具有一定的空间分布特点^[9]，因此某一种空间结构特征在一定采样尺度下才能表现出来^[10]，并且由于成土环境、条件不同，大多数土壤理化性质都存在较大差异^[11-12]。现有研究表明，土壤有机质的空间变异性与土壤类型、取样尺度、取样方法^[13]有关。土壤有机质的空间异质性研究方法最初主要是传统统计学方法，该方法只能描述土壤全局性的变化特征^[14]，不能精准描述土壤空间变异的特性^[15]。近几年，地统计学方法广泛应用于土壤学^[16]、气象学^[17]和生态学^[18]等领域，目前最常用的是GIS与地统计学相结合的方法^[14]，但是已有的研究侧重于用单一的地统计学方法研究土壤有机质表层的空间分布特征，对研究区不同土壤剖面的研究相对较少^[19]。此外，各种方法的研究结果迥异，单一的地统计学方法并不能准确描述土壤理化性质的变化规律。因此，本研究拟运用4种地统计学方法揭示玛依湖区不同土壤剖面有机质的空间分布规律，其中，趋势分析法能够根据采样点土壤有机质含量的高低，直观反映出湖区土壤有机质在东西、南北方向的空间变化规律；反距离权重插值法根据已有采样点数据，插值出湖区土壤有机质整体的分布状况；空间自相关法可以反映出湖区土壤有机质在空间上的聚集程度和空间差异性；半变异函数分析法能够定量描述土壤有机质的空间变异特征。

玛依湖是新疆克拉玛依市重要的生态屏障，湖区的生态环境关乎克拉玛依市主城区的生态安全，近年来，由于上游补水量变化较大，对玛依湖区湿地植被群落产生极大影响。此外，湖区周围存在大量的油气资源，油气资源的合理开发与生态屏障保护之间的协调发展变得十分重要。本文以新疆克拉玛依市玛依湖区为研究区域，以湖区周围的典型土壤为研究对象，运用趋势法、空间插值法、空间自相关法和半变异函数分析法4种地统计学方法分析上、中、下湖区和不同土层深度土壤有机质的分布状况，探讨影响土壤有机质空间异质性的因素，以便为湖区生态恢复、防护林选种、育种以及合理的资源开发利用提供科学依据和理论支撑。

1.   材料与方法

1.1   研究区概况

玛依湖区(85°7′~85°26′E，45°14′~45°35′N)，位于克拉玛依市东南部、玛纳斯河下游古河道区，东部濒临古尔班通古特沙漠^[20]。该区域位于中纬度内陆地区，属于典型的温带大陆性气候。多年平均气温为8.6 ℃，1月为最冷月，多年月平均气温−15.4 ℃，极端低温−40.5 ℃，7月为最热月，多年月平均气温27.9 ℃，极端高温46.2 ℃；多年平均降水量108.9 mm，年潜在蒸发量高达1 492.0 mm，是同期降水量的13.7倍；全年日照2 704.0 h，7月最多，达302.5 h，12月最少，仅99.8 h，无霜期232.3 d，年平均风速为3.2 m·s⁻¹。湖区内梭梭Haloxylon ammodendron、柽柳Tamarix chinensis、西伯利亚白刺Nitraria sibirica等荒漠植物生长茂盛，土壤类型以灰漠土、草甸土和沼泽土为主。

1.2   样品采集与处理

2017年11月，围绕整个玛依湖区进行随机取样，采样距离、采样范围适用于所使用的地统计学方法。利用GPS定位，在选定点挖取近百个土壤剖面，并最终筛选出60个具有代表性的土壤剖面进行采样，并按照土壤层次的取样原则，在0~20、20~40、40~60和60~80 cm土层对每个采样点进行采样，共计240个土壤样品，将土壤样品带回实验室，用重铬酸钾容量−外加热法测定土壤有机质。

1.3   地统计学方法

趋势分析法和反距离权重插值(Inverse distance weight，IDW)法可以反映研究区土壤理化性质的分布状况和变化规律，均在Arcgis 10.2相应模块中进行。

空间自相关法可以反映某一属性在空间上的聚集程度和空间差异性^[21-22]，其计算公式如下：

式中，I为莫兰(Moran)指数，n为空间数据的个数， $ {x}_{i} $ 和 $ {x}_{j} $ 分别为i区和j区的属性值， $\overline x $ 为所有要素的平均值，W_ij为空间权重矩阵。I>0表示空间正相关，I=0表示不相关，I<0表示空间负相关，I值越大，空间相关性越明显^[22]。

半变异函数是地统计学中研究土壤变异性的关键函数，它能够定量描述土壤属性的空间变异规律^[23]，其计算公式如下：

式中，γ(h)是样点间距的半方差，该值随着h的增加而增加， $ Z({x}_{i}) $ 、 $ Z({x}_{i}+h) $ 表示在 $ {x}_{i} $ 和 $ {x}_{i}+h $ 时的土壤属性观测值， $ N\left(h\right) $ 是间距为h的计算对数。半变异函数分析的结果主要由块金值（C₀）、基台值(C₀+C)和块金效应[C₀/(C₀+C)]体现，一般C₀/(C₀+C)越大，表明相应物理量空间分布越分散，相关性越弱。

1.4   数据分析

在SPSS 21中对数据进行描述性统计分析、正态分布检验和多重比较分析；在Arcgis 10.2中对整个流域4个土层的土壤有机质进行趋势分析和IDW值计算；在Geoda软件中对土壤有机质进行空间自相关分析；在GS+9.0中构建土壤有机质半变异函数。

2.   结果与分析

2.1   不同土层土壤有机质含量描述性统计

通过SPSS 21中Kolmogorov-Smirnow对玛依湖区由浅到深4个土层的有机质含量进行非参数检验，P值分别为0.168、0.341、0.837和0.646，均大于0.05，因此可以认为，整个流域不同土层的土壤有机质含量均符合正态分布，可进行描述性统计分析，但是其只能描述土壤有机质全局性特征，不能精准描述土壤有机质空间变化规律，故需要地统计学方法进行相应分析。由表1可以看出，玛依湖区不同土层土壤有机质平均含量以0~20 cm土层最多，其平均质量分数为8.72 g·kg⁻¹，其次为20~40 cm土层，平均质量分数为6.36 g·kg⁻¹，40~60 cm土层的土壤有机质含量最少，平均质量分数为5.12 g·kg⁻¹。当土壤有机质含量变异系数(CV)<10%时，为弱变异；当10%≤CV<100%时，为中等变异；当CV≥100%时为强变异^[24]，而玛依湖湖区各个土层的有机质含量变异系数分别为72.39%、77.84%、62.11%和63.93%，均呈中等程度变异。

表  1  不同土层土壤有机质含量统计特征值

Table  1.  Statistical characteristic values of soil organic matter content at different soil layers

土层/cm Soil layer	w(有机质)/ (g·kg−1) Organic matter content					统计特征值 Statistical characteristic value
	最小值 Min.	最大值 Max.	平均值 Average	标准差 Standard deviation		方差 Variance	变异系数/% Coefficient of variation	偏度 Skewness	峰度 Kurtosis	P
0~20	0.41	35.15	8.72	6.31		39.87	72.39	1.50	3.88	0.168
20~40	0.28	25.68	6.36	4.95		24.51	77.84	1.54	3.71	0.341
40~60	0.19	13.20	5.12	3.18		10.10	62.11	0.58	−0.17	0.837
60~80	0.33	15.23	5.30	3.39		11.48	63.93	0.77	0.73	0.646

下载: 导出CSV 

| 显示表格

2.2   不同土层土壤有机质含量趋势法分析

通过Arcgis10.2中的趋势分析模块，可以得到不同土层土壤有机质含量的空间分布特征及其异向性分布特征参数，能为变异函数模型的获取和插值图的绘制提供基础，进而更准确地描述研究区土壤有机质含量的空间变异规律。

图1是所有采样点不同土层在南北、东西方向上的土壤有机质含量变化趋势。从图1中可以看出，0~20 cm土层中，土壤有机质含量在南北方向是逐渐增加的，而在东西方向是凸线，呈现先增加后减少的趋势；20~40 cm土层中，南北方向和东西方向都是凸线，南北方向的土壤有机质含量的增加速度逐渐缓慢，东西方向呈现先缓慢增加后快速减少的趋势；40~60 cm土层中，南北方向的土壤有机质含量呈现直线式快速增加，东西方向呈现先缓慢增加，后急速减少的趋势；60~80 cm土层中，南北方向的土壤有机质含量呈先快速增加后缓慢增加的趋势，东西方向呈快速减少的趋势，由趋势分析法可以看出，玛依湖湖区土壤有机质含量在0~20、20~40、40~60和60~80 cm土层变化速率虽然存在差异，但是所有土层土壤有机质含量均表现为：南北方向呈增加趋势，东西方向呈减少趋势。

图  1  不同土层土壤的有机质含量变化趋势

Figure  1.  The variation trend of soil organic matter content at different soil layers

下载: 全尺寸图片 幻灯片

2.3   不同土层土壤有机质含量反距离权重插值法分析

在Arcgis 10.2的空间插值模块中对土壤有机质含量进行反距离权重插值分析，从图2中可以看出，在0~20 cm土层中，土壤有机质含量反距离权重极小值出现在上湖，有机质含量在上湖较少，中湖和下湖均存在极大值，有机质含量整体水平较上湖有明显的增加；在20~40 cm土层中，有机质含量上湖较低，极大值出现在中湖和下湖，中湖的整体有机质含量高于上湖和下湖；在40~60 cm土层中，有机质含量极小值主要出现在上湖，极大值主要出现在中湖和下湖，中湖和下湖的有机质含量整体高于上湖；在60~80 cm土层中：有机质含量极小值主要出现在上湖，极大值主要出现在中湖和下湖，中湖和下湖的有机质含量整体高于上湖。此外，从图2中还可以看出，各土层中，土壤有机质含量在南北方向是逐渐增加的，而在东西方向呈现减少的趋势。此结果与趋势法分析结果基本一致。

图  2  不同土层土壤有机质含量的反距离权重插值(IDW)

Figure  2.  Inverse distance weight(IDW) interpolation method of soil organic matter content at different soil layers

下载: 全尺寸图片 幻灯片

2.4   不同土层土壤有机质含量的空间自相关法分析

在GeoDa软件中对玛依湖区各土层的土壤有机质含量进行空间自相关分析，结果如图3。由图3可以看出：整个湖区由浅到深各土层土壤有机质含量的Moran指数(I)均大于0，分别为0.1643、0.1236、0.1955和0.2461，均在空间上呈现正相关；随后对4个土层的Moran指数采用蒙特卡洛模拟的方法进行检验，在α=0.05水平下，4个土层对应的P值分别为0.007、0.016、0.003和0.001，表明4个土层土壤有机质含量的空间自相关是显著的，具有良好的空间结构。图3中I_e为I的期望值，通过I_e值计算出4个土层的Z值(Z值能反映空间正相关或者负相关)分别为3.1510、2.5934、3.5903和4.6355，说明了60~80 cm土层土壤的有机质含量空间正相关最显著，空间聚集程度最高，40~60 cm土层次之，20~40 cm土层的Z值最小。分析图3可以看出，在0~20、20~40和40~60 cm土层，采样点主要分布在第1、第2和第3象限，表明土壤有机质含量的高值区域分布较广，研究区土壤有机质含量的空间分布受高值分布影响较大；60~80 cm土层的采样点主要分布在第1、第2和第4象限，表明土壤有机质含量的高低值区域分布较为均匀，土壤有机质含量的空间自相关性受高值和低值共同影响^[22]。

图  3  不同土层土壤有机质含量的空间自相关分析

A1、B1、C1和D1显示各土层莫兰指数(I)；A2、B2、C2和D2显示各土层莫兰指数经置换后对应的期望值(I_e)

Figure  3.  Spatial autocorrelation analysis of soil organic matter content at different soil layers

A1, B1, C1和D1 show Moran index (I) values at different soil layers; A2, B2, C2和D2 show the expected values (I_e) of Moran index after iteration at different soil layers

下载: 全尺寸图片 幻灯片

2.5   不同土层土壤有机质含量的半变异函数法分析

利用GS+9.0软件对玛依湖区各个土层土壤的有机质含量构建半变异函数并进行相关分析，得到半变异函数模型(图4)及模型相关参数(表2)，可以看出，0~20 和20~40 cm土层土壤有机质含量的变异函数模型均为球状模型，40~60和60~80 cm土层土壤有机质含量的变异函数模型分别为高斯模型和指数模型，且40~60 cm土层土壤有机质含量高斯模型的决定系数最高，R²为0.661，拟合效果较好，0~20、20~40和60~80 cm土层土壤有机质含量模型的拟合效果较差。块金效应[C₀/(C₀+C)]反映了土壤有机质含量的空间依赖性，块金效应的大小反映了土壤物理变量空间相关性的强弱，一般C₀/(C₀+C)越大，表明相应物理量空间分布越分散，相关性越弱。从表2可以看出，C₀/(C₀+C)随着土层深度的增加迅速增加，0~20和20~40 cm的土层中，C₀/(C₀+C)分别为0.033和0.045，表明土壤有机质含量的空间相关性较强，40~60和60~80 cm土层中，C₀/(C₀+C)分别为0.427和0.420，表明土壤有机质含量具有一定的空间相关性，但会受到随机因素的影响；4个土层的块金值(C₀)分别为1.30、1.08、4.23和5.32，说明存在由采样和试验引起的误差，并且40~60和60~80 cm土层中误差较大；4个土层的基台值(C₀+C)分别为38.96、24.22、9.91和12.68，表明40~60 cm土层系统内随机性变异和结构性变异的总和较小，而0~20、20~40和60~80 cm土层中，随机性变异和结构性变异的总和较大；由浅到深4个土层的变程分别为1510.00、1250.00、7170.69 和39660.00 m，表明土壤有机质含量在该范围内具有空间相关性，60~80 cm土层的变程最大，20~40 cm土层的变程最小，表明研究区4个土层的空间变异尺度具有一定的差异。

表  2  各土层土壤有机质含量变异函数模型参数

Table  2.  Variogram model parameters of soil organic matter content at different soil layers

土层/cm Soil layer	变异函数模型类型 Variogram model	块金值(C0) Nugget	基台值(C0+C) Sill	变程/m Range	C0/(C0+C) Nugget/Sill	决定系数(R2) Coefficient of determination
0~20	球状模型 Spherical model	1.30	38.96	1 510.00	0.033	0.177
20~40	球状模型 Spherical model	1.08	24.22	1 250.00	0.045	0.061
40~60	高斯模型 Gaussian model	4.23	9.91	7 170.69	0.427	0.661
60~80	指数模型 Exponential model	5.32	12.68	39 660.00	0.420	0.413

下载: 导出CSV 

| 显示表格

图  4  不同土层土壤有机质含量的变异函数模拟模型

Figure  4.  Variogram model of soil organic matter content at different soil layers

下载: 全尺寸图片 幻灯片

2.6   土壤有机质含量的垂直分布特征

从图5可以看出，不同土层土壤有机质含量的垂直分布无明显规律，只有在极个别点才有明显的垂直分布特征，比如8、12、21、33、35、52和58号采样点，均是0~20 cm土层的土壤有机质含量较高，质量分数分别为21.23、18.82、16.83、19.21、17.14、14.77和17.69 g·kg⁻¹，20~40、40~60和60~80 cm土层的土壤有机质含量较低；46号采样点的0~20 cm土层的土壤有机质含量是一个极大值，质量分数高达35.15 g·kg⁻¹，而该采样点20~40、40~60和60~80 cm土层的土壤有机质含量较低，质量分数仅为6.44、3.21和3.71 g·kg⁻¹；28、29、30和31号采样点，所有土层的土壤有机质含量均较低。由此可见，不同土层土壤有机质不存在明显的垂直分布规律，但是不同采样点的不同土层土壤有机质含量存在明显差异。因此可以得出，玛依湖区不同土层土壤有机质含量的空间分布以水平差异性为主，在水平差异性存在的情况下，局部地区土壤有机质含量存在明显的垂直分布特征。

图  5  不同土层土壤有机质含量的垂直分布特征

Figure  5.  Vertical distribution characteristics of soil organic matter content at different soil layers

下载: 全尺寸图片 幻灯片

为了探究玛依湖区土壤有机质垂直分布特征，进一步对上湖、中湖和下湖不同土层的土壤有机质含量做了显著性分析(图6)，结果表明(在0.05的显著性水平下)：上湖各土层的土壤有机质含量无显著性差异，有机质含量垂直分布特征不明显；中湖0~20 cm与20~40 cm土层的有机质含量无显著差异，40~60 cm与60~80 cm土层土壤有机质含量无显著差异，0~20 cm与40~60和60~80 cm土层的土壤有机质含量有显著差异；下湖中0~20 cm与20~40、40~60和60~80 cm土层的土壤有机质含量有显著差异，20~40、40~60和60~80 cm土层的土壤有机质含量无显著差异。

图  6  不同湖区不同土层土壤有机质含量多重比较分析

各图中，柱子上方的不同小写字母表示差异显著(P<0.05，Duncan’s法)

Figure  6.  Multiple comparative analysis of soil organic matter content at different soil layers of different lake areas

In each figure, different lowercase letters on bars indicate significant differences(P<0.05, Duncan’s test)

下载: 全尺寸图片 幻灯片

3.   讨论与结论

玛依湖位于克拉玛依市的东南部，是克拉玛依市东部重要的生态屏障，近几年由于强烈的人为活动和自然因素的干扰，使得湖区生态环境发生了退化，进而改变了土壤属性，而目前关于玛依湖区土壤属性分布特征、影响因素的研究基本空白。周倩倩等^[25]研究发现有机质的分布受土地利用类型、土壤颗粒组成、土壤质地的影响。此外，张海鸥等^[26]研究发现，不同种植年限下土壤有机质与全氮含量呈显著正相关。张枝枝等^[3]研究表明，高程因子与土壤有机质表现为显著相关。袁晓宇等^[1]研究了土壤有机质含量对抗风蚀的作用。这些研究都为玛依湖区土壤属性的进一步研究提供了思路。

土壤有机质含量、空间分布受诸多因素的影响。鲍丽然等^[27]研究显示，渝西北的土壤有机质空间变异主要受地形坡度、土壤类型、土地利用类型和耕作措施的影响，玛依湖区地势平坦，大部分区域的坡度在0~2°之间，土壤有机质空间变异受坡度的影响较小，土壤类型主要有草甸土、灰漠土、风沙土、灰棕漠土、盐土和沼泽土，土壤有机质的空间分布受土壤类型的影响较大。杜佩颖等^[28]研究显示，平原丘陵过渡区土壤有机质的空间分布与黏粒含量密切相关，黏粒含量越高，土壤对有机质的吸附作用越强，有利于有机质的积累，砂粒含量越高，有机质矿化速率越快，不利于有机质的积累，玛依湖区不同土层土壤粒径分布差异较大，土壤主要由粉粒组成，黏粒和粉粒含量自上湖到下湖不断增加，砂粒含量不断减少，而有机质含量自上湖到下湖逐渐增加，有机质的分布受黏粒、砂粒含量影响，这与杜佩颖等^[28]的研究相一致。董莉丽等^[29]研究显示，陕西省土壤有机质空间分布与地形、植被类型、气候和土壤其他属性有关，并且发现坡向不同的地方，植被生长环境不同，植被覆盖也不同，导致表层土壤养分含量和分布状况的异质性，而玛依湖区由于地势平坦，坡向对有机质的影响较小，但是湖区外围不同的植被类型，如落叶阔叶林、荒漠、灌丛、草甸、沼泽等接受到的太阳辐射不同，地下生物种类、数量不同，地下生物量不同，对有机质的分布特征产生的影响也不同；此外，玛依湖湖水受上游补水量的影响较大，湖区面积在不断变化，湖区面积的扩张与收缩对有机质的空间分布也会产生影响。

不同的研究方法由于假设前提、参数、模型结构不同而有细微差别。趋势分析和反距离权重插值分析能够直观地反映出土壤有机质含量的空间变异规律，但是无法准确定量描述差异的大小，趋势分析法和反距离权重插值可以相互验证，在本研究中，这2种方法分析的玛依湖区不同土层土壤有机质含量的变化趋势基本保持一致，均表现为南北方向的土壤有机质含量呈增加趋势，东西方向呈减少趋势，趋势法能够反映出土壤有机质变化趋势的快慢(通过趋势线的斜率)，反距离权重插值能够更加直观地反映整个研究区土壤有机质含量的变化趋势，2种方法相辅相成；半变异函数受到采样点分布、研究区大小的影响，存在随机误差；空间自相关性分析中，空间权重矩阵的构建受到人为因素的影响较大，半变异函数和空间自相关法可以相互验证，在本研究中，空间自相关法分析出0~20、20~40、40~60和60~80 cm土层的Z值分别为3.1510、2.5934、3.5903和4.6355，半变异函数分析出的块金效应分别为0.033、0.045、0.427和0.420，表明研究区不同土层土壤有机质含量在表层(0~20和20~40 cm)空间相关性较弱，空间聚集程度较低，在底层(40~60和60~80 cm)空间相关性较强，空间聚集程度较高。此外，土壤属性垂直分布规律没有更好的方法去描述。本文采用不同的方法、多角度分析了克拉玛依市玛依湖流域不同土层土壤有机质的空间分布规律，多种方法相互印证，揭示了研究区不同土层土壤有机质的空间分布规律。

本文运用4种地统计学方法揭示了玛依湖区土壤有机质的空间异质性，趋势分析法和反距离权重插值法表明，玛依湖区不同土层土壤有机质水平分布差异较大，局部地区土壤有机质含量存在明显的垂直分布特征，土壤有机质含量南北方向呈增加趋势，东西方向呈减少趋势。空间自相关法和半变异函数分析法表明，表层(0~20和20~40 cm土层)土壤有机质空间相关性不显著，空间聚集程度较低，底层(40~60和60~80 cm土层)土壤有机质空间正相关较显著，空间聚集程度最高。玛依湖区土壤有机质含量的空间变化特征受土壤类型、土壤质地、外围植被类型以及湖区面积变化的影响较大，此外，土壤有机质是土壤养分的重要组成部分，是各种植物生长的营养来源，在湖区外围生态屏障建设时，防护林树种、种植深度、种植密度的选择应当结合土壤有机质的空间变化规律进行，如在南部和西部有机质含量较低的区域，种植对有机质含量要求较低的乡土树种，并适当减小种植密度，增加种植深度。