李建盟，殷 勇，于慧春，袁云霞，李 迎

（河南科技大學食品與生物工程學院，河南 洛陽 471023）

番茄是我國主要種植蔬菜作物之一，以其豐富的營養(yǎng)價值和良好的風味深受人們喜愛[1-2]。番茄中較高的水分含量且表面附著的微生物易導致其腐敗變質(zhì)，使貨架期縮短[3-4]；若不能及早預警腐敗會造成大量經(jīng)濟損失。因此，貯藏期間對番茄品質(zhì)的監(jiān)控及早期預警尤為重要。

目前，在番茄的品質(zhì)監(jiān)控方面，主要通過硬度、可溶性固形物含量、pH值、總糖含量、番茄紅素含量等指標表征番茄品質(zhì)[5-7]。近些年，近紅外光譜、高光譜等技術(shù)結(jié)合理化指標也逐漸被用于表征番茄品質(zhì)[8-10]。但是這些方法由于取樣的隨機性和樣本間差異性不能較準確地表征貯藏期間番茄整體質(zhì)量變化狀況。3D熒光技術(shù)作為近年來興起的檢測技術(shù)，可以實現(xiàn)微量級甚至痕量級的定性及定量分析，并以其靈敏度高、特異性強而廣泛應用于眾多領(lǐng)域，如醫(yī)學[11-12]、環(huán)境[13-14]、化工[15-16]、食品[17-18]等。3D熒光光譜由激發(fā)波長、發(fā)射波長及熒光強度組成，其包含了被測對象豐富的信息，有利于獲得可靠的分析結(jié)果。

不同品種及產(chǎn)地番茄的品質(zhì)存在差異，但番茄在貯藏期間通過呼吸作用交換到貯藏室氣氛中的物質(zhì)（如醇、醛、酮類等揮發(fā)性物質(zhì)[19-20]）以及表面逸散到氣氛中的微生物與其品質(zhì)變化往往是對應的。微生物內(nèi)具有色氨酸、酪氨酸和苯丙氨酸等可產(chǎn)生熒光的物質(zhì)[21-23]，且醇、醛、酮類等揮發(fā)性物質(zhì)中的苯環(huán)結(jié)構(gòu)也可產(chǎn)生熒光信息。因此，利用貯藏室氣氛信息可實現(xiàn)番茄貯藏過程中品質(zhì)變化的監(jiān)控。貯藏室氣氛是一個隨機的均勻整體，將貯藏室氣氛作為研究對象不僅可以消除取樣隨機性和樣本間差異性的影響，而且能表征番茄貯藏期間整體品質(zhì)的狀況。此外，由于生理鹽水不產(chǎn)生熒光，且微生物可在生理鹽水中存活，番茄產(chǎn)生的具有羥基和羧基的揮發(fā)性物質(zhì)也可溶于生理鹽水，因此生理鹽水可作為收集載體，進行熒光檢測和分析，以實現(xiàn)番茄品質(zhì)的監(jiān)控及早期預警。

番茄屬于呼吸躍變型蔬菜[24]，品質(zhì)變化會出現(xiàn)拐點?；诖?，本實驗借助貯藏室氣氛的熒光信息，運用小波包分解方法選擇表征番茄品質(zhì)變化的3D熒光光譜特征波長，嘗試根據(jù)特征波長確定番茄的腐敗基準，最后以特征波長作為不同貯藏日的光譜表征信息，并運用馬氏距離實現(xiàn)番茄腐敗的早期預警。

1 材料與方法

1.1 材料與試劑

‘左粉’番茄于2020年9月13日購于洛陽大張超市，總質(zhì)量為350 kg。

NaCl（分析純） 天津市北辰方正試劑廠。

1.2 儀器與設備

HKM-II型空氣浮游微生物采樣器 廣東環(huán)凱微生物科技有限公司；TG16A臺式高速離心機 上海盧湘儀離心儀器有限公司；紫外增強型三維熒光光譜儀（光源為150 W紫外增強型連續(xù)氙燈，波長精準度為±1 nm，信噪比大于800∶1）由實驗室自制。

