苹果样本于2020年10月13日购买于当地丹尼斯超市, 产地为山东烟台, 品种为烟台红富士, 共350 kg。 由于研究的目的是构建苹果腐败预警模型, 并考虑到贮藏条件不会对研究结果产生影响, 只会加速或减缓样本的腐败进程^[1]。 因此, 为加快试验进程, 样本贮藏室温度设置为常温。 试验过程中, 每日固定时间采集样本高光谱图像信息以及色差值和失水率数据。 贮藏室温度在18.8~26.2 ℃之间, 相对湿度在53.0%~64.9%RH之间。

每日随机从贮藏室不同位置取出18~19个苹果, 每个苹果在不同的部位切下2~3个测试样本, 共制备50个样本。 切割苹果时, 将苹果果柄向上, 果蒂向下放置, 并沿果柄轴线方向自上而下切下厚度为15 mm且包含果皮的切片作为测试样本。 然后将每一个测试样本的切割面向下, 分别放在培养皿中, 最终得到50张高光谱图片。

图1是样本的高光谱信息获取系统, 该系统主要由4部分构成: (1)高光谱相机(光谱仪)为德国的Inno-Spec公司制造, 型号是IST50-3810型, 该相机可在371.05~1 023.82 nm范围内获取光谱信息, 每两个波长之间间隔0.51 nm; (2)高光谱数据处理软件(SICap-STVR V1.0.x); (3)用于系统照明的卤素灯, 该灯共有4个, 每个功率为500 W, 由德国的Esylux公司制造, 型号是RK90000420108型; (4)送样平台, 其速度可在0~10.00 mm· s^-1之间调节。 试验中所使用高光谱系统由本实验室自主设计。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 高光谱系统Fig.1 Hyperspectral image acquisition system

为确保高光谱图像采集质量, 通过预试验可知, 送样平台的输送装置速度设为1.75 mm· s^-1, 高光谱相机物距设为260 mm, 每个测试样本共获得1 288个波长的高光谱图像信息。

研究发现, 针对苹果样本, 仅采用平均反射光谱信息来表征高光谱信息时所构建的预警模型不够稳健、 准确。 分析认为, 平均反射值光谱信息仅反映了整体光谱信息, 属于静态信息; 而光谱信息的动态变化情况没有得到反映, 致使高光谱信息的表征不够充分。 因高光谱图像灰度值的方差可以表征光谱信息的波动情况, 体现了光谱信息的动态性。 因此, 拟采用高光谱图像灰度值的平均值联合其方差来共同表征苹果样本的高光谱信息, 以构建稳健的预警模型。

1.3.1 基于均值表征的光谱曲线校正

由于试验过程中其他光源以及暗电流的影响, 需要对高光谱图像灰度值的均值和方差分别进行黑白标定^[1]。 均值表征的光谱信息校正公式如式(1)

R=(Rsample-Rdark)(Rwhite-Rdark)(1)

式(1)中: R为黑白校正后样本灰度值平均值; R_sample为样本原始灰度值平均值; R_dark为全黑标定灰度值平均值; R_white为全白标定灰度值平均值。

1.3.2 基于方差表征的光谱曲线校正

基于方差表征的高光谱信息的校正方法是, 先分别对高光谱图像n个像素点在b_i(i=1, 2, …, 1 288)波长下的灰度值a_k(k=1, 2, …, m, m≤ n× i)进行黑白校正, 然后再计算校正后像素灰度值的方差值。 简述如下:

由于每个波长的高光谱图像黑白校正原理相同, 因此以第b_i波长下像素点灰度值a_j(j=1, 2, …, n)为例进行推导。 又因黑白标定是全黑或全白标定高光谱图像, 所以第b_i波长的高光谱图像所有像素点的黑板和白板灰度值分别相等, 记为R_darki和R_whitei。 校正时, 首先根据式(1)可得出第b_i波长的每个像素点灰度值的校正值A_j; 然后, 根据A_j得出第b_i波长的像素点灰度值的方差V_i; 最后, 根据ENVI5.1软件获得的样本像素点在第b_i波段灰度值均值 a̅i和灰度值标准差σ _i。 经推导可得到全波段灰度值方差的校正公式, 如式(2)所示。

V=σ2(Rwhite-Rdark)2(2)

式(2)中: V为经黑白校正后的样本灰度值方差; σ 为样本全波段下灰度值标准差; R_dark和R_white含义与式(1)一致。

