小麥加工鏈中重金屬鎘含量的深度網(wǎng)絡(luò)預(yù)測模型

金學(xué)波，張佳帥，郭天洋，王小藝，蘇婷立，賴燕群，孔建磊，白玉廷

（1.北京工商大學(xué)人工智能學(xué)院，北京 100048；2.北京工商大學(xué) 中國輕工業(yè)工業(yè)互聯(lián)網(wǎng)與大數(shù)據(jù)重點實驗室，北京 100048；3.北京工商大學(xué)食品與健康學(xué)院，北京 100048；4.北京服裝學(xué)院，北京 100105；5.中糧糧油工業(yè)（荊州）有限公司，湖北 荊州 434300）

小麥?zhǔn)且环N在世界各地廣泛種植的谷類作物，在我國小麥的生長范圍廣泛，是人們喜愛的主食之一，也是牲畜的主要飼料用料之一。據(jù)《國家統(tǒng)計局關(guān)于2021年糧食產(chǎn)量數(shù)據(jù)的公告》，2021年全國小麥播種面積23 568.4千 hm，總產(chǎn)量為12 226萬 t。

在小麥的生產(chǎn)中，重金屬污染問題引起了人們的廣泛關(guān)注。若食用被重金屬污染的小麥，當(dāng)污染物生物累積性達(dá)到一定濃度后會引起慢性中毒甚至癌癥、生理畸形等病癥。如重金屬鎘在攝入后可在人體的肝、腎等器官組織中蓄積，造成各器官組織的損傷，當(dāng)達(dá)到一定濃度后會破壞骨骼中鈣的正常補給，嚴(yán)重者會造成骨軟化癥等惡性疾病。因此，通過預(yù)測方法掌握食品加工鏈重金屬污染發(fā)展規(guī)律并及時加以預(yù)防控制具有極大的實際應(yīng)用意義。

目前用于預(yù)測建模的方法主要是多元線性回歸方法。例如，牛耀星等采用一級動力學(xué)模型結(jié)合Arrhenius方程建立了基于金針菇子實體品質(zhì)指標(biāo)的貨架期預(yù)測模型，并對模型的預(yù)測精確度進(jìn)行了驗證及評價。李婕使用Huang模型、Baranyi模型、修正的Gompertz模型對雙孢蘑菇表面背景菌群和假單胞菌的生長情況進(jìn)行建模，模型適合描述溫度對雙孢蘑菇樣品表面假單胞菌生長速率的影響。宋波等采用逐步多元回歸分析建立不同鎘濃度等級下水稻籽粒重金屬鎘含量預(yù)測模型，模型能夠預(yù)測稻米中鎘的累積量，為廣西鎘異常區(qū)內(nèi)其他水稻產(chǎn)地的安全生產(chǎn)提供了參考。王燁鋒采用多元逐步線性回歸分析方法對土壤-水稻體系中鎘的積累情況建立預(yù)測模型，決定系數(shù)（）均達(dá)0.7以上。通過實測值反推模型，3 組模型擬合散點值均勻分布在＝線兩側(cè)，并通過采集田間樣品對模型進(jìn)行了驗證。也有一些學(xué)者進(jìn)行了神經(jīng)網(wǎng)絡(luò)建模的嘗試，取得了良好的效果。例如，于灝等針對小麥中鎘的含量分別建立了多元回歸模型與淺層神經(jīng)網(wǎng)絡(luò)模型，結(jié)果表明基于神經(jīng)網(wǎng)絡(luò)的小麥鎘含量預(yù)測模型的精度總體優(yōu)于多元線性回歸模型。

綜上可知，目前對于原糧產(chǎn)地小麥中鎘含量分析的研究較多，卻鮮有研究涉及小麥加工過程中鎘含量的變化。而在加工過程中，鎘含量的變化與小麥粉的品質(zhì)直接相關(guān)。本研究基于小麥原糧（即原麥）的鎘含量數(shù)據(jù)預(yù)測小麥加工成品——小麥粉中的鎘含量，以此為依據(jù)指導(dǎo)加工環(huán)節(jié)中的重金屬危害管控，從而避免小麥粉的重金屬含量超標(biāo)，保證食品安全。

