基于无人机遥感的玉米水分利用效率与生物量监测

韩文霆 汤建栋 张立元 牛亚晓 王彤华

(1.西北农林科技大学水土保持研究所，陕西杨凌 712100; 2.西北农林科技大学水利与建筑工程学院，陕西杨凌 712100；3.西北农林科技大学机械与电子工程学院，陕西杨凌 712100)

0 引言

青贮玉米是指在籽粒蜡熟期后可全株收获制作成青贮饲料的玉米品种，是目前用于饲料的主要玉米品种，在我国农业畜牧业发展和生态建设中具有重要作用[1]。地上生物量是指植物组织地面以上部分的质量，由植物进行光合作用的干物质积累形成，是反映青贮玉米长势、预测产量和收益[2]以及评价农作区生态状况的重要指标。内蒙古地区是我国青贮玉米的主要种植地，该地区处于典型的干旱半干旱气候区，其农业发展在很大程度上依赖于水资源的开发利用，存在水资源短缺、利用率低和农业水资源管理落后等问题[3]，严重制约了当地农业的可持续发展。对干旱半干旱地区青贮玉米生物量及其水分利用状况进行快速、准确的评估可以为农业生产提供管理调控的科学依据，对精准农业生产和水资源高效利用具有重要意义。

传统的生物量和作物水分利用效率估测方法有田间采样法和作物模型模拟法。田间采样法需要破坏性取样，该方法耗时、耗力，且很难获得详细的作物生长时空分布状况[4]；作物模型模拟法则需要通过大量的野外观测数据对模型进行参数校正，如土壤特性、作物参数和气象数据等[5]。随着科技的发展，卫星遥感估测技术已成为大面积作物长势和农情监测重要的研究手段。其中，将卫星遥感数据与作物生长模型相结合来估测作物生物量和水分利用状况是一种主要方法[6-7]。文献[8]结合卫星遥感数据与双作物系数模型估算了玉米和甜菜的蒸散发(Evapotranspiration，ET)，再通过实测的产量数据估算作物水分利用效率(Crop water use efficiency，CWUE)。虽然卫星遥感数据与辐射传输、作物生长和土壤水分平衡模型的集成已被广泛用于区域和田间尺度的生物量和作物水分利用效率制图[9]，但卫星遥感存在重访周期较长、时间分辨率较低的缺陷，而作物模型对遥感数据的观测时间、空间密度以及反演参数的精度要求较高，这些因素会直接影响田块尺度上作物生长模拟的精度[10]。无人机遥感平台具有高频、迅捷、低成本、高空间分辨率、受天气影响小等特点，在农田信息精准监测方面具有良好的应用前景[11-14]。

研究表明，利用无人机遥感技术能够很好地反演多种作物生物量[15-17]，但目前主要还是通过回归分析建立经验模型来估测作物生物量[18-19]。这些研究均以数据为驱动，而无法从作物的生长发育机制来解释[20]和分析作物整个生长过程的水分利用状况。目前，研究重点主要是评估相关的作物生长性状，即通过光谱指数分析不同生长条件下的植物水分状况。比如，文献[21]使用无人机平台结合能量平衡模型与光能利用效率模型估测了柳树的瞬时净辐射、蒸散发、作物总初级生产力、光能利用效率和水分利用效率。但此研究描述的是在单一时间点上的作物水分利用状态，而作物的水分消耗是一个随时变化的过程[9]。因此，将基于作物生长过程的模拟模型与反映作物生长状况的光谱特征相结合会更好地表达作物生长和水分消耗的本质[22]。

目前，将无人机遥感平台与作物模型相结合估算田间作物生物量及水分利用效率的研究较少。通常，用于常规无人机系统监测生物量和水分利用效率的模型应该易于操作、所需参数有限且容易获取。AquaCrop模型[23]是一个不断改进更新、以水分为驱动的作物生长模型，文献[24]将卫星遥感数据输入基于Aquacrop简化的水分利用效率模型中，用来估算作物水分生产力和生物量。通过无人机遥感平台与作物生长模型的集成可以更好地发挥两者优势，对提高田间尺度上作物生长参数的模拟精度、促进无人机遥感平台在精准农业中的应用具有重要意义[25]。

