环境科学  2023, Vol. 44 Issue (12): 6909-6920   PDF    
引用本文
赵明松, 陈宣强, 徐少杰, 邱士其, 王世航. 基于MGWR的土壤pH值空间建模及其影响因素分析[J]. 环境科学, 2023, 44(12): 6909-6920.
ZHAO Ming-song, CHEN Xuan-qiang, XU Shao-jie, QIU Shi-qi, WANG Shi-hang. Spatial Prediction Modeling for Soil pH Based on Multiscale Geographical Weighted Regression(MGWR) and Its Influencing Factors[J]. Environmental Science, 2023, 44(12): 6909-6920.
基于MGWR的土壤pH值空间建模及其影响因素分析
赵明松1,2,3, 陈宣强1,2,3, 徐少杰1, 邱士其1, 王世航1,2,3     
1. 安徽理工大学空间信息与测绘工程学院, 淮南 232001;
2. 矿山采动灾害空天地协同监测与预警安徽省教育厅重点实验室, 淮南 232001;
3. 矿区环境与灾害协同监测煤炭行业工程研究中心, 淮南 232001
收稿日期: 2022-12-04; 修订日期: 2023-02-14
基金项目: 安徽省自然科学基金项目（2208085MD88）；国家自然科学基金项目（41501226）；安徽理工大学人才引进项目（ZY020）
作者简介: 赵明松(1983~), 男, 博士, 副教授, 主要研究方向为数字土壤制图和土壤空间变异, E-mail: zhaomingsonggis@163.com
摘要: 以安徽、河南、江苏和山东省为研究区，利用599个土壤样点数据，从地形、气候和生物等方面选取与土壤pH相关的9个环境因子，采用多尺度地理加权回归（MGWR）、混合地理加权回归（Mixed GWR）、地理加权回归（GWR）和多元线性回归（MLR）这4种模型对研究区土壤pH空间分布进行建模，并结合MGWR与分位数回归揭示环境因子对土壤pH作用的空间差异性.结果表明：①研究区土壤pH在不同空间距离上呈不同程度的显著全局和局部空间自相关性，聚集特征明显.② 4种模型中MGWR模型最优，MGWR、Mixed GWR、GWR和MLR的建模集Radj2为0.64、0.62、0.59和0.48.MGWR的残差独立分布性最强，其空间自相关性最弱，全局Moran's I仅为0.07.③ 3种GWR预测结果显示，研究区土壤pH值空间分布总体由北至南逐渐降低，河南北部最高，安徽南部最低.④ MGWR回归结果表明年均降雨量（MAP）、多尺度谷底平坦度（MRVBF）和海拔对土壤pH的影响较强，且存在较强的空间异质性.在江苏北部和山东大部分地区，MAP对土壤pH的影响较强；在江苏北部和山东西部，MRVBF对土壤pH的正向作用较强；在江苏北部和中部，海拔对土壤pH的负向作用最强.⑤ MAP对不同分位数水平上的土壤pH均呈显著负作用，作用强度随分位数水平增加呈减弱趋势；MRVBF对低分位数水平（θ为0.1~0.4）上的土壤pH呈显著负作用，对高分位数水平（θ为0.5~0.9）的土壤pH作用不显著.研究结果可为利用MGWR开展大区域土壤属性影响因素分析及预测制图提供参考.
关键词: 土壤pH      多尺度地理加权回归(MGWR)      混合地理加权回归      分位数回归      数字土壤制图     
Spatial Prediction Modeling for Soil pH Based on Multiscale Geographical Weighted Regression(MGWR) and Its Influencing Factors
ZHAO Ming-song1,2,3 , CHEN Xuan-qiang1,2,3 , XU Shao-jie1 , QIU Shi-qi1 , WANG Shi-hang1,2,3     
1. School of Geomatics, Anhui University of Science and Technology, Huainan 232001, China;
2. Key Laboratory of Aviation-aerospace-ground Cooperative Monitoring and Early Warning of Coal Mining-induced Disasters of Anhui Higher Education Institutes, Huainan 232001, China;
3. Coal Industry Engineering Research Center of Collaborative Monitoring of Mining Area's Environment and Disasters, Huainan 232001, China
Abstract: Anhui, Henan, Jiangsu, and Shandong provinces were selected as the study area. A total of 599 soil samples and nine environmental factors of soil pH were collected. The spatial distribution of soil pH was modeled based on multi-scale geographically weighted regression(MGWR), mixed geographically weighted regression(Mixed GWR), geographically weighted regression(GWR), and multiple linear regression(MLR) models. Then, the spatial difference in the effect of environmental factors on soil pH was revealed using MGWR and quantile regression models. The results showed that: ① soil pH showed significant global and local spatial autocorrelation at different spatial distances, and the clustering characteristics were obvious. ② The MGWR model was the best among the four models, and the Radj2 of MGWR, Mixed GWR, GWR, and MLR were 0.64, 0.62, 0.59, and 0.48, respectively. The residual of MGWR had the strongest independent distribution and the weakest spatial autocorrelation with a global Moran's I of 0.07. ③ Three types of GWR predictions showed that the spatial distribution of soil pH decreased gradually from north to south in the study area, with the highest in northern Henan and the lowest in southern Anhui. ④ MGWR modeling results showed that there was strong spatial heterogeneity of mean annual precipitation(MAP), multi-resolution valley bottom flatness(MRVBF), and elevation affecting soil pH. MAP had a stronger effect on soil pH in northern Jiangsu and most parts of Shandong. The positive effect of MRVBF on soil pH was stronger in northern Jiangsu and western Shandong. The negative effect of elevation on soil pH was stronger in northern and central Jiangsu. ⑤ The quantile regression analysis showed that the mean annual precipitation had a significant negative effect on soil pH at different quantile levels of soil pH, and influence intensity decreased with the increase in pH quantile level. MRVBF had a significant negative effect on soil pH at a low quantile level(θ=0.1 to 0.4) but had no significant effect on soil pH at a high quantile level(θ=0.5 to 0.9). These results can provide an important reference for mapping soil properties and analyzing its influence factors based on the MGWR model in large regions.
Key words: soil pH      multi-scale geographically weighted regression(MGWR)      mixed geographically weighted regression      quantile regression      digital soil mapping     

