基于光谱成像的猪肉新鲜度空间分布预测评价方法

赵茂程 吴泽本 汪希伟,2 邢晓阳 陈加新 唐于维一

(1.南京林业大学机械电子工程学院， 南京 210037；2.南京林业大学机电产品包装生物质材料国家地方联合工程研究中心， 南京 210037)

0 引言

猪肉新鲜度的评价指标包括感官、挥发性盐基氮(Total volatile basic nitrogen，TVB-N)含量和微生物含量等指标[1-2]，本文选择TVB-N含量客观评价猪肉新鲜度。

高光谱成像技术在农林产品质量安全检测等方面获得了较多应用[3-7]。光谱图像中每个像素点均包含特定信息，因此可以实现指标空间分布预测的可视化[8-17]，即基于兴趣区域的均值光谱建立化学计量学预测模型，并将其应用于像素光谱中，从而得到指标空间分布预测图像。这一方法已被应用于农林产品品质检测工作中[18]。但像素位置缺乏微观指标参考值，无法对可视化图像进行直接评价，因此目前的研究仅停留在建立可视化图像层面，未对指标空间分布预测的效果进行评价。

然而，将众多准度相近的化学计量学模型应用到像素光谱进行指标空间分布预测时，像素预测结果存在明显差异，甚至出现预测结果不具备统计意义的情况[19-20]。究其原因，其差异主要表现在两方面：一是准度，即各像素位置微观预测值的统计均值与理化检测值的偏差程度；二是精度，即根据指标理论允许范围，异常点在兴趣区域内所占比值。因此，本文以猪肉新鲜度指标空间分布预测为例，建立基于不同光谱预处理及化学计量学模型所得新鲜度指标空间分布预测图像，从准度和精度两方面，对指标空间分布预测评价的方法进行研究。

1 材料与方法

1.1 实验材料

实验所用猪肉样本购买于江苏省南京市麦德龙下关商场，分割出5条猪肉背长肌作为研究对象，用内置冰袋的保温箱运送至南京林业大学逸夫科技实验楼。将每条猪肉背长肌分割成18块厚度约为10 mm的肉块，共获取90个实验样本。分割完成后，将实验样本放置于4℃的冰箱中冷藏保存，以维持猪肉内部的自然组织结构状态。

1.2 高光谱图像采集及处理

1.2.1高光谱成像系统

高光谱成像系统由光谱成像单元、照明系统、计算机和辅助支架组成，其中光谱成像单元包括工作在可见-近红外波段(550～1 000 nm)的CVA-200型声光可调谐滤波器(BRIMEROSE公司，美国)、ORCA-R2型可见-近红外相机(HAMAMATSU公司，日本)和可变焦镜头(18～200 mm，尼康公司，日本)；照明系统包括C3K型UPS稳压电源(山特电子(深圳)有限公司)输出的12个50 W卤素灯和1个半球形漫反射穹顶。由计算机(Yangtian A4600t， Windows XP系统)通过基于Visual Studio 2008平台编写的采集程序实现高光谱图像采集控制。

1.2.2高光谱图像采集

每隔24 h从冰箱中取出10个实验样本，用保温箱运送至高光谱成像系统处采集光谱数据。在550～1 000 nm波段范围内，区域灰度均值范围设定为1 600～2 400，采用变曝光时间采集方式，相机曝光时间如图1所示。采集步长为3 nm，曝光物距为860 mm，图像分辨率为516像素×672像素，单次采集实验样本的高光谱图像。

图1 曝光时间Fig.1 Exposure time

由于970～1 000 nm波段内图像信噪比过低，故舍弃该部分图像，仅保留550～970 nm共141个波段的光谱图像。

在猪肉样本高光谱图像采集结束后，先采集标称值75%的标准反射率标定板(SRT-75-100型，Labsphere公司，美国)的光谱图像；后盖上镜头盖，采集此时的光谱图像，即暗噪图像。将以上两种光谱图像作为当日采集猪肉样本光谱图像相对反射率校正的基准。

1.2.3高光谱图像预处理

(1)光谱均值滤波

