基于高光谱数据与针叶电阻率的松材线虫病早期诊断研究

刘 飞,谈家金*,谢谷艾,苏 峻,叶建仁

1. 南京林业大学林学院,南方现代林业协同创新中心,江苏 南京 210037 2. 江西省林业科学院林业有害生物防治研究所,江西 南昌 330013 3. 南京林业大学理学院,江苏 南京 210037

引 言

松材线虫(Bursaphelenchusxylophilus)是目前对我国森林危害最严重的林业检疫性有害生物,由其引起的松材线虫病是一种松树毁灭性病害[1-3]。我国自1982年发现该病以来,病情经过41年不断发展,目前全国共有19个省(自治区、直辖市)701个县级行政区发生疫情,松材线虫病的准确诊断对该病的有效防控至关重要。该病现行的诊断技术包括林间症状诊断法、病原线虫鉴定法、流胶法等,这些方法受到条件或技术的限制而不完善;为了能够有效发现和控制松材线虫病,其爆发前的早期诊断更是成为当前研究的热点[4]。如高瑞贺等发现健康马尾松感染松材线虫病后,其针叶的蒸腾速率和净光合速率变化显著,可作为松材线虫病早期诊断的指标[5];Gaspar等在接种松材线虫后,发现松树3-蒈烯挥发量相比健康松树显著提高,可用于松材线虫病早期诊断的指标[6]。

大量研究表明,植物在受生物胁迫或非生物胁迫过程中,其自身色素、含水量等生理生化指标会发生一定的变化[7-8]。而植物对电磁辐射的吸收与反射特性会随着这些生理生化参数的改变而改变,在遥感上表现为不同的光谱反射率。因而遥感在研究植物胁迫方面应用广泛,可见光区和近红外区是指示植物色素含量变化和健康情况的关键波段[9-11]。健康马尾松在感染松材线虫病后,其叶片内部结构会被破坏,光合色素和水分含量也会发生变化[12];伴随病害不断的发展,马尾松的光谱反射率呈现出较大变化[13]。松材线虫病在病情发展的早期阶段,感病松树自身不会发生明显症状,且叶绿素和水分含量变化不明显,但当肉眼可以观察到针叶叶色明显变色时,感病松树会在短期内萎蔫死亡。目前研究中,无论是基于光谱本身的直接诊断方法,或是基于生理生化组分的间接诊断方法在病害的早期阶段效果并不理想[14-15]。

为了实现对感病植株的早期诊断,可以考虑寻找一种基于光谱的新方法用于判断植株是否感病,这将对松材线虫病的防治具有重要的实践意义与应用价值。马尾松感染松材线虫病后除生理生化指标发生变化外,其生物电也必然会发生一定的变化[16]。本研究拟从电学角度出发,在研究感病松树针叶电阻率变化的基础上,将基于马尾松针叶反射光谱数据建立针叶电阻率的预测模型,以实现在感病后的早期阶段对马尾松松材线虫病的快速判别,为利用无人机、卫星遥感对松林进行松材线虫病监测提供理论依据。

1 实验部分

1.1 松材线虫接种、病情观察及样品采集

实验地设在安徽省休宁县商山镇五方水库的一处马尾松林内(118°21′63″E,29°69′23″N;156.7 m)。2021年8月17号采用注干法对8～9年生马尾松接种松材线虫,实验共设置处理和对照两组,处理组和对照组的马尾松分别设置6棵和4棵。在松树主干距离地面1.5 m处,使用直径6 mm的钻头(钻头与主干夹角约45°)钻入木质部4 cm,在钻孔内放入适量的无菌棉球;使用1 mL的枪头(枪头前端已被修剪为注射器针头状,增加接触面积)插入钻孔,用微量进液器注入线虫悬浮液,处理组接种线虫量为30 000条·株-1(浓度30 000条·mL-1),对照组接等量无菌水。接种后每4～7天对马尾松发病状况进行观察并记录感病指数。马尾松病情分级与感病指数按谈家金等的方法[17],略加修改后计算:0级正常,针叶绿色,代表值为0;1级为1/2以下针叶褪绿、1/4以下叶发黄,代表值为1;2级为1/2以上叶褪绿、1/4～3/4叶发黄,代表值为2;3级为3/4以上叶发黄、1/2以下叶变红,代表值为3;4级为1/2以上叶变红、植株濒死或死亡,代表值为4。计算公式:感病指数=Σ(各病级株数×病级代表值)×100/(总株数×最高病级的代表值)。接种当天以及接种后7、11、15、19和26 d,采集外观表现健康的针叶进行光谱和电阻的测定。样品采集时间设置在上午10:00—12:00,选择晴朗无风无云或少云天气,从马尾松向阳面3个位置(即冠层上、中、下3个位置各取一个点)采集当年生成熟针叶,将采集好的样本针叶用密封袋封口,放置在冰盒内带回室内进行光谱和电阻的测定。