土壤pH是重要的土壤属性, 它与土壤的肥力水平、微生物与动物群的活动和腐殖质的形成息息相关, 并直接影响土壤养分存在的形态和有效性[1~3].同时土壤pH还能够对土壤中重金属元素的存在状态、有效性和迁移转化特性进行控制, 从而对区域生态环境质量产生重要影响[3~5].因此, 预测土壤pH的空间分布, 探究不同环境因子对土壤pH空间分布的影响, 对于合理规划利用土壤资源、提升土壤环境质量有重要意义[4, 6].

土壤形成发育过程中环境因子与土壤属性的关系存在空间差异性, 随着空间位置的变化, 环境因子对土壤属性的作用强度不同, 在大区域中这种关系更为明显, 传统的线性或全局模型难以处理它们间的复杂关系[7, 8].如有研究表明, 中国农田土壤pH和有机质在旱地和水旱轮作的土壤上呈显著负相关, 而在水田上则呈显著正相关[8].此外, 不同的环境因子对土壤属性的作用具有尺度依赖性[9~11].有研究表明, 北京顺义区土壤镉的空间分布存在多尺度空间变异, 在小尺度上主要受土地利用强度的影响, 大尺度上主要受土壤类型的影响[11].不同空间尺度下, 土壤类型、成土母质、地形、土地利用以及耕作施肥等因素均会对土壤pH值的空间变异产生不同程度的影响, 即使在同一地区, 不同时期的土壤pH值空间变异程度和主控因素也不同[12~14].

近年来地理加权回归(geographically weighted regression, GWR)模型在土壤属性的空间建模中应用较多, 且建模精度较高[15~18].GWR模型中自变量的回归系数随空间位置而变化, 能够直观地探测不同空间位置上环境因子对目标变量的影响程度[19].普通GWR用一个固定的带宽来概括所有自变量对因变量作用的空间尺度, 实际上它们作用的空间尺度并不一致.在此基础上, 学者们改进了普通GWR模型, 如混合地理加权回归(mixed geographically weighted regression, Mixed GWR)[20]和多尺度地理加权回归(multi-scale geographically weighted regression, MGWR)[21].Mixed GWR模型中自变量分为全局项和局部项, 即采用两种带宽拟合模型, 兼具了自变量对因变量作用的全局性和局部性.与普通GWR模型相比, Mixed GWR的建模精度有所提高[20, 22].MGWR模型考虑了不同环境因子的差异化作用尺度, 能识别不同空间尺度下各自变量影响因变量的地理过程[23].近年来MGWR模型逐渐应用于多个领域, 如土壤有机碳多尺度建模[24]、空气质量影响因素的空间异质性[25, 26]、房价影响因素的空间异质性[27]和工业发展效率驱动因素的尺度差异[28]等.MGWR可以更加精确地建模, 定量识别自变量对因变量作用的空间异质性和尺度差异性, 在空间数据关系异质性建模方面将发挥更大的作用.

安徽、河南、江苏和山东地处黄淮海平原和长江中下游平原, 该区域气候、成土母质和地形变化多样, 土壤类型丰富, 是我国重要的粮食生产基地.本研究以这4个省份的表层土壤pH值为研究对象, 利用MGWR、Mixed GWR、GWR和多元线性回归(multiple linear regression, MLR)进行空间预测建模, 并对比4种建模结果差异, 同时通过MGWR和分位数回归结果探究环境因子对土壤pH影响的空间异质性, 以期为利用MGWR模型在大区域开展土壤属性空间变异和预测制图提供参考.

1 材料与方法 1.1 研究区概况

研究区介于110°21′~122°42′E, 29°41′~38°24′N, 包括安徽、河南、江苏和山东, 总面积5.676×105 km2.研究区地处亚热带向温带的过渡带, 年均降雨量550~1 800 mm, 年均气温11~16.4℃.地形西高东低, 海拔0~2 237 m, 地势高的区域集中在河南西部、安徽的皖南、皖西山区和山东的中部, 其他区域地形平坦, 海拔在100 m以下(图 1).土地利用以耕地类型为主, 主要分布在地形平坦区域; 林地主要分布在安徽皖西和皖南丘陵山地区, 河南西部山地, 江苏宁镇扬丘陵区和沂沐丘陵地区, 山东中部地区[29~32].

图 1 研究区位置与样点分布示意 Fig. 1 Location of study area and soil samples

1.2 数据来源

(1) 土壤样点数据来源于《中国土系志·安徽卷》[29]、《中国土系志·河南卷》[30]、《中国土系志·江苏卷》[31]和《中国土系志·山东卷》[32], 采样时间为2010~2011年6~11月, 样点采集按照优势土壤类型(面积较大)、地形-母质-土地利用等景观要素组合, 共计599个土壤剖面.由于土壤剖面按照发生层采样, 表层土壤的厚度不相同, 本研究按照深度加权后取0~20 cm土层的土壤pH作为数字土壤制图对象.

