基于Logistic模型的盐渍农田向日葵根长密度分布模拟

任志鹏，何良伟

（1.黑龙江省农垦科学院，哈尔滨 150000；2.武汉大学水利水电学院，武汉 430072）

0 引 言

土壤盐渍化是制约我国耕地生产力的重要因素[1]。以内蒙古河套灌区为例，盐渍化面积约33.33万hm2，占总耕地面积的63.8%，是限制当地耕地资源扩展和生态文明建设的重要因素[2-4]。根系作为作物从土壤汲取水分和养分的重要器官，是农田水肥管理、作物生长模拟关注的重点。然而，由于根系藏于地下，难以与土壤分离，相比于作物地上部分，对根系形态的观测和定量研究十分有限[5,6]。

根系是由主根和不同级别的侧根组成的多级分叉结构，在空间上呈非均匀分布。为量化根系空间分布，学者们尝试利用根系深度的数学函数对特定深度的根系指标（RD模型）或者累积根系指标（CRD模型）进行表征。RAATS[7]、VILIAM[8]和KAGE[9]基于根系随土层深度迅速减少的事实，分别提出根表面积、根重量或根长密度（RLD）的指数RD模型。GALE、GRIGAL[10]和SCHENK、JACKSON[11]分别利用指数函数和Logistic函数提出根长或根重的CRD模型。而不论是RD或是CRD模型，研究表明模型的参数除了与土壤水盐环境相关，还受生育期、作物遗传信息和气候的影响[10-14]。多数研究虽提出模型与环境变量的定性关系，但很少有研究考虑模型和生育期内环境动态变化的定量关系。考虑到作物生育期内根系分布的动态变化，WU[15]改进了RD模型，假定作物不同生长阶段的根长密度分布具有较好相似性，提出归一化根长密度分布模型（NRLD）用以表征作物不同生长阶段下根长密度的相对分布。ZUO[16]发现小麦的归一化根长密度分布不受气候条件、土壤条件和品种的影响，具有统一的变化趋势，以此建立了小麦的NRLD模型。NING[6]提出了改进的多项式NRLD模型，验证在棉花、小麦、玉米和水稻等作物上具有良好的适用性。

然而，研究发现在不同土壤水盐环境下作物根系相对分布呈动态变化。一方面土壤盐浓度增加，抑制作物根系生长，导致根生物量、根表面积、根长等指标降低[17-20]。另一方面作物在面对盐胁迫时，根系生长过程会做出响应，以抵御盐胁迫。研究发现盐胁迫下，作物会优先发育中、细根，改变粗细根的比例和根系分布形状适应胁迫环境[21,22]，这导致盐渍农田下的作物根系生长和分布与正常环境下有所不同。以上因素都可能使作物在不同生长的根长密度分布的相似性变差，影响NRLD模型的适用性。因此，针对盐渍农田中作物根系的生长分布规律需要进一步研究，同时，NRLD模型在盐渍环境中对作物生育期根系分布模拟的适用性，还需要进一步分析。

基于此，本文选取具有一定耐盐性的经济作物向日葵为研究对象，比较分析了2016-2017年不同盐分田块向日葵根长密度的动态变化规律，提出了Logistic根系动态分布模型（LRDG），并建立了模型参数与土壤水盐和积温之间的定量关系，实现水盐胁迫下向日葵根长密度分布的动态模拟。

1 材料和方法

1.1 试验区概况

试验于2016-2017年在内蒙古河套灌区义长灌域永联试验站开展。河套灌区地处东经106°20′-109°19′，北纬40°19′-41°18′，温带大陆性季风气候，具有沙漠型气候特点，日照充足，干燥多风，降雨量少。年平均气温6.1 ℃，无霜期117～136 d，封冻期达180 d，年平均太阳辐射的总量为6.5 MJ/m2，年平均日照时数3 200 h。该地区年降水量139～222 mm，多集中在夏季，年蒸发量高达2 200～2 400 mm。灌区主要的土壤质地类型有粉壤土、粉土和黏壤土等，土壤田间持水率在23.67%～30.24%，耕地含盐量约0.08%～0.65%。灌区的主要作物包含向日葵、玉米、小麦和西瓜等，试验选用中等耐盐耐旱的经济作物向日葵作为试验对象。

1.2 试验设计与观测内容

