孔慶明，谷俊濤,2，高 睿，李澤東，馬 錚，蘇中濱*

(1.東北農(nóng)業(yè)大學 電氣與信息學院，黑龍江 哈爾濱 150030;2.黑龍江省網(wǎng)絡空間研究中心，黑龍江 哈爾濱 150090)

玉米種植占我國糧食種植面積的42%，2019年我國玉米產(chǎn)量達2.57億噸，與此同時也遺留下2.41億噸秸稈，如何有效利用玉米秸稈成為行業(yè)難題[1]。玉米秸稈約占整株玉米生物量的50%，目前多數(shù)國家將秸稈作為反芻動物的重要粗飼料來源，秸稈資源飼料化成為現(xiàn)代農(nóng)業(yè)的重要發(fā)展方向[2]。但玉米秸稈的粗蛋白含量低，纖維含量高，適口性差，限制了其有效利用[3-4]。目前秸稈的氨、堿化處理結合超聲技術可有效打破纖維素與木質(zhì)素之間的連接，使玉米秸稈更易為微生物附著，從而加速秸稈消化過程[5]，達到提高粗蛋白含量，改善適口性，增加采食量，提高消化率的最終目的。但在進行最優(yōu)氨堿超聲條件判定過程中，秸稈中粗蛋白、纖維含量的測定基于化學檢驗，方法速度慢，成效性差，不適合大批量樣品的測定和篩選。目前針對秸稈成分的快速檢測方法較多，從非接觸式方法來看光譜技術成為首選，高光譜及近紅外檢測方法應用較多[6-11]。本文結合近紅外光譜分析技術對氨、堿化玉米秸稈粗蛋白含量快速檢測方法開展研究，并對模型構建過程中數(shù)據(jù)冗余問題進行探究，采用組合間隔偏最小二乘(SIPLS)變量優(yōu)選方法進行特征波長選擇。通過構建快速檢測模型實現(xiàn)氨堿最優(yōu)條件的判定，解決秸稈低能量、多纖維、難消化的營養(yǎng)特性問題，對充分發(fā)揮秸稈潛在營養(yǎng)優(yōu)勢、解決人畜爭糧問題、推動低碳農(nóng)業(yè)發(fā)展具有重要意義。

1 實驗部分

1.1 樣品采集、制備及標定

1.1.1 樣品采集選取東北農(nóng)業(yè)大學高標準試驗田的優(yōu)質(zhì)玉米秸稈(先玉335，中單909-第一積溫帶上限主栽品種)共計107株為實驗樣品，樣品制備儀器采用微型植物粉碎機、超聲儀和50 mL離心管；化學試劑為5%尿素溶液和4%氫氧化鈉溶液。將收獲后的秸稈選取根上40 cm處切割成6組3 cm長度樣品(每3組做1個平行)，分別對平行樣本進行氨、堿加工處理。

堿化處理：取秸稈干物質(zhì)，將4%氫氧化鈉溶液溶于水，均勻噴灑在秸稈上，實驗室條件下保持含水量為45%左右，密封保存，經(jīng)堿化處理30 d后，打開密封袋，放置排氨72 h，待排氨后放置于65 ℃烘箱中烘干48 h至恒重，再研磨粉碎過200目篩[12]。

無處理空白樣：晾曬制風干樣秸稈，未經(jīng)其他處理，稱重后將秸稈置于65 ℃烘箱中烘干48 h至恒重，再研磨粉碎過200目篩。氨、堿化處理與無處理空白樣處理同時進行。

1.1.3 超聲處理取過200目篩的秸稈粉末在水浴條件下進行超聲波處理。超聲功率分別為60、75、90、105、120 W，超聲時間分別為5、10、15、20、25 min，固液比分別為1∶5、1∶10、1∶15、1∶20、1∶25，容器內(nèi)聲功率密度分別為1.2、1.5、1.8、2.1、2.4 W·mL。處理完成后烘干并嚴格密封，保證良好的厭氧環(huán)境，在室溫靜置穩(wěn)定后測定秸稈中水分、粗蛋白、半纖維素、纖維素和木質(zhì)素的含量。