鉴于以上研究现状和存在问题，本文利用无人机遥感平台获取不同生育期大田玉米冠层高分辨率多光谱遥感影像，结合当地气象数据和土壤水分数据，将其整合输入简单的作物水分利用效率模型中，并获取大田玉米水分利用效率与生物量的空间分布，研究分析本文方法在不同水分胁迫下的响应，以期为干旱和半干旱地区农田作物水分利用和灌溉管理提供理论依据和技术支持。

1 材料与方法

1.1 试验地概况

试验于2018年在内蒙古自治区鄂尔多斯市达拉特旗昭君镇(40°26′0.29″N，109°36′25.99″E，海拔1 010 m)开展，该区域属于典型温带大陆性气候。试验地土壤为砂壤土，0～120 cm平均田间持水率为13%(体积含水率)，凋萎系数为5.6%，土壤容重为1.56 cm3[26]。试验玉米品种为“钧凯918”，试验田为半径约60 m的圆形地块(图1)，出苗日期5月18日，收获日期为9月8日。其中玉米播种行距为0.58 m，株距0.25 m，行向由东至西，采用中心轴式喷灌机进行灌溉。

试验区设置5个水分处理区，每个水分处理区下设置3个6 m×6 m的采样区，也是3个重复区，共15个小区(图1)。对扇形区域TRT1进行充分灌溉(田间持水率的95%)，对其他4个区分别以TRT1处理区的不同百分比进行不同梯度的水分胁迫处理，具体水分处理及生育期内降雨及灌溉情况见表1(如TRT3的3个生育期灌溉处理(40/100/80)分别为TRT1的40%、100%、80%)。其中玉米所处的生育期阶段是按照“叶领法”[26]进行划分，分别为营养生长期(Vegetative stage，V期)、生殖生长期(Reproductive stage，R期)和成熟期(Maturation stage，M期)，营养生长期是根据发育完全的叶片的数量来确定的，生殖生长期(R期)是从玉米抽穗至谷粒成熟，R4-R6也为玉米的成熟期。本试验进行水分胁迫处理时间为：V7期开始至整个生殖生长期和成熟期。实际灌溉量和降雨量分别通过安装在喷灌机上的流量计(MIK2000H型)和试验地相邻的标准气象站采集。各区其他管理措施除灌溉量外完全相同。

表1 不同玉米生育期各个处理区灌溉处理和降雨量情况

1.2 无人机系统与影像数据采集

数据采集采用自主研发的无人机多光谱影像采集系统，如图2所示。该系统采用开源飞控PIXHAWK，经纬M600型机架，搭载RedEdge型多光谱相机(MicaSense，USA)，该相机具有5个波长的光谱采集通道：475 nm(蓝光，B)、560 nm(绿光，G)、668 nm(红光，R)、717 nm(红边，RE)、840 nm(近红外，NIR)，且配有光强传感器和固定反射率校正板(Group Ⅷ)，无人机系统参数可参照文献[25]。

在玉米生育期(V-R6)进行影像采集，间隔为5～7 d。采集时间选定当天11：00—13：00，此时天气晴朗、阳光强度稳定、无云且风力较小。无人机飞行高度为70 m，多光谱相机镜头垂直向下，地面分辨率为4.77 cm/像素。每次都采用固定航线飞行，航向和旁向重叠度均为85%。每次采集后的无人机影像采用Pix4DMapper软件(Pix4DInc.，Switzerland，https:∥www.pix4d.com/)进行拼接。主要拼接过程为：导入图像和地面控制点坐标，生成点云，辐射校正，生成数字表面模型、正射影像及归一化植被指数[27](Normalized differential vegetation index，NDVI)分布图。使用15个样区的shp文件利用R语言程序(R-3.5.1，https:∥www.r-project.org/)对植被指数分布图进行批量化裁剪，提取各样区的植被指数。

1.3 地面数据获取

(1)气象站数据采集：农业气象站由河北省清易电子科技有限公司组装和调试，位于距离试验区1 km处，其下垫面为苜蓿，高度保持在12 cm左右。监测参数包括空气温度(Air temperature，Ta)、相对湿度(Relative humidity，RH)、2 m处风速(u2)、太阳净辐射(Netsolar radiation，Rn)和降雨量(Precipitation，P)，除降雨量外，每间隔30 min采集一次数据。