作为示例, 图2为1个测试样本在653.65 nm波长下的高光谱图像灰度图, 其大小为1 032 pixel× 900 pixel。 由于测试样本其自身近似球型的几何结构, 会使光线在其局部表面发生散射, 这将严重影响高光谱测试数据的精确性和有效性^[2]。 因此, 在避开散射区域的基础上, 以样本上部空间顶点为中心选择如图2红色区域所示的矩形区域作为感兴趣区域(ROI), 该区域大小为900 pixel× 650 pixel。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 苹果样本图像感兴趣区域提取Fig.2 Extraction of interest region in apple sample image

由于检测器异常或电流不稳定, 采集的光谱会出现尖刺等噪声信息。 为提高测试数据的准确性, 采用Savitzky-Golary(SG)平滑^[3] 、 标准正态变量变换(SNV) ^[4]、 多元散射校正(MSC)^[4]、 二阶导数(Der2)^[5] 、 一阶导数(Der1)^[5]、 局部加权散点平滑(LOWESS)^[6]共6种方法, 分别对高光谱图像灰度值均值和方差表征光谱曲线进行预处理, 并对比预处理结果选出最优方法。

苹果表皮出现的颜色转变现象, 是其在贮藏过程中最直观的品质变化。 因此, 选用苹果表皮色差变化作为表征苹果贮藏过程中品质变化的一种理化指标, 该指标也是水果成熟度判断、 等级分类时常用的指标^[7]。 已有的相关文献中均采用一种理化指标, 如可溶性固形物含量(SSC)^[8]、 失水率^[9]来表征样本贮藏中的品质变化。 分析认为, 尽管色差指标属于一种综合指标, 但也仅是从一种角度表征样本的综合属性, 仅采用色差一种理化指标, 表征信息可能不够充分; 同时, 考虑到失水率也是苹果在贮藏过程中变化比较直观的一个理化指标, 因此, 本工作对苹果样本失水率变化也进行了测定。

1.6.1 色差

对固定的6个苹果样本测定色差值, 每日固定时间使用美国X-rite公司生产的色差仪, 测定样本表皮的色差值, 色差仪型号为X-Rite Color i5型。 每次测量色差前, 都要对仪器进行黑白校正。 在每个样本表皮的不同位置, 重复3次测量, 取其平均值即为样本色差值。 根据色差值数据中的a^* 值、 b^* 值, 通过CIELAB颜色体系^[7]中的公式可得到苹果样本色调值H^* (也被称为色调角)。 表1给出了样本色调角的测试数据。 色调角越小, 说明苹果样品表皮颜色越红, 越大说明表皮颜色越黄。 色调角计算公式如式(3)

H* =tan-1b* a*(3)

表1

Table 1

表1(Table 1)

表1 苹果样本色调角变化(° ) Table 1 Hue angle change of apple samples 贮藏日数/d H* 值 贮藏日数/d H* 值 1 14.43 15 22.23 2 15.19 16 21.75 3 15.02 17 22.18 6 18.68 20 22.88 7 16.63 21 24.45 8 16.90 22 24.41 9 17.34 23 24.18 10 18.68 24 26.17 13 20.30 27 25.99 14 20.85

表1 苹果样本色调角变化(° ) Table 1 Hue angle change of apple samples

式(3)中: a^* 值代表红度, 其数值越大表明颜色越红, 越小颜色越绿; b^* 为黄度, 其数值越大表明颜色越黄, 越小颜色越蓝。

1.6.2 失水率

对固定的3个苹果样本测定失水率; 每日固定时间使用佳衡日平仪器有限公司制造的电子天平, 天平型号是DT-500B型。 失水率计算公式如式(4)所示。 表2给出苹果样本失水率测试结果。

表2

Table 2

表2(Table 2)

表2 苹果失水率变化(%) Table 2 Water loss rate change of apple samples 贮藏 日数/d 失水率 贮藏 日数/d 失水率 贮藏 日数/d 失水率 1 0.00 10 0.45 19 1.01 2 0.04 11 0.52 20 1.11 3 0.08 12 0.59 21 1.18 4 0.13 13 0.64 22 1.28 5 0.17 14 0.70 23 1.36 6 0.22 15 0.76 24 1.42 7 0.26 16 0.82 25 1.52 8 0.31 17 0.88 26 1.61 9 0.38 18 0.94 27 1.71