1.3 方法

1.3.1 熒光光譜采集

采用常溫貯藏加速番茄貯藏品質(zhì)變化。貯藏期間，每天記錄貯藏環(huán)境溫、濕度，貯藏溫度波動范圍24.5～27.1 ℃，相對濕度波動范圍64.0%～66.1%。每天固定時間采集3D熒光數(shù)據(jù)，為確保腐敗基準確定的可靠性，實驗截止時間從感官出現(xiàn)腐?。ǖ?個貯藏日）后延續(xù)至番茄出現(xiàn)大面積腐爛，使實驗數(shù)據(jù)能更好地覆蓋腐敗基準。實驗共持續(xù)8 d。

以0.9 g/100 mL NaCl溶液經(jīng)121 ℃、20 min滅菌后作為采集貯藏室氣氛樣本的載體。將采樣器置于貯藏庫中央，放入裝有0.9 g/100 mL NaCl溶液的培養(yǎng)皿采樣，氣體采樣量為1 000 L。采樣結(jié)束后將樣品移至離心管，并用0.9 g/100 mL NaCl溶液定容到4 mL，10 000 r/min離心1 min，棄去2 mL上清液，將剩余樣本振蕩均勻后進行3D熒光數(shù)據(jù)采集。激發(fā)波長范圍200～450 nm，步長5 nm；發(fā)射波長范圍260～800 nm，步長0.4 nm，積分時間設置為1 s。

1.3.2 光譜預處理

由于采集的3D熒光光譜數(shù)據(jù)存在瑞利散射[25]，其數(shù)值遠超光譜區(qū)域數(shù)值，會掩飾光譜的真實信息，因此采用Matlab2014a軟件中eemscat軟件包去除瑞利散射。去除瑞利散射后再對光譜進行Savitzky-Golay（S-G）平滑處理[26]。

1.3.3 光譜小波包分解及表征

3D熒光光譜中相鄰或相近的激發(fā)波長具有相似性和漸變性，相鄰激發(fā)波長變化不明顯。由于小波包分解具有信號“顯微放大功能”，因此采用小波包分解將各激發(fā)波長對應的發(fā)射光譜分解后進行能量值計算，以分析不同激發(fā)波長間的差異。以3 層小波包分解為例進行說明，如圖1所示，原始信號S經(jīng)第1層分解得到兩個頻段，低頻段A1和高頻段D1；第2層分解時，由A1和D1分解得到AA2、DA2、AD2、DD2這4 個頻段，AA2為最低頻段，其他頻段對應的頻率范圍依次遞增；同理進行第3層分解可得到8 個頻段，AAA3為最低頻段。小波包分解后得到的每個頻段為一個系數(shù)集，如式（1）所示。

式中：X為頻段對應系數(shù)集；m為頻段內(nèi)的系數(shù)個數(shù)。

各頻段信號能量可以用小波包能量表征，各頻段小波包能量計算如式（2）所示。

式中：Ej為第j頻段小波包能量。

小波包總能量（E）為各個頻段能量的總和，如式（3）所示。

式中：n為當前小波包分解層數(shù)下頻段數(shù)量。

各頻段小波包能量占比（P）按式（4）計算。

圖1 3 層小波包分解示意圖Fig. 1 Schematic diagram of three-layer wavelet packet decomposition

小波包分解中，基函數(shù)和分解層數(shù)的選擇對于光譜信息的表征具有一定影響，根據(jù)已有研究，選擇4 階Symlet小波作為小波基函數(shù)[27]。另外，小波包分解層數(shù)太少不利于表征光譜特征信息，分解層數(shù)過多會導致計算量過大和有用信息丟失，因此選擇3 層分解。

小波包能量可表征光譜信息，而光譜能量主要集中在低頻段，高頻段一般為噪聲和干擾；低頻段小波包能量占比越高表明對應的激發(fā)波長更能反映光譜的真實信息。經(jīng)計算發(fā)現(xiàn)，小波包分解后最低頻段的小波包能量占比最高，因此，用最低頻段的能量值作為激發(fā)波長的表征量，通過分析最低頻段的小波包能量表征熒光光譜中各激發(fā)波長下光譜的變化狀況。