本文采用光谱均值滤波提高像素光谱质量。设置固定宽度的滤波器窗口，沿原始光谱滑动，不同窗口的跨度相互重叠，取窗口内k个波段的平均值代替第1个波段的数据。滤波后第i个波段的光谱数据计算方法为

(1)

式中Rout(i)——光谱均值滤波后第i个波段的光谱图像

Rin(h)——第h个波段的原始光谱图像

m0——原始光谱图像的波段数量

本文采用的声光可调谐滤波器每个波段的带宽最大为20 nm。为使光谱均值滤波带宽与光学滤镜带宽相当，本文设置滤波器宽度为6、18、30、42、54 nm，对应波段数量分别为2、6、10、14、18。

(2)相对反射率校正

为进一步提高信噪比，消除图像采集过程中暗电流、背景光强度及光源分布不均匀等产生的噪声影响，需要对采集的高光谱图像进行相对反射率校正[21]。利用当日系统暗噪图像和标称值为75%的标准反射率标定板的光谱图像，对采集的样本光谱进行相对反射率校正，计算公式为

(2)

式中RT——相对反射率校正后的光谱图像

Ro——相对反射率校正前的光谱图像

Rb——当日系统暗噪图像

Rw——标称值为75%的标准反射率标定板光谱图像

w——根据标准反射率，用于Rw的标准反射率标定板各波段实际反射率

(3)肌肉有效兴趣区域创建

TVB-N含量指标主要反映肌肉中蛋白质的分解程度，因此将兴趣区域定义为肌肉有效兴趣区域(Eligible muscle region of interest，EMROI)，即仅保留肌肉区域，排除脂肪、筋膜及肉皮等组织成分的干扰。猪肉中肥肉和肌肉部分与背景在图像特定通道上有不同的反射率，因此，首先采用固定阈值法在727 nm处人工设置反射率阈值为0.3，对高光谱数据集内所有图像进行分割，实现猪肉与背景的分离；其次，由于部分猪肉样本的侧面也被采集至图像中，用半径为8像素的圆形模板对所有分离后猪肉图像进行边缘腐蚀操作，以消除对预测模型的影响，仅保留待检的上表面区域；最后在805 nm处人工设置反射率阈值0.45对所有腐蚀后猪肉图像进行分割，以剔除肥肉部分，实现肌肉的精确分离。

(4)特征光谱提取

通过求取不同波长下EMROI内所有空间像素位置处的光谱数据均值，提取得到猪肉样本不同预处理方法下的特征光谱{si}，即均值光谱。计算公式为

(3)

式中si——特征光谱(均值光谱)在第i个波段的反射率

aij——EMROI内第j个像素点在第i个波段的反射率

n——EMROI内像素点数量

m——光谱图像的波段数量

(5)特征波段选取

本文采集的高光谱图像数据量较大，包含550～970 nm范围内的141个波段信息，其中包含很多冗余信息，降低了模型运算效率。因此，需要采用合适的方法选取相关性较高的特征波段。采用连续投影算法(Successive projection algorithm，SPA)[22]对基于EMROI的全波段光谱进行特征波段提取，该算法是一种前向循环的使矢量空间共线性最小化特征变量筛选方法，可以筛选出有效信息，降低数据之间共线性影响，增加模型的鲁棒性和泛化性。

为便于工业现场应用，实现样本的快速检测，常需要组建专用的多光谱成像系统。基于扫描式光谱成像，光谱波段数量一般不超过6个；基于快照式光谱成像，相机在一次曝光周期内同时采集图像的波段数量通常为20个左右。因此，本文分别优选出6个和20个特征波段，用于建立猪肉新鲜度预测模型。

1.3 挥发性盐基氮含量测定

参照对鲜冻肉的评价标准[1]，TVB-N含量(质量比)不超过15 mg/(100 g)可认定为新鲜肉。将采集完高光谱图像的样本用保温箱运送至理化检测实验室，以测定TVB-N含量。参照对挥发性盐基氮测定的国标方法[23]，取样本表面5 mm的纯肌肉部分，剁成肉粉，采用半微量定氮法进行测定。在同等条件下，对获得的两次独立测定结果计算平均值。若两次结果的绝对差值不超过平均值的10%，以该平均值作为该实验样本的TVB-N含量理化检测值。