另一方面，加工過程中的危害物含量具有強非線性、強隨機性的特征，因此線性回歸、淺層神經(jīng)網(wǎng)絡(luò)等模型建模精度不高，預(yù)測性能有待提升。將深度學(xué)習(xí)網(wǎng)絡(luò)應(yīng)用于強非線性、強隨機性數(shù)據(jù)預(yù)測領(lǐng)域的研究有很多。Haider等建立了一種基于長短期記憶（long and short term memory，LSTM）網(wǎng)絡(luò)的小麥產(chǎn)量預(yù)測模型。牛哲文等建立了基于門控循環(huán)單元（gated recurrent unit，GRU）的預(yù)測模型，使用具有波動性和不確定性的風(fēng)速、風(fēng)向等歷史天氣數(shù)據(jù)對風(fēng)力發(fā)電的功率進(jìn)行了預(yù)測。但是，由于這些深度網(wǎng)絡(luò)模型對噪聲非常敏感，并不適用于小麥加工過程中鎘含量預(yù)測模型的建立。

為了提升小麥加工過程中鎘含量預(yù)測的準(zhǔn)確度，本研究基于小麥加工鏈鎘含量數(shù)據(jù)進(jìn)行深度學(xué)習(xí)網(wǎng)絡(luò)建模，并使用正則化方法優(yōu)化其損失函數(shù)，通過加入噪聲懲罰項防止訓(xùn)練時模型對噪聲的過度擬合，并使用貝葉斯超參數(shù)優(yōu)化方法對模型超參數(shù)進(jìn)行優(yōu)化。

1 材料與方法

1.1 材料

小麥為某企業(yè)收購原糧。

1.2 設(shè)備

Y7000P便攜式電腦（配備16 GB運行內(nèi)存和RTX3050Ti顯卡） 聯(lián)想控股股份有限公司。

1.3 方法

1.3.1 小麥加工鏈及關(guān)鍵節(jié)點

小麥加工鏈?zhǔn)侵笍脑湹叫←湻鄣乃猩a(chǎn)環(huán)節(jié)，簡化示意圖如圖1所示。該加工鏈共有清理、潤麥、皮磨等10 個加工環(huán)節(jié)。原糧委托某公司按圖1進(jìn)行加工，并對各加工環(huán)節(jié)進(jìn)行抽檢。抽檢方法為：每隔2 h對全加工鏈抽檢1 次，共抽取16 次，每次在各環(huán)節(jié)分別抽取3～4 份（每份1 kg）平行樣品。為便于表示各環(huán)節(jié)的樣本，設(shè)定原糧、清理、潤麥、皮磨一（1M芯）、皮磨二（2M芯）、皮磨三（3M芯）、皮磨四（4M芯）、皮磨五（5M芯）、皮磨六（6M芯）、包裝、存放環(huán)節(jié)分別為～。

圖1 小麥加工鏈簡化示意圖Fig. 1 Flow chart of wheat processing chain

1.3.2 鎘含量的測定

取1.3.1節(jié)樣品，委托某公司參照GB/T 20380.4—2006《淀粉及其制品 重金屬含量 第4部分：電熱原子吸收光譜法測定鎘含量》測定鎘含量（以鎘元素計），單位為mg/kg。

1.4 數(shù)據(jù)處理與分析

1.4.1 鎘含量數(shù)據(jù)集構(gòu)建

取樣時設(shè)置3～4 個平行，各環(huán)節(jié)每次取樣的鎘含量結(jié)果以平行樣品鎘含量的平均值表示，共獲得16 條全加工鏈的鎘含量數(shù)據(jù)。由于這些數(shù)據(jù)數(shù)量不足以訓(xùn)練神經(jīng)網(wǎng)絡(luò)，利用等差插值的方法對數(shù)據(jù)進(jìn)行增廣處理。

1.4.2 基于正則化損失函數(shù)的GRU深度網(wǎng)絡(luò)模型分析

首先用門控循環(huán)單元建立深度學(xué)習(xí)模型，其次，設(shè)計基于正則化方法優(yōu)化的損失函數(shù)，該損失函數(shù)在模型原損失函數(shù)的基礎(chǔ)上添加訓(xùn)練數(shù)據(jù)中噪聲的懲罰項，達(dá)到在訓(xùn)練過程中減小對噪聲擬合程度的目的。因此訓(xùn)練過程中損失函數(shù)會在一個穩(wěn)定值附近波動，該穩(wěn)定值用來量化數(shù)據(jù)中的噪聲。

