基于高光譜成像的大米中蛋白質含量的可視化研究

王朝輝，趙層，趙倩，王艷輝，賴翰卿，王曉東，王靖會

（1.吉林農業大學食品科學與工程學院，吉林長春130118；2.廣東地球土壤研究院，廣東廣州510385；3.吉林省長春市凈月開發區福祉街道辦事處，吉林長春130122；4.吉林農業大學信息技術學院，吉林長春130118）

中國稻米的產量占世界稻米產量的30%以上[1]，而吉林省“梅河大米”是中國粳米地理標志產品，其產地位于世界黃金糧食生產帶（北緯45°）。在實際生活中梅河大米的種類繁多，通常采用凱氏定氮法和分光光度法等化學類方法測定不同品種中大米的蛋白質含量，但是這些傳統的化學方法不但對樣品本身具有破壞性，而且步驟繁瑣，檢測周期過長。紅外光譜技術作為一種快速無損的檢測手段，已經廣泛的應用于大米主要成分（蛋白質[2]、脂肪[3]、淀粉[4]、水分[5]）的檢測，但是其只能根據光譜信息得到組分的含量，無法實現更加直觀的表達即含量的可視化。

高光譜是包括圖像信息和光譜信息[6]的三維立方體數據，得到的高光譜圖像既包含了大米的內部信息（內部物理結構、化學成分的信息）也包含了大米的外部信息（粒型、缺陷等），能彌補近紅外不能快速識別某一物質的空間分布情況[7-8]即圖像的缺失。目前，高光譜成像技術在農產品方面的主要應用有瘀傷識別分類[9]、籽粒的識別[10]、品種的檢測[11]、含水量預測[12]等，而高光譜成像技術在大米中的主要應用有堊白檢測[13]、水分含量檢測[14]、摻假大米的檢測[15]等。雖然高光譜檢測技術在組分的檢測方面已有諸多應用，但關于高光譜成像技術應用于大米中蛋白質含量分布可視化的研究還未見相關報道。

本文選取了吉林省梅河市4 個產區的3 個品種（稻花香、秋田小町、吉粳60）的大米為研究對象，利用高光譜成像技術對采集的大米進行檢測，求取大米感興趣區域的平均光譜，為了降低光譜的信噪比和得到相對穩健的模型，對光譜進行了卷積平滑（savitzkygolay，S-G）、均值中心化（mean centering）、多元散射校正 （multiplicative scatter correction，MSC）3 種算法的預處理，建立了大米蛋白質含量的偏最小二乘回歸（partial least-squares regression，PLSR）、主成分回歸（principal component regression，PCR）、誤差反向傳播神經網絡（back propagation neural networks，BP-NN）3 種預測模型。采用SPA 選取特征波長，建立特征波長模型，并將大米高光譜圖像轉變為蛋白質含量分布圖，實現不同品種產地大米蛋白質含量的可視化。

1 材料與方法

1.1 材料與儀器

以吉林省梅河市種植的稻米品種（稻花香、秋田小町、吉粳60）為研究對象，通過空間分層隨機布點的原則采集了包括稻花香、秋田小町、吉粳60 3 個品種在內的120 個稻米樣品。試驗前先將采集的稻米進行脫殼、去糙，之后剔除霉變、未成熟和粒型不完整的籽粒，選擇大小均勻、飽滿的籽粒進行高光譜掃描和進行大米蛋白質含量的化學測定。采集的120 個大米樣品中86 個用于校正集，34 個用于檢驗集。

ImSpector V10E-QE 型高光譜圖像系統：芬蘭Spectral Imaging 有限責任公司；9589（EKE-ER）全光譜鹵素燈：荷蘭Philips 公司；PSA200-11-XZolix 電控位移臺：北京Zolix 公司；C8484-05G 相機：日本Hamam atsu Photonics；GZ02DS20 可升降樣品臺：北京廣正儀器有限公司；JLGJ4.5 型檢驗礱谷機、JNMJ3 型檢驗碾米機：臺州市糧儀廠；JXFM110 錘式旋風磨：上海嘉定糧油儀器有限公司；K1100 全自動凱氏定氮儀、SH220 石墨消解儀：濟南海能儀器股份有限公司；PHS-3C 型酸度計：上海儀電科學儀器股份有限公司。