1.4 基于均值光谱的猪肉新鲜度预测模型建立与评价

偏最小二乘回归(Partial least squares regression，PLSR)[24]是一种广泛应用于光谱分析的线性回归建模方法，是典型的相关分析和主成分分析的集成和发展，可以同时实现提取变量特征、分析变量间相关性和回归建模，是光谱分析领域应用最广泛的方法之一[25]。其思路是从自变量集合中选取主成分，然后建立主成分与自变量的回归方程，公式为

(4)

式中Y——预测值

β0——增益常数

βi——第i个波段的增益系数

Xi——第i个波段的平均反射率

本文基于Matlab工具箱中偏最小二乘回归函数plsregress建立预测模型，得到回归参数向量β，作为化学计量学模型的波段增益。其中，βi值越大，对预测结果的贡献度越高，相应波段的放大倍率也越高。

1.5 基于像素光谱的猪肉新鲜度空间分布预测评价方法

1.5.1猪肉新鲜度空间分布预测

将基于均值光谱的化学计量学模型{β0，β1,…，βm}应用到像素光谱进行新鲜度指标空间分布预测，计算出相应像素点的TVB-N含量预测值，得到灰度图像。各像素位置微观预测值的计算公式为

(5)

式中Pj——EMROI内第j个像素点的微观预测值

根据数值大小结合伪色彩处理技术，绘制得到猪肉新鲜度可视化图像。肉样表面预测值由高(腐败)至低(新鲜)分别以红色至蓝色显示，背景设置成黑色，肥肉部分设置成白色，样本原始尺寸(即包含样本侧面的未腐蚀图像尺寸)在肉样外侧用白色线条表示，色带范围设置为0～25 mg/(100g)。

1.5.2空间分布预测评价方法

由于无法对各像素位置的微观新鲜度TVB-N含量指标进行直接评价，因此从准度及精度两方面，对基于不同光谱预处理及化学计量学模型所得肉样检测区域新鲜度指标空间分布预测结果进行评价。

(1)准度评价

(2)精度评价

求取不同波长下EMROI内所有空间像素位置处光谱数据均值得到均值光谱的过程，相当于经过一个尺寸为EMROI大小的空间均值滤波器降噪，可有效提高信噪比。在将基于均值光谱所建化学计量学模型应用到像素光谱时，由于未经过该空间均值滤波，像素光谱质量明显低于区域均值光谱，导致模型的预测精度下降，从而造成猪肉新鲜度TVB-N含量预测值出现负值。根据新鲜度理化指标的理论允许范围，以TVB-N含量微观预测值小于零的像素点在兴趣区域内所占比值作为精度评价指标，即异常预测点占比越小，空间预测精度越高。

1.5.3空间分布预测精度影响因素

(1)像素光谱质量

通过理想均值滤波对光谱图像进行6～54 nm的5种不同带宽的光谱滤波降噪，可以有效提高像素光谱质量。将暗噪图像各波段EMROI内所有像素点亮度值的标准差，记为暗噪标准差，以传感器输出值(Digital number，DN)作为计量单位，从而量化图像随机噪声抑制效果，以评价像素光谱质量。通常来讲，暗噪标准差越小，图像的噪声水平越低，像素光谱质量越高。

(2)模型波段增益

对于一般的化学计量学模型，波段增益是光谱图像每一个波段的放大倍率。本文应用PLSR建立预测模型，m个波段的增益系数分别为PLSR回归系数{β1，β2，…，βm}。

1.6 数据处理软件

光谱图像预处理、PLSR预测模型建立、空间分布预测、数据分析及图形绘制均采用Matlab R2017b软件；数据汇总整理由Excel 2016软件完成。

2 结果与分析

2.1 样本集划分

