无人机光谱分析在水稻产量预测中的应用

隋丽娜,房 建,郭立峰

(河北民族师范学院 数学与计算机科学学院,河北 承德 067000)

0 引言

水稻在我国广泛种植,是我国主要的粮食作物,水稻的稳产高产直接关系到国家粮食安全[1]。由于种植条件的特殊性,在东北等主要水稻产地,水稻主要为大规模的集中式种植[2]。该种植模式有助于管理,但当某个生长时期出现病灾和虫灾时,会对水稻产量产生严重影响[3]。因此,急需建立一套水稻全生长周期的产量预测模型,当一个生长时期出现问题时,通过其他生长周期进行调整,保证水稻稳产。目前,水稻产量模型主要包括单株图像采集法[4]、飞行器图像分析[5]和卫星光谱分析[6]。单株图像分析对一株水稻不同生长时期进行图像处理,分析其结节数、植株高度和植株宽度等参数,然后对其产量进行建模分析。该方法建模精度高、针对性强,但受限于样本数量及样品选择,对于大范围的水稻产量预测容易出现偏差[7]。飞行器图像分析方法对水稻冠层图像进行分析,背景噪声处理难度较大,且受到稻田环境因素影响,算法精度有待提高[8]。卫星光谱分析采样范围过大,光谱信息庞杂,建模难度大,且获得光谱数据成本高[9]。无人机技术的发展,使获得低成本光谱数据成为可能。以无人机搭载光谱仪进行40m高度均速光谱采集,通过分析水稻不同生长时期的光谱,得到水稻亩产预测模型。同时,进行多组“留一法”实验,结果表明:该方法具有良好的预测精度,性能稳定,预测精度高。该模型数据获得方便,预测精度较高,适于大规模一线推广。

1 系统组成

系统包括无人机光谱系统、水稻产量模型系统、模型建立依据和产量模型建立4部分,如图1所示。

图1 系统结构图Fig.1 Structure of system。

无人机光谱系统保持无人40m高度匀速飞行,采用机载反射光谱传感器检测水稻冠层叶绿素光谱。模型建立依据,叶绿素可以表征植物生长状态,冠层叶绿素含量可以预测水稻产量。通过无人机机载光谱仪检测水稻冠层叶绿素光谱,由于叶绿素在红边区域、近红外区域均有特征吸收波长,因此采用差值植被指数综合考虑光谱包含信息。计算水稻全生命周期中的差值植被指数及其对应的亩产量,对其进行共线性分析,确定模型建立因子,并采用逐步线性拟合的方式建立模型。同时,采用“留一法”进行亩产预测,对模型进行线性拟合评价和统计指标评价,测试模型可靠性。

2 建模原理

太阳与地球间距离遥远,照射到地球上的太阳光是理想的白色平行光,当白色平行光照射到水稻上时,叶绿素会吸收特定波长的光,然后将白光反射到定高飞行的无人机光谱传观器上。由于特定波长的光被吸收,在光谱传感器输出光谱曲线上形成波谷,通过分析光谱即可确定叶绿素含量。叶绿素吸收光能,将水和CO2转化为有机物和氧气,实现光合作用,可以作为衡量水稻长势及产量的标志。因此,通过无人机光谱遥感计算叶绿素含量可以预测水稻最终产量。

2.1 叶绿素含量与水稻产量

水稻不同生长周期叶绿素SPAD值与亩产关系如图2所示。

图2 水稻不同生长周期叶绿素SPAD值与亩产关系Fig.2 Relationship between chlorophyll SPAD and rice yield in different growth period。

叶绿素吸收光能并离子化,能量被储存在三磷酸腺苷中,最终将水和CO2转化为有机物和氧气,实现光合作用。由于无人机反射光谱检测重点针对冠层叶绿素含量,现采用SPAD-502型叶绿素检测仪对不同生长时期的水稻冠层叶片进行检测,每叶片3点取样,计算水稻冠层叶片SPAD值,即

(1)

其中,i为单株水稻冠层叶片数量,测试结果如图2(a)所示。

叶绿素SPAD值从拔节期开始显著增长,拔节期到孕穗期增长速率最大,到抽穗期达到最大,随后含量开始下降。这是由于在抽穗期之前光合作用主要集中在植物生长上;从抽穗期开始,转向生殖生长,叶片中的氮开始向稻穗集中转移,叶绿素含量逐渐降低,此时水稻处于生殖生长阶段。设冠层叶片叶绿素含量A为

(2)

其中,LAI为冠层叶面指数。不同生长时期冠层叶片叶绿素含量A与亩产拟合直线,相关系数如图2(b)所示。其在孕穗期达到最高,乳熟期次之。这是由于该时期叶片发育茂盛,叶绿素含量能较好反映水稻生长状况,抽穗期与乳熟期主要进行稻穗生殖生长,叶片变黄,叶绿素含量降低。总体上叶绿素含量与亩产相关系数均达到0.72以上,拟合精度较高,冠层叶绿素含量可以作为水稻产量的评价对象。

