盐碱土壤Kostiakov入渗模型参数的BP预报模型

李昊哲，樊贵盛

(太原理工大学水利科学与工程学院，太原 030024)

0 引 言

盐碱地指的是土壤中盐分含量过高，以致其直接影响了各类作物的正常生长。据统计，中国约有9 913 万hm2的盐碱地，占世界盐碱地的1/10，其中有1 333.33 万hm2盐碱地被认为具有农业利用潜力[1]。中国盐碱地不仅资源量大且分布广泛，在东北、华北、西北内陆地区以及长江以北沿海地带都分布着大量盐碱土壤。对盐碱地的合理开发和有效治理，已经成为解决目前人口增加和耕地逐渐减少矛盾的重要突破口。自然盐碱地的形成主要是由于地下水和底层土壤中的盐分，随着土壤毛细管作用随地下水运移至地表，并在水分蒸发后留在土壤表层[2]。对此，目前人们利用大量淡水灌溉使这些盐分逆向运移回地下，这是改良盐碱地最重要的水利措施。因此，盐碱地土壤水分入渗特征的研究对盐碱地的开发和改良有着重要的意义。

目前土壤盐碱化已经上升为全球问题，并受到学术界的重点关注，众多学者围绕土壤盐碱化的问题展开了众多研究和探讨，相对而言，对盐碱地土壤水分入渗的研究还比较少。胡顺军[3]针对盐碱地的冲洗定额做了大量研究和探索，找到了理论上更严密，操作上更简便的计算方法。褚琳琳[4]研究了喷灌对盐碱地土壤湿润锋的运移、水分再分布特征及盐分淋洗的影响，得出了喷灌中对湿润锋运移影响最大的是喷灌强度和黏粒含量。王一民[5]通过室内的模拟实验，模拟膜下滴灌技术对盐碱地进行改良，得出土壤初始含水率和滴头流量对盐碱土水盐运移的影响规律。张金龙[6]分析了盐碱土灌溉冲洗与排水改良技术参数，并结合天津滨海区的实际情况，提出了适应当地自然环境的灌排改良工程技术参数的估算方法。

土壤传输函数法[7]是指利用易获得的土壤基本理化参数，如土壤含水率、土壤容重、土壤有机质、土壤含盐量等，预测土壤不易获得的其他物理参数的方法。黄元仿[8]成功利用土壤理化参数估算了不同含水率或基质势下的土壤导水率。王志强[9]通过土壤基本物化参数，准确估算饱和含水率、田间持水率等水力特征参数。也有人在土壤水分入渗参数预测方面应用了土壤传输函数，如郭华[10]研究了冻融条件下的水分入渗，利用非线性方法对Kostiakov二参数入渗模型参数做出了预测。但有关于盐碱地入渗参数的预测，还未曾有过报道。本文以盐碱地土壤为研究对象，利用大田盐碱地土壤水分入渗资料，建立了BP神经网络训练样本，试图利用BP神经网络的方法，实现对土壤水分入渗Kostiakov模型参数的预测，为盐碱地的开发和改良提供技术支持。

1 材料与方法

1.1 试验区土壤条件

盐碱荒地土壤入渗试验选在山西省应县盐碱荒滩。这里位于山西省北部，朔州市东部，气候寒冷，年平均温度7 ℃左右，属北温带大陆性气候，年降雨量360 mm。试验在盐碱地大面积分布的4个试验点进行，分别位于应县的4个乡镇：臧寨、杏寨、大黄崴、大临河。所选试验点全部为原生盐碱荒地土壤，土壤盐分约在2.1～5.1 g/kg，pH在7.1～8.8，属于中度盐碱地，其上只生长着稀疏的盐碱植物。地下水埋深在1.0～2.0 m，土壤有机质含量为0.9～1.4 g/kg左右，4个试验区的土壤类型都有较多的砂粒含量，占60%左右，黏粒含量最少约占7%，均属于砂壤土。

1.2 试验方案与方法

本次盐碱荒地土壤入渗试验主要利用双套环入渗仪进行，其中外环直径60 cm，内环直径26.0 cm，内外环高度均为25 cm，下环深度约20 cm。试验时为使内外环在同样的水位入渗，采用水位平衡装置保证内外环水位齐平。一般来说，土壤入渗试验达到稳渗状态的时间一般在60 min左右，为确保试验成功，本试验水分稳定入渗的时间取90 min。