1.3.4 特征激發(fā)波長的選擇

1.3.4.1 特征激發(fā)波長初選

采集的3D熒光光譜中含有冗余信息，需提取可表征番茄品質(zhì)變化的特征波長。為篩選有用的熒光區(qū)域激發(fā)波長，將最低頻段小波包能量占比超過99%所對應的激發(fā)波長作為有用熒光區(qū)域信號，而低于此值舍去。利用最低頻段小波包能量占比實現(xiàn)激發(fā)波長的初選。

1.3.4.2 特征激發(fā)波長終選

以初選后各激發(fā)波長對應發(fā)射光譜小波包分解的最低頻段能量值作為表征量，繪制不同激發(fā)波長-對應能量值的變化曲線。能量曲線中峰值一側(cè)能量值波動大，另一側(cè)則變化平穩(wěn)、波動小。能量波動大的一側(cè)所反映的光譜變化規(guī)律不明顯；而能量波動小的一側(cè)，光譜變化規(guī)律性較強。因此，選擇曲線峰值波動小的一側(cè)所對應激發(fā)波長作為特征激發(fā)波長范圍。

1.3.5 特征發(fā)射波長的選擇

數(shù)據(jù)采集時發(fā)射波長步長為0.4 nm，剔除無熒光信息背景區(qū)域后仍有726 條發(fā)射波長，發(fā)射波長數(shù)量多不利于分析。在化工[28]、環(huán)境[29]、食品[30]等領(lǐng)域研究中，3D熒光發(fā)射波長步長一般設置為1、2、5 nm。而本研究中發(fā)射波長步長為0.4 nm，雖分辨率高可細化光譜，被測對象可被更精細描述，但考慮到相鄰波長的相似性和漸變性，可將相鄰波長進行歸組，對特征激發(fā)波長下的發(fā)射光譜進行波段劃分，通過對各波段的分析選取特征發(fā)射波長。采用取整原則，由于發(fā)射波長2 nm尺度中僅包含5 個步長，數(shù)量過少不足以表征熒光信息，因此選擇尺度為4、8、12 nm進行波段劃分，通過對比分析確定波段尺度。

將每個波段中發(fā)射波長按對應激發(fā)波長從小到大順序首尾相連構(gòu)成一個新的數(shù)據(jù)鏈，分別對各波段相應數(shù)據(jù)進行3 層4 階Symlet小波包分解，以其小波包分解最低頻段的能量值作為該波段的表征量，特征激發(fā)波長下發(fā)射光譜可表征為各波段對應的小波包表征能量變化曲線。小波包能量值大表明其對應波段信號強，信息更為豐富，因此，選取不同貯藏日表征能量變化曲線峰值所對應的波段進行特征發(fā)射波長的選擇。

1.3.6 腐敗基準的確定

貯藏期間番茄品質(zhì)不斷變化，當番茄品質(zhì)達到腐敗狀態(tài)時，腐敗前后熒光光譜信息不同，表現(xiàn)為熒光光譜值出現(xiàn)異?；蛲蛔儭Ｋ刑卣骷ぐl(fā)波長對應能量值的總和即為該貯藏日樣本光譜檢測的表征量。當能量總和變化出現(xiàn)突變或拐點時即表明熒光光譜信息出現(xiàn)突變，說明番茄開始出現(xiàn)腐敗，由此確定腐敗基準。

此外，由于番茄腐敗前后熒光光譜信息不同，采用歐氏距離對不同貯藏日特征激發(fā)波長下發(fā)射光譜的表征能量值進行聚類分析。將每日所采集番茄熒光光譜信息作為單獨一類，計算類間歐氏距離，并將距離最小的兩類合并為一個新類，然后計算新類和其他類的距離再進行合并，重復此過程至僅剩兩類。最后這兩類即為腐敗前、后的光譜信息表征。