GRU在時序預(yù)測任務(wù)中有很好的效果，此前已經(jīng)有大量的研究中都使用GRU作為預(yù)測模型并且證明了其性能。GRU模型中只有兩個門，分別是更新門和重置門。更新門控制前一時刻的狀態(tài)信息帶入到當(dāng)前狀態(tài)的程度，該門的數(shù)值越大就說明前一時刻的信息帶入越多。重置門用于控制前一時刻狀態(tài)的遺忘程度，數(shù)據(jù)越大說明忽略越少。

GRU模型結(jié)構(gòu)如圖2所示，根據(jù)該結(jié)構(gòu)進(jìn)行前向傳播時每個單元中參數(shù)計算，如式（1）～（4）所示。

式中：x為時刻的輸入；z為更新門衰減系數(shù)；r為重置門衰減系數(shù)；h是-1時刻的狀態(tài)向量；h為時刻狀態(tài)向量；W、U為更新門權(quán)重；W、U為重置門權(quán)重；W、U為待選a的權(quán)重；b、b、b為偏置向量；為Sigmoid激活函數(shù)；*表示逐元素乘法；tanh為tanh激活函數(shù)。

圖2 GRU模型結(jié)構(gòu)圖Fig. 2 Structure diagram of GRU model

但GRU模型在訓(xùn)練時易對鎘含量數(shù)據(jù)中含有的非線性噪聲過度擬合，降低模型的泛化能力，從而影響模型的預(yù)測性能，因此需要采取措施防止模型過度擬合。正則化的作用在于限制模型中的參數(shù)，讓模型的參數(shù)不會過大，從而減少模型過度擬合的可能。在此步驟中，為減少模型因過度學(xué)習(xí)加工鏈數(shù)據(jù)中的隨機噪聲而導(dǎo)致模型過度擬合，設(shè)計噪聲平滑損失函數(shù)，如式（5）所示。

噪聲平滑損失函數(shù)的計算公式中包括兩個部分，第一部分衡量數(shù)據(jù)擬合程度，表示將實際序列與擬合序列之間的殘差平方和最小化的目標(biāo)；第二部分衡量數(shù)據(jù)平滑程度，表示訓(xùn)練過程中對序列平滑度的需求。其中，擬合程度用平均絕對誤差表示，平滑程度通過計算P的范數(shù)來實現(xiàn)。是正則化項，將其視為在擬合與平滑兩個目標(biāo)之間的權(quán)衡取舍。

P的計算是實現(xiàn)噪聲平滑損失函數(shù)的關(guān)鍵，計算時首先需要定義矩陣，該矩陣對數(shù)據(jù)進(jìn)行平滑度懲罰。矩陣如式（6）所示。

矩陣的維度由輸入數(shù)據(jù)決定，假設(shè)輸入數(shù)據(jù)長度為，則的維度為（，），然后通過與輸入數(shù)據(jù)組成矩陣的乘法計算得到矩陣P，最后通過矩陣范數(shù)的計算得到可以表示數(shù)據(jù)每3 個點間平滑程度的范數(shù)。

在訓(xùn)練過程中，本實驗所設(shè)計的損失函數(shù)可以有效地考慮噪聲對模型的影響從而達(dá)到減少過度擬合的作用。當(dāng)公式（5）中為0時，網(wǎng)絡(luò)學(xué)習(xí)原始數(shù)據(jù)，隨著不斷增大，擬合線趨向平滑。所以，本實驗所設(shè)計的損失函數(shù)可以使模型對噪聲具有更強的魯棒性，有效減小噪聲在訓(xùn)練過程中的影響，從而提高模型的預(yù)測準(zhǔn)確度。

1.4.3 應(yīng)用網(wǎng)絡(luò)的貝葉斯優(yōu)化方法分析

深度學(xué)習(xí)模型的超參數(shù)選擇直接決定了模型的性能。本研究基于序列模型的優(yōu)化（sequence model based optimization，SMBO）方法進(jìn)行模型超參數(shù)優(yōu)化，優(yōu)化的超參數(shù)主要包括GRU中的神經(jīng)元數(shù)、Dropout率、訓(xùn)練次數(shù)、批處理大小、優(yōu)化器。