土壤的常规理化性质参数包括土壤含水率、土壤干容重、土壤质地、土壤有机质含量、土壤含盐量、土壤pH值等。土壤含水率测定用传统烘干称重法测量获得；土壤干容重通过腊封法进行测定得到；土壤质地通过筛分+比重计法得到筛分曲线，然后分析土壤的颗粒级配，进而确定土壤质地；试验点土壤有机质含量的测定是利用化学的方法通过重铬酸钾容量法来测定；试验点土壤总盐量的测定利用烘干质量法；土壤pH值由pH计测定。

1.3 Kostiakov二参数入渗模型

土壤水分入渗指的是降雨或灌溉过程中，水分从地表进入土壤内部的过程。随着土壤状况和外界条件的不同，土壤水分入渗的过程也会发生变化。对于非盐碱土壤的入渗过程，常用Kostiakov二参数模型和Kostiakov-Lewis三参数模型。而对于盐碱地，选用Kostiakov二参数模型更加简便，且精度较高[11]。因此，本文选Kostiakov二参数模型：

I=ktα

(1)

式中：I为累积入渗量，cm；k为经验入渗系数，一般指入渗开始后第一个单位时段末的累积入渗量，在数值上等于第一个单位时段内的平均入渗率，cm/min；α为经验入渗指数，反映土壤入渗能力的衰减速度。

1.4 BP神经网络方法

BP神经网络指的是(Back Propagation)误差反传，即误差反向传播算法。这种算法在本质上是一种神经网络学习的数学模型，所以有时也称为BP模型。它是利用最速下降法，通过反向传播来不断调整网络的权值和阈值，使网络的误差最小。BP神经网络模型包括输入层、隐含层和输出层。本文中的输入层为土壤常规理化参数，包括含水率、容重、质地、有机质、pH值等，输出层为kostiakov二参数入渗模型中的入渗系数k与入渗指数α，其结构示意图如图1所示。

图1 BP神经网络结构图

由于BP神经网络神经元的响应函数为Sigmoid函数，要求输出输入变量的值在(-1，1)之间，这就要求在将样本数据导入到神经网络模型之前，需要对数据样本进行归一化处理，本文选取如下公式进行：

(2)

同时，根据数据样本的特点，本文隐含层的激活函数选tansig函数，输出层的激活函数选择purelin函数。BP神经网络训练的学习率一般取值0.01～0.8；BP神经网络训练的最大学习迭代次数一般取值1 500～3 000；BP神经网络训练要求精度一般取值0.001～0.000 1。

2 BP神经网络预报模型的建立

2.1 入渗参数k、α的影响因素

影响入渗的因素有很多，其中以表征土壤基本理化性质的含水率、容重、质地和有机质等为主。土壤体积含水率越高，其水吸力越小，土壤水分下渗的速度越慢，入渗参数k和α越小。土壤干容重越大，土壤越密实，土壤孔隙率越小，连通性越差，单位面积入渗通量就越小，k和α越小。土壤质地是土壤固相物质各粒级土粒的配合比例，一般用土壤中黏粒、粉粒及砂粒的百分含量来表征土壤质地，其中黏粒含量越高，土壤越重，土壤颗粒越细小，吸附能力越强，导致同样时间内水分通量减小，土壤水分的入渗速度减小，k和α的值也会随之变小。土壤有机质有助于土壤团聚体的形成，能有效增强土壤的稳定性，使土壤结构稳定、孔隙均匀，因此有机质含量较低的土壤在入渗过程中，由于其存在较大空隙，会导致初始入渗速率会很大，即k的值很大，但随着水分的下渗，会导致土壤内部发生挤压和坍塌，堵塞部分大孔隙，减缓土壤入渗的速率，使土壤入渗速率的衰减速度变快，α变小。

相对于普通非盐碱土壤，盐碱土壤的含盐量较高，因而其pH值也较大，而土壤的含盐量和pH值也和土壤水分的入渗有很大关系[12]。土壤含盐量指的是土壤中常见的八大离子(Na+、K+、Ca2+、Mg+2+、Cl-、HCO-3、CO2-3、SO2-4)含量之和，盐碱土壤中的大量盐分会直接改变土壤结构，特别是Na+含量的升高会引起土壤颗粒的膨胀与分散，从而使土壤的透水性变差，影响土壤水的渗透。土壤pH值反映了土壤溶液的酸碱性，盐碱土壤的碱度较高，而高碱度条件下，土壤胶体中的Ca2+、Mg2+离子易被土壤溶液中的Na+离子置换，并生成沉淀，阻塞土壤孔隙，影响水分通过，从而影响水分的入渗过程。