運用特征激發(fā)波長對應能量值總和的變化與特征發(fā)射波長光譜表征能量值的聚類分析兩種方法可最終確定腐敗基準；并且這兩種方法所得結(jié)果的一致性可充分說明腐敗基準的準確性和合理性。

1.3.7 根據(jù)腐敗基準反選特征發(fā)射波長

根據(jù)1.3.5節(jié)所述，以各貯藏日能量變化曲線峰值對應波段進行特征發(fā)射波長的選擇，而不同貯藏日能量變化曲線峰值對應波段不同，因此，為進一步縮小范圍，按照1.3.6節(jié)方法，采用歐氏距離對不同貯藏日特征激發(fā)波長下發(fā)射光譜的表征能量值進行聚類分析，選擇腐敗前、后兩類距離最大波段中的發(fā)射波長作為終選特征發(fā)射波長，以此實現(xiàn)運用腐敗基準反選特征發(fā)射波長，并最終確定特征發(fā)射波長。

1.3.8 早期腐敗預警

腐敗基準確定后，可用腐敗基準日的特征發(fā)射波長光譜表征能量值作為基準信息向量，計算腐敗基準日前各貯藏日特征發(fā)射波長光譜表征能量值到基準信息向量的馬氏距離。馬氏距離越大，表明番茄品質(zhì)越好，遠未達到腐敗狀態(tài)；馬氏距離越小，表明番茄品質(zhì)越差，越接近腐敗狀態(tài)。根據(jù)馬氏距離實現(xiàn)番茄的品質(zhì)監(jiān)控及早期預警。馬氏距離按式（5）計算。

式中：D為第i個貯藏日特征發(fā)射波長光譜表征能量值到基準信息向量的馬氏距離；xi（x1,x2, …,xm）表示第i個貯藏日特征發(fā)射波長光譜表征能量值；y為腐敗基準日特征熒光光譜信息向量；S為腐敗基準日特征熒光光譜信息向量的協(xié)方差。

綜合上述研究方法，番茄貯藏過程中，基于貯藏氣氛的3D熒光特征波長選擇方法及實現(xiàn)腐敗預警的研究流程如圖2所示。

圖2 基于貯藏氣氛的3D熒光特征波長選擇方法及實現(xiàn)腐敗預警的研究流程Fig. 2 Flow chart of tomato spoilage early warning based on 3D fluorescence feature wavelength selection of storage room gas

2 結(jié)果與分析

2.1 光譜預處理結(jié)果

如圖3A所示，原始3D熒光光譜僅可觀察到一道散射墻，其熒光強度遠超周圍熒光區(qū)域，無法獲取有用信息。對原始光譜數(shù)據(jù)進行去除瑞利散射和S-G平滑處理，圖3B顯示出熒光信號區(qū)域及熒光峰的位置，有利于后續(xù)研究。

圖3 原始（A）和去除瑞利散射、S-G平滑處理后（B）熒光光譜Fig. 3 Original fluorescence spectra (A), and fluorescence spectra subjected to Rayleigh scattering removal and S-G smoothing (B)

2.2 特征激發(fā)波長選擇結(jié)果

2.2.1 特征激發(fā)波長初選結(jié)果

對3D熒光光譜進行分析時發(fā)現(xiàn)，發(fā)射波長550～800 nm為背景區(qū)域，不存在熒光信號，因此去除。將各激發(fā)波長對應的發(fā)射光譜進行小波包分解，小波包分解產(chǎn)生頻段的能量占比如圖4所示。熒光區(qū)域300 nm激發(fā)波長對應發(fā)射光譜最低頻段的能量占比超過99%，其他頻段能量占比極小；而非熒光區(qū)域215 nm激發(fā)波長對應發(fā)射光譜雖然最低頻段能量占比最高，但是與熒光區(qū)域差異較大。最低頻段能量占比顯示了不同激發(fā)波長間差異性，以最低頻段能量占比超過99%為準則，選取激發(fā)波長范圍260～405 nm中30 個激發(fā)波長作為初選激發(fā)波長。