(2)土壤含水率数据采集：采用干燥法测量试验区玉米根系层的土壤含水率。从V7期开始，大约每隔7 d在每一个采样区中心位置采集一次土样，并在灌溉前后和有效降雨后进行加测。取土深度分别为10、20、30、45、60、90、120 cm，将7个不同深度土层的土壤含水率取平均值，代表6 m×6 m采样区土壤含水率的均值，并采用线性插值法对土壤含水率进行处理，得到各样区土壤含水率的逐日数据序列。

(3)地上部玉米干生物量数据采集：在无人机光谱数据采集完后，在采样区内随机选取一个0.5 m×0.5 m的区域进行地上玉米生物量的采集。将玉米从茎基部剪下，获得完整的冠部，编号分装好后放入干燥箱内，在105℃下杀青2 h，恒温(80℃)干燥至恒质量后称量。称量时采用精度为0.01 g的电子天平，获取采样区的地上部玉米干生物量。

1.4 基于无人机遥感的水分利用效率与生物量估算方法

1.4.1基于生物量的水分利用效率及标准化水分利用效率模型

植物生长和生物量积累是光合作用过程中CO2同化和水分通过植物气孔的蒸腾交换的结果，其中所需的能量由太阳辐射提供[28]。考虑到这一过程的生理基础，生长模拟的经典方法包括将蒸腾的水分转换为生物量的效率，这种经典的模型被称为水利用效率(Water use efficiency，WUE)模型，是许多作物模型的基础[29]。WUE表示作物单位水分消耗获得的干生物量，生物量是实际蒸腾量Tc,adj时间积分与WUE乘积的结果[30]，即

(1)

式中mBio——t0～t时间段内增长的地上部干生物量，g/m2

t0——出苗日期

t——生物量获取日期

WUE——作物水分利用效率，kg/m3

由于WUE容易受不同的区域、年份和管理措施等的影响而产生差异，文献[31]指出，针对不同气候环境下的蒸发需求，需要进行标准化处理。已有许多研究使用饱和气压差(Vapor pressure deficit，VPD)和参考蒸发蒸腾量(Reference evapotranspiration，ET0)来进行标准化处理，在联合国粮农组织发布的AquaCrop模型中认为针对气候差异使用ET0是最佳水分利用效率标准化方法[23,32]。据此提出了标准化水分利用效率(Normalization of water productivity，WP*)的计算模型[33]。文献[34]研究表明，在内蒙古地区多变的气候条件下，WP*是一个比WUE更为稳定评估作物生产力和生长速率的参数。

(2)

其中

Tc,adj=ktkstkswET0

(3)

式中WP*——标准化水分利用效率，g/m2

kt——蒸腾系数

kst——温度胁迫系数

ksw——水分胁迫系数

1.4.2无人机遥感平台与水分利用效率模型的结合

水分利用效率与生物量关系的核心参数是随时间变化的蒸腾量Tc,adj。Tc,adj又由蒸腾系数kt决定，而基于传统的方法很难获取大面积农田上连续的作物蒸腾系数kt，因无人机遥感平台可以获取高时空分辨率的作物冠层信息，为获取kt信息提供新的思路。研究表明，通过基于植被指数经验模型可以很好地估测作物蒸腾系数，提供田间尺度上详细连续的kt分布[35]。文献[8]提出了可以适用于不同气候条件下的、基于归一化植被指数NDVI的蒸腾系数计算模型

(4)

式中kt,max——作物达到完全覆盖时的最大蒸腾系数

NDVImax——作物达到全覆盖时最大NDVI

NDVImin——裸土状态下NDVI[27]

α——模型经验系数

为实现上述计算方法，本研究实测的NDVImax、NDVImin分别为0.82、0.14；根据FAO-56作物需水量计算手册[36]中双作物系数计算kt,max，目前在中国利用双作物系数法估算和区分灌溉条件下农田蒸腾和蒸发量的研究较多[37]，通过FAO-56提供的计算公式，使用气象与作物参数数据，可以计算得到kt,max为1.2；本研究采用简单线性模型，因此α取1。代入式(4)得