(2) 环境因子数据土壤pH作为土壤制图的目标变量.选择地形、气候和生物因子作为环境因子用于土壤pH值的建模.环境因子均为栅格形式的空间数据, 其中所有地形因子、生物因子和气候因子的原始空间分辨率分别为90 m、250 m和1 000 m.利用ArcGIS软件提取所有采样点的环境因子值, 在IBM SPSS Statistics 26进行多重共线性检测后, 最终筛选5个地形因子, 2个气候因子和2个生物因子, 所选环境因子数据描述见表 1.

表 1 环境因子数据描述 Table 1 Description of environment factor data

地形因子主要利用SAGA v.2.1.4[33]软件从SRTM数字高程模型(digital elevation model, DEM)提取.其中地形湿度指数(TWI)反映每个像元土壤水分的干湿程度, 其值越大代表该地集水能力越强.多尺度谷底平坦指数(MRVBF)识别多分辨率下的平坦和低地形进而识别谷底, 多尺度脊顶平坦指数(MRRTF)识别高平坦区域来识别山脊.当多尺度谷底(脊顶)平坦指数的值为0时, 表示该地区为沉积(侵蚀)地区; 值等于或大于1时, 表示该地区是逐渐增大的平坦区域[34].

所有环境因子重采样为1 000 m空间分辨率, 土壤pH预测结果的空间分辨率为1 000 m.

1.3 研究方法 1.3.1 空间自相关

GWR建模要求因变量存在显著的全局和局部空间自相关性[35].采用全局莫兰指数和局部莫兰指数分析土壤pH值的全局和局部空间自相关性.