1.1.4 粗蛋白含量標準測定依據(jù)國標GB 5009.3-2010《食品中蛋白質(zhì)的測定標準》[13]中的凱氏定氮法對氨堿及超聲處理后的玉米秸稈固體粉末進行粗蛋白測算，得107個樣本的粗蛋白含量分布區(qū)間為2.481 8%～6.519 3%，平均值為3.724 6%。

1.2 光譜采集

經(jīng)氨、堿及超聲處理后的秸稈樣本為固態(tài)粉末，采用Antaris Ⅱ光譜儀對秸稈固體樣品在4 000～12 000 cm-1波段范圍內(nèi)進行漫反射方式掃描，掃描分辨率為4 cm-1，采用積分球掃面方式，背景掃描及樣品掃描均設定為64次，共計獲得2 075個特征波段點。

圖1 玉米秸稈近紅外光譜圖像Fig.1 Near infrared spectra image of corn straw

2 結果與討論

2.1 模型構建前處理

2.1.1 消除異常樣本107個玉米秸稈光譜掃描結果如圖1所示，部分樣品明顯偏離樣本群體中心，因此采用基于3D視圖的雙殘差融合杠桿值法及單殘差法進行異常樣本判定[14]，綜合二者判斷結果進行異常樣本剔除，對粗蛋白樣本構建交互驗證(Cross validation)定量分析模型，統(tǒng)計樣本X殘差(代表光譜)、樣本Y殘差(代表粗蛋白含量)及杠桿值(Leverage表示差異性)。由圖可見，45、99樣品具有明顯的異常分布(距離群體中心較遠，圖2A)且具有較高的殘差(圖2B)，均方根誤差(RMSE)分別為1.690 3、1.443 2，因此將其剔除，對剩余105個樣本進行正態(tài)分析，統(tǒng)計樣本選擇及分布，判定其是否具有代表性。

2.1.2 光譜去噪及樣品分類首先采用小波變換方法進行噪聲處理[15]，在小波降噪過程中閾值的選取的方法直接影響降噪質(zhì)量，而閾值方法中軟、硬閾值各有優(yōu)缺點，其中經(jīng)軟閾值去噪后信號較為平滑，但會丟失信號的部分特征，即相似性欠佳；硬閾值雖可保留信號特征，但平滑性不足。通常軟閾值去噪相似性上的誤差在允許范圍內(nèi)，所以較硬閾值法應用更廣泛，本文選取軟閾值法進行光譜去噪。采用基于DaubechiesN(DBN)系列中DB2 4層分解對比缺省閾值、Bridge-Massart策略及Penalty閾值3種數(shù)學模型的信號重構結果與原信號的信號噪聲標準差作為評價標準，結果表明：上述3種模型下信號標準差分別為0.012 36、0.063 67、0.031 07。將基于缺省閾值重構的光譜數(shù)據(jù)與粗蛋白含量數(shù)據(jù)進行交互驗證，從表1中可以看出,經(jīng)異常樣本剔除及小波去噪后粗蛋白模型的決定系數(shù)(R2CV)從0.788 9提升至0.920 8，交互驗證均方根誤差(RMSECV)從0.475 2降至0.329 1，表明以上方法可有效提升模型準確性。圖3為原始光譜及3種閾值方法信號重構后的光譜圖及105個樣本的正態(tài)分析結果，可看出該樣品集分布均勻，具有較好的代表性。

表1 玉米秸稈光譜預處理后模型交互驗證結果Table 1 Corn straw model validation results after pretreatment

對105個樣品進行分類，依據(jù)Kennardstone進行樣品集分類，最終選取70 個樣品為校正集，35 個為驗證集，其統(tǒng)計分析結果如表2所示。

