王晓珂 刘婷婷 许桂玲 冯跃华,3* 彭金凤 李杰 罗强鑫 韩志丽 卢苇 PHONENASAY Somsana
(1贵州大学农学院,贵阳550025;2黔西南州农业农村局,贵州 兴义562400;3贵州大学/山地植物资源保护与种质创新教育部重点实验室,贵阳550025;第一作者:411282a22he.cdb@sina.cn;*通讯作者:fengyuehua2006@126.com)
氮素是作物生长发育过程中需求量最大的矿质元素之一[1],与光合作用中光合色素的合成密切相关[2],能显著影响作物产量和品质[3]。氮素的丰缺会对环境和作物本身造成不同程度的影响[4]。因此,快速、准确、实时监测作物氮素营养状况,是氮肥合理使用的前提[5]。
随着高光谱遥感的发展,不同学者进行了一些关于作物氮素状况反演模型的研究,主要集中在作物叶片氮素含量、叶片氮素积累量和地上部氮素积累量方面。RAJEEV等[6]通过高光谱遥感构建小麦孕穗期叶片氮素含量和地上部氮素积累量的估算模型,结果表明,植被指数对叶片氮素含量的预测精度要优于地上部氮素积累量;FITZGERALD等[7]利用冠层叶绿素含量指数(CCCI)预测小麦冠层氮素营养状况时认为,采用高光谱遥感方法估测氮素时,应使用地上部氮素积累量作为度量单位;XUE等[8]研究表明,与水稻叶片氮素含量相比,叶片氮素积累量与冠层反射率在全生育期内相关性一致,叶片氮素积累量估算模型的预测能力较好;ZHAO等[9]研究表明,水稻叶片氮素积累量对冠层光谱参数的响应比叶片氮素含量更敏感。由此可见,关于表征作物氮素营养状况的最佳指标,不同的学者研究结果不尽相同。为此,本研究设置不同施氮素水平和不同水稻品种的大田试验,来构建水稻氮素营养诊断模型,通过比较得到最佳的诊断模型以及最适的氮素营养指标。
试验于2019年在贵州省黄平县旧州镇寨碧村进行。试验田耕层土壤理化指标:pH 5.02,有机质18.38 g/kg,速效氮209.50 mg/kg,速效钾65.73 mg/kg,速效磷4.56 mg/kg。
供试水稻品种为Q优6号和宜香优2115;供试氮肥为尿素(含N质量分数为46.2%)、磷肥为过磷酸钙(含P2O5质量分数为16.0%)、钾肥为氯化钾(含K2O质量分数为60.0%)。
试验采用裂区设计,主区处理为水稻品种(V),设置2个水平,分别为V1(Q优6号)、V2(宜香优2115);副区处理为施氮素量(N),设置5个水平,分别为N0(0 kg/hm2)、N1(75 kg/hm2)、N2(150 kg/hm2)、N3(225 kg/hm2)、N4(300 kg/hm2),副区面积25.9 m2;每种施氮量处理采用分次施肥法,基肥、分蘖肥、促花肥、保花肥的施氮量分别占总施氮量的35%、20%、30%、15%;磷肥和钾肥的用量分别为P2O596 kg/hm2、K2O 135 kg/hm2,磷肥作基肥一次施入,钾肥作基肥、促花肥各施50%。每个处理重复3次。4月18日播种,5月28日移栽,行株距30 cm×20 cm,每丛插1苗。田间管理按照高产栽培管理措施进行。
1.3.1 冠层光谱
水稻冠层光谱测量采用美国Analytical Spectral Device(ASD)公司生产的FieldSpec@4 Standard-Res型地物光谱仪。波长范围350~2 500 nm。冠层光谱测定选择在天气晴朗风力小于3级时进行,测定时间范围为北京时间10∶00—15∶00,以保证有较高的太阳高度角。传感器探头(视场角为25°)垂直向下,距冠层顶垂直高度约0.75 m,地面视场范围直径为0.33 m。分别在水稻的孕穗期(7月23日)、抽穗期(8月4日)和成熟期(9月19日)进行冠层光谱测定,每个小区记录4个观测点,每个观测点测量5条光谱曲线,以其平均值作为该小区的光谱数据。测量过程中,每个小区测量前对仪器进行标准白板校正,以消除环境变化所带来的影响。
1.3.2 水稻地上部干物质量
水稻地上部生物量样品的采集与冠层光谱测定同时进行。采用破坏性取样的方法,根据各个小区的平均茎蘖数,每小区选择生长一致且具有代表性的水稻植株4株,按叶片、茎鞘、穗等器官进行分类,在105℃下杀青0.5 h后于80℃下烘干至恒质量,分器官称其干质量,再计算单位面积地上部各器官的干物质量(kg/hm2)。
1.3.3 叶片氮素含量(Leaf N content,LNC)、叶片氮素积累量(Leaf N accumulation,LNA)和地上部氮素积累量(Aerial part N accumulation,APNA)
将各器官的烘干样品粉碎,采用H2SO4-H2O2消化-凯氏定氮法测定样品的氮素含量(%),之后计算各器官氮素积累量(kg/hm2),地上部氮素积累量(kg/hm2)为各器官氮素积累量之和。
1.4.1 植被指数的计算
光谱数据用光谱仪自带的ViewSpec Pro软件进行处理并导出,剔除不稳定的波段,即350~400 nm、1 350~1 480 nm、1 780~1 990 nm和2 400~2 500 nm[10];选择4种常用的植被指数进行计算,分别是RVI(比值植被指数)、DVI(差值植被指数)、NDVI(归一化植被指数)和PVI(垂直植被指数),其计算公式见表1。
表1 植被指数计算公式
1.4.2 氮素营养模型建立及模型评价
构建任意两波段组合的4个植被指数,通过相关性分析,利用决定系数(r2)确定每个氮素营养模型(LNC反演模型、LNA反演模型和APNA反演模型)的最优波段组合。
以最优波段组合构建的植被指数为自变量,以氮素营养指标为因变量进行回归分析构建反演模型,包括线性函数模型、指数函数模型、对数函数模型、幂函数模型、抛物线函数模型。利用含量梯度法[14]将数据集按2∶1的比例分为建模集和测试集,每个氮素营养指标获得构建模型和测试模型的数据分别为60个和30个。利用相关指数(R2)、均方根误差(RMSE)和相对误差(RE)来评价模型的精度。
式中,yi为实测值,为预测值,为实测值的平均值,n为数据集的样本量。
图1 水稻冠层原始光谱反射率与LNC的相关系数图(n=90)
水稻冠层原始光谱反射率与LNA的相关性变化趋势(图2)和水稻冠层原始光谱反射率与LNC的相关性变化趋势一致,但在不同波段的相关程度略有差异。在401~1 332 nm、1 481~1 540 nm和1 991~2 400 nm波段范围内,水稻冠层原始光谱反射率与LNA的相关程度整体上低于水稻冠层原始光谱反射率与LNC的相关程度;在1 333~1 349 nm和1 541~1 779 nm,水稻冠层原始光谱反射率与LNA的相关程度整体上高于水稻冠层原始光谱反射率与LNC的相关程度。水稻冠层原始光谱反射率与LNA呈负相关关系的波段为401~731 nm、1 333~1 349 nm、1 481~1 779 nm和1 991~2 399 nm,负相关性达到极显著水平的波段为401~726 nm、1 481~1 640 nm、1 712~1 779 nm、1 991~2 399 nm,负相关程度最大的波段为674 nm,相关系数为-0.835。水稻冠层原始光谱反射率与LNA呈正相关关系的波段为732~1 331 nm,正相关性达到极显著水平的波段为737~1 148 nm,正相关程度最大的波段为780 nm,相关系数为0.772。
图2 水稻冠层原始光谱反射率与LNA的相关系数图(n=90)
水稻冠层原始光谱反射率与APNA的相关性(图3)整体上呈现正相关关系,相关性的变化趋势与前述变化趋势有明显差异。水稻冠层原始光谱反射率与APNA的相关程度均较差,整体上未达到极显著水平,相关程度最大的波段为401 nm,对应的相关系数为0.310。
图3 水稻冠层原始光谱反射率与APNA的相关系数图(n=90)
将350~2 500 nm(剔除不稳定的波段)高光谱波段两两组合构成的4个高光谱植被指数RVI、PVI、DVI和NDVI与水稻氮素营养指标进行皮尔逊相关性分析,得到决定系数矩阵图(图4~图6),根据植被指数的计算公式和决定系数矩阵图的对称性,其中子图(a)RVI和(b)PVI为全矩阵图,而子图(c)DVI和(d)NDVI为三角矩阵图。图中颜色越深的区域表示植被指数与水稻氮素营养指标的相关性越高。
图4 4个植被指数和水稻叶片氮素含量的决定系数矩阵图(n=90)
图5 4个植被指数和水稻叶片氮素积累量的决定系数矩阵图(n=90)
图6 4个植被指数和水稻地上部氮素积累量的决定系数矩阵图(n=90)
图4为4个植被指数和水稻叶片氮素含量(LNC)的决定系数矩阵图。由图4可知,4个植被指数和水稻叶片氮素含量相关性较高(r2>0.7)的波段组合分布区域不同。对于RVI,相关性较高的波段组合集中分布在400~770 nm和750~1 350 nm范围内,RVI和LNC决定系数最大值为0.886,其对应的波段组合为694 nm和763 nm;对于PVI,相关性较高的波段组合集中分布在700~1 350 nm和400~700 nm范围内,PVI和LNC决定系数最大值为0.869,其对应的波段组合为864 nm和483 nm;对于DVI,相关性较高的波段组合较为分散,主要是可见光和近红外波段的部分波段组合,DVI和LNC决定系数最大值为0.883,其对应的波段组合1 292 nm和1 258 nm;对于NDVI,相关性较高的波段组合集中分布在750~1 350 nm和400~770 nm范围内,NDVI和LNC决定系数最大值为0.881,其对应的波段组合为1 296 nm和1 220 nm。
图5为4个植被指数和水稻叶片氮素积累量(LNA)的决定系数矩阵图。比较图4和图5可知,图5中4个植被指数和LNA的相关性比图4有所降低,但是4个植被指数和LNA的决定系数较高的波段组合分布区域与图4一致。RVI和LNA决定系数最大值为0.798,其对应的波段组合为692 nm和775 nm;PVI和LNA决定系数最大值为0.794,其对应的波段组合为864 nm和455 nm;DVI和LNA决定系数 最大值为0.819,其对应的波段组合为1 294 nm和1 258 nm;NDVI和LNA决定系数最大值为0.791,其对应的波段组合为780 nm和660 nm。
图6为4个植被指数和水稻地上部氮素积累量(APNA)的决定系数矩阵图。由图6可知,4个植被指数和水稻地上部氮素积累量的相关性均较差。RVI和APNA决定系数最大值为0.201,其对应的波段组合为1 196 nm和1 193 nm;PVI和APNA决定系数最大值为0.211,其对应的波段组合为700 nm和401 nm;DVI和APNA决定系数最大值为0.234,其对应的波段组合为1 261 nm和1 172 nm;NDVI和APNA决定系数最大值为0.200,其对应的波段组合为1 196 nm和1 193 nm。虽然4个植被指数与APNA的决定系数最大值均达到极显著水平(r2=0.073,n=90,p<0.01),但是对APNA的解释程度均较低,解释能力仅仅在20%左右。
以上述决定系数最大值对应的波段组合构建的植被指数为自变量,以氮素营养指标为因变量进行回归分析构建反演模型,其结果如表2~表4所示。
表2 叶片氮素含量的估算模型
2.3.1 水稻叶片氮素含量诊断模型
我们没有见过面,但彼此都有耳闻。高文鹏过来时,主动和我打招呼。他说看到我身上的工衣,就像见到了故人。我们就像故人一样坐下来,要了两个菜,喝啤酒。初次见面,好像也没什么不自在的。高文鹏给我的印象就是高。他说,一米八二。我开个玩笑,太高了,景花厂那座小庙,哪能容得下你?高文鹏笑笑,个子高顶什么用,我的水平也不高,不如你呀,话又说回来,我对阿花可是忠心耿耿啊,像我这样的人也很难找啊。可阿花居然找到了,这不,你来帮她了,我就放心了。我说,你说你忠心耿耿,为何又离她而去呢?高文鹏一声喟叹,说我离开也是情非得已啊。
由表2可知,4个植被指数的线性和非线性模型均可以很好地估算水稻叶片氮素含量,建模集R2在0.458~0.904之间,测试集R2在0.521~0.895之间,模型预测值和实测值之间的R2均达到极显著水平。根据测试集R2值最大、RMSE和RE值最小的原则确定最优的估算模型。对于RVI和PVI,最优的反演模型均为抛物线函数,模型表达式分别为y=-1.67x2-4.14x+3.26和y=-87.71x2-64.78x-4.46;对于DVI和NDVI,最优的反演模型均为线性函数,模型表达式分别为y=-193.88x+1.33和y=-89.35x+3.48。估算水稻叶片氮素含量最优的模型为NDVI的线性模型,其建模集R2、RMSE和RE分别为0.873、0.366%和16.664%,其测试集R2、RMSE和RE分别为0.895、0.331%和15.110%。
2.3.2 水稻叶片氮素积累量诊断模型
表3为叶片氮素积累量的估算模型,与叶片氮素含量的估算模型结果一致,模型预测值和实测值之间的R2均达到极显著水平,除RVI的幂函数模型外,4个植被指数的线性和非线性模型均可以很好估算水稻叶片氮素积累量,建模集R2在0.656~0.798之间,测试集R2在0.361~0.880之间。最优估算模型的确定方法同上,对于RVI、DVI和NDVI,最优的估算模型为抛物线函数,模型表达式分别为y=-146.59x2-105.61x+92.93、y=3E+04x2-6E+03x+24.39和y=-233.53x2+474.56x-153.01;对于PVI,最优的估算模型为线性函数,模型表达式为y=2E+03x+1E+03。估算水稻叶片氮素积累量最优的模型为DVI的抛物线函数,其建模集R2、RMSE和RE分 别 为0.790、16.887 kg/hm2和27.736%,其测试集R2、RMSE和RE分别为0.880、12.810 kg/hm2和21.002%。
表3 叶片氮素积累量的估算模型
2.3.3 地上部氮素积累量诊断模型
从表4可见,地上部氮素积累量诊断模型的建模集和测试集R2整体上均达到显著水平,说明所选植被指数对地上部氮素积累量的估测有一定的效果,但是R2的值均较低,建模集R2均在0.3以下,测试集R2的最大值在0.5左右。
表4 地上部氮素积累量的估算模型
2.3.4 几种诊断模型的比较
为进一步比较上述氮素诊断模型的预测能力,对用于评价模型精度的指标(R2、RMSE和RE)进行单因素完全随机设计试验资料的方差分析。使用SSR法进行多重比较。
首先对氮素营养指标进行单因素完全随机设计试验资料的方差分析,因3种氮素营养指标对应的RMSE单位不同,所以选择R2和RE进行比较,结果如表5所示。由表5可知,3种氮素营养指标建模集和测试集R2整体上表现为LNC>LNA>APNA,且建模集R2两两之间均达到显著差异,测试集R2在LNC和LNA之间未达到显著差异,但LNC和LNA极显著高于APNA。建模集RE整体上表现为LNC<APNA<LNA,LNC和APNA之间差异不显著,但两者均极显著低于LNA;测试集RE整体上表现为LNC和APNA极显著低于LNA,其中LNC和APNA差异不显著。综合比较,3种氮素营养指标诊断模型的预测能力LNC>LNA>APNA,选用LNC诊断模型有较好的预测效果。
表5 氮素营养指标对氮素营养诊断模型的影响
其次对LNC的4个植被指数进行单因素完全随机设计试验资料的方差分析。结果(表6)表明,建模集和测试集的R2、RMSE和RE在4个植被指数之间差异不显著。建模集R2整体上表现为DVI>PVI>NDVI>RVI,而建模集RMSE和RE趋势正相反,表现为DVI<PVI<NDVI<RVI;测试集R2整体上表现为DVI>NDVI>PVI>RVI,而测试集RMSE和RE趋势正相反,表现为DVI<NDVI<PVI<RVI。整体上DVI表现较好,但由于DVI无法拟合对数函数和幂函数,模型不健全,所以综合分析预测能力最好的植被指数为NDVI。
表6 植被指数对氮素营养诊断模型的影响
最后,将NDVI 5种回归方程的建模集和测试集评价指标整体进行单因素完全随机设计试验资料的方差分析。如表7所示,建模集和测试集R2整体上表现为抛物线>线性>指数>对数>幂,而RMSE和RE表现为抛物线<线性<指数<对数<幂。综合比较,预测能力最好的回归模型为抛物线函数。
表7 回归方程对氮素营养诊断模型的影响
氮素是水稻植株内多种有机物质的组成成分,氮素的丰缺直接影响与光合作用有关色素的生物合成,从而影响水稻的生长发育和产量的形成[15],所以,利用遥感技术实时监测水稻氮素状况具有重要意义[16]。
前人研究指出,水稻冠层氮素敏感波段为520~550 nm、580~690 nm和740~1 070 nm[17]。本研究结果表明,与LNC和LNA负相关性最大的波段为674 nm,与LNC和LNA正相关性最大的波段分别为779 nm和780 nm,与前人的研究结果相似。
前人研究表明,蓝紫光波段和红光波段是叶绿素主要的吸收波段[18],叶绿素含量随着氮素含量增加而增加[19],说明蓝紫光和红光波段反射率的变化主要是通过氮素影响叶绿素含量表现出来的;近红外波段的光谱特征主要受叶片内部构造控制[20],氮素能够改变叶片结构[21],从而造成近红外波段光谱反射率的差异。短波红外波段在大气窗口中透过率较高(超过90%),植被的反射信号较强,一定程度上包含了氮素的光谱反射特征[22]。本研究发现,在所选择的4个植被指数中,最佳波段组合包含蓝紫光波段、红光波段、近红外波段以及短波红外波段,与绿色植物的生理结构及光谱波段的特性相吻合。
关于利用植被指数反演氮素营养状况的最佳波段组合,TIAN等[23]研究表明,绿光波段的553 nm和537 nm波段组合构建的比值植被指数能够很好地反演不同栽培条件下水稻冠层叶片氮素含量,而本研究结果表明,比值植被指数最优的波段组合为红光波段(694 nm)和近红外波段(763 nm)的波段组合,与TIAN等[23]的研究结果不同,这可能是不同水稻品种的光谱特征不同以及氮素敏感波段存在差异导致的[24]。CHU等[25]研究表明,770 nm和752 nm波段组合构建的比值植被指数能够很好地反演水稻冠层叶片氮素积累量,与本研究筛选的比值植被指数最优波段组合692 nm和775 nm较为一致,波段均分布在红光和近红外波段范围内。
从整体上看,水稻冠层氮素营养状况估算模型效果以叶片氮素含量估算模型最好,叶片氮素积累量模型次之,地上部氮素积累量模型最差,这与RANJIAN等[6]的研究结果相似。究其原因,可能是由于水稻地上部器官与空气的接触面积大小导致的,与冠层顶端空气接触面积最大的器官是叶片,所以光谱仪传感器所接收的反射光主要是由叶片产生,从而通过高光谱估测叶片氮素营养状况要比估测植株的氮素营养状况好;叶片氮素积累量综合了叶片氮素含量、比叶重和叶面积指数等因素,受品种、冠层结构等多重因素的影响,相对于叶片氮素含量其估测难度较大[26],因而叶片氮素含量估算模型要优于叶片氮素积累量估算模型。
本研究得出最优的氮素诊断模型为叶片氮素含量诊断模型,其模型表达式为LNC=1E+03NDVI2-132.55NDVI+3.72,建 模 集R2、RMSE和RE分 别 为0.879、0.357%和16.267%,测试集R2、RMSE和RE分别为0.895、0.331%和15.136%。