全局莫兰指数(global Moran's I)用于描述整个区域内地理实体间空间自相关程度, 计算公式如下：

(1)

式中, n为样本总数, wi, j为样本i和j之间的空间权重, S0为所有样本之间空间权重的和, xi为要素i的属性值, X为对应属性的平均值.全局Moran's I的值介于-1~1之间, 正值表明土壤pH整体上以空间聚类特征为主, 负值则表明以离散分布为主, 且绝对值越大表明空间相关程度越高; 越趋于0则表明越偏向于随机分布.

局部莫兰指数(local Moran's I)用于描述局部区域内地理实体间的空间自相关程度, 计算公式如式(2)、式(3).

(2)
(3)

式中, n、wi, j、xi和X含义同公式(1).局部Moran's I的值为正, Z(统计检验值)得分为正时, 为“高-高(H-H)”聚类, 即土壤pH的高值样点聚集分布; Z得分为负时, 为“低-低(L-L)”聚类, 即土壤pH的低值聚集分布.局部Moran's I的值为负, Z得分为正时, 为“高-低(H-L)”聚类, 即土壤pH的高值点被低值点包围; Z得分为负时, 为“低-高(L-H)”聚类, 即土壤pH的低值点被高值点包围.具体含义、计算和显著性检验过程参见文献[12, 36].

1.3.2 分位数回归

分位数回归采用加权残差绝对值之和的方法估计参数.当自变量对不同分位数上的因变量的分布产生不同的影响时, 例如出现左偏或右偏的情况时, 分位数回归能更加全面地刻画分布特征, 描述自变量对因变量概率分布的影响.公式如下[37]:

(4)

式中, S为不同分位点估计值; θ为分位点; Hθ和Hij分别为第i个样点下第j个自变量在不同分位点θ的预测值与实际值; dij为第i个样点第j个自变量的数值.

1.3.3 地理加权回归(GWR)

MLR假定自变量对因变量的作用在全局上是同质的, 其回归系数在空间上一致, 忽视了空间位置对作用程度的影响.GWR是典型的局部回归模型, 其公式如下:

(5)

式中, (ui, vi)为第i个样点的空间坐标; Yi和Xik为因变量Y和自变量集Xk在空间位置处的实测值; β0(ui, vi)为空间位置(ui, vi)处的常数项; βk(ui, vi)为连续函数βk(u, v)在i点处的值; εi为随机误差.选用自适应bi-square空间权函数来拟合回归点的回归系数; 选用改进的赤池信息量准则(corrected akaike information criterion, AICc)来确定最优带宽[38].

1.3.4 混合地理加权回归(Mixed GWR)

Mixed GWR将地理变异性较强的自变量作为局部变量, 其余作为全局变量, 其公式如下:

(6)

式中, zil为具有固定回归系数γl的自变量, 即全局变量, 与传统回归模型的系数一样在全局范围保持不变其余参数含义与公式(5)中相同.Mixed GWR认为部分自变量对因变量的作用程度在整体区域内保持恒定, 应将其作为全局变量, 而剩余自变量对因变量的作用程度随着空间位置而变化, 则作为局部变量.

1.3.5 多尺度地理加权回归(MGWR)

MGWR区分各自变量对因变量作用的空间尺度, 允许不同的地理过程在各自最佳的带宽上运行, 增加了整体模型的可靠性, MGWR的公式如下:

(7)

式中, (ui, vi)为第i个样点的空间坐标; yi和xij为因变量和自变量集在i处的值; β0(ui, vi)为空间位置(ui, vi)处的常数项; εi为随机误差项, βbwj(ui, vi)为使用特定带宽bwj拟合的空间位置i处第j个自变量的回归系数.

1.3.6 精度评价

模型精度评价选用以下指标: Lin氏一致性协调相关系数(Lin's concordance correlation coefficient, LCCC)、平均绝对误差(mean absolute error, MAE)、均方根误差(root mean square error, RMSE)以及调整后决定系数(Radj2).LCCC是比较实际值和预测值之间一致性的方法; MAE用来衡量实际值与预测值之间的偏差; Radj2反映因变量的全部变异被自变量解释的比例.

(8)
(9)
(10)
(11)

式中, n为样本总数, ρ为实际值和预测值的相关系数, Oi和Pi为样点i的实际值和预测值; O和P分别为实际值与预测值的均值; σo和σp分别为实际值和预测值的标准差.

空间模型的优劣也可以根据模型的残差趋势进行评估[19], 若模型残差空间自相关性下降可以间接证明模型的改善, 残差存在较强的空间相关性则说明模型还存在一些未捕捉到的可预测信息.除上述4个精度评价指标外, 本研究还通过对比土壤pH值建模残差的空间自相关性的强弱, 来比较模型的优劣.

1.4 数据处理与分析

① 首先根据土壤pH值进行随机抽样, 将599个样点按照4∶1比例分为建模集(479个)和验证集(120个); ②对所选的环境因子进行了容差和方差膨胀因子检验, 剔除出环境因子间的多重共线性, 在IBM SPSS Statistics 26软件中完成; ③分析不同空间距离上土壤pH全局和局部空间自相关性, 在ArcGIS 10.7软件中完成; ④利用筛选出的环境因子进行MGWR、Mixed GWR、GWR与MLR建模, 在MGWR 2.2和GWR 4.0软件中完成, 比较建模精度; ⑤结合分位数回归和MGWR回归系数分析环境因子对土壤pH影响的空间异质性, 分位数回归在IBM SPSS Statistics 26软件中完成.

2 结果与讨论 2.1 土壤pH值描述统计特征

研究区土壤pH的平均值为7.10, 变化范围在4.20~9.18, 变异系数为15.35%, 属于中等变异(表 2).训练集和验证集的平均值、中位数、标准差和变异系数都非常接近.研究区碱性土壤样本(pH>7.5)占42.57%, 主要分布在北部; 中性样点(6.5≤pH≤7.5)占28.55%; 酸性土壤样点(pH＜6.5)占28.88%, 主要分布在安徽和江苏南部(图 2).从各省土壤pH均值来看, 大小排序为: 河南>山东>江苏>安徽, 变异系数大小排序为: 安徽>江苏>河南>山东, 表明安徽的土壤pH值变异程度最大, 山东的土壤pH值变异程度最小.

表 2 土壤pH值统计特征 Table 2 Statistical characteristics of soil pH value

图 2 研究区样点分布 Fig. 2 Distribution of soil samples in the study area

2.2 土壤pH值空间自相关分析

对总样本的土壤pH进行全局空间自相关检验, 在默认空间距离阈值为56.4 km下存在显著的全局自相关性, 全局Moran's I为0.70(P＜0.01), 总体上土壤pH以聚集分布为主.在100、200、300、400、500和600 km的空间距离阈值下, 全局Moran's I值分别为0.65(P＜0.01)、0.58(P＜0.01)、0.52(P＜0.01)、0.47(P＜0.01)、0.43(P＜0.01)和0.39(P＜0.01).整体上在不同空间距离阈值上, 土壤pH均呈现显著的空间自相关性, 且随着空间距离阈值的增大, 全局Moran's I值逐渐下降, 即全局空间自相关性逐渐减弱.

在不同空间距离阈值上, 研究区土壤pH也呈现出相似的局部空间聚类特征.本文仅展示默认空间距离阈值下的局部空间自相关聚类图(图 3).研究区中有298个样点存在显著局部空间聚类, 其中“高-高”聚类点共156个(土壤pH值较高且空间差异较小的样点聚集分布)主要分布于河南北部、安徽北部; “低-低”聚类点共142个(土壤pH值较低且空间差异较小的样点聚集分布)主要分布在江苏和安徽南部.

图 3 土壤pH局部聚类 Fig. 3 Local spatial cluster of soil pH

2.3 基于多种GWR的土壤pH空间建模

综合多种精度指标比较, 建模精度由高到低依次为: MGWR>Mixed GWR>GWR>MLR(表 3).3种GWR模型的精度高于MLR模型, GWR区分了不同空间位置上环境因子对土壤pH的作用程度, 模型解释能力得到很大的改善, 这与大量的研究结果一致[15~17, 39].GWR模型较MLR模型, 其AICc值降低了74.31, 一般认为AICc值减小2以上时, 认为改进后模型明显更优[19]; 残差平方和(residual sum of squares, RSS)降低了88.42, Radj2增加了0.11, 说明GWR模型对预测变量进行局部回归的方式有效改进了MLR全局回归模型, 有效呈现了在不同空间位置上环境因子对土壤pH作用程度的差异性, 较大地提升了模型的预测精度.

表 3 各模型的建模结果 Table 3 Modeling results of each model

图 4为验证结果, 综合所有指标来看验证精度由高到低依次为: MGWR>Mixed GWR>GWR>MLR.MGWR在验证集的Radj2为0.57, 对pH变异的解释能力达到了57%, MAE和RMSE分别为0.53和0.66, LCCC值为0.76, 与Mixed GWR和普通GWR的精度较接近.从散点图的分布来看MGWR的验证结果在1 ∶1线附近分布得更均匀集中.因此, 总体上MGWR模型的精度最优, 这与多数研究结果一致[27, 28]. 3种GWR模型训练集和验证集的Radj2变化在0.05~0.07以内, 表明GWR及其改进模型均有较高的稳定性.

图 4 各模型验证结果散点图 Fig. 4 Scatter plots of validation of each model

对建模集的土壤pH残差进行全局和局部空间自相关检验, 结果表明所有模型的残差都存在一定的空间自相关性(图 5), 其空间自相关性由强到弱依次为: MLR>GWR>Mixed GWR>MGWR.MGWR的残差空间自相关性最弱, 残差独立分布性最强, 全局Moran's I仅为0.07(P＜0.05), 在空间上存在高-高聚类点19个, 低-低聚类点4个.其次是Mixed GWR模型, 其全局Moran's I为0.10(P＜0.01), 在空间上存在高-高聚类点20个, 低-低聚类点22个.GWR模型的全局Moran's I为0.17(P＜0.01), 在空间上存在高-高聚类点40个, 低-低聚类点30个.MLR模型的全局Moran's I为0.33(P＜0.01), 在空间上存在高-高聚类点85个, 低-低聚类点54个.MGWR不仅能反映空间对象局部效应, 还考虑了不同影响因子的差异化作用尺度, 有效降低了模型的残差自相关性, 模型更优.

图 5 各模型残差局部聚类 Fig. 5 Local clustering maps of residual of each model

2.4 基于MGWR的土壤pH预测制图

基于4种模型和环境因子栅格图预测了研究区土壤pH空间分布(图 6).研究区土壤总体呈南酸北碱的格局, 安徽南部地区土壤pH最低, 河南地区最高.MGWR和Mixed GWR模型的预测制图结果的空间格局相似, 结果较为平滑.普通GWR和MLR预测结果呈碎斑化, 即出现了较多低值被高值包围和高值被低值包围的斑块区域.土壤pH的空间分布图汇总统计结果显示(表 4), 4种模型预测制图的pH均值介于7.13~7.20, 均略高于采样点土壤pH均值7.10, 表明预测制图结果均存在不同程度的高估现象.土壤样点pH实际值的标准差为1.09, 预测制图结果的标准差介于0.67~0.77, 表明制图结果存在一定程度的平滑效应.MGWR预测制图结果的均值、标准差以及各pH值分级比例与土壤pH实际值最接近, 制图结果最优.

图 6 土壤pH制图结果 Fig. 6 Results of soil pH mapping

表 4 各模型制图结果pH统计特征 Table 4 Statistical characteristics of mapping results of each model

2.5 环境因子对土壤pH影响的空间异质性 2.5.1 多尺度地理加权回归系数分析

MGWR的标准化回归系数统计结果见表 5.全局带宽为样本数479, MRRTF、坡度、EVI和MAT接近全局尺度, 其回归系数值变化也较小, 表明这些环境因子趋向于在全局范围上影响土壤pH值.TWI、MRVBF、海拔和MAP的带宽较小, 表明其对土壤pH的影响存在较强的空间异质性, 回归系数在研究区上变化差异也较大.从影响因子的标准化回归系数的平均值和中值的绝对值来看, MAP、海拔、MRVBF和MRRTF对研究区土壤pH空间分布的影响程度较高.下文以MAP、海拔、MRVBF和MRRTF的标准化回归系数分布(图 7), 阐述环境因子对土壤pH影响的空间差异性.

表 5 MGWR标准化回归系数1) Table 5 Standardized regression coefficients of MGWR

图 7 MGWR的标准化回归系数分布 Fig. 7 Standardized regression coefficient distribution in MGWR model

MAP的回归系数值介于-0.679~-0.175之间[图 7(a)], 均值和中值为-0.435和-0.494(表 5), 说明MAP对土壤pH呈负向作用, 其它环境因子不变的情况下, 随着MAP的增大, 土壤pH逐渐下降.我国主要农田酸化的原因是氮沉降和盐基离子淋失[40], 而这两者都与降雨密切相关, 当降雨量远大于蒸发量时, 淋溶作用的发生让土壤酸化.在江苏北部和山东大部分地区, MAP的回归系数绝对值较高, 表明这些区域的MAP对土壤pH的负向作用较强; 低值集中在安徽南部和江苏南部地区, 在这些区域MAP对土壤pH的影响较弱.

海拔对土壤pH的影响强度存在强烈的地区差异, 回归系数介于-1.233~0.068之间, 回归系数均值和中值为-0.362和-0.174(表 5), 以负向作用为主[图 7(b)], 即其它环境因子不变的情况下, 随海拔增加土壤pH呈降低趋势.这与陈清霞等[12]和吴正祥等[14]研究的结果相似.研究区中的安徽南部和皖西地区属山地丘陵区, 海拔较高300~1 700 m.随着海拔的升高, 降雨量不断增大, 雨水冲刷和淋失作用不断增强, 土壤盐基不饱和度变大, 致使土壤不断向酸性方向发展[14], 该地区土壤pH预测值低于6.5, 主要呈弱酸性和强酸性(图 6), 这体现了MAP和海拔对土壤pH的负向作用.回归系数正值集中在河南西部和山东北部, 其余地区海拔对土壤pH呈负作用, 其中在江苏北部海拔对土壤pH的负向作用较强.

MRVBF的回归系数值介于-0.073~0.325之间[图 7(c)], 均值和中值为0.143和0.186(表 5), 表明了MRVBF对土壤pH的正向作用.MRVBF越小表示该位置越接近山谷地貌, 山谷地区蓄水能力较强, 淋溶作用更加剧烈, 土壤pH相对较小.回归系数正值的高值集中在山东西部和江苏北部, 在这些区域MRVBF对土壤pH的正向作用较强, 低值集中在山东东部、江苏中部和河南东北部地区, 这些区域MRVBF对土壤pH的正向作用较弱.

MRRTF的回归系数值介于-0.171~-0.089[图 7(d)], 均值和中值为-0.140和-0.151(表 5), 为负向作用.MRRTF越小即该位置越接近山脊地貌, 山脊地区相对于平坦区域和山谷地区的蓄水能力更差, 淋溶作用相对更弱, 所以土壤pH值相对较大.回归系数绝对值呈自西向东逐渐增大的趋势, 绝对值高值集中在山东半岛、安徽东部和江苏, 这些区域MRRTF对土壤pH的影响较强; 低值集中在河南西北部地区, 表明在这些区域MRRTF对土壤pH的影响较弱.

2.5.2 分位数回归结果

以土壤pH值的0.1~0.9共9个分位点为例, 分析不同pH水平下, 环境因子对其作用的差异, 建模结果及参数检验结果见表 6和图 8.当土壤pH值处于不同水平上时, 各因素对其影响有较大差异.TWI、坡度、海拔和MAT的两种回归系数的置信区间重合较多, 分位数回归结果与MLR结果没有显著差异.

表 6 土壤pH分位数回归系数1) Table 6 Soil pH quantile regression coefficient

水平红色实直线表示MLR回归参数估计值, 虚线为MLR回归参数95%的置信区间上下限; 黑色实折线表示分位数回归参数估计值, 蓝色区域表示分位数回归各分位点95%的置信区域; 当两个模型置信区间重叠时, 可以认为两个结果差异没有统计学意义 图 8 土壤pH分位数回归系数变化 Fig. 8 Variation in regression coefficient of soil pH quantile

在MLR模型中, MRVBF对土壤pH呈显著负作用.在分位数回归中, MRVBF对低分位数(θ为0.1~0.4)水平上的土壤pH呈显著负作用, 且回归系数的绝对值大于MLR回归中的系数绝对值, 表明MRVBF作用强度大于MLR模型; 在土壤pH的高分位数水平时(θ为0.5~0.9), MRVBF的回归系数不显著, 且其对土壤pH作用强度小于MLR模型中的作用强度(表 6).在MLR中坡度对土壤pH呈显著负作用, 在分位数回归中坡度整体上对土壤pH作用显著, 其中当土壤pH处在0.8分位数水平时其作用强度最高, 在其它分位数水平下坡度对土壤pH作用强度保持稳定.

在MLR模型和分位数回归模型中, EVI和NDVI对土壤pH的作用均不显著.分位数回归结果可以看出NVDI在0.1~0.3和0.9的分位数水平上对土壤pH呈显著负作用, 其它分位数水平上NDVI对土壤pH呈负作用, 但不显著.

在MLR模型中, MAP对土壤pH呈显著负作用; 分位数回归结果中, MAP对土壤pH的作用强度随着pH分位数水平的增加呈现逐渐减弱的趋势.在0.1~0.7分位数水平上, MAP的回归系数绝对值较大, 表明其对土壤pH的作用强度大于MLR模型; θ >0.7分位数水平时MAP的分位数回归系数略小于MLR中的回归系数, 表明其对pH作用强度小于MLR模型中的作用强度.

上述结果分析表明, 环境因子对不同的土壤pH等级的影响强度及作用关系(正向或负向作用)存在不同程度的异质性.研究区土壤pH空间分布规律明显(图 6), 由北向南土壤pH值降低.不同土壤pH等级分布在不同的地理区域, 因此不同分位数上的回归结果, 可以视为具有一定空间作用尺度的局部模型.与MGWR模型的作用类似, 本研究中分位数回归可以评价环境因子对pH作用的空间异质性.

3 结论

(1) 研究区土壤pH均值为7.10, 变化范围在4.20~9.18.在不同空间尺度上, 土壤pH存在不同程度的全局和局部空间自相关性, 以空间聚集特征为主.

(2) 3种GWR模型可以有效地预测土壤pH.MGWR优于其它模型, 可以解释建模集土壤pH的64%的变异, 验证集57%的变异, 其残差分布独立性相较其他模型更强.MGWR的土壤pH预测制图结果表明, 研究区土壤pH值空间分布呈由北至南逐渐降低的趋势.

(3) MGWR结果表明, MAP、MRVBF、海拔和MRRTF对研究区土壤pH影响较强, 且存在较强的空间差异性.在江苏北部和山东大部分地区MAP和MRVBF对土壤pH的作用较强, 江苏北部海拔对土壤pH的负向作用强.

(4) MAP和MRVBF等对不同分位数上的土壤pH的作用强度不同.MAP对土壤pH呈显著负作用, 其作用强度随着pH分位数的增加逐渐减弱; MRVBF对低分位数上的土壤pH呈显著负作用.

参考文献
[1] Chen S C, Liang Z Z, Webster R, et al. A high-resolution map of soil pH in China made by hybrid modelling of sparse soil data and environmental covariates and its implications for pollution[J]. Science of the Total Environment, 2019, 655: 273-283. DOI:10.1016/j.scitotenv.2018.11.230
[2] Li Q Q, Li S, Xiao Y, et al. Soil acidification and its influencing factors in the purple hilly area of southwest China from 1981 to 2012[J]. CATENA, 2019, 175: 278-285. DOI:10.1016/j.catena.2018.12.025
[3] Zeng F R, Ali S, Zhang H T, et al. The influence of pH and organic matter content in paddy soil on heavy metal availability and their uptake by rice plants[J]. Environmental Pollution, 2011, 159(1): 84-91. DOI:10.1016/j.envpol.2010.09.019
[4] 陈怀满. 土壤中化学物质的行为与环境质量[M]. 北京: 科学出版社, 2002.
[5] 杨忠平, 卢文喜, 辛欣, 等. 长春市城区土壤中Pb、Cu、Zn、Cd及Cr化学形态特征[J]. 土壤通报, 2011, 42(5): 1247-1255.
Yang Z P, Lu W X, Xin X, et al. Chemical morphological features of Pb, Cu, Zn, Cd and Cr in urban soil of Changchun City, China[J]. Chinese Journal of Soil Science, 2011, 42(5): 1247-1255.
[6] 杨忠芳, 陈岳龙, 钱鑂, 等. 土壤pH对镉存在形态影响的模拟实验研究[J]. 地学前缘, 2005, 12(1): 252-260.
Yang Z F, Chen Y L, Qian X, et al. A study of the effect of soil pH on chemical species of cadmium by simulated experiments[J]. Earth Science Frontiers, 2005, 12(1): 252-260.
[7] 赵明松, 张甘霖, 李德成, 等. 江苏省土壤有机质变异及其主要影响因素[J]. 生态学报, 2013, 33(16): 5058-5066.
Zhao M S, Zhang G L, Li D C, et al. Variability of soil organic matter and its main factors in Jiangsu province[J]. Acta Ecologica Sinica, 2013, 33(16): 5058-5066.
[8] 韩天富, 柳开楼, 黄晶, 等. 近30年中国主要农田土壤pH时空演变及其驱动因素[J]. 植物营养与肥料学报, 2020, 26(12): 2137-2149.
Han T F, Liu K L, Huang J, et al. Spatio-temporal evolution of soil pH and its driving factors in the main Chinese farmland during past 30 years[J]. Journal of Plant Nutrition and Fertilizers, 2020, 26(12): 2137-2149.
[9] 徐慧秋, 黄银华, 吴志峰, 等. 广州市农业土壤As和Cd污染及其对景观异质性的多尺度响应[J]. 应用生态学报, 2016, 27(10): 3283-3289.
Xu H Q, Huang Y H, Wu Z F, et al. Agricultural soil contamination from As and Cd and its responses to landscape heterogeneity at multiple scales in Guangzhou, China[J]. Chinese Journal of Applied Ecology, 2016, 27(10): 3283-3289.
[10] Liu Y, Lv J S, Zhang B, et al. Spatial multi-scale variability of soil nutrients in relation to environmental factors in a typical agricultural region, eastern China[J]. Science of the Total Environment, 2013, 450-451: 108-119. DOI:10.1016/j.scitotenv.2013.01.083
[11] 刘伟, 郜允兵, 周艳兵, 等. 农田土壤重金属空间变异多尺度分析——以北京顺义土壤Cd为例[J]. 农业环境科学学报, 2019, 38(1): 87-94.
Liu W, Gao Y B, Zhou Y B, et al. Multi scale analysis of spatial variability of heavy metals in farmland soils: case study of soil Cd in Shunyi District of Beijing, China[J]. Journal of Agro-Environment Science, 2019, 38(1): 87-94.
[12] 陈清霞, 陆晓辉, 涂成龙. 安顺市土壤pH空间变异及影响因素分析[J]. 环境科学, 2022, 43(4): 2124-2132.
Chen Q X, Lu X H, Tu C L. Spatial variation and influencing factors of soil pH in Anshun city[J]. Environmental Science, 2022, 43(4): 2124-2132.
[13] 李伟峰, 叶英聪, 朱安繁, 等. 近30 a江西省农田土壤pH时空变化及其与酸雨和施肥量间关系[J]. 自然资源学报, 2017, 32(11): 1942-1953.
Li W F, Ye Y C, Zhu A F, et al. Spatio-temporal variation of pH in cropland of Jiangxi province in the past 30 years and its relationship with acid rain and fertilizer application[J]. Journal of Natural Resources, 2017, 32(11): 1942-1953.
[14] 吴正祥, 周勇, 木合塔尔·艾买提, 等. 鄂西北山区耕层土壤pH值空间变异特征及其影响因素研究[J]. 长江流域资源与环境, 2020, 29(2): 488-498.
Wu Z X, Zhou Y, Muhtar·Amat, et al. Spatial variability of soil pH value and its influencing factors in the soil layer of northwestern Hubei Province[J]. Resources and Environment in the Yangtze Basin, 2020, 29(2): 488-498.
[15] 赵明松, 刘斌寅, 卢宏亮, 等. 基于地理加权回归的地形平缓区土壤有机质空间建模[J]. 农业工程学报, 2019, 35(20): 102-110.
Zhao M S, Liu B Y, Lu H L, et al. Spatial modeling of soil organic matter over low relief areas based on geographically weighted regression[J]. Transactions of the Chinese Society of Agricultural Engineering, 2019, 35(20): 102-110.
[16] Zhang C S, Tang Y, Xu X L, et al. Towards spatial geochemical modelling: use of geographically weighted regression for mapping soil organic carbon contents in Ireland[J]. Applied Geochemistry, 2011, 26(7): 1239-1248.
[17] Wang D, Li X X, Zou D F, et al. Modeling soil organic carbon spatial distribution for a complex terrain based on geographically weighted regression in the eastern Qinghai-Tibetan Plateau[J]. CATENA, 2020, 187. DOI:10.1016/j.catena.2019.104399
[18] Song X D, Brus D J, Liu F, et al. Mapping soil organic carbon content by geographically weighted regression: a case study in the Heihe River Basin, China[J]. Geoderma, 2016, 261: 11-22.
[19] Fotheringham A S, Brunsdon C, Charlton M. Geographically weighted regression — the analysis of spatially varying relationships[M]. Chichester: John Wiley & Sons Inc, 2002.
[20] Brunsdon C, Fotheringham A S, Charlton M. Some notes on parametric significance tests for geographically weighted regression[J]. Journal of Regional Science, 1999, 39(3): 497-524.
[21] Fotheringham A S, Yang W B, Kang W. Multiscale geographically weighted regression(MGWR)[J]. Annals of the American Association of Geographers, 2017, 107(6): 1247-1265.
[22] Yang S H, Liu F, Song X D, et al. Mapping topsoil electrical conductivity by a mixed geographically weighted regression kriging: a case study in the Heihe River Basin, Northwest China[J]. Ecological Indicators, 2019, 102: 252-264.
[23] 卢宾宾, 葛咏, 秦昆, 等. 地理加权回归分析技术综述[J]. 武汉大学学报(信息科学版), 2020, 45(9): 1356-1366.
Lu B B, Ge Y, Qin K, et al. A review on geographically weighted regression[J]. Geomatics and Information Science of Wuhan University, 2020, 45(9): 1356-1366.
[24] 南富森, 李宗省, 张小平, 等. 基于MGWR的黄河流域土壤有机碳密度和影响因素分析[J]. 环境科学, 2023, 44(2): 912-923.
Nan F S, Li Z S, Zhang X P, et al. Spatial distribution characteristics and influencing factors of soil organic carbon density in Yellow River Basin based on MGWR model[J]. Environmental Science, 2023, 44(2): 912-923.
[25] 周梦鸽, 杨依, 孙媛, 等. 2016~2020年山东省空气质量时空分布特征及影响因素分析[J]. 环境科学, 2022, 43(6): 2937-2946.
Zhou M G, Yang Y, Sun Y, et al. Spatio-temporal characteristics of air quality and influencing factors in Shandong province from 2016 to 2020[J]. Environmental Science, 2022, 43(6): 2937-2946.
[26] 周志凌, 程先富. 基于MGWR模型的中国城市PM2.5影响因素空间异质性[J]. 中国环境科学, 2021, 41(6): 2552-2561.
Zhou Z L, Cheng X F. Spatial heterogeneity of influencing factors of PM2.5 in Chinese cities based on MGWR model[J]. China Environmental Science, 2021, 41(6): 2552-2561.
[27] 沈体雁, 于瀚辰, 周麟, 等. 北京市二手住宅价格影响机制——基于多尺度地理加权回归模型(MGWR)的研究[J]. 经济地理, 2020, 40(3): 75-83.
Shen T Y, Yu H C, Zhou L, et al. On hedonic price of second- hand houses in Beijing based on multi-scale geographically weighted regression: Scale law of spatial heterogeneity[J]. Economic Geography, 20, 40(3): 75-83.
[28] Liu K, Qiao Y R, Zhou Q. Analysis of China's industrial green development efficiency and driving factors: research based on MGWR[J]. International Journal of Environmental Research and Public Health, 2021, 18(8). DOI:10.3390/ijerph18083960
[29] 李德成, 张甘霖, 王华. 中国土系志·安徽卷[M]. 北京: 科学出版社, 2017.
Li D C, Zhang G L, Wang H. Soil series of China, Anhui[M]. Beijing: Science Press, 2017.
[30] 吴克宁, 李玲, 鞠兵, 等. 中国土系志·河南卷[M]. 北京: 科学出版社, 2019.
Wu K N, Li L, Ju B, et al. Soil series of China, Henan[M]. Beijing: Science Press, 2019.
[31] 黄标, 潘剑君. 中国土系志·江苏卷[M]. 北京: 科学出版社, 2017.
Huang B, Pan J J. Soil series of China, Jiangsu[M]. Beijing: Science Press, 2017.
[32] 赵玉国, 宋付朋. 中国土系志·山东卷[M]. 北京: 科学出版社, 2020.
Zhao Y G, Song F P. Soil series of China, Shandong[M]. Beijing: Science Press, 2020.
[33] Conrad O, Bechtel B, Bock M, et al. System for Automated Geoscientific Analyses(SAGA) v. 2.1.4[J]. Geoscientific Model Development, 2015, 8(7): 1991-2007.
[34] Gallant J C, Dowling T I. A multiresolution index of valley bottom flatness for mapping depositional areas[J]. Water Resources Research, 2003, 39(12). DOI:10.1029/2002WR001426
[35] 刘琼峰, 李明德, 段建南, 等. 农田土壤铅、镉含量影响因素地理加权回归模型分析[J]. 农业工程学报, 2013, 29(3): 225-234.
Liu Q F, Li M D, Duan J N, et al. Analysis on influence factors of soil Pb and Cd in agricultural soil of Changsha suburb based on geographically weighted regression model[J]. Transactions of the Chinese Society of Agricultural Engineering, 2013, 29(3): 225-234.
[36] 陈彦光. 基于Moran统计量的空间自相关理论发展和方法改进[J]. 地理研究, 2009, 28(6): 1459-1463.
Chen Y G. Reconstructing the mathematical process of spatial autocorrelation based on Moran's statistics[J]. Geographical Research, 2009, 28(6): 1459-1463.
[37] 陈建宝, 丁军军. 分位数回归技术综述[J]. 统计与信息论坛, 2008, 23(3): 89-96.
Chen J B, Ding J J. A review of Technologies on quantile regression[J]. Statistics & Information Forum, 2008, 23(3): 89-96.
[38] 张连发. 基于流数据的地理加权回归建模方法的研究[D]. 武汉: 武汉大学, 2019.
Zhang L F. Research on Geographically weighted regression modeling approach based on flow data[D]. Wuhan: Wuhan University, 2019.
[39] 郭龙, 张海涛, 陈家赢, 等. 基于协同克里格插值和地理加权回归模型的土壤属性空间预测比较[J]. 土壤学报, 2012, 49(5): 1037-1042.
Guo L, Zhang H T, Chen J Y, et al. Comparison between Co-Kriging model and geographically weighted regression model in spatial prediction of soil attributes[J]. Acta Pedologica Sinica, 2012, 49(5): 1037-1042.
[40] Guo J H, Liu X J, Zhang Y, et al. Significant acidification in Major Chinese croplands[J]. Science, 2010, 327(5968): 1008-1010.