表中：ρ為土體密度，ω為天然含水率，ωp為土體塑限，ωL為土體液限，c為粘聚力，φ為內(nèi)摩擦角，k為土體的滲透系數(shù)。

表2 玉米秸稈粗蛋白定量分析樣品集分類Table 2 Classification of corn straw crude protein sample set

2.2 特征波段選取

玉米秸稈光譜含有大量冗余信息，降低了模型解析速度，需對其進行特征波段優(yōu)選，常用方法有主成分分析法(PCA)、相關系數(shù)法(CC)、模擬退火法(SA)、連續(xù)投影算法(SPA)、遺傳算法(GA)、間隔偏最小二乘(IPLS)、后向區(qū)間間隔偏最小二乘(BIPLS)、組合間隔偏最小二乘(SIPLS)、競爭性自適應重加權采樣法(CARS)等[16-19]。其中IPLS是將整個光譜區(qū)域劃分為等寬的若干個子區(qū)間，然后在每個子區(qū)間用PLS方法建立回歸模型，并用交互驗證法確定最佳主因子數(shù)，以交互驗證均方根(RMSECV)作為局部模型的精度衡量標準，取精度較高的局部模型所在的子區(qū)間作為特征子間。SIPLS是在IPLS基礎上聯(lián)合同一次區(qū)間劃分中精度較高的幾個局部模型所在的子區(qū)間作為一個整體建立模型[20]，從中選擇使RMSECV值最小的組合區(qū)間即為最優(yōu)的特征子區(qū)間。組合區(qū)間的數(shù)量一般取2～4個，但隨著組合區(qū)間的增加，其計算量及計算時間也隨之呈指數(shù)級增長。BIPLS可彌補IPLS單一區(qū)間的問題，也可彌補SIPLS計算時間較長的問題，其原理是將光譜區(qū)間分為n個等寬的子區(qū)間，每次預留出1個子區(qū)間，把剩余的n-1個子區(qū)間作為一個大區(qū)間來建立偏最小二乘回歸模型，計算模型RMSE值，使剩余的n-1個子區(qū)間中具有最小RMSECV值的子區(qū)間即為第一個排除區(qū)間，直到計算只剩一個區(qū)間為止。本研究將分別論證IPLS、BIPLS以及SIPLS的特征優(yōu)選實現(xiàn)過程及結果。

2.2.1 IPLS將整個譜圖區(qū)間依據(jù)間隔數(shù)10、20、30、40進行區(qū)間劃分并分別構建交互驗證模型，各間隔區(qū)間的RMSECV分布見圖4，從圖中可以看出，RMSECV最小值明顯集中于9 600～10 400 cm-1。圖5為10～40間隔數(shù)下所選取最優(yōu)區(qū)間及最優(yōu)區(qū)間下校正集及驗證集樣本分布，其數(shù)據(jù)結果見表3。結合圖5和表3可見，10間隔數(shù)下第8間隔在主因子數(shù)為7的情況下驗證效果最好，驗證集相關系數(shù)(rp)為0.966 6，驗征集決定系數(shù)(R2P)為0.934 3，驗征集均方根誤差(RMSEP)為0.277 8。由此可見，IPLS的精準間隔劃分可為特征波段的選取提供精準定位的判定依據(jù)。