表2 苹果失水率变化(%) Table 2 Water loss rate change of apple samples

α=W-wW×100%(4)

式(4)中: W为苹果样本原始重量; w为第一日之后每日的样本重量; α 为样本的失水率(%)。

1.7.1 色调角

苹果样本在贮藏过程中其表皮会发生一系列变化, 表皮颜色会逐渐由红色转变为褐色或淡黄色, 表皮天然蜡质会发生分解并使表皮变得油腻^[10]; 同时, 表皮开始逐渐出现黑色斑点, 直至苹果发生腐败。 对样本观察发现, 在第21日苹果开始出现腐败, 因此, 初步将第21日确立为苹果样本的腐败基准日。 为保证后续研究数据的覆盖性、 完整性, 试验在苹果完全腐败时(第27日)结束。

根据苹果样本的a^* 值、 b^* 值以及式(3), 可得到样本色调角随贮藏时间变化的折线图。 通过对折线图和样本实际观察情况的综合分析, 将图中变化趋势出现转折的数据点所在贮藏日确定为腐败基准日。

1.7.2 叶绿素特征吸收波长

苹果在贮藏过程中果皮发生的变化, 与果皮中花青苷以及叶绿素二者含量的变化密切相关^[11]。 因此, 根据测得的叶绿素特征吸收波长(675 nm左右)^[12]处光谱反射值, 作出其均值随不同贮藏天数的变化趋势图。 将图中数值变化趋势出现转折的天数与色调角变化折线图进行对比分析, 以进一步确立腐败基准日。

为消除高光谱图像灰度值的均值和方差表征光谱信息的共线和冗余信息, 采用连续投影算法(SPA)^[13]并结合样本色调角和失水率, 提取预处理后的光谱信息的特征波长。 色调角和失水率分别是从不同的角度反映苹果样本品质的综合属性, 为尽可能使模型具有一定的鲁棒性和泛化能力, 同时缩小两种理化指标之间的差异性对模型产生的不利影响, 选取可表征两个理化指标的共有波长为特征波长。

巴氏距离(BD)算法可以准确计算分类差异^[14], 有利于区分不同天数的光谱信息之间的差异性, 因此可根据这种差异性表征不同天数苹果样本的品质变化。 模型构建时, 首先确定样本的腐败基准日, 并确定能够表征该基准日样本品质属性的特征波长; 然后, 计算不同天数特征波长与该基准特征波长之间的巴氏距离值; 最后, 根据对不同天数与基准日之间距离值的大小进行比较, 即可有效判断贮藏期间任意天数样本接近腐败的程度。 不同贮藏天数样本与基准日样本之间的距离值越小, 说明其品质越差。 巴氏距离计算公式如式(5)

Bst=(ut-us)24(σs2+σt2)+12lnσs2+σt22σsσt(5)

式(5)中: B_st(s=1, 2, …, 21, t=1, 2, …, 21)为第s个贮藏日样本特征波长与第t个贮藏日样本特征波长之间的巴氏距离; u_s和u_t为第s和t个贮藏日的样本高光谱图像灰度值均值; σ _s和σ _t为第s和t个贮藏日样本高光谱图像灰度值标准差。

通过分析比较6种数据预处理结果, 发现SG平滑可显著消除数据噪声。 图3、 图4分别给出苹果样本基于均值和方差原始光谱图及SG处理后光谱图。 图中均给出第1, 6, 11, 16和第21天样本的光谱信息。 由于1~99波段光谱反射值较高, 考虑到坐标刻度对图像视觉效果的影响, 图中仅给出100~1 288波段的光谱数据。 对比分析图3(a, b)和图4(a, b)表明, SG处理后光谱曲线尤其是曲线两端的毛刺显著减少, 光谱噪声得到了有效处理。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 苹果样本基于均值表征的光谱值SG处理前后结果 (a): SG处理前的光谱; (b): SG处理后的光谱Fig.3 The spectral data of apple samples based on the mean representation before and after SG treatment (a): Before SG smoothing; (b): After SG smoothing

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 苹果样本基于方差表征的光谱值SG处理前后结果 (a): SG处理前的光谱; (b): SG处理后的光谱Fig.4 The spectral data of apple samples based on the variance representation before and after SG treatment (a): Before SG smoothing; (b): After SG smoothing

2.2.1 色调角的基准日