kt=1.8NDVI-0.25

(5)

通过2个水分利用效率模型和无人机遥感植被指数NDVI之间关系的描述，本研究基于无人机遥感结合水分利用效率模型估算大田青贮玉米水分利用效率及生物量的流程如图3所示。具体流程为：①通过无人机遥感平台获取玉米各个生育期不同日期整个地块的NDVI分布。②通过式(4)、(5)获取整个地块的蒸腾系数kt，通过时间线性插值获取整个生育期内每天的kt分布图。③结合每天的气象数据和土壤水分数据得到校正系数kst、ksw和ET0，即可通过式(3)获得整个生育期玉米不同时期各个处理区∑ktkswkst和∑Tc,adj分布。④使用实测的生物量数据通过式(1)即可计算每个处理下的WUE和WP*分布。虽然已有很多研究计算过玉米的WUE，但是WUE随环境的易变性和WP*数据的稀缺性，需要在当地对其进行校正[35]。⑤基于以上步骤，同样可以计算各小样区的WUE和WP*经验值，然后采用WUE和WP*经验值进行生物量的估算并验证，从而获取整个处理区的生物量分布图。

1.4.3环境胁迫因子kst、ksw计算

以上kt表示的是标准状况下没有胁迫状况的作物蒸腾系数，而实际的作物生长会在不同的环境、管理措施下偏离正常的生长。本文根据FAO-66手册，将温度胁迫系数kst和水分胁迫系数ksw计入模型中。

温度胁迫系数kst在0～1之间，0表示温度太低不能满足作物生长，大于0小于1表示由于受作物生长所需热量不能完全满足，作物蒸腾量只能部分转换为作物生物量，1表示作物生长所需热量能完全满足，作物蒸腾量可全部转换为作物生物量[38]。本文采用文献[39]提出的计算方法计算kst。

(6)

式中SrelT——相对水分胁迫水平

T0——玉米最适宜生长的温度,℃

Tb——不能满足玉米生长的临界温度,℃

Tm——实测的日平均气温,℃

根据文献[39]：T0、Tb分别取30、8℃。当Tm>T0时，kst=1；当Tm≤T0时，kst计算式为

(7)

式中Sn——kst上限，取1

Sx——kst下限，取0.001

r——速率因子，取15[40]

由于本试验设置了不同的水分胁迫处理，同时存在未及时灌溉的情况，玉米会产生水分胁迫。土壤根系层低于作物适宜生长的含水率时会产生水分胁迫，影响作物生物量累积。水分胁迫系数计算式为[41]

(8)

式中θ——作物根系层土壤含水率，%

θfc——田间持水率(体积含水率)，取13%[26]

θwp——凋萎系数(体积含水率)，取5.6%[26]

θj——适宜土壤含水率，取8.9%

p——发生水分胁迫之前根系中所消耗水量与土壤总有效水量的比值，取0.55

1.5 精度评价及方法验证

为验证不同水分处理下本文方法的适用性，将TRT1、TRT2和TRT4每个处理区中的2个小样区(2小样区各参数平均值)作为水分利用效率模型建立区，将TRT1、TRT2和TRT4处理区剩余的一个小样区和TRT3和TRT5区平均值分别作为验证数据集。采用决定系数R2、均方根误差(Normalized root mean square error，RMSE)、一致性指数d(Index of agreement)3个指标进行精度评价。一致性指数越大表示监测值与预测值的一致性程度越高，其计算式为

(9)

Pi——预测值n——样本数量

2 结果分析

2.1 玉米不同生育期的基础气象数据与作物参数变化

2018年玉米生育期参考蒸发蒸腾量ET0、降雨量和灌溉量的变化过程见图4。ET0整体呈下降趋势，其中5月末至6月初呈波动增大趋势，而在7—9月初呈减小趋势，青贮玉米整个生育期ET0为554 mm。青贮玉米主要降雨量集中在7月玉米营养生长期的后期，在降雨量不足时进行补充灌溉。灌溉量主要集中于玉米生育期前期，不同处理区灌溉量在V7期之前一致，从V7-R6期按TRT1的不同百分比进行灌溉(表1)。