1.2 针叶光谱反射率测定

采用Ocean Optics USB2000+光谱仪对针叶光谱反射率进行测定,该光谱仪由美国Ocean Insight公司生产,其作用波段范围为350～1 000 nm,光谱分辨率为1.5 nm,光谱采样间隔为0.44 nm。光谱仪探头视场角25.4°,探头距离样本表面距离8 cm,采用10 W卤素灯为光源,方位角为45°(光源与样本表面的夹角)。每次测定前,使用白色面板进行光谱校正。测定仪器各配件如图1所示,开关电源(a)、探头(b)、样本槽(c)和光源(d)均放置在内部密闭黑色的箱子内,如图1(A)所示。测定过程中设定光谱平均次数为10次,积分时间100 ms。对针叶样本进行测定时,使用黑色双面胶将针叶中间部分整齐平铺在黑色样本板上,并剪去多余部分,针叶之间不能有缝隙,如图1(B)所示。将冠层上、中、下3个位置光谱反射率平均值作为该植株的光谱反射率。

图1 针叶光谱反射率测定仪器各配件(A)黑色箱子内部各组件:(a)开关电源;(b)探头;(c)样本槽;(d)光源;(B)样品板与待测针叶样本Fig.1 Accessories of needle spectral reflectance measurement instrument(A) Internal components of the black box:(a) Switching power supply;(b) The probe;(c) The sample slot;(d) The light source;(B) Sample plates and sample of needles to be tested

1.3 针叶电阻率测定

采用M4070 LCR测试仪测定针叶电阻,该仪器是由中国晶研仪器科技有限公司生产,其量程范围为0.00 Ω～10.000 MΩ,基本测量准确度0.5%。使用3D打印机打印黑色薄板,将两个电极插入对应位置,在半球形凹槽内放入导电膏(济宁高新区金诺特医用凝胶厂生产),方便针叶电阻的测量,如图2所示。将针叶横截面近似看成半椭圆形,剪取针叶中部4 cm,使用游标卡尺(分辨率0.01 mm)对针叶的宽度和厚度进行测量。将针叶与M4070 LCR测试仪连接后3～5 s的读数作为其电阻值,从而计算出其电阻率。冠层向阳面上、中、下3个位置各取5根针叶进行电阻的测定,将5根针叶电阻率平均值(电阻率异常值的样本除外)作为对应冠层位置针叶电阻率,将上、中、下3个位置电阻率平均值作为该植株的针叶电阻率。

图2 针叶电阻测定过程Fig.2 Needle resistance measurement process

1.4 数据分析

1.4.1 数据预处理

地面光谱与针叶电阻同步采样共6次,将获取的数据去除明显异常值后,共获得52组数据。其中,处理组含34组数据,对照组含18组数据。使用Excel数据分析工具对光谱数据和电阻数据处理汇总;使用SPSS软件对电阻率数据进行分析,并在Origin 2018中绘图。

在MATLAB R2021b中调用Savitzky-Golay函数对光谱数据进行平滑处理,既可以消除数据中的噪声部分,又能将数据中重要的信息保留,提高了数据的平滑性。将平滑后的光谱定义成原始光谱(OR)。光谱数据自身具有复杂性和高维数特点,为了提取更多有用的信息,需要对光谱数据进行一些复杂的数学变换,为预测模型的创建提供充足的信息源。目前通常的光谱变换包括:一阶导数(FD)、二阶导数(SD)、对数变换(LOG)、倒数变换(1/R)和连续统去除法(CR)。