2.2.2 BIPLS采用BIPLS依據(jù)間隔數(shù)10、20、30、40進行區(qū)間劃分，每次預留出1個子區(qū)間，計算其余所有區(qū)間的均方根誤差(RMSE)，并不斷排除多余區(qū)間，由于運算數(shù)據(jù)較多，運算量較大，因此本文僅展示間隔數(shù)40的波段選擇過程，其特征選擇的預留間隔、子區(qū)間、RMSE及變量數(shù)的變化結果見表4。從表中可看出RMSE隨變量的不斷減少而逐漸降低，直至預留間隔數(shù)9時RMSE達最小值0.249 0，此后RMSE隨變量的繼續(xù)減少發(fā)生波動并呈上升趨勢，因此判定該間隔為其最優(yōu)區(qū)間解。選擇預留間隔數(shù)為9時的468個變量(實際間隔24～29)作為間隔數(shù)40的最優(yōu)選擇，其余間隔數(shù)10、20、30情況下依此類推，統(tǒng)計結果見表5。表中數(shù)據(jù)可見BIPLS在間隔數(shù)為10時選擇7、8間隔(8 808～10 404 cm-1)具有最優(yōu)模型驗證精度，rp為0.978 0，R2P為0.956 4，RMSEP為0.229 1。其中BPLS間隔數(shù)為10時最優(yōu)定量分析模型校正及驗證結果見圖6。

表3 IPLS優(yōu)選后定量分析模型驗證結果Table 3 Validation results after IPLS selection

表4 間隔40條件下RMSE及變量數(shù)量變化情況Table 4 Change of RMSE and variable number under 40 interval

表5 BIPLS優(yōu)選后定量分析模型驗證結果Table 5 Validation results after BIPLS selection

圖6 BIPLS特征間隔選取及定量分析模型構建結果Fig.6 Optimal feature interval selection and model building results of BIPLSA:optimal feature interval selection(最優(yōu)特征區(qū)間選取結果);B:predicted results of calibration set(校正集預測結果);C:predicted results of validation set(驗證集預測結果)

2.2.3 SIPLSSIPLS以組合波段為主，其運算過程較繁雜，以間隔數(shù)40為例，其預留間隔、RMSE及變量數(shù)量變化見表6。選取2～4個變量組合并提取最優(yōu)驗證結果，在采用SIPLS方法特征選取的過程中發(fā)現(xiàn)，在間隔數(shù)10～30情況下選取2～4個變量組合，每次交互驗證運算時間為0.5～1.0 min；而隨著變量增大其模型驗證時間呈指數(shù)級增長，在間隔數(shù)30時4變量組合執(zhí)行時間明顯增大，間隔數(shù)40時3變量組合模型運算次數(shù)達27 405次，4組合(間隔數(shù)40)時模型運算次數(shù)達91 300次，運算時間激增至24 min，因此本實驗不計入間隔數(shù)40條件下的變量優(yōu)選結果。另外，即使選擇不同間隔數(shù)，其特征波段區(qū)間依然較為明顯和集中，表6數(shù)據(jù)可見在執(zhí)行間隔數(shù)30條件下2變量組合24、28具有最優(yōu)驗證結果，由此可見，不同間隔數(shù)及不同變量組合條件下其特征波長的選取均具有高度一致性，基本集中于9 800～10 400 cm-1，但執(zhí)行間隔數(shù)40的多變量組合方式的運算時間過長。SIPLS方法30間隔數(shù)條件下2變量組合(24、28)的特征波長區(qū)間選取及定量分析模型校正及驗證結果見圖7。

表6 SIPLS優(yōu)選后定量分析模型驗證結果Table 6 Validation results after SIPLS selection

圖7 SIPLS特征間隔選取及定量分析模型構建結果Fig.7 Optimal feature interval selection and model building results of SIPLSA:optimal feature interval selection(最優(yōu)特征區(qū)間選取結果);B:predicted results of calibration set(校正集預測結果);C:predicted results of validation set(驗證集預測結果)

2.3 分析討論

為進一步驗證玉米秸稈粗蛋白含量快速檢測中SIPLS特征優(yōu)選方法的優(yōu)勢，對多種特征波長選取方法及結果進行綜合比較，PCA、相關系數(shù)法、CARS、GA、MWPLS等特征優(yōu)選結果見圖8。驗證結果統(tǒng)計見表7。結果顯示：采用SIPLS選取波段10 128～10 398 cm-1及11 196～11 462 cm-1時具有最優(yōu)模型驗證結果，rp為0.978 4(R2P為0.957 2)，RMSEP為0.221 1。