2.2 冠层叶绿素无人机光谱

水稻处于不同生长时期的无人机光谱如图3所示。其中,总体上拔节期、孕穗期、抽穗期和乳熟期光谱曲线走势一致,叶绿素对于绿光波段吸收能力很弱,故在530～600nm形成波峰,对于蓝光和红光吸收能力较强,在400～500nm形成较小波谷;在620～680nm形成较深波谷,同时在700～800nm红光波段与近红外线波段交接区域,反射率急剧增长,形成红边区域;在900～960nm近红外波段形成波谷。由图3可知:叶绿素具有3个特征吸收波段,其中红边区域与近红外波段吸收效果显著,因此采用采用红边区域与近红外波段组合参数模式,可提高对于叶绿素含量检测的稳定性与敏感性。

图3 不同生长时期冠层叶绿素光谱Fig.3 The rice canopy chlorophyll spectral in different growth period。

由于水稻产量与水稻4个生长时期的冠层叶绿素含量拟合,相关系数达到0.7以上,通过冠层叶绿素含量可以有效预测水稻产量;同时,冠层叶绿素含量的无人机光谱数据具有良好的特征吸收波长,通过红边区域与近红外波段组合参数模式可以有效预测其含量,因此通过分析无人机光谱可以有效的预测水稻产量。

3 建模建立

叶绿素具有3个特征吸收波段,在红边区域和近红外波段吸收效果显著,采用单一波段不能概括光谱曲线全部信息,因此采用多波段组合进行组合计算,得到水稻植被组合因子。在水稻生长全生命周期中,拔节期到孕穗期过程中,水稻冠层叶片处于高速生长状态;当处于抽穗期和乳熟期时,冠层叶片由绿色向黄色发展,要充分考量冠层叶绿素含量的变化,因此选用对于背景变化敏感的差值植被指数[10]作为模型建立因子,计算公式如式(3)所示。其中,DNNIR为近红外光谱反射率;DNR为红边区域光谱反射率。

DVI(NIR,R)=DNNIR-DNR

(3)

由于水稻生长分为拔节期、孕穗期、抽穗期和乳熟期,为了综合不同时期对最终亩产的影响,将4个时期的差值植被指数作为模型建立的因子。

3.1 不同生长时期差值植被指数与亩产关系

水稻不同生长时期差值植被指数与亩产之间的相关系数R灰度云图,如图4所示。

图4 不同时期水稻DVI(NIR,R)与产量相关系数云图Fig. 4 Cloud chart of linear correlation coefficient between DVI(NIR,R)and rice yield in different growth period。

图4中,红边波长为横坐标,近红外特征波长为纵坐标,每个坐标相关系数R数值用该点的灰度表示,即坐标点越亮,该点对于的相关系数R越大。具体数值参照右侧灰度强度条。

图4(a)为水稻拔节期差值植被指数与亩产相关系数云图,该时期拟合相关系数R普遍在0.4以下,最大值为0.631 6。水稻处于拔穗期时,冠层叶片发育不成熟,叶片发育程度较低,土壤背景杂波处于主要因素。图4(b)为孕穗期相关系数R云图,总体上较拔节期有所提高,大部分处于0.5以上,最大值为0.771 3。这是由于冠层叶片在孕穗期显著生长,叶绿素含量显著增加,红边区域吸收能力增强。图4(c)为抽穗期相关系数R云图,大部分处于0.6以下,最大值为0.634 1,出现在1 080nm和1 089nm处。由于两特征吸收波段相近,造成红边向量与近红外向量线性相关较大,因此拟合关系系数R较低。图4(d)为乳熟期相关系数R云图,最大值出现在900nm及1 100nm处,该波段为淀粉特征吸收波段,红光和蓝光波段相关系数强度下降,叶绿素含量降低。在乳熟期水稻生长以生殖生长为主,稻穗中淀粉开始生成。不同生长时期差值植被指数与亩产相关系数最大值及出现波长如表1所示。

表1 不同生长时期DVI(NIR,R)特征波长及与亩产相关系数Table 1 Characteristic absorption wavelength and correlation coefficient between DVI(NIR,R)and rice yield in different growth period。

3.2 共线性判定

水稻在生长过程中分为4个时期,现对4个时期的差值植被指数与最终亩产量进行线性拟合,得到水稻产量模型,如式(4)所示。其中,X1～X4分别为拔节期、孕穗期、抽穗期和乳熟期差值植被指数DVI(NIR,R)向量。

Y=aX1+bX2+cX3+dX4+m