大田试验在低（S1）、中（S2）和高（S3）3个盐分水平展开，2016年播种前各田块0～100 cm土壤饱和电导率（ECe）依次为1.87～2.71、3.51～8.28和7.54～10.32 dS/m；2017年灌水间隔时间和取样时间都较短，播种前各田块0～100 cmECe相对较小，依次为1.07～2.84、1.30～3.28和1.05～5.95 dS/m。在大田试验开始前进行机械化的表层松土工作，松土厚度10 cm，之后会在各区田的长边方向铺置3条宽度为70 cm的地膜。地膜铺设完毕后，需要对区田进行充分灌溉，在土壤水分达到适宜播种的土壤含水量后，与当地农民同步播种。播种方式选用膜上穴播，一穴3粒，出苗后约1周间苗，保留长势最好的一株，每条地膜播种2行，播种的行距和株距为50 cm×50 cm。其他田间管理措施一致。

大田试验的观测指标主要包括：气象、地下水位、作物指标和土壤理化指标等。气象数据由永联试验基地内的自动气象站观测得到，观测的指标包含空气温度、太阳辐射、空气湿度、大气压力、降雨量和蒸发量等指标。地下水位由每个田块旁设置的地下水位观测井测得。土壤含水率采用称重法测量，具体而言利用土钻取得不同深度的新鲜土壤，放置于铝盒内密封，称重后置于烘箱，105 ℃烘干8 h，取出再次称重，确定土壤的质量含水率，由对应深度的土壤容重确定土壤体积含水率（θv）。土壤含盐量由ECe表征，ECe由土水比为1∶5的土壤混合溶液的上层清液电导率（EC1∶5）转换确定[23]。

作物根系由ROOT-700根系监测系统进行观测。2016年的观测时间分别为播种后的29、36、43、50、57、64和105 d，2017年的观测时间为播种后的28、38、48、58、68、78和110 d。在植株播种前，先在各区田的地膜中轴位置平行地安装4根间距为100 cm，长度200 cm，内径7.1 cm的微根管，埋入位置大约在地下1.5 m处，露出地表0.5 m，露出的部分要用黑色的胶带和棉布包裹住，避免微根管透光。微根管埋设角度为与地表呈45°的夹角，垂直埋深在100 cm左右。微根管在测定时，把长约35 cm的扫描探头深入微根管中，以每20 cm作为一个观测深度控制逐层地深入微根管中，分层扫描，扫描之后重新将微根管重新用黑色胶带密封，以防止水汽地进入，最后用扫描得到的图像利用WinRHIZO和WinRHIZOTron根系分析软件对向日葵的根系进行识别，然后对根长等根系指标进行提取。

1.3 归一化根长密度模型

为探究NRLD模型在盐渍农田向日葵根长密度分布中的适用性，对向日葵归一化根系深度zr处的标准根长密度NRLD（zr）进行了计算：

式中：Lr是最大根深。

在计算NRLD（zr）的基础上，选取ZUO[24]、NOVAK[8]和WU[15]建立的幂（Power）、指数（Exponential）和多项式（Polynomial）模型对作物不同深度的NRLD（zr）进行了拟合。

Power：

Exponential：

Polynomial：

式中：a为幂函数模型参数；b和c为指数模型参数；d、e、f和g为多项式模型参数。

1.4 Logistic累积根长密度模型及参数量化

Logistic模型参数物理意义清晰，并且能够表征复杂的根系分布状况。基于SCHENK等[11]和FAN等[13]的研究，本文尝试利用Logistic模型对向日葵根系分布进行拟合：

式中：z为土壤深度，cm；Y(z)为由地表至z深度处的根长累积占比，无量纲；参数α决定了表层根系占比，参数α越小，表示表层根系占比越大；参数β决定了根系聚集生长的土层位置，β越小，表示根系聚集生长所在的土层位置越深。

Logistic输出为根长累积占比（CRD），为了与其他NRLD模型比较，应将其输出转化为标准化根长密度：式中：n为土层层数，本试验中n=7；zi为第i个土层的深度，cm（i=1，2，…，n），z0=0；Lr为作物最大根深，取100 cm。

考虑到土壤水盐以及作物生长对根长密度分布的影响，尝试建立Logistic模型参数与上述影响因素之间的定量关系，从而实现不同土壤环境以及作物生长阶段根长密度分布的动态表征。其中，作物生长过程与积温密切相关[24,25]，因此选用积温（AT）对作物生长阶段进行定量表征。考虑到作物根区的水盐非均匀分布[26]，分别尝试利用表层（0～10 cm）土壤水盐分浓度CML1、下层（占据最大根深50%的土壤深度）土壤水盐分浓度CML2和平均（整体根区）土壤水盐分浓度CML3。其中，土壤水盐分浓度CML（单位：mg/cm3）的计算方法如下：