图5为TRT1、TRT2和TRT4处理区青贮玉米土壤含水率变化曲线。本试验从V7期开始(图5中绿色虚线处)进行水分胁迫处理，在水分胁迫处理前各试验区的降雨灌溉量一致，但由于存在不同处理区间的土壤异质性，TRT2在进行水分处理前和进行水分处理后期某些时间段含水率大于TRT1和TRT4处理区，但在水分胁迫处理后大部分生育期各处理区的土壤含水率与灌溉梯度一致，从大到小依次为TRT1、TRT2、TRT4。

图6为玉米生育期内的地上部干生物量的变化趋势。本试验从2018年6月19日，V6初期开始采集生物量数据，随着时间推移而逐渐增大，大约至R5(蜡熟期)末期，玉米的地上部干生物量基本不再增加。因为在蜡熟期后期，干物质累积过程结束，干生物量达到最大值。不同水分处理下的干生物量增长趋势也不同，总体上TRT1区的生物量高于其余2个处理区。在不同的水分胁迫处理下，相比于TRT1处理区，TRT2和TRT4两处理区最终的累积干生物量分别减少20%和22%。

图7为玉米生育期内冠层植被指数NDVI的变化曲线。在所有的16次测量中，3个处理区的NDVI都呈现出先缓慢增长(V2-V5)，快速增长(V6-VT)，平缓不再增长(R1-R4)到逐渐减少(R5-R6)的变化趋势。NDVI主要与作物的覆盖度有关，由于在V6-VT期，玉米的叶片发育迅速，因此玉米的覆盖度和冠层NDVI快速增长。VT期的玉米达到全覆盖，NDVI达到最大值(约为0.82)且基本不再增长。在R5-R6期玉米的叶片开始部分枯萎发黄和脱落，导致覆盖度和NDVI减小。且从图7可知NDVI对水分胁迫响应明显，3个处理区的NDVI在不同的水分胁迫状况下具有明显的梯度。

2.2 水分胁迫系数和温度胁迫系数

图8a为TRT1、TRT2和TRT4处理区通过对土壤含水率插值的水分胁迫系数ksw计算结果，水分处理从V7期开始，之前生育期都为充分灌溉，且灌溉频率高，ksw为1。由于青贮玉米在VT-R6生育期内不同处理区间存在不同水分处理，具有不同的胁迫梯度，具体表现为：在V7-VT期，TRT1、TRT2和TRT4处理区平均ksw分别为1.00、0.97和1.00；在R1-R3期分别为0.93、0.90和0.83；在R4-R6期分别为0.94、0.72和0.69。虽然TRT1区是按照100%进行灌溉，但可能由于采样点土质等其他原因也出现了较长时期的水分胁迫。同时，由于V7-VT期灌溉降雨量较多(图4)，3个处理区ksw差异不大，都基本不存在水分胁迫状况。

图8b为玉米生育期每天的平均气温与相应温度胁迫系数kst，在温度影响下kst整体表现为生育期前期和末期小，6、7、8月这3个月值大，有利于玉米生长，也是玉米的主要生长期。由于5月下旬和9月上旬昼夜温差大且日间温度相对较低导致一天平均温度过低，导致kst很小，一部分时间甚至接近于0，整个生育期的kst平均值为0.78，说明当地气候对于玉米生物量累积具有较大的影响。

2.3 水分利用效率(WUE)

图9为各区根据实测的各生育期地上部干生物量与相对应日期基于光谱数据计算的玉米生育期内实际蒸腾量(∑Tc,adj)拟合结果；结果表明不同水分处理下生物量与蒸腾量(∑Tc,adj)之间具有极显著的正相关关系(P<0.001)，决定系数R2分别为0.92、0.95和0.93。但从图9中可以看出，在玉米的生育期后期(R5、R6)，相关性明显较低。在不同的水分处理下，TRT1、TRT2和TRT4 3个处理区整个生育期的蒸腾量分别为257、253、183 mm。在水分效率分析中线性关系的截距与0没有区别，因此将斜率作为WUE[24,35]，TRT1、TRT2和TRT4处理区的水分利用效率分别为7.4、6.9、7.9 kg/m3。不同水分处理下的WUE没有明显的趋势，TRT4处理区虽然具有较高的WUE，但其生物量没有TRT1和TRT2处理区高，TRT1处理区的生物量最高且有较高的WUE。