貝葉斯優(yōu)化方法在確定模型參數(shù)時，使用代理模型來擬合真實的目標(biāo)函數(shù)，并根據(jù)擬合結(jié)果主動選擇最“潛在”的評估點。所以需要定義一個目標(biāo)函數(shù)（）和優(yōu)化的超參數(shù)空間。目標(biāo)函數(shù)（式（7））表示通過貝葉斯優(yōu)化需要達(dá)到的最小化的目標(biāo)，通過設(shè)定模型的目標(biāo)函數(shù)，找到此度量上產(chǎn)生最佳得分的模型超參數(shù)。

貝葉斯優(yōu)化的代理函數(shù)如式（8）所示。

式中：*為貝葉斯優(yōu)化確定的最優(yōu)參數(shù)；為輸入的一組超參數(shù)；為多維超參數(shù)的參數(shù)空間；（）為目標(biāo)函數(shù)。

貝葉斯優(yōu)化主要由兩個步驟組成：首先通過第+1步估計和更新高斯過程，然后通過最大化代理函數(shù)來指導(dǎo)超參數(shù)的采樣。在高斯過程中，設(shè)置目標(biāo)函數(shù)（）服從高斯分布（式（9））。

式中：（）為目標(biāo)函數(shù)；（）為目標(biāo)函數(shù)（）的均值；（，’）為目標(biāo)函數(shù)（）的協(xié)方差矩陣；表示高斯分布。初始的（，’）如式（10）所示。

式中：為超參數(shù)的數(shù)量；表示協(xié)方差矩陣；（，）、（，w）、（w，）、（w，w）為矩陣內(nèi)的不同超參數(shù)的組合。

在進(jìn)行貝葉斯優(yōu)化時，高斯過程的協(xié)方差矩陣會隨著迭代過程改變，假設(shè)在第+1步輸入的一組參數(shù)為w，則此時協(xié)方差矩陣如式（11）所示。

式中：’為+1時刻的協(xié)方差矩陣；為初始的協(xié)方差矩陣；w為+1時刻的超參數(shù)值；為超參數(shù)的組合矩陣，為的轉(zhuǎn)置，的組合公式為＝[（w，），（w，），…，（w，w）]，此時可以得到目標(biāo)函數(shù)的后驗概率（式（12））。

得到后驗概率后，通過超參數(shù)搜索方法來尋找最優(yōu)超參數(shù)，本實驗使用上置信界（upper confidence bound，UCB）采集函數(shù)（式（13））來完成超參數(shù)搜索。

1.4.4 預(yù)測效果評價標(biāo)準(zhǔn)與模型對比驗證

以80%的鎘含量數(shù)據(jù)用作訓(xùn)練集，對模型的性能進(jìn)行訓(xùn)練優(yōu)化；以20%的鎘含量數(shù)據(jù)數(shù)據(jù)作為驗證集，驗證模型的性能。采用均方根誤差（root mean square error，RMSE）、平均絕對誤差（mean absolute error，MAE）、皮爾遜相關(guān)系數(shù)評估預(yù)測的性能，計算分別如式（14）～（16）所示。三者均反映模型預(yù)測值與真實值之間的關(guān)系。RMSE、MAE越小，說明模型預(yù)測結(jié)果與真實情況越接近。越大，說明模型預(yù)測結(jié)果與真實情況越接近。

采用線性回歸、循環(huán)神經(jīng)網(wǎng)絡(luò)（recurrent neural network，RNN）、LSTM和GRU模型與本研究所得模型進(jìn)行對比。通過比較各模型基于小麥加工鏈鎘含量數(shù)據(jù)的預(yù)測效果評價指標(biāo)（RMSE、MAE、）驗證模型的預(yù)測性能，驗證本研究所得模型在實際小麥加工鏈鎘含量預(yù)測中的可靠性。

2 結(jié)果與分析

2.1 小麥加工鏈中小麥的鎘含量數(shù)據(jù)集構(gòu)建

小麥加工鏈各環(huán)節(jié)樣品的鎘含量數(shù)據(jù)經(jīng)增廣處理后最終獲得2 057 條全加工鏈的鎘含量數(shù)據(jù)。以其中3 條全加工鏈的鎘含量變化為例，如圖3所示，在清理環(huán)節(jié)、潤麥環(huán)節(jié)中鎘含量略有增加，之后各環(huán)節(jié)鎘含量的變化相對平穩(wěn)。

圖3 各加工環(huán)節(jié)不同樣本中鎘含量Fig. 3 Cadmium contents in different samples measured during processing

2.2 不同模型預(yù)測的性能對比