1.2 高光譜圖像的采集與校正

測量的光譜范圍是408.360 0 nm～1 007.220 0nm，包含477 個波段，系統采集大米樣品前通過對比多次預實驗結果來調整曝光時間、物距等各項參數。采集圖像時物距為13.5 cm，曝光時間為15 ms，位移臺移動速度為1.62 mm/s。設定好儀器參數后將放有大米樣品的黑色底板平鋪在載物臺上（黑色底板的反射率接近于0，大米樣本按照網格單粒隨機的方式擺放至載物臺黑板上），進行大米高光譜數據的采集。高光譜成像系統的結構示意圖如圖1 所示。

圖1 高光譜成像系統示意圖Fig.1 Schematic diagram of hyper-spectral imaging system

在高光譜圖像采集過程中由于會受到外界（自然光光源分布不均）和系統（暗電流噪聲）等因素的干擾[16]，導致得到的高光譜圖像包含大量的無效信息，為了減少由于上述因素導致的噪聲增強現象，因此在采集大米高光譜圖像之前，必須先進行原始大米圖像的黑白校正。校正公式如下

式中：Ra為原始大米的圖像；Rd為黑板圖像；Rb為白板圖像；Rn校正后大米圖像。

1.3 蛋白質含量的測定

將采集完光譜數據的大米樣品進行磨粉過篩，采用GB5009.5-2016《食品安全國家標準食品中蛋白質的測定》測量大米中的蛋白質含量，之后把稱量好的大米樣品和催化劑放入準備好的消化管中（大米樣品重量 0.5g 左右，CuSO4和 K2SO4的重量分別是 0.2 g 和3 g），準備完成后在消化管中倒入10 mL 的濃硫酸進行消化，并做空白管試驗，當消化管內溶液的顏色變為綠色透明時證明其消化完全，則可通過全自動凱氏定氮儀對其進行蛋白質含量的測定，每個樣品測量3 次取其平均值作為該大米樣品蛋白質含量的化學值。

1.4 數據處理

1.4.1 感興趣區域與平均光譜的提取

在ENVI5.0 軟件中對獲得的高光譜圖像進行感興趣區域的選取，選取3×15 網格內所有單粒米樣中間部位11×11 像素的正方形區域作為1 個感興趣區域（region of interest，ROI），ROI 全部像素點的平均光譜作為建立模型的光譜數據。

1.4.2 光譜預處理

為了獲得相對穩定和準確的模型，從原始光譜中摒棄噪音等過多無用信息的干擾，從而放大差異，提高分辨率。采用卷積平滑（savitzky-golay，S-G）、均值中心化（mean centering）、多元散射校正（multiplicative scatter correction，MSC）3 種光譜預處理方法。通過3種預處理算法的對比，以期找到適合本次試驗的光譜預處理方法。

1.4.3 模型的建立

采用偏最小二乘回歸（partial least-squares regression，PLSR）、主成分回歸（principal component regression，PCR）、誤差反向傳播神經網絡（back propagation neural networks，BP-NN）3 種常用的算法建立蛋白質的預測模型，采用不同建模算法，最終呈現的建模效果也各不相同。本文以大米為試驗對象，以校正集和檢驗集樣品的相關系數（Rc2、Rp2）和均方根誤差（RMSEC、RMSEP）為評價模型性能的指標，通過對比最終得到穩健性好的大米蛋白質預測模型。

1.4.4 偏最小二乘回歸模型

偏最小二乘回歸 （partial least square regression，PLSR）是多因變量對多自變量的回歸模型，是常用的一種光譜建模分析方法[17]。PLSR 通過對多自變量和多因變量進行壓縮分解，分解時得到兩者最強對應關系的主成分，加強模型的預測能力。

1.4.5 主成分回歸模型

主成分回歸（principal component regression，PCR）是一種多元線性回歸分析方法，利用主成分分析將樣品光譜矩陣進行分解，把得到的變量組合成相互無關的新變量[18]，用較少的變量代表原有的變量，能有效解決預測變量高度相關和共線的問題。

1.4.6 誤差反向傳播神經網絡模型

誤差反向傳播神經網絡（back propagation neural networks，BP-NN）是一種按誤差逆傳播算法訓練的多層前饋網絡[19]，是目前應用最廣泛的神經網絡模型之一。本文選用有導師學習算法來調整神經元間的聯接權，使得網絡輸出更符合實際。