(4)

在进行多元线性拟合前必须确保X1-X4线性无关,即4向量不共线。判定依据为容忍度T和主成分分析特征值。其中,容忍度T计算公式为

(4)

其中,Rij为Xj与Xi的线性相关系数。容忍度越小,该向量与其他向量的线性先关性越高。当T<0.1时,认为线性相关性非常严重;主成分分析特征值E越小,线性相关性越高。现计算拔节期、孕穗期、抽穗期和乳熟期差值植被指数DVI(NIR,R)向量的容忍度和主成分特征值,如图5所示。其中,拔节期容忍度为0.097,主成分特征值为0.14;抽穗期容忍度为0.051,主成分特征值为0.09。差值植被指数DVI(NIR,R)向量在拔节期与抽穗期高度共线,应给予舍弃。因此,以孕穗期和乳熟期差值植被指数DVI(NIR,R)向量作为模型建立因子。

图5 共线性判定Fig.5 The collinear test。

3.3 逐步回归分析

由共线性分析可知:拔节期与抽穗期容忍度和主成分特征值很小,线性相关性高,在模型中剔除。因此,水稻产量模型简化为

Y=aX1+bX2+m

(5)

现采用逐步线性拟合的方法进行求解,主要步骤如下:

1)将孕穗期差值指数、乳熟期差值指数与亩产量进行线性拟合,对回归系数βi预测值与实际值进行F检验。

y=Xiβi+εi(i=1,2)

(6)

孕穗期与乳熟期对应的检验结果为F1、F2,取F(1)=max(F1,F2),查表的置信区间为1-a情况下,临界值为F1,当F(1)>F1时,引入第1个变量。

2)引入第2个量,建立亩产量y与X1和X2的二元线性回归方程,对回归系数预测值与实际值进行F检验,计算结果为F(1)。查表的置信区间为1-a情况下,临界值为F2,当F(2)>F2引入第2个变量,完成产量预测模型。

现取8组样本光谱数据,采用“留一法”进行交叉预测,以其中7组样品光谱作为建模样本,另一组作为预测样本,进行8次模型预测,结果表2所示。表2中,8组预测值和实际值平均相对误差值为7.102%,表明模型具有较高的预测精度。

表2 水稻产量模型系数及其预测值Table 2 Coefficient and predictive value of production prediction of rice。

4 模型评价

4.1 线性拟合评价

多元线性拟合的模型评价标准为线性决定系数R2和均方根差RMSE。其中,线性决定系数R2表示模型预测值与实际值的符合程度,计算公式为

(7)

均方根差RMSE为预测值与真实值的误差平方根的均值,计算公式为

(8)

现计算采用采用“留一法”进行交叉预测,计算每一组模型R2和RMSE值,如图6所示。

图6 模型线性性能评价Fig. 6 Linear fitting evaluation。

其中,R2值分布在0.607～0.648之间,均值为0.630;RMSE值分布在21.612～25.145之间,均值为24.300。这说明,模型拟合线性性能良好。

4.2 统计指标评价

统计学中常将总体相对误差Rs和预估精度P作为模型预测精度的评价标准,两标准的计算公式为

(9)

采用式(9)计算得到8组样品的Rs为0.0142,预估精度P为93.71。测试结果表明,该模型在统计学上具有良好的预测精度。

5 结论

为了实现大范围的水稻产量预测,基于无人机水稻反射光谱分析建立了水稻产量模型。由于冠层叶绿素浓度可以有效表征植物生长情况,预估水稻产量。通过分析无人机测量的水稻冠层叶绿素光谱发现:在620～680nm形成较深波谷;在700～800nm红光波段与近红外线波段交接区域,反射率急剧增长,形成红边区域;在900～960nm近红外波段形成波谷,包含丰富信息。因此,采取多波段组合进行综合计算,以完整地提取光谱信息;选用包含红边区域与近红外区域的水稻差值植被指数,作为建模因子。同时,分析差值植被指数在拔节期、孕穗期、抽穗期和乳熟期的全生长周期中的特征吸收波长,将其和亩产进行线性拟合,得到线性相关系数。对其在不同生长时期进行共线性判定,发现差值植被指数在拔节期与抽穗期的容忍度和主成分特征值均很低,存在严重的线性相关,因此给予舍弃,并对孕穗期与乳熟期差值植被指数与水稻亩产进行逐次线性拟合。对8组样本进行“留一法”预测表明:模型的线性决定系数R2值分布在0.607～0.648之间,均值为0.630;RMSE值分布在21.612～25.145之间,均值为24.300。由此表明,模型精度较高。计算统计量总体相对误差Rs=0.0142和预估精度P=93.71,表明该模型在统计学上具有较高精度。本方法数据采集方便,光谱处理简单,适用于大规模推广。