式中：λ为转化系数，由试验拟合得到，本研究选用0.432 4；ρb为土壤干容重，g/cm3；θv为土壤体积含水率。

考虑到模型参数α与AT、CML之间的相互关系，采用xex函数和幂函数或指数函数相结合的方程分别对参数α进行拟合：

式中：l、m和k为拟合系数，其中，l和m为考虑水盐胁迫的增大系数，k为决定参数α达到峰值的系数，与向日葵营养生长期到生殖生长期的转变节点有关，转变节点为1/k。

考虑到模型参数β与AT、CML之间的相互关系，本文采用线性函数和幂函数或指数函数相结合的方程分别对参数β进行拟合：

式中：n、p为拟合参数，其中，n表示参数β随AT下降的缩小系数，p为考虑了水盐胁迫对参数β的增大系数。

2 结 果

2.1 土壤水盐动态变化

作物生育期内各田块的土壤剖面θv整体呈下降趋势（见图1）。2016年的播种后天数DAS=1～27 d和2017年DAS=1～30 d，根区的θv下降显著，尤其2016年S1田块0～25 cm土层的θv由0.35 cm3/cm3骤降至0.25 cm3/cm3，而对于S2和S3田块，0～15 cm土层的θv下降0.07～0.10 cm3/cm3，其余土层θv略微下降。2016年DAS=51～69 d和2017年DAS=60～70 d，S1田块70～90 cm土层的θv在0.20～0.25 cm3/cm3，而S2和S3田块除了0～25 cm土层，其余土层仍大致保持0.25 cm3/cm3以上。2016年DAS=69～105 d和2017年DAS=80～112 d，θv下降速率减缓，而盐分处理越高，θv相对越高，S1田块保持在0.2 cm3/cm3左右，S2和S3田块基本处于0.25～0.30 cm3/cm3。

2016年DAS=1～51 d和2017年DAS=1～60 d，S1田块的表层和深层土壤的ECe相差不显著，分布均匀，整个土壤剖面的ECe小于3 dS/m而S2和S3田块随着深度增加，土层的ECe增长明显（见图2）。2016年DAS=51～69 d和2017年DAS=60～80 d，S1的土壤盐分分布与初始状态接近，整体盐分含量变化不大。但在S2和S3田块的表层ECe逐渐增加，尤其在S3田块中，表层土壤的ECe迅速增加，0～15 cm土层ECe能到10 dS/m以上。

图2 2016-2017年不同盐分田块（S1、S2和S3）土壤电导率（ECe）动态变化Fig.2 Soil conductivity dynamics in plots with different salinity (S1, S2 and S3) from 2016-2017

2.2 根长密度分布

各盐分田块作物不同深度RLD及其占比见图3、图4。可以看到2016年DAS=29 d和2017年DAS=28 d，各田块RLD相近，主要集中在0～30 cm土层，70%～90%的根系主要集中在0～45 cm土壤。2016年DAS=36～50 d和2017年DAS=38～48 d，根系下扎迅速，S1和S2田块的向日葵根系在此阶段达到100 cm，而S3田块在2016年DAS=64 d和2017年DAS=68 d才达到100 cm土层。2016年DAS=50 d和2017DAS=58 d以后，各田块中部土层出现根系聚集生长现象：S1和S2田块，0～30 cm土层根系占比在40%左右，30～75 cm土层根系占比在40%以上；相比田块S1和S2田块，S3田块作物表层根系占比骤减，中部土层聚集现象更为显著，如2016年S3田块峰值在5 cm/cm3以上。2016年DAS=64～105 d和2017年DAS=78～110 d，S1和S2田块的0～30 cm土层根系占比有所回升，S3田块0～30 cm根系占比趋于稳定，保持在20%左右。

图3 2016-2017年不同盐分田块（S1、S2和S3）各土层作物根长密度（RLD）动态变化Fig.3 Dynamics of crop root length density of each soil layer in different salinity plots (S1, S2 and S3) from 2016-2017

图4 2016-2017年不同盐分田块（S1、S2和S3）各土层作物根长占比变化Fig.4 Change in the percentage of crop root length of each soil layer in different salinity plots (S1, S2 and S3) from 2016-2017

2.3 静态根长密度分布模拟