本文共采集90个样本。图像采集时，第21号样本20%表面肌肉区域被不慎掉落的脂肪组织覆盖，造成图像光谱数据与理化检测值不对应，被判为无效样本排除，剩余有效样本共89个。选择合适的定标集方法划分样本可以有效提高模型的预测能力。本文按照TVB-N含量排序，采用留出法[26]确定样本集，其中67个样本为训练集，22个样本为预测集。如表1所示，训练集、预测集和总样本的TVB-N含量变化范围为3.49～26.19 mg/(100 g)，样本标签值跨度范围在猪肉TVB-N含量正常范围内，并且数据集样本之间差异明显，具有广泛性；平均值和标准差相近，保证训练集与预测集处于相同的分布状况。数据集具有较强的泛化性和一致性，为建立可靠的猪肉新鲜度预测模型提供了数据基础。

表1 猪肉样本训练集与预测集TVB-N含量统计结果Tab.1 Statistical results of TVB-N value in training set and prediction set of pork samplesmg/(100 g)

2.2 光谱曲线分析

对550～970 nm内141个波段求取EMROI内所有空间像素位置处的光谱数据均值，可获得猪肉样本不同预处理方法下的光谱反射率曲线。图2分别表示猪肉样本原始光谱反射率曲线及滤波器宽度为6、18、30、42、54 nm的滤波后光谱反射率曲线。由图2a可以看出，不同猪肉样本之间的光谱曲线具有相似的趋势。575 nm附近的吸收峰是脱氧肌红蛋白与氧合肌红蛋白的合频，肌红蛋白的存在状态决定了鲜肉的肉色[27]。

图2 89个猪肉样本EMROI内光谱反射率曲线Fig.2 Spectral reflectance curves of 89 pork samples in EMROI

对比图2a～2f发现，随着滤波器宽度的增加，光谱均值滤波有效降低了图像噪声，抑制了原始光谱曲线的细小波动，使光谱曲线变得更加光滑，预测模型的稳健性得以提升；但曲线峰谷的数量及波动幅度出现下降，特别是575 nm附近的吸收峰，当滤波器宽度为54 nm时该吸收峰已经消失，说明噪声抑制能力增强的同时可能会导致部分细节信息的丢失。

2.3 猪肉新鲜度预测模型的建立

采用偏最小二乘回归法，分别基于原始光谱和滤波后5种不同带宽的光谱对全波段、精选出的20、6个特征波段建立新鲜度预测模型。为了得到可靠稳定的模型，采用交叉验证法对训练集进行验证，以减小训练过程中的过拟合。设置交叉验证次数为5，分别建立PLSR猪肉新鲜度预测模型，模型性能如表2所示。

表2 PLSR猪肉新鲜度预测模型准度Tab.2 Accuracy of PLSR pork freshness model

2.4 猪肉新鲜度空间分布预测评价

2.4.1猪肉新鲜度空间分布预测准度评价

图3 均值光谱预测值与像素光谱预测值对照结果Fig.3 Comparison results of mean spectral prediction and pixel spectral prediction

从原理上对样本基于均值光谱的预测值和各像素位置微观预测值的统计均值相等的原因进行探讨。偏最小二乘回归建模过程中，利用线性化学计量学算法，样本基于均值光谱的预测值公式为

(6)

将基于均值光谱所建模型应用到像素光谱进行新鲜度指标空间分布预测时，利用线性化学计量学算法，各像素位置微观预测值的统计均值公式为

(7)

式中P2——各像素位置微观预测值的统计均值

将P2化简，可得

(8)

由式(8)可知，对于线性化学计量学算法，样本基于均值光谱的预测值恒等于各像素位置微观预测值的统计均值。

文献[19]利用类似光谱成像系统进行过新鲜度分布预测。与本文不同，其所得样本表面新鲜度像素级定量可视化预测均值与有效兴趣区域特征光谱预测值之间存在明显差异，而该差异是由在光谱预处理环节采用的标准正态化处理(Standard normal variate，SNV)算法所引入的非线性所致。

2.4.2猪肉新鲜度空间分布预测精度评价

(1)各种新鲜度分布可视化图像感官对比