1.4.2 特征波段选择

高光谱数据具有高维、复杂的特点,为了筛选出与针叶电阻率相关性高的波段,采用随机森林进行特征波段筛选。随机森林(random forest,RF)是一种非线性模型,以决策树作为基学习器,使用装袋(Bagging)原理来处理训练数据集[18]。由于其计算过程相对简单,对样本和变量的选择过程是随机的,因而能够在多种领域应用并且均可以得到很好的分类效果[19]。这些特点使得RF特别适用于光谱数据的研究。本次实验使用RF进行特征波段的筛选,RF中树的数量nTrees设置300,输入的特征数量10。

1.4.3 最小二乘支持向量机回归模型

最小二乘支持向量机(least squares support vector machine,LSSVM)算法是基于支持向量机算法的一种改进,极大提高了求解问题的速度和收敛精度[20],且不再受到样本数量的影响,具有很强的泛化能力。

采用LSSVM算法对针叶光谱反射率与针叶电阻率进行回归分析,从而对马尾松针叶电阻率进行反演。以光谱变换下的特征波段为自变量,针叶电阻率为因变量,基于LSSVM算法构建马尾松针叶电阻率预测模型。为了使得模型测试的结果更客观,80%的样本选作建模集,20%的样本用作验证集。主要用以下三个指标评估马尾松针叶电阻率预测模型的精度:决定系数(coefficient of determination,R2)、均方根误差(root mean square error,RMSE)、平均绝对误差(mean absolute error,MAE)。当预测模型的R2趋近1、RMSE和MAE越小时,表明其精度越高。

2 结果与讨论

2.1 接种后马尾松病情发展

接种后11 d,发现有1棵接种线虫的马尾松1个2级侧枝针叶开始褪绿,此时,病情指数为4;接种后15 d,该侧枝针叶变为黄褐色,此时,病情指数为4;接种后19 d,该侧枝针叶变成红褐色,针叶已枯死,此时,已有3棵接种线虫的马尾松针叶变色,病情指数为12.5。接种后26 d,已有4棵接种线虫的马尾松极少数针叶变色,病情指数为16.7。接种后67 d,接种线虫的6棵马尾松少数针叶均已变色,病情指数为25。

2.2 接种后马尾松针叶电阻率的变化

接种实验结果显示,马尾松显症早期,针叶电阻率显著升高。由图3可知,接种后15 d,此时仅1棵接虫马尾松的1个侧枝针叶变色,接虫株和对照株的电阻率开始达到极显著差异,这种差异一直维持到接种后26 d(p<0.01)。整个实验期间,相比接虫株的针叶电阻率,对照株的针叶电阻率波动较小。因此,从电阻率的大小可以初步判断当马尾松针叶电阻率大于323 Ω·m时,该松树已经感病。

图3 接种松材线虫后马尾松针叶电阻率的变化**表示在0.01水平差异显著Fig.3 Changes in needle resistivity of Pinus massoniana inoculated with B. xylophilus** Shows significant correlation difference at 0.01 level

2.3 接种后马尾松针叶光谱数据特征分析

由于部分波段区域出现大量噪声,主要对490～980 nm范围内的马尾松针叶反射光谱进行分析。为了丰富基于光谱建立针叶电阻率模型的信息,使用5种数学变换方式对所有有效光谱数据进行光谱预处理,仅对数变换没有改变光谱曲线的趋势[图4(a—f)]。

图4 不同预处理方法的光谱变化(a):OR;(b):FD;(c):SD;(d):LOG;(e):1/R;(f):CRFig.4 Spectral changes of different pretreatment methods(a):OR;(b):FD;(c):SD;(d):LOG;(e):1/R;(f):CR