圖4 熒光區(qū)域300 nm激發(fā)波長（A）和非熒光區(qū)域215 nm激發(fā)波長（B）對應發(fā)射光譜的小波包能量占比Fig. 4 Wavelet packet energy ratio of emission spectra corresponding to excitation wavelength of 300 nm in fluorescent region (A) and 215 nm in non-fluorescent region (B)

2.2.2 特征激發(fā)波長終選結(jié)果

如圖5所示，不同貯藏日的能量曲線峰值對應的激發(fā)波長各不相同，圖中右側(cè)能量曲線變化趨勢平穩(wěn)，所以選擇靠近右側(cè)峰值對應的激發(fā)波長335 nm作為選取激發(fā)波長的基點，共選取基點右側(cè)包括335 nm在內(nèi)的15 個激發(fā)波長。

圖5 不同貯藏日特征激發(fā)波長-對應能量值的變化曲線Fig. 5 Variation curves of wavelet packet energy versus feature excitation wavelength at different storage times

2.3 特征發(fā)射波長選擇結(jié)果

以12 nm尺度劃分的波段為例，對該波段相應數(shù)據(jù)進行小波包分解，不同貯藏日特征激發(fā)波長下發(fā)射光譜該波段對應的小波包能量變化曲線如圖6所示。貯藏過程中不同貯藏日波段對應的能量變化曲線峰值在第10～14個波段之間，因此選擇第10～14連續(xù)5 個波段中的發(fā)射波長作為特征發(fā)射波長。

圖6 不同貯藏日特征激發(fā)波長下發(fā)射光譜波段對應的小波包能量變化曲線（劃分尺度12 nm）Fig. 6 Variation curves of wavelet packet energy versus emission spectral bands corresponding to different feature excitation wavelengths at different storage times (wavelength intervals 12 nm)

以尺度4、8、12 nm對特征激發(fā)波長下發(fā)射光譜進行波段劃分，各波段對應發(fā)射光譜進行小波包分解，以最低頻段能量作為波段表征量，特征激發(fā)波長下發(fā)射光譜可表征為各波段對應的能量變化曲線。能量變化曲線峰值對應波段的發(fā)射波長范圍如表1所示。整個貯藏過程中，尺度4、8、12 nm劃分波段能量變化曲線峰值對應發(fā)射波長范圍分別為376～416、372～428、368～428 nm。因尺度4 nm和8 nm選取的發(fā)射波長范圍包含在尺度12 nm所確定的發(fā)射波長范圍內(nèi)，故取12 nm作為最終波段劃分尺度。尺度12 nm所確定發(fā)射波長范圍368～428 nm作為特征發(fā)射波長。

表1 不同貯藏日特征激發(fā)波長下發(fā)射光譜能量變化曲線峰值對應波段發(fā)射波長Table 1 Emission wavelength corresponding to peak values of emission spectral energy change curves at different feature excitation wavelengths at different storage times nm

2.4 腐敗基準確定結(jié)果

2.4.1 能量總和變化確定的腐敗基準

所有特征激發(fā)波長下發(fā)射光譜小波包分解最低頻段能量值的總和即為該貯藏日樣本光譜檢測的表征量。由圖7可知，隨著貯藏時間延長，能量總和呈下降趨勢。第5個貯藏日開始曲線變化趨于平穩(wěn)，說明熒光光譜信息在第5個貯藏日出現(xiàn)拐點，表明番茄品質(zhì)變化達到腐敗狀態(tài)。因此，確定第5個貯藏日為腐敗基準日。

圖7 不同貯藏日特征激發(fā)波長下發(fā)射光譜小波包分解最低頻段能量總和變化Fig. 7 Changes in the lowest wavelet packet frequency band energy of emission spectra corresponding to different feature excitation wavelengths as a function of storage time

2.4.2 腐敗基準確定結(jié)果

采用歐氏距離對不同貯藏日特征激發(fā)波長下發(fā)射光譜表征能量值進行聚類分析。由圖8可知，第1～4個貯藏日的熒光光譜特征信息聚為一類，第5～8個貯藏日聚為一類，說明前4 個貯藏日的熒光光譜特征信息相近，而在第5個貯藏日番茄品質(zhì)發(fā)生變化，因此確定腐敗基準日為第5個貯藏日。