将2016-2017年向日葵大田试验的根长密度数据归一化处理后，分别对各生育阶段用Power、Exponential、Polynomial和Logistic模型进行拟合，分析各模型在不同盐分水平下对向日葵静态根长密度分布模拟的适用性。根据拟合结果（见图5和表1），Power和Exponential模型在盐渍农田下，对向日葵各生育阶段的静态根长密度分布模拟较差，Polynomial和Logistic模型对向日葵静态根长密度分布模拟效果更好。

表1 不同根长密度分布模型在2016-2017年各田块的模拟结果Tab.1 Simulation results of different root length density distribution models in S1、S2、S3 plots in 2016-2017

图5 各根长密度分布模型模拟结果比较Fig.5 Comparison of simulation results of each root length density distribution model

Power、Exponential、Polynomial和Logistic模型拟合的决定系数R2分别为0.548、0.444、0.869和0.757，RRMSE分别为0.460、0.471、0.225和0.338，Polynomial和Logistic模型对向日葵静态根系的整体拟合效果明显好于Power和Exponential模型。从不同盐分水平上看，除2017年S1田块，Power、Exponential和Polynomial模型的拟合精度随盐分增大呈下降的趋势，在2016年S1田块，RRMSE分别为0.334、0.393和0.093，而在2016年S3田块，RRMSE分别为0.571、0.645和0.300。Logistic模型除2017年S1和S3田块外，其他田块拟合RRMSE皆小于0.300。

2.4 NRLD模型对根系动态分布模拟

NRLD模型中的Polynomial模型能够拟合盐渍农田向日葵各生育阶段的静态根长密度分布状况，但能否用Polynomial模型来模拟整个生育期向日葵动态根长密度分布，还要考察模型参数随生育阶段推进以及在不同盐分条件下的变化。以下统计各盐分田块在所有生育阶段拟合的模型参数（见表2）。

表2 Polynomial模型参数在2016-2017年各田块中整个生育期内的变化情况Tab.2 Polynomial model parameters throughout the growth period of sunflower in the S1、S2、S3 plots in 2016-2017

从不同盐分田块看，比较各田块参数的平均值、最大值和最小值，Polynomial模型的参数在各田块的变化范围在不同区间，同时参数浮动范围大，说明模型参数受盐分水平的影响。其中，参数d影响表层标准化根长密度，2016年S2和S3田块的参数d小于0是不合实际情况的。从生育期内的参数变化上看，Polynomial模型的参数随生育阶段变化幅度大，差异系数较高，尤其在2016年的S2田块参数d差异系数达到-1.44，在2016的S1田块，参数e、f和g差异系数分别为-7.92、-5.57和3.11，说明模型参数也受生育阶段推进的影响。

2.5 LRDG模型对根系动态分布模拟

2016年各盐分田块Logistic模型参数α和β变化情况见图6。不同盐分田块参数α在幼苗期缓慢提高，在蕾期达到峰值，且田块盐分水平越高，α的增幅越明显，其中S2、S3田块α整体增幅约为11和25；花期至收获期，α参数值逐渐下降。不同田块参数β在整个生育期总体呈现逐渐降低趋势，盐分水平越高，参数β整体越高。花期后，参数β趋于平缓，S1、S2和S3田块在收获期分别为0.038、0.061和0.076。

图6 Logistic模型α和β参数在2016年不同盐分田块（S1、S2和S3）的动态变化Fig.6 Dynamics of logistic model α and β parameters in different salinity plots (S1, S2 and S3) in 2016

结合参数变化规律，利用2016年试验数据，构建参数α和β与AT、CML的定量关系；采用xex函数和幂函数或指数函数相结合的方程分别对参数进行拟合，拟合效果见表3。

表3 参数α的量化模型及量化效果Tab.3 Quantification model of parameter α and quantification effect

采用线性函数和幂函数或指数函数相结合的方程分别对参数β进行拟合，拟合效果见表4。

由表3和表4可知，采用表层土壤水盐作为解释变量对参数α和β的拟合效果最佳。参数α以AT、CML1为解释变量的lCMLmATe-kAT函数形式的拟合效果最好，R2和RRMSE分别为0.884和0.211。参数β以AT和CML1为解释变量的e-nAT+pCML函数形式的拟合结果最好，R2和RRMSE分别为0.831和0.161。由此，得到可描述盐渍农田向日葵根系动态分布的LRDG模型：

表 4 参数β的量化模型及量化效果Tab.4 Quantification model of parameter β and quantification effect

根据构建的LRDG模型，计算2017年大田试验的作物根系在各生育阶段的土层分布，计算结果见图7。结果显示，LRDG模型计算得到的参数与Logistic模型与实测根长密度分布直接拟合的参数接近，参数α的R2和RRMSE分别为0.880和0.425，参数β的R2和RRMSE分别为0.760和0.297。