图5给出样本色调角在贮藏期间变化的折线图。 贮藏期间样本色调角随天数增加呈上升趋势, 表明样本表皮颜色由红色逐渐转变为黄色, 与实际观察结果相符。 但第6日样本色调角度数值异常升高, 这可能是贮藏前期, 样本为应对叶绿素分解, 自身形成了可以合成更多叶绿素的相关调节机制^[11]。 随着前期叶绿素的迅速合成, 使得其短期内的合成速度已经超过了其分解速度, 导致第6日叶绿素的累积含量达到最高, 从而导致在样本b^* 值变化不明显的前提下, a^* 值大幅度降低。 此外, 发现第4、 5和6日叶绿素特征吸收波长光谱平均反射值逐渐降至最低点的现象(见图6), 同样表明该时间段内样本表皮叶绿素总量在逐渐增加, 吻合对图5进行分析的结果。 由于色调角在第21日之后出现较大转变, 并且连续3天数值稳定在24.40° 左右, 因此, 拟定第21日为苹果样本的腐败基准日。 第24和第27日色调角升至26° 左右, 说明苹果样本腐败加剧。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 样本色调角随贮藏天数的变化Fig.5 The change of hue angle of apple samples with the storage days

2.2.2 叶绿素特征吸收波长的基准日

图6给出了675.11 nm处苹果表皮叶绿素的平均光谱反射值变化趋势图。 由图6可知, 贮藏前期, 平均光谱反射值变化趋势出现先迅速上升然后迅速下降, 之后再快速上升再缓慢下降的现象。 这可能是由于苹果为应对其表皮叶绿素分解而自发形成了相关的调节机制, 在该机制的调节下, 贮藏过程中样本叶绿素的合成与分解同时发生, 且前期叶绿素的合成量高于分解量, 所以光谱信息出现波动。 特别是第6日, 由于样本表皮叶绿素累积量达到最高, 导致光谱反射值降到最低点, 这也与图5中第6日色调角异常偏高的试验现象相吻合。 由于后期调节机制逐渐失灵, 苹果品质开始恶化, 其叶绿素的合成逐渐减少, 表皮叶绿素总量逐渐降至低点, 所以样品对655.11 nm的波长吸收减弱, 反射升高。 从图6可以看出, 反射值在第21日出现明显变化, 21和22日保持在较高值0.76左右, 故, 将第21日设定为腐败基准日。 第23和24日光谱反射值下降, 25至27日平稳保持在0.75左右, 这可能是随着样品的腐败, 其表皮逐渐增多的黑色菌斑吸收了部分655.11 nm的波长。 因此, 根据叶绿素特征吸收波长平均光谱反射值的变化也可确定腐败基准, 并与依据色调角变化所确定的腐败基准是一致的。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 675.11 nm处的基于灰度值均值的苹果样本光谱值Fig.6 The spectral data of apple samples at 675.11 nm based on mean representation

为更加全面的获得表征苹果样本品质属性的特征波长, 采用全波段光谱信息进行特征波长提取。 将经SG处理后的全波段灰度值均值光谱信息结合样本色调角和失水率, 采用SPA提取特征波长。 样本灰度值均值结合样本色调角提取了17个波长, 结合样本失水率提取了12个; 最终确定7个二者共同波长为特征波长。 同理, 苹果样本灰度值方差结合色调角提取了15个波长, 结合失水率提取了13个波长; 最终确定8个二者共同波长为特征波长。 表3给出特征波长选取结果。 从表3可以看出, 基于灰度值均值和方差所获取的特征波长基本都集中在600~700 nm范围内。 这可能是因为该波段内原始光谱曲线波动性较大, 可有效体现不同样本之间的区别, 这也符合图3、 图4中光谱信息均在600~700 nm范围内出现较大波动的现象。

表3

Table 3

表3(Table 3)

表3 特征波长 Table 3 Feature wavelengths 表征类型 特征波长/nm 均值 428.04 635.76 652.12 676.14 689.41 705.73 960.97 方差 440.09 479.52 502.58 624.50 631.15 649.05 663.88 689.92

表3 特征波长 Table 3 Feature wavelengths

2.4.1 均值表征的建模