1.4.7 連續投影算法（successive projection algorithm，SPA）特征波長的選取

連續投影算法（successive projection algorithm，SPA）是一種使矢量空間共線性最小化的前向變量選擇算法[20]，采用SPA 對原始高光譜數據進行特征波長的選取，解決了高光譜數據自身維度高，變量多等問題，由于后期模型的建立使用消除了大量冗余信息的高光譜數據，因此可以根據挑選的特征波長（保留光譜信息中最優的變量組）建立簡化高效的模型。

1.4.8 大米中蛋白質含量的可視化

通過對比得到最適合本研究的大米蛋白質含量的預測模型，將像素點的光譜數據帶入預測模型，可得到各像素點位置的蛋白質含量。通過各像素點位置蛋白質含量的高低，使用MATLAB 對灰度圖像進行各像素點顏色的設定。光譜預處理、預測模型的建立和大米蛋白質可視化的實現均在MATLAB2016a 軟件中完成。

2 結果與分析

2.1 大米的原始光譜曲線

圖2 是選取ROI 全部像素點的平均值得到的反射率光譜（408.360 0 nm～1 007.220 0 nm），120 個大米全波長光譜的趨勢大致相同，說明大米在全波長的情況下具有相同的反射特性，但在反射強度上存在差異。

圖2 120 個大米樣本光譜圖Fig.2 Spectral map of 120 rice samples

2.2 蛋白質含量統計分析

對凱氏定氮法測量的大米蛋白質的化學值進行統計分析，得到校正集、檢驗集和全部樣品的最小值、最大值、平均值和標準偏差，結果如表1 所示。

表1 校正集和檢驗集樣品蛋白質含量Table 1 Calibration set and test set sample protein content

從表1 中可以看出校正集蛋白質含量范圍大于檢驗集的蛋白質含量范圍，說明校正集和檢驗集的劃分合理。

對120 個大米樣品的蛋白質含量進行統計分析得到如圖3 所示的柱形圖。

圖3 不同品種和產地的蛋白含量柱狀圖Fig.3 Bar graph of protein content of different varieties and regions

如圖所示3 個品種的大米在吉樂鄉產區的蛋白質含量明顯高于其他產區，不同品種大米的蛋白質含量在產區間的分布規律也不相同，稻花香從曙光到黑山頭蛋白質含量表現為依次遞增，而秋田小町和吉粳60從曙光、灣龍再到黑山頭蛋白質含量表現為先增加后減少。對品種和產地進行雙因素方差分析，得到關于品種因素、產地因素和品種產地的交互對應的P 值，因為所得三者的P 值均遠遠小于0.01，因此品種和產地對大米蛋白質含量均有顯著影響。

2.3 光譜預處理

通過對比不同預處理算法建立的PLSR 模型的相關系數和均方根誤差，得到最優的光譜預處理算法。表2 是光譜預處理的結果，通過對比可以發現MC 預處理后的光譜所建立模型最優，其校正集和檢驗集的相關系數為0.927 5 和0.907 9，SG 和MSC 預處理后所得模型檢驗集相關系數下降為0.807 0 和0.897 1。

表2 基于不同預處理方法的PLSR 模型結果Table 2 Results of PLSR models based on different preprocessing methods

2.4 全波段模型的建立

本文以477 個波段的光譜數據為自變量，以測量的蛋白質含量為因變量，建立PLSR、PCR、BP 神經網絡預測模型，基于全波段不同建模算法的結果見表3。

表3 基于全波段不同建模算法的結果Table 3 Results based on different modeling algorithms for the full band

如表3 所示，建立PLSR 模型時確定的最優潛在變量個數是13，其校正集的相關系數和均方根誤差RMSEC 分別為0.927 5 和0.034 6，其檢驗集的相關系數和均方根誤差RMSEP 分別為0.907 9 和0.050 8；建立PCR 模型時確定的最優主成分個數是17，其校正集的相關系數和均方根誤差RMSEC 分別為0.920 0和0.037 9，其檢驗集的相關系數和均方根誤差RMSEP 分別為0.897 1 和0.055 7；建立BP 神經網絡模型時，神經網絡隱含層神經元數20，其校正集的相關系數和均方根誤差RMSEC 分別為0.968 5 和0.015 7，其檢驗集的相關系數和均方根誤差RMSEP 分別為0.862 5 和0.074 5。通過對3 種建模算法的結果對比得到最優的建模算法為PLSR。