入渗系数k在数值上接近于第一个单位时段末(1 min)的累积入渗量，故k应只与地表表层的土壤理化指标有关，这里取0～20 cm的土壤理化指标作为输入变量，即0～20 cm土壤体积含水率θ0、容重γ0、有机质含量G、黏粒含量w1、粉粒含量w2、全盐量ε以及pH值δ。而对于入渗指数α反映的是土壤入渗能力的衰减程度，不仅与耕作层有关还与犁底层有关，即取20～40 cm土壤体积含水率θ1、土壤容重γ1、黏粒含量w3、粉粒含量w4作为入渗指数α的输入变量。

本文所选入渗模型为Kostiakov二参数入渗模型，故选取入渗系数k、入渗指数α为本次试验的输出变量。

2.2 预测样本与检验样本的确定

根据大田原生盐碱荒地土壤的系列入渗试验得到累积入渗量I与入渗历时t的数据样本，结合Kostiakov二参数入渗公式，借助MATLAB7.0拟合得到入渗参数k和α，建立了120组累积入渗量I90与入渗参数间的对应关系如表1所示。

表1 Kostiakov二参数入渗模型参数值

再根据试验中配套测定的土壤基本理化参数，得到上述120组入渗参数所对应的土壤理化参数如表2所示。

表2 土壤理化参数

根据表1和表2所示入渗模型参数与土壤基本理化参数之间的一一对应关系建立本次试验的数据样本。最后通过筛选，在120组数据样本中，选取110组作为建模训练样本，6组作为预测模型的检验样本。

2.3 BP神经网络模型的建立与检验

2.3.1 BP神经网络模型的建立

采用Matlab7.0建立BP模型结构如下：

net=newff(minmax(traininput),

[m,1],{'tansig','purelin'},'trainlm')

式中：newff为建立前馈神经网络的函数；traininput为模型的输入样本，minmax(traininput)为样本的范围；[m,1]为隐含层的神经元节点个数，本文在预测入渗系数k时，试算确定的节点数为[20，1]，预测入渗指数α时，试算确定的节点数为[15，1]；tansig为隐含层的激活函数；purelin为输出层的激活函数；trainlm为BP神经网络的反向传播训练函数。

将样本输入网络进行训练，神经网络学习训练经过迭代后收敛，并且网络完全准确地识别了学习样本，建立了输入参数与输出参数之间复杂的非线性映射关系，训练结果如下：