根据2.2已确定样本的腐败基准日为第21日, 应用前21日(含第21日)灰度值均值表征的光谱信息来表示各贮藏日的光谱信息。 建模时, 为使样本数据尽可能反映苹果样本贮藏品质变化, 提高模型可靠性, 每天应用50个测试样本的特征光谱信息的平均值来表征该天的光谱信息测量结果; 那么7个特征波长下每个贮藏日光谱信息表征结果为7维向量。 图7给出各贮藏日巴氏距离的计算结果, 横坐标为贮藏天数, 纵坐标为巴氏距离值。 根据图7可知, 前期巴氏距离值出现较小程度的波动, 这可能是样本中含有个别受到损伤的苹果, 从而使均值表征光谱信息出现一定程度的波动。 中期距离值又出现较小程度的波动, 这可能是因为此时部分样本表皮因蜡质分解而逐渐变得油腻, 使得样本对光线的吸收与反射受到干扰。 后期距离值变化基本趋于平缓, 这可能是因为样本中表皮变得油腻的苹果数量在逐渐增多, 从而增加了对建模结果的干扰强度。 对模型整体分析发现, 随着贮藏时间的增长, 巴氏距离值呈下降趋势, 说明苹果样本越来越接近腐败, 符合苹果样本在贮藏过程中的实际变化情况。 但是, 模型整体变化趋势具有一定的波动性, 容易造成误判。 这可能是提取的基于均值的高光谱表征信息不能充分反映苹果贮藏过程中的品质变化情况, 仅反映样本变化的整体信息, 而细节变化的动态信息没有被反映出来。 所以, 需补充表征信息以更充分地反映苹果在贮藏过程中的品质变化。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 基于均值表征的高光谱信息的巴氏距离建模Fig.7 The Bhattacharyya distance results of hyperspectral information based on mean representation

2.4.2 方差表征的建模

为获得苹果贮藏过程中品质变化的动态信息, 采用灰度值方差表征的光谱信息来建立预警模型。 8个特征波长每个贮藏日光谱信息表征结果为8维向量。 图8给出基于灰度值方差的各贮藏日巴氏距离计算结果。 从图8可以看出, 前6日巴氏距离值波动性较大, 可能是试验前期的样本中含有个别受到机械损伤的苹果, 从而造成方差表征的信息出现较大波动(均值表征信息受此影响不显著)。 随着贮藏时间的延长, 样本巴氏距离值呈下降趋势, 且逐渐接近腐败基准, 与实际贮藏过程相符。 说明基于方差的光谱信息也可以表征苹果品质的变化, 从动态信息方面予以体现。 对比图7和图8, 可以发现图7中在贮藏的中后期巴氏距离值偏高, 这可能是试验中后期, 样本表面的蜡质成分开始变得油腻, 并逐渐加重, 使测试样本更容易发生散射现象, 导致样本光谱反射率增强, 从而引起巴氏距离值偏高, 这也充分说明均值表征的信息属整体信息, 因此体现出了与方差表征信息建立巴氏距离模型所不同的变化趋势。

图8

Fig.8

	Figure Option ViewDownloadNew Window
	图8 基于方差表征的高光谱信息的巴氏距离建模Fig.8 Bhattacharyya distance of hyperspectral information based on variance representation

2.4.3 均值融合方差联合表征的建模

为获得更全面的反映苹果贮藏中品质变化的光谱信息, 将均值表征融合方差表征的光谱信息, 采用与2.4.1小节相同的方法, 计算了各贮藏日的巴氏距离。 联合表征下的特征波长为15个, 那么每个贮藏日光谱信息表征结果为15维向量。 图9给出均值融合方差二维表征光谱信息的巴氏距离计算结果。 从图中可以看出, 各贮藏日苹果样本的巴氏距离整体呈下降趋势, 符合苹果样本的实际腐败进程。 综合对比图7、 图8和图9发现, 图9中模型前期的波动性得到有效减弱, 说明均值融合方差联合表征的方法可有效消除个别样本损伤对模型准确性的影响。 但是, 模型中期仍存在一定程度的波动性, 这说明表皮油腻的样本数量增加的现象, 会对模型的建立产生影响。 对模型整体分析发现, 巴氏距离整体变化趋势的波动性减弱。 这说明均值融合方差表征的光谱信息, 可以有效减少样本预警结果的误判, 充分证明了运用二维表征光谱信息进行预警的有效性和稳健性。

图9

Fig.9

	Figure Option ViewDownloadNew Window
	图9 均值融合方差二维表征光谱信息的巴氏距离结果Fig.9 Bhattacharyya distance results based on 2D hyperspectral information representation with mean and variance