根据各像素点预测值结合伪色彩处理技术，绘制得到猪肉新鲜度可视化图像。同一肉样经不同光谱预处理及化学计量学预测模型所得的新鲜度指标空间分布可视化图像存在明显差异。整体来看，随滤波器宽度的增加，大部分化学计量学模型的可视化效果呈现出逐渐提升的趋势；但存在个别模型游离在趋势之外，可视化效果不升反降。选取一块肥瘦分明的典型猪肉样本，利用不同的化学计量学模型对新鲜度空间分布预测，在同一维度上，即特征波段数量相等或同为全波段建模时，各选择3个典型模型建立可视化图像，如图4所示。图4a～4c中，图4c的可视化效果最好；图4d～4f中，图4f的可视化效果最好；图4g～4i中，图4h的可视化效果最好。肉样表面预测值由高(腐败)至低(新鲜)分别以红色至蓝色显示，背景设置成黑色，肥肉部分设置成白色，样本原始尺寸(即包含样本侧面的未腐蚀图像尺寸)在肉样外侧用白色线条表示，色带范围设置为0～25 mg/(100g)。颜色越接近红色，预测值越高，猪肉越腐败；颜色越接近蓝色，预测值越低，猪肉越新鲜。可视化图像右侧上方附以该模型的各波段暗噪标准差和化学计量学模型增益，暗噪标准差用蓝色线条表示，各波段增益值用红色线条表示；右侧下方附以肉样表面TVB-N含量预测值分布直方图，预测值等于零用红色线条标出，预测值小于零像素占比用红色数字标出。

图4 猪肉新鲜度预测模型空间可视化图像Fig.4 Spatial visualization images of pork freshness prediction model

(2)新鲜度空间分布预测精度对比

将准度相当的化学计量学模型应用到像素光谱，指标空间分布预测精度会存在明显差异，图5为猪肉新鲜度18种空间分布预测的精度与像素光谱质量及预测模型波段增益的关系图。图5a、5c、5e、5g、5i、5k、5m、5p、5r分别与图4a～4i所示模型对应。从中可以看出，大部分模型随滤波器宽度的增加，预测精度逐渐提升；但当滤波器宽度增加到一定程度时，预测精度开始下降。图5a～5f、图5g～5l、图5m～5r分别表示基于6个特征波段、20个特征波段以及全波段建模时，原始光谱及滤波器宽度依次为6、18、30、42、54 nm滤波后光谱的空间分布预测精度与像素光谱质量及预测模型波段增益的关系。椭圆的横轴半径表示该模型各波段增益标准差，纵轴半径表示暗噪标准差；箱线图上下边缘分别表示预测值小于零像素占比的最大值和最小值，箱子上下两端边位置分别表示预测值小于零像素占比的上下四分位数，箱子中间红色线条表示预测值小于零像素占比的中值，极端异常值用红色“+”标出。可以发现，在同一维度上，即特征波段数量相等或同为全波段建模时，椭圆的横轴半径和纵轴半径越小，箱子面积就越小，箱子位置越低，模型的预测精度越高。

图5 猪肉新鲜度18种空间分布预测的精度与像素光谱质量及预测模型波段增益的关系Fig.5 Relationships between precision of 18 spatial distribution prediction of pork freshness and spectral quality of pixels and band gain of prediction model

基于6个特征波段建模时，随滤波器宽度的增加，椭圆纵轴半径逐渐减小，暗噪标准差逐渐降低，由20.58 DN下降至7.91 DN，像素光谱质量逐渐提升；图5a～5e模型椭圆横轴半径无明显差异，各波段增益标准差均处于107.17～149.45的范围内。此时箱子面积逐步减小，箱子位置逐渐降低，预测精度逐步提升。而图5f模型相较于图5e模型，椭圆横轴半径有所增大，各波段增益标准差达到205.59，增长幅度为37.56%；在暗噪标准差减小幅度仅为4.24%的情况下，箱子面积增大，箱子位置升高，预测精度略有下降。

基于20个特征波段建模时，相同的预测精度变化趋势再次出现。随滤波器宽度的增加，椭圆纵轴半径逐渐减小，暗噪标准差逐渐降低，由20.51 DN下降至7.86 DN，像素光谱质量逐渐提升；图5g～5k模型椭圆横轴半径虽略有起伏，但各波段增益标准差均处于45.66～90.97的范围内。此时箱子面积逐步减小，箱子位置逐渐降低，预测精度逐步提升。而图5l模型相较于图5k模型，椭圆横轴半径有所增大，各波段增益标准差达到99.55，增长幅度为31.51%；在暗噪标准差减小幅度仅为4.97%的情况下，箱子面积增大，箱子位置升高，预测精度略有下降。