图10为各区根据实测的收获时地上部干生物量与基于光谱数据的玉米生育期总的实际蒸腾量∑Tc,adj通过式(1)计算获取的不同水分处理区(TRT1、TRT2和TRT4)的WUE分布图，具体处理步骤见1.4.2节流程图及其说明。由表2可知，3个处理区TRT1、TRT2、TRT4 WUE均值分别为7.6、6.7、8.6 kg/m3，其中TRT1、TRT2 WUE平均值与图9中拟合的结果相近，而水分胁迫更严重的TRT4较之于TRT1、TRT2表现出了很大的差异。WUE平均值随着水分胁迫并没有明确的变化趋势，但是随着水分胁迫的增加，WUE空间变异性呈现出增长的趋势。如表2所示，TRT1、TRT2和TRT4处理区WUE的变异系数分别为10.2%、12.6%和17.2%。

表2 处理区TRT1、TRT2和TRT4 WUE分布统计参数

2.4 标准化水分利用效率(WP*)

图11为采用1.4.2节的方法对水分利用效率进行标准化处理获得的WP*，生物量与∑ktkswkst之间同样具有极显著的正相关关系(P<0.001)，且R2均不小于0.93。不同水分胁迫下的WP*分别为31.6、29.4、33.7 g/m2，都与FAO-66手册推荐值33.7 g/m2相近。

图12为根据实测的收获时地上部干生物量与基于光谱数据的玉米生育期总的∑ktkswkst通过式(2)计算获取的不同水分处理区(TRT1、TRT2和TRT4)的WP*分布图，具体处理步骤见1.4.2节流程图及其说明。由表3可知3个处理区TRT1、TRT2和TRT4的均值分别为37.9、31.1、39.5 g/m2。其中TRT2的WP*与图11中拟合的结果相近，TRT1、TRT4表现出了相对较大的差异，可能是由于受土壤背景的影响和使用收获时单个生物量计算使得相对于图12中整个生育期拟合结果偏大。但从图12中可以看出，大部分WP*还是分布在30～35 g/m2区间内。WP*平均值随着水分胁迫并没有明确的变化趋势。但是随着水分胁迫的增加，与WUE类似，WP*分布的标准差和变异系数呈现出增长的趋势。如表3所示，TRT1、TRT2和TRT4处理区WP*的变异系数分别为8.6%、9.7%和13.5%。

表3 处理区TRT1、TRT2和TRT4 WP*分布统计参数

2.5 模型检验及生物量估算

为验证该模型的适用性和可重复性，本文使用2.3节和2.4节模型拟合的WUE与WP*对验证集的生物量进行估算，并通过实测的地上部生物量进行验证。数据验证集包括TRT1、TRT2和TRT4每个处理区下预留的一个采样区和与TRT1、TRT2、TRT4水分胁迫模式不同的TRT3和TRT5处理区。图13a、13b分别为使用2.3节和2.4节中模拟的结果(TRT1、TRT2和TRT4的WUE和WP*)进行生物量估测的验证结果。结果显示在整个生育期内基于标准化水分利用效率WP*的验证精度(RMSE为256 g/m2，d=0.96)略高于WUE，2种方法误差主要来自于玉米地上部干生物量大于1 500 g/m2的生育期后期(R5-R6)。V-R4期采用2种方法进行验证时，RMSE分别降低为91.7 g/m2和126 g/m2，一致性指数d上升为0.98。

由于TRT3和TRT5处理区的水分胁迫状况不同于TRT1、TRT2和TRT4处理区，所以分别采用TRT1、TRT2和TRT4区的WUE平均值(7.4 kg/m2)与WP*平均值(31.5 g/m2)对TRT3和TRT5处理区的生物量估算验证，验证不同水分状况下WUE与WP*估测生物量的适用性。图13c、13d为TRT3和TRT5处理区的验证结果，总体上WP*验证精度(RMSE为195 g/m2，d=0.97)高于WUE(RMSE为306 g/m2，d=0.93)。在使用WUE对不同水分胁迫下的TRT3和TRT5处理区进行验证时，在全生育期都会呈现出高估的现象。而WP*的验证精度远高于WUE，但使用WP*进行估测时发现，当实测生物量小于1 000 g/m2时表现出高估而大于1 000 g/m2时低估的现象。综上，无论在相同水分胁迫还是不同水分胁迫条件下，在估测生物量时采用WP*效果比WUE更好。