蛋白質PLSR 模型的實測值和預測值的散點圖見圖4。

圖4 蛋白質PLSR 模型的實測值和預測值的散點圖Fig.4 Scatter plot of measured and predicted values of protein PLSR model

2.5 SPA特征波長的選取

特征波長的計算過程在MATLAB 中實現，設定最小和最大的選定波長數為5 到50，最終得到的結果如圖5 所示。

圖5 SPA 篩選特征波長Fig.5 SPA screening characteristic wavelength

通過SPA 篩選出27 個特征波長，分別為：416.690 0、429.840 0、426.250 0、459.910 0、469.590 0、519.590 0、861.890 0、881.170 0、892.750 0、910.760 0、927.480 0、939.060 0、959.650 0、966.080 0、971.230 0、975.080 0、978.940 0、982.800 0、990.510 0、993.090 0、995.660 0、1 000.800 0、1 003.370 0、948.070 0、747.900 0、487.810 0、1 005.940 0 nm。

PLSR 通過計算分別得到合適的X 和Y 的潛變量P 和Q，通過迭代法調整得分矩陣T 和U，同時保持殘差矩陣E 和F 的值接近0，求得T 和U 的內在聯系，間接得到X 和Y 的映射關系[21]。

基于SPA 篩選的27 個特征波長的反射值作為自變量X，全自動凱氏定氮儀測定的大米蛋白質含量為因變量Y 建立PLSR 模型，得到27 個特征波長的反射值和蛋白質含量的內在關系，并用此關系在已知波譜反射值的情況下對蛋白質含量值進行預測。圖6 所示SPA-PLSR 模型的校正集和檢驗集的相關系數分別為0.923 2 和 0.903 9。

圖6 蛋白質SPA-PLSR 模型的實測值和預測值的散點圖Fig.6 Scatter plot of measured and predicted values of the protein SPA-PLSR model

2.6 大米中蛋白質含量分布的可視化

選取稻花香、秋田小町和吉粳60 的大米高光譜圖像進行成像分析，通過ENVI5.0 從大米樣本的高光譜圖像中提取27 個特征波長下的高光譜圖像，然后提取特征圖像中每個像素點的光譜反射值，將提取的所有像素點的光譜反射值代入2.5 已建立的PLSR 模型，最終預測到每個像素點的蛋白質含量。將高光譜灰度圖像轉入MATLAB2016a 中，高光譜灰度圖像中所有像素點的灰度值用得到的對應像素點的蛋白質含量值代替，生成如圖7 所示的大米蛋白質含量分布圖。

圖7 不同品種產地蛋白質含量分布可視化圖Fig.7 Visualization of protein content distribution in different varieties of production areas

如圖7 所示3 個品種在吉樂鄉產區的蛋白質含量值最高，在曙光、灣龍和黑山頭產區蛋白質含量無明顯差異，因為吉樂鄉采集的大米地處山區與其他3 個產區相比距離最遠，地貌差異最大，因此形成的大米蛋白質含量差異明顯，表明產地對大米蛋白質的含量是有影響的，與圖3 所示的規律基本相同。

從全部大米蛋白質含量分布圖可以看出米粒內部大米蛋白質的分布是不均勻的，但是總分布趨勢基本相同，胚內多于胚乳，大米蛋白質含量的高低主要表現為胚乳內蛋白質含量分布情況，這與Little and Dawson（1960）公布的大米蛋白質染色圖得到的結論基本相同。

3 結論

采用高光譜成像技術對大米中蛋白質含量分布的可視化進行了可行性研究，通過MC 光譜預處理方法和SPA 特征波段的選取，得到了簡化高效的PLSR蛋白質含量預測模型，基于該定量模型的基礎上對不同品種和不同產地大米中蛋白質含量的分布進行可視化研究。由于同一品種間大米形狀相近，難以通過普通RGB 圖像區分，通過對蛋白質含量分布成像可以為大米產地的識別提供思路，而對比不同品種間大米的蛋白質含量分布圖，可以為后期選育大米品種提供依據。