基于全波段建模时，随滤波器宽度的增加，椭圆纵轴半径逐渐减小，暗噪标准差逐渐降低，由20.41 DN下降至7.83 DN，像素光谱质量逐渐提升；图5m～5p模型椭圆横轴半径相差无几，各波段增益标准差最大仅为43.85。此时，箱子面积逐步减小，箱子位置逐渐降低，预测精度逐步提升。而图5q模型相较于图5p模型，椭圆横轴半径有所增大，各波段增益标准差达到114.85，增长幅度为1 036.00%，远超过7.23%的暗噪标准差减小幅度；图5r模型椭圆横轴半径进一步增大，各波段增益标准差达到287.54，增长幅度为150.36%，而暗噪标准差减小幅度仅为4.63%；与此同时，图5q和图5r模型的箱子面积明显增大，箱子位置升高，预测精度显著下降。

总之，预测精度与像素光谱质量的变化趋势非常接近，即像素光谱质量越高，预测精度越高；而异常模型的出现则是由于各波段增益标准差明显增大，增长幅度远大于暗噪标准差减小幅度所造成的。

像素光谱质量与空间分布预测精度的相关度散点图如图6所示。暗噪标准差与预测值小于零像素占比中值的拟合方程为Y1=1.295X-6.672，大部分蓝色圆点分布趋势贴近拟合直线，其相关系数为0.72，表明两者具有良好的一致性，即像素光谱质量与空间分布预测精度相关度较高。而左上方界外点的出现，则是受到各波段增益的影响，说明单纯考虑化学计量学模型的准度指标，忽视其波段增益数值偏大，可能导致其应用到像素光谱时精度极差。暗噪标准差与像素预测值标准差中值的拟合方程为Y2=0.591X+0.456，橙色三角形分布非常贴近拟合直线，且没有界外点出现，其相关系数为0.76，说明两者同样具有良好的一致性，即随像素光谱质量的提升，暗噪标准差减小，此时像素预测值标准差随之减小，兴趣区域内各像素点预测值分布范围逐渐减小，分布更加精确集中，并不是整体向预测值增大的方向移动。

图6 像素光谱质量与空间分布预测精度相关度Fig.6 Pearson correlation coefficient between pixel spectral quality and spatial distribution prediction precision

综上所述，新鲜度指标空间分布预测的精度明显受到像素光谱质量及化学计量学模型波段增益值的共同影响，其中像素光谱质量占主导作用。因此，在应用线性化学计量学算法，将基于兴趣区域均值光谱所建化学计量学模型应用到像素光谱进行新鲜度指标空间分布预测时，应充分提高像素光谱信噪比，同时注意限制化学计量学的波段增益值，以提高空间分布预测的精度。

3 结论

(1)虽然受限于缺乏像素位置微观指标参考值无法进行直接评价，但通过对像素位置微观预测值的统计均值进行准度评价及根据指标理论正确值的有效范围进行精度评价的方法，能够对基于光谱成像的化学计量学指标空间分布预测质量进行综合评价。

(2)在应用线性化学计量学算法，将基于兴趣区域均值光谱所建化学计量学模型应用到像素光谱进行新鲜度指标空间分布预测时，其准度不会下降。

(3)新鲜度指标空间分布预测的精度明显受到像素光谱质量及化学计量学模型波段增益值的共同影响，其中像素光谱质量占主导作用(R=0.72)。

(4)在应用线性化学计量学算法，将基于兴趣区域均值光谱所建化学计量学模型应用到像素光谱进行新鲜度指标空间分布预测时，高准度化学计量学模型不一定适用于像素光谱进行空间分布预测，区域预测效果与空间分布预测效果并无必然联系；可以通过提高像素光谱信噪比和限制模型波段增益提高预测的精度。