圖8 不同貯藏日特征激發(fā)波長下發(fā)射光譜表征能量值的聚類分析Fig. 8 Cluster analysis of energy values of emission spectra at different feature excitation wavelengths as a function of storage time

2.5 腐敗基準反選特征發(fā)射波長結(jié)果

根據(jù)能量變化曲線峰值選取的特征發(fā)射波長在第10～14個波段（劃分尺度12 nm）中（圖6），對各波段對應發(fā)射波長光譜表征能量值進行聚類分析，發(fā)現(xiàn)第11個波段對應發(fā)射波長光譜表征能量值聚類分析腐敗前、后的距離最遠，說明此波段的發(fā)射波長能夠最優(yōu)表征貯藏過程中番茄品質(zhì)的變化。第11個波段對應發(fā)射波長光譜表征能量值聚類分析結(jié)果如圖9所示，第1～4個貯藏日為腐敗前光譜信息，第5～8個貯藏日為腐敗后光譜信息，與能量總和變化確定的腐敗基準結(jié)果相印證。因此，選擇第11個波段中的發(fā)射波長作為熒光光譜的特征發(fā)射波長。

圖9 最優(yōu)波段對應發(fā)射波長光譜表征能量值的聚類分析Fig. 9 Cluster analysis of spectral energy values of emission wavelengths corresponding to optimal band

2.6 馬氏距離對番茄早期腐敗的預警結(jié)果

由圖10可知，隨著貯藏時間延長，貯藏日特征發(fā)射波長光譜表征能量值到基準信息向量的馬氏距離呈下降趨勢，表明番茄逐漸趨于腐敗，當馬氏距離接近零時說明番茄開始腐敗。雖然隨著貯藏的進行，腐敗基準距離出現(xiàn)波動，但整體呈下降趨勢，因此應用馬氏距離的變化進行番茄貯藏過程中品質(zhì)的監(jiān)控及早期預警可行且有效。

圖10 腐敗前各貯藏日與腐敗基準間的馬氏距離變化趨勢Fig. 10 Trend of Mahalanobis distance between storage days before spoilage and spoilage benchmark as a function of storage time

3 結(jié) 論

通過采集不同貯藏時間番茄貯藏室氣氛3D熒光光譜數(shù)據(jù)，經(jīng)原始光譜數(shù)據(jù)去除瑞利散射、平滑和剔除無信息熒光背景區(qū)域處理，對每個激發(fā)波長對應發(fā)射光譜進行小波包分解，以最低頻段小波包能量占比進行特征激發(fā)波長初步篩選，再根據(jù)激發(fā)波長小波包能量變化優(yōu)選出15 個特征激發(fā)波長。然后，以尺度12 nm對特征激發(fā)波長下發(fā)射光譜進行波段劃分，在獲取波段小波包能量變化曲線峰值的前提下，選取5 個峰值波段中的發(fā)射波長作為特征發(fā)射波長。在此基礎(chǔ)上，根據(jù)不同貯藏日，運用特征激發(fā)波長對應能量值總和的變化與特征發(fā)射波長光譜表征能量值的聚類分析均可確定第5個貯藏日為番茄腐敗基準日。另外，根據(jù)腐敗基準前、后信息歐氏距離反選特征發(fā)射波長，進一步篩選出最優(yōu)波段對應發(fā)射波長作為終選特征發(fā)射波長。最后，將腐敗前各貯藏日與腐敗基準間的馬氏距離用于監(jiān)控貯藏過程中番茄品質(zhì)變化及腐敗預警。結(jié)果表明，運用小波包分解可選擇特征激發(fā)和發(fā)射波長，根據(jù)特征激發(fā)波長下發(fā)射光譜的總能量變化及其聚類分析兩種方法均可確定準確、可靠的腐敗基準，且所建立的馬氏距離預警方法具有可行性。以貯藏室氣氛作為研究對象克服了取樣隨機性和樣本間差異性的影響，使研究結(jié)果更具有魯棒性。