由表1可知,通过随机森林算法对各预处理方法进行特征波段筛选后发现存在较大不同。原始光谱(OR)的特征波段集中在近红外区域,一阶导数和连续统去除法光谱的特征波段集中分布在绿光和近红外区域,对数变换和倒数变换筛选出的波段主要分布在红光和近红外区域。综合比较发现,特征波段主要分布在近红外区域,说明病害发生后的早期针叶内细胞结构可能发生了变化,这与高瑞贺等研究的感染松材线虫病引起马尾松针叶光合作用相关机构改变的结果相一致[5]。大量研究表明,在可见光和近红外区域的反射光谱在病虫害监测中具有很高的应用价值,它们通常与病害的发展过程密切相关[9-10]。

表1 基于随机森林算法筛选的特征波段Table 1 Feature bands screening based on random forest algorithm

2.4 马尾松针叶电阻率预测模型分析

基于原始光谱进行五种光谱变换,使用随机森林对原始光谱和各预处理后的光谱数据进行筛选,从而得到与针叶电阻率相关特征波段,最终利用LSSVM算法建立的针叶电阻率预测模型。该模型的预测效果汇总于表2,发现不同预处理间预测效果的综合表现具有差异,其中二阶导数模型预测效果综合表现最好,建模集与验证集的平均R2、RMSE和MAE分别为0.848、32.332和7.067。

表2 基于不同光谱变换筛选输入变量的模型表现Table 2 Performance of model based on the different spectral transformation input variable selection

各预测模型间除了倒数变换和连续统去除法变换外,其他光谱变换后的预测模型比原始光谱预测模型验证集的R2均有所提升,RMS和MAE均有所降低。一阶导数预测模型的改善效果较小,建模集与验证集的平均R2仅提升1.5%,但RMSE和MAE却分别升高28%和1.7%。而二阶导数模型预测效果提升较为明显,R2提升了4.2%,RMSE和MAE分别降低2.5%和18.9%。对数模型预测效果也有所提升,R2提升了3.1%,但RMSE和MAE分别升高了57.3%和11.3%。综合比较,基于二阶导数选择的特征变量在最小二乘支持向量机算法预测模型中综合表现最优,特征波段为594.986、646.107、646.451、782.896、784.841、839.164、863.890、902.021、947.901、962.315;其次为LOG、FD、1/R和CR,这表明对光谱预处理的必要性,处理后的光谱数据增强了样本的光谱反应特征,进一步提升反演模型的精度。

图5为基于二阶导数(SD)变换的LSSVM建模方法反演结果散点图,从图5中可以清晰看出,该模型的实测值样点和预测值样点多数分布在y=x线附近,能够较好地对马尾松针叶电阻率进行估测。因此,认为SD-RF-LSSVM可作为马尾松针叶电阻率大小的最优估测模型。相较于传统的林间症状诊断法、病原线虫鉴定法以及流胶法[4],本文提出的基于SD-RF-LSSVM模型的光谱检测方法不仅可以实现对该病害的快速无损诊断;此外,该方法能够在病树极少数针叶变色时就完成对病害的监测,相比之下,传统的遥感监测方法需要等到病树针叶大部分已经明显变色才能发现病害[16]。因此,该方法可以在病树出现极少数针叶变色时进行早期预警和监测,及时发现病害的发生和变化情况,从而更好地制定病害防治措施。

图5 基于二阶导数数据变换的最小二乘支持向量机模型检验结果Fig.5 Test results of LSSVM model based on the SD data transformation

3 结 论

多年来,我国在松材线虫病的致病机理及防治方法上取得许多重要进展,然而在该病的早期诊断方面一直未有突破,而利用遥感技术对该病的早期监测和预警长期处在探索阶段。本研究表明,使用针叶反射光谱反演针叶电阻率是可行的,且SD-RF-LSSVM建立的预测模型精度最高,可用于针叶电阻率的快速估测,为实现基于遥感的松材线虫病早期诊断与监测提供了思路与方法。

实验仅对针叶叶片层面的光谱特征进行分析,虽可为将来冠层尺度研究提供参考,但冠层尺度的松树冠层监测更复杂,光谱识别与分析受到冠层形状、土壤、植被等因素影响。需要在冠层尺度比较健康和感病马尾松的光谱特征差异,从而为使用遥感技术实现病害的快速有效识别提供方法。此外样本数量较少,模型需要进一步检验和验证。同时,环境条件对结果的影响需要进一步研究,增强模型的稳定性和适用性是今后研究的重点。