綜合對比IPLS/BIPLS及SIPLS等其他特征波長選擇方法總結如下：

(1) 主成分分析(PCA)：主成分分析法在進行特征選取中依據(jù)計算協(xié)方差等統(tǒng)計特征代表性實現(xiàn)數(shù)據(jù)降維，選取7個有效主成分即可代表所有特征波段(圖8A)，實現(xiàn)了數(shù)據(jù)重構及數(shù)據(jù)降維的特征優(yōu)選。但相較于偏最小二乘建模方法，其特征選擇僅針對光譜，不考慮樣本特性與光譜的關聯(lián)性，使得模型驗證精度降低，且較低的重構主成分數(shù)無法實現(xiàn)對大量光譜信息的完整表達，無法完全描述非線性問題，表7中實驗結果也印證了PCA在玉米秸稈粗蛋白模型構建中特征優(yōu)選效果不夠理想。

(2) 相關系數(shù)法：相關系數(shù)法的思想提出較早，是早期光譜模型構建的特征判定方法，分為自相關性(光譜-光譜)及互相關性(光譜-樣品特性)，通常自相關較低、互相關較高會具有較好的模型精度。但從相關研究可以看出[21]，由于近紅外吸收特性較低，多種樣品間自相關與互相關性通常呈矛盾狀態(tài)，如4 000～8 000 cm-1波數(shù)范圍內(nèi)相關系數(shù)普遍低于0.3，粗蛋白與光譜呈現(xiàn)較差的弱相關性，r大于0.7的強相關性區(qū)間為11 352～11 997 cm-1，相比于其他特征波段選擇方法，具有明顯的差異性(圖8B)。

表7 不同特征波段優(yōu)選定量分析模型及驗證結果Table 7 Quantitative analysis model and validation results of different characteristic wavenumber

(3) CARS與GA：二者基本原理類似，均為基于進化理論而衍生的算法，其中CARS依據(jù)“適者生存”原理以指數(shù)衰減函數(shù)變化為判定依據(jù)，篩選出每次循環(huán)構建的回歸系數(shù)絕對值最大的變量，如圖8C所示，共選取146 個特征變量。GA通過不斷種群迭代形成最優(yōu)子代，通過信息的選擇、變異、交叉來實現(xiàn)，本研究設定適應度函數(shù)為RMSE，通過不斷迭代變異找到適應性最強即RMSE最小的特征變量，可見9 883、10 256 cm-1等8個波段構成的特征波段組合具有最小的RMSEP(圖8D)。但CARS和GA兩種方法均需借助先前經(jīng)驗，影響選擇結果的因素較多，需設計與調(diào)整的參數(shù)(如選擇方式、交叉重組方式等)較多，且選擇結果的隨機性較大，需多次驗證(依據(jù)經(jīng)驗10 次平行)方可篩選出最優(yōu)特征波段點，因此所選取結果不具有唯一性。

(4) 移動窗口最小二乘法MWPLS：移動窗口偏最小二乘基于窗口沿著光譜軸連續(xù)移動，每移動一個波長點，采用交互驗證方式建立一個模型，得到系列不同窗口(移動波長點)和主因子數(shù)對應的殘差平方和(PRESS或SSR)，但移動窗口大小的寬度需不斷嘗試(本文通過多次試驗選取31窗口，圖8E)。從試驗過程來看，窗口大小設置過大容易忽略和錯過特征波長點，設置過小則容易陷入局部循環(huán)。