2.6 生物量空间分布

由2.5节结果可知，使用标准化水分利用效率WP*模型对生物量有良好的估算效果，本文采用式(2)对TRT3和TRT5处理区进行生物量估算，获取不同水分胁迫下的生物量空间分布状况。虽然表1中的TRT3和TRT5处理区水分处理相差不多，但是由于土壤异质性导致TRT3实际的土壤含水率大于TRT5：水分处理后的TRT3和TRT5ksw分别为0.91和0.75。由图14可知，在不同水分胁迫下的TRT3和TRT5处理区的生物量同样差异明显，TRT3和TRT5处理区生物量均值分别为1 788、1 401 g/m2。

3 讨论

生物量与作物水分利用效率是反映青贮玉米产量以及指导灌溉的重要参数，而植被不同时期的生长状况可以通过其冠层光谱特征来反映，再结合作物生长模型，可从作物生长发育机制的角度更为精准地估算作物生物量及水分利用状况。本文利用无人机平台获取玉米全生育期内高时空分辨率的多光谱图像，并将气象数据与土壤水分数据一起输入到水分利用效率模型中，较为精准地估算了作物生物量及水分利用效率，并获取了高空间分辨率的大田青贮玉米水分利用效率分布图。结果表明使用基于无人机遥感数据的∑Tc,adj及∑ktkswkst与实测地上部生物量拟合的水分利用效率(WUE)和标准化水分利用效率(WP*)具有极显著的相关性，决定系数R2都在0.92以上，这与文献[24]利用卫星数据构建土壤调节植被指数SAVI计算WP*精度相一致。并且较文献[35,40]利用多个遥感平台估算水分利用效率，本研究通过单一无人机平台就能获取更加细致的作物信息，减少融合多个卫星遥感平台数据带来的误差，同时获取田块尺度上更为详细的WUE和WP*的空间分布，有利于评估田块尺度内大田玉米水分利用效率的空间变异性(图10和图12)。本研究在田块尺度中水分利用效率空间变异性一部分由土壤异质性带来的，本试验采样区虽然是均匀地布置在试验区内，但由于土壤的空间异质性会导致作物生长差异较大，对数据采集分析具有一定的影响。同时，由表2和表3可以看出，不同的水分胁迫状况(即干旱)同样可以导致作物水分利用效率产生较大的空间变异性，这与文献[42]分析相一致。从另一方面也可以看出基于高分辨率无人机多光谱遥感平台相对于卫星遥感的优势[43]。不同土壤状况下的作物长势(生物量、水分利用效率)的高度变化信息的获取，意味着使用无人机遥感结合作物模型有利于更加精细的灌溉管理。

目前国内外使用无人机遥感平台估算作物生物量主要还是在特定的生育期寻找特定的光谱指数与特定日期的生物量建立经验模型，如文献[44]利用棉花主要生育时期的无人机近红外影像数据，提取4种不同的植被指数，通过与棉花地上生物量的实测值建立拟合关系，分析了不同植被指数在棉花各生育时期的估算效果并对其进行了验证，最高决定系数R2为0.86。文献[45]利用平均绿色归一化植被指数(GNDVI)估算了干豆开花和补荚中期的生物量，相关系数最高为0.73。本文通过植被指数NDVI计算作物关键参数Tc,adj输入到简单的作物生长模型中，最终拟合结果决定系数R2都在0.92以上，较之前的研究相关性更高。STEDUTO等[29,46]认为大田作物绝大部分生育期的WP*是不变的，因此，本文基于无人机遥感数据拟合出单一的WUE或WP*可以较为准确地同时获取多个生育期的生物量，在这方面本研究极大地提高了基于无人机遥感平台估算作物生物量的效率和准确度。