k=purelin(iw2×(tansig(iw1×F+b1))+b2

α=purelin(iw2×(tansig(iw1×E+b1))+b2

其中：

F=[θ0，γ0，w1，w2，G，ε，δ]

E=[θ0，θ1，γ0，γ1，w1，w2，w3，w4，G，ε，δ]

式中：iw2为模型隐含层到输出层的权值；iw1为模型输入层到隐含层的权值；b1为模型输入层到隐含层的阈值；b2为模型隐含层到输出层的阈值。

2.3.2 入渗模型参数的预测

(1)入渗系数k的预测。本文对入渗系数k进行预测时，BP神经网络参数值设定如下：训练的学习率设定为0.01；最大学习迭代次数设定为1 500；要求精度设定为0.000 1；隐含层节点数最后经过试算确定为20，最后当计算步数为873时，训练精度达到要求，训练结束，训练结果如表3所示。

(2)入渗指数α的预测。本文对入渗指数α进行预测时，BP神经网络参数值设定如下：训练的学习率设定为0.01；最大学习迭代次数设定为1 500；要求精度设定为0.000 1；隐含层节点数最后经过试算确定为15，最后当计算步数为231时，训练精度达到要求，训练结束，训练值如表4所示。

(3)入渗参数的误差分析。根据上述对Kostiakov入渗参数的预测得到入渗系数k和入渗指数α的预测值，与实测值进行比较，计算其绝对误差与相对误差，计算结果如表5所示。

表3 入渗系数k的预测矩阵表

表4 入渗指数α的预测矩阵表

表5 k、α误差分析表

由表5可以看出：入渗系数k预测值与实测值相比绝对误差在0～0.061 9 cm/min之间，平均值为0.003 7 cm/min；相对误差在0～3.31%之间，平均值为0.29%；入渗指数α预测值与实测值相比绝对误差在0～0.044 5之间，平均值为0.005 2；相对误差在0.01%～9.40%之间，平均值为1.28%，两入渗参数预测的相对误差平均值较小，预测精度较高。

(4)90 min累积入渗量的综合误差分析。根据Kostiakov二参数入渗公式，结合上述入渗参数的预测结果，计算出累积入渗量I90的预测值，与样本实测值相比较，误差分析结果如表6所示。

由表6可以看出：90 min累积入渗量I90的绝对误差平均值为0.239 9 cm，最大值为2.209 5 cm，最小值为0.000 3 cm；相对误差的平均值为2.37%，最大值为9.04%，最小值为0，联合误差完全在可接受的范围之内。由此表明所建立的BP神经网络预报模型具有较高的精度。

表6 I90误差分析表

2.3.2 预报模型的检验

用试验预留的6组数据检验模型的预测精度，结果如表7所示。

表7 检验结果表

根据检验结果可以看出：6组检验样本中，入渗系数k的相对误差最大值为2.07%，最小值为0.29%，平均值为0.90%；入渗指数α的相对误差最大值为1.68%，最小值为0.09%，平均值为0.91%；90 min累积入渗量I90的相对误差最大值为5.59%，最小值为0.02%，平均相对误差为1.89%。实例检验结果表明Kostiakov入渗模型的两入渗参数以及90 min累积入渗量的相对误差能控制在6%以下，精度较高，由此说明应用所建立的BP预测模型能获得较好的效果。

3 结 语

本文建立的基于BP神经网络模型的土壤水分入渗参数预测模型可以较好地反映Kostiakov二参数模型中入渗参数与多种影响因素之间的复杂的非线性关系，入渗系数k、入渗指数α和90 min累积入渗量I90的相对误差平均值分别为0.29%、1.28%和2.37%。同时检验结果的相对误差平均值分别为0.90%、0.91%和1.89%，误差都在理想精度范围内。说明利用盐碱土壤基本理化参数对入渗参数的预测是可行的，所建模型能较好预测出土壤水分入渗参数，能为干旱地区盐碱土壤的改良提供技术支持，同时在一定程度上丰富了土壤传输函数的理论成果。

本文成功利用土壤传输函数法对盐碱土壤的水分入渗参数做出预报，所选择的输入变量获得了满意的结果。经过反复试算，确定用土壤含水率、容重、质地、有机质、全盐量以及pH等参数对Kostiakov二参数土壤入渗模型参数k，α进行预测是可行的。

□

[1] 王遵亲.中国盐渍土[M].北京:科学出版社,1993.

[2] Doran J C, Turnbull J W. Australian trees and shrubs: species for land rehabilitation and farm planting in the tropics[R].Canberra: Australian Centre for International Agricultural research,1997.

[3] 胡顺军,田长彦,宋郁东.基于土壤饱和入渗理论计算盐碱地冲洗定额[J].土壤学报,2010,47(3):563-567.

[4] 褚琳琳,康跃虎,陈秀龙,等.喷灌强度对滨海盐碱地土壤水盐运移特征的影响[J].农业工程学报,2013,29(7):76-82.

[5] 王一民,虎胆·吐马尔白,弋鹏飞,等.盐碱地膜下滴灌水盐运移规律试验研究[J].中国农村水利水电,2010,(10):13-17.

[6] 张金龙,张 清,王振宇.天津滨海盐碱土灌排改良工程技术参数估算方法[J].农业工程学报,2011,27(8):52-55.

[7] 朱安宁,张佳宝,陈效民,等.封丘地区土壤传递函数的研究[J].土壤学报,2003,40(1):54-59

[8] 黄元仿,李韵珠.土壤水力性质的估算——土壤转换函数[J].土壤学报,2002,39(4):517-522.

[9] 王志强,刘廷玺,等.寒冷干旱区土壤水力特性参数的模拟估算[J].沈阳农业大学学报,2004,10,35(5):426-428.

[10] 郭 华,樊贵盛.冻融土壤Kostiakov入渗模型参数的非线性预报模型[J].节水灌溉, 2015,(11):1-4.

[11] 樊贵盛,李 尧,苏冬阳,等.大田原生盐碱荒地入渗特性的试验[J].农业工程学报,2012,28(19):63-70.

[12] 王 雪,樊贵盛.改善原始盐碱荒地入渗能力措施的实验研究[J].灌溉排水学报,2009,28(6):46-49.