(5) SIPLS在原有IPLS基礎上實現(xiàn)間隔波段自由組合，特征組合更為自由化，在實現(xiàn)特征間隔變化的基礎上可實現(xiàn)特征波長的精確定位，不會出現(xiàn)GA、CARS等單獨波段點提取情況，同時由于多組合、可變特征間隔相比于MWPLS的固定窗口大小更具靈活性的特性，不會破壞光譜數(shù)據(jù)的連貫性和特征吸收的漸變性，使得整個光譜區(qū)間分為多個模塊進行分層次交互驗證，因此在多組合、多間隔情況下模型表現(xiàn)較好，表7中最優(yōu)模型(SIPLS 30 間隔，2組合)校正結果也論證了該方法的特性。但同時該法也存在弊端，主要表現(xiàn)為多間隔、多波段組合條件下特征選取運算量較大，系統(tǒng)響應較慢，大量特征組合變量形式構建交互驗證模型消耗的時間隨間隔增加呈指數(shù)級增長(40間隔4組合運算次數(shù)91 300次，耗時24 min)。目前，基于SIPLS可實現(xiàn)全譜波段4 000～12 000 cm-1下低于30間隔數(shù)4變量組合的特征優(yōu)選，并呈現(xiàn)較好的準確性與實時性，但間隔數(shù)超40且變量組合數(shù)量大于4時現(xiàn)有服務器難以實現(xiàn)特征波段的快速、準確選取，特征選取響應時間較長。因此下一步工作會重點結合深度學習構建多層解析網(wǎng)絡，整合優(yōu)化非線性算法來實現(xiàn)多變量組合的快速選取，解決特征波段快速、精準定位的難題。

3 結 論

本文構建了玉米秸稈中粗蛋白近紅外光譜定量分析模型，對樣本進行異常剔除及光譜去噪后，對特征波段選取方法進行了探討，重點對IPLS及其改進型方法BIPLS、SIPLS的波長選取原理、選取過程進行闡述，實現(xiàn)了定量分析模型的快速準確定標，為秸稈氨堿化處理最優(yōu)條件判定提供了數(shù)據(jù)支持，研究結果如下：

(1) 數(shù)據(jù)預處理(剔除異常樣本、光譜去噪)可有效提高模型樣本的分布均勻性和代表性。本研究對107個樣本進行殘差計算，剔除2個具有明顯較大殘差的樣本，交互驗證模型決定系數(shù)R2CV從0.788 9提升至0.894 8，RMSECV從0.475 2下降至0.343 4。使用小波變換去除光譜噪聲后，R2CV從0.894 8提升至0.920 8，RMSECV從0.343 4下降至0.329 1。

(2) 特征波段選取可有效提取代表性數(shù)據(jù)，降低數(shù)據(jù)冗余，提升模型解析速度并提高模型驗證精度。本文采用IPLS、BIPLS、SIPLS對光譜區(qū)間重新排列，依據(jù)交互驗證誤差根結果進行多種特征波段選取(PCA、CARS、相關系數(shù)法、GA、MWPLS等)，結果顯示除PCA、GA及相關系數(shù)法外，其他特征優(yōu)選驗證結果都符合預期，可有效準確找到特征波段區(qū)間，其中在SIPLS 30 間隔情況下選取24、28間隔(共計138個特征波段點)時定量分析模型具有最優(yōu)驗證結果，rp為0.978 4，R2P為0.957 2，RMSEP為0.221 1。由此可見，SIPLS在間隔30以內(nèi)可有效提取出2 組變量組合最優(yōu)樣品解，但在40 間隔以上擬合時間較長，變量組合超過4，運算時間呈現(xiàn)指數(shù)級增長，無法應用于實驗室檢測。要完成在線檢測(大量樣本、快速分析)則有待與其他方法如深度學習等進行非線性組合優(yōu)化及融合。

本文構建了氨、堿處理后玉米秸稈粗蛋白含量的快速定量分析模型，發(fā)現(xiàn)基于SIPLS的特征優(yōu)選方法可有效、快速地測定粗蛋白含量，模型精度及誤差滿足要求，該方法可為氨堿化秸稈處理提供重要的快速判定方法和檢測手段。