3 讨 论

比较不同盐分田块作物RLD发现，S3田块RLD总体高于S1和S2田块，这与作物抗逆生长策略有关。研究表明部分作物为降低盐胁迫对地上部分的生长抑制作用，会将更多的资源发育地下部分，增强根系吸收养分和水分的能力，提高对逆境胁迫的适应性[27-29]。并且，相比于S1和S2田块，S3田块最大根深变化速度更为缓慢，这是由于作物为适应盐胁迫会优先发育侧根，侧根上的细根是地下部分最具活力和生理活性的部分，在水分和养分的提取中起重要作用[21,22]。

在生育初期，作物最大根深变化慢，根系在土壤表层聚集，盐分胁迫一定程度上抑制了向日葵根系向深层土壤生长，因此S3田块根系在表层占比最高。在蕾期，作物根系下扎迅速，根系在中部土层横向延伸发展，表层根系占比减少，中部土层根系占比增多，甚至大于表层根系，根系集中分布于30～75 cm土层，盐分处理水平提升，刺激根系的侧根延伸生长，导致中部土层根系集中生长的现象越明显。根系最适宜生长的土层为30～75 cm，这是因为积盐的作用使得土壤表层的盐胁迫愈发显著，而深层土壤的含水率较高，通气孔隙较少，这些限制条件使得中部土层更适宜根系生长[30]；在生育末期，整个土壤剖面的根系衰老，表层土壤根系占比有所回升，但高盐田块表层根系占比回升不显著。可以看出在作物整个生育期，受到土壤盐分的影响，根系的相对分布形状发生了改变，由生育前期的上层分布多下层少，到生育后期的上层分布少中下层分布多。

不同模型的静态根长密度分布模拟结果显示，Power、Exponential等NRLD模型不能很好地表征盐渍农田的向日葵根长密度分布。这主要由于Power和Exponential模型仅能模拟根长密度随土层深度减少的情形，而难以反映生育中期向日葵根系在中间土层集中生长的现象。而参数较多的Polynomial模型可以较好地反映向日葵根系集中生长现象，但模型参数所表达的物理意义并不明确。并且，在盐渍农田向日葵根系在生育期内的标准化根长密度变化幅度较大，尤其在部分土层，标准差达到0.630以上，盐渍农田的向日葵标准根长密度分布相似性差，不符合NRLD模型的假定。Logistic模型能够模拟盐渍农田的向日葵静态根长密度分布，各参数也有明确的物理意义，其中参数α表示表层根系含量，参数β与根系聚集生长的土层位置相关，由此参数能够与环境因素建立定量关系。建立的LRDG模型能够考虑根系在整个生育期对土壤水盐运移变化的响应，适用于盐渍农田向日葵根系动态分布的定量表征。

4 结 论

本文基于2 a的野外大田试验，分析了不同盐分田块土壤水盐运移和向日葵RLD的变化，提出了适用于表征盐渍农田向日葵RLD动态变化的Logistic模型，并对模型参数进行了量化，得到如下主要结论。

（1）盐分处理水平越高，土壤ECe上升越显著，但θv下降越缓。在作物收获时，S2和S3田块在0.25 cm3/cm3以上，而S1田块在0.2 cm3/cm3左右，S1、S2和S3田块ECe分别能达到3.0、6.5和9.0 dS/m。

（2）S3田块的RLD整体明显高于S1和S2田块，尤其在成熟期，RLD在土壤剖面上的分布峰值在5 cm/cm3以上；而S3田块的最大根深变化则较为缓慢，晚于S1和S2田块10～20 d后达到最大根深。并且，盐分胁迫影响根系在生育期的相对分布，幼苗期根系占比随土层深度减少，0～45 cm土层根系占比达到70%以上，蕾期后，根系在30～75 cm土层聚集生长，根系占比达到40%以上。

（3）Polynomial和Logistic模型对盐渍农田的向日葵静态根长密度分布模拟效果较好，RRMSE分别为0.225和0.338。同时，盐渍农田的向日葵的标准根长密度分布相似性较差，部分标准化根深处的NRLD（zr）标准差在0.630以上，不符合NRLD模型的假定。本文利用CML1和AT对Logistic模型参数α和β进行了定量表征，并以2017年大田试验验证，参数α和β的量化效果RRMSE分别为0.425和0.297，以此得到可描述盐渍农田向日葵根系动态分布的LRDG模型。