基于2.1節(jié)小麥加工鏈10 個環(huán)節(jié)中2 057 條全加工鏈的鎘含量數(shù)據(jù)進(jìn)行實驗。將本研究所得模型與線性回歸、RNN、LSTM、GRU模型進(jìn)行對比，采用RMSE、MAE、評估模型的性能，結(jié)果如表1所示。本研究所得模型的RMSE（0.259 5）最小、MAE（0.161 6）最小且最大，比GRU模型的預(yù)測性能更好（RMSE降低了46.37%）。

表1 不同預(yù)測方法基于鎘數(shù)據(jù)的性能對比Table 1 Performance comparison of different prediction methods based on cadmium data

通過驗證集的驗證與對比分析可知，本研究所得模型在預(yù)測任務(wù)中表現(xiàn)更為出色，擁有更強的含噪時序數(shù)據(jù)分析能力。由此可見，本研究所得模型可以有效地分析含噪時間序列數(shù)據(jù)，對于小麥加工鏈危害物等含噪數(shù)據(jù)有良好的預(yù)測性能。

2.3 模型預(yù)測結(jié)果

基于2.1節(jié)中20%的原糧鎘含量數(shù)據(jù)，對后續(xù)清理、潤麥、皮磨等10 個環(huán)節(jié)鎘含量進(jìn)行預(yù)測，其中6 個樣本預(yù)測結(jié)果如圖4所示。其中，清理環(huán)節(jié)和潤麥環(huán)節(jié)樣品的鎘含量略有上升，尤其是清理環(huán)節(jié)，由于鎘元素在小麥表皮的富集，這兩個環(huán)節(jié)的浸潤對鎘元素有一定的析出作用。皮磨一及皮磨二環(huán)節(jié)后鎘含量持續(xù)下降，說明皮磨對危害物鎘有較好的去除作用。預(yù)測結(jié)果還表明，～環(huán)節(jié)鎘含量基本趨于平穩(wěn)。

GB 2762—2017《食品安全國家標(biāo)準(zhǔn) 食品中污染物限量》規(guī)定小麥原糧及成品糧的鎘含量不得高于0.1 mg/kg。由本研究所得模型可知，當(dāng)小麥原糧中鎘含量低于GB 2762—2017規(guī)定的0.1 mg/kg時，加工鏈中清理環(huán)節(jié)和潤麥環(huán)節(jié)的鎘含量可能會使其暫時超標(biāo)，但最終成品小麥粉符合GB 2762—2017要求；而當(dāng)原麥中鎘含量高于0.1 mg/kg時，成品小麥粉中鎘含量將高于0.1 mg/kg。因此考慮檢測數(shù)據(jù)偏差的情況下，為確保小麥粉的食用安全，原糧中危害物鎘含量高于0.1 mg/kg時不能進(jìn)行后續(xù)加工流程。

圖4 模型預(yù)測小麥加工鏈中樣品鎘含量的變化規(guī)律Fig. 4 Trend of cadmium content in wheat processing chain predicted by the model

3 結(jié) 論

本實驗基于小麥加工鏈各個環(huán)節(jié)小麥中的鎘含量。采用深度學(xué)習(xí)的建模方法對加工鏈各環(huán)節(jié)樣品中的鎘含量進(jìn)行建模預(yù)測。所建模型首先使用正則化方法修改深度學(xué)習(xí)模型的損失函數(shù)，通過加入噪聲懲罰項來減小噪聲項對模型預(yù)測性能的影響。同時使用貝葉斯優(yōu)化算法優(yōu)化模型超參數(shù)。經(jīng)過對比驗證，本研究所建模型對于實際的小麥加工鏈鎘含量數(shù)據(jù)預(yù)測準(zhǔn)確度較高。

應(yīng)用該預(yù)測模型可以準(zhǔn)確地預(yù)測出加工鏈中從原麥到小麥粉的各個關(guān)鍵環(huán)節(jié)鎘含量的變化，能夠?qū)π←溂庸ゆ湹逆k危害物防控實踐進(jìn)行指導(dǎo)。預(yù)測結(jié)果表明，通過加工，原麥中的鎘含量呈整體下降趨勢。但是，在清理和潤麥的環(huán)節(jié)，鎘含量會略有上升。由預(yù)測模型可知，如果原麥的鎘含量滿足GB 2762—2017要求，即小麥中的鎘含量不高于0.1 mg/kg，則經(jīng)過加工的成品小麥粉基本滿足GB 2762—2017的要求。