结果显示，本文方法估计生物量在大部分生育期虽然较为准确，但是在高值区，即生育期后期(R5-R6)误差明显，导致整体精度并不高。同一水分胁迫下，分别使用TRT1、TRT2和TRT4处理区WUE和WP*进行各区的生物量验证，RMSE只有259 g/m2和256 g/m2(图13a、13b)，低于预期。这与文献[24]的研究结果相似，在R5期，其研究中雨养条件下的玉米生物量误差最高能达到800 g/m2左右。造成这种高值区误差较大的原因可能是因为在模型中∑Tc,adj和∑ktkswkst为单调递增，且本研究玉米生育后期的水分胁迫较严重(图8)，所以模拟生物量的增长率会随着ksw急剧减小而快速下降，但可能由于玉米的抗旱性，实际的玉米生长过程生物量增长率并不会如此快速地下降，因此出现低估现象。针对这一问题，文献[47]采用Logistic和Sigmoid方程提出了一种非线性的WP*模型，并将其替换了原来Aquacrop模型中WP*的算法，提高了生物量的模拟精度。基于此，后续研究可以使用非线性的WUE或WP*模型与无人机遥感平台相结合来提高此方法对后期高值区玉米生物量的估算精度。尽管如此，水分利用效率线性模型在同样水分处理条件下，基于拟合的WUE、WP*能够精确地估计大部分生育期的大田玉米生物量，具有很强的应用价值。

以往基于该模型计算水分利用效率和生物量都是在充分灌溉条件下或是没有完整分析水分胁迫处理对该模型的影响[35,40]。本文通过设置不同的水分胁迫处理分析了模型的适用性。如前所述，在同种水分胁迫处理下，使用拟合的WUE与WP*在V-R4期有很好的精确性，RMSE分别为126、91.7 g/m2，在R5-R6生育期的误差较大。而使用TRT1、TRT2和TRT4的WUE、WP*平均值验证TRT3、TRT5处理区时，WUE(RMSE为306 g/m2)的估算精度远不如WP*(RMSE为195 g/m2)(图13c、13d)，说明水分胁迫对WUE的影响较大，因此使用WUE估算不同水分状况下的玉米生物量会造成较大的误差。同样，在玉米生长期V-R4内，使用TRT1、TRT2和TRT4处理区WP*平均值估测不同水分状况下的TRT3和TRT5处理区时的精度(图13d)不如同种水分胁迫下的估算精度(图13b)且会出现高估的现象，原因是TRT3、TRT5处理区的V7-VT期的灌溉水量仅为充分灌溉的40%，低于TRT1、TRT2和TRT4处理区的灌溉水量(100%、65%、65%)。低灌溉降雨量下的TRT3、TRT5处理区在生育期前期因水分胁迫比其它处理区严重，导致这时期的水分生产力低于TRT1、TRT2和TRT4处理区，所以才会出现高估现象。综上水分胁迫对于该模型的WUE生物量验证精度具有较大的影响，但是使用WP*对不同水分胁迫下的生物量进行估算具有可以接受的精度。

4 结论

(1)在不同水分胁迫下，将基于无人机多光谱遥感植被指数NDVI估算的蒸腾系数kt及气象、土壤水分数据输入到作物模型中，模拟得出的作物水分利用效率WUE与标准化水分利用效率WP*与实测生物量具有极显著的相关性，WUE和WP*决定系数都在0.92以上。在估测WUE和WP*的基础上，获取了大田玉米WUE和WP*的高空间分辨率分布图，同时分析土壤水分胁迫对WUE和WP*空间异质性的影响，表明土壤水分胁迫增加会导致其空间异质性增大。

(2)在同一水分胁迫下，使用拟合的WUE和WP*对生物量的估测精度几乎相同，在玉米V-R4生育期估测精度较高(RMSE分别为126、91.7 g/m2，d均为0.98)，但在R5-R6生育期内精度不高。尽管如此，基于同一水分胁迫下拟合的WUE和WP*仍具有很好的应用性。

(3)在不同水分胁迫下，使用WUE和WP*估测生物量时，WUE容易受到水分胁迫影响，精度不理想，而WP*的精度较高(RMSE为195 g/m2，d=0.97)。因此，WP*更适合用于估算不同水分胁迫下的生物量。