謝秋菊 馬超凡 王圣超 包 軍 劉洪貴 于海明

(1.東北農(nóng)業(yè)大學(xué)電氣與信息學(xué)院, 哈爾濱 150030; 2.東北農(nóng)業(yè)大學(xué)動(dòng)物科技學(xué)院, 哈爾濱 150030;3.農(nóng)業(yè)農(nóng)村部生豬養(yǎng)殖設(shè)施工程重點(diǎn)實(shí)驗(yàn)室, 哈爾濱 150030;4.教育部北方寒區(qū)智能化繁育與養(yǎng)殖工程研究中心, 哈爾濱 150030)

0 引言

豬舍空氣質(zhì)量是影響生豬呼吸道健康的關(guān)鍵因素,是豬舍環(huán)境控制關(guān)注的熱點(diǎn)之一。在生豬養(yǎng)殖過程中,豬呼吸產(chǎn)生的CO2、豬糞尿等排泄物和飼料殘?jiān)纸猱a(chǎn)生的NH3、H2S、CH4等有害氣體對(duì)豬的生長、繁殖、健康造成極大危害。其中NH3與CO2的氣體濃度備受關(guān)注。因?yàn)镹H3是一種無色具有刺激氣味的有害氣體,不同濃度的NH3會(huì)對(duì)豬的呼吸系統(tǒng)、黏膜等造成不同程度的損傷,并且NH3排放還會(huì)造成周圍環(huán)境污染[1];豬舍中的CO2濃度是衡量豬舍內(nèi)空氣質(zhì)量狀況的重要標(biāo)志,它主要來自于豬的呼吸,舍內(nèi)CO2濃度過高不僅會(huì)造成豬頭昏和生產(chǎn)性能下降,更有可能對(duì)豬造成生命威脅[2]。因此,對(duì)豬舍中NH3和CO2這兩種氣體濃度進(jìn)行監(jiān)測與控制,是實(shí)現(xiàn)豬舍環(huán)境控制的重要任務(wù)。通過大量的現(xiàn)場實(shí)驗(yàn)得出,在豬舍內(nèi)進(jìn)行一系列控制,例如開啟風(fēng)機(jī)1 h,舍內(nèi)的溫濕度、氣體濃度會(huì)有較明顯的變化。通過對(duì)舍內(nèi)氣體濃度變化進(jìn)行預(yù)測,可以實(shí)現(xiàn)豬舍環(huán)境控制設(shè)備運(yùn)行狀態(tài)的提前調(diào)控,從而彌補(bǔ)在豬舍環(huán)境實(shí)際控制過程中的環(huán)境控制效果滯后的問題。

目前,畜禽舍內(nèi)有害氣體濃度變化模型主要有兩種,一種是基于理化統(tǒng)計(jì)的模型,這類模型通過對(duì)氣體的產(chǎn)生機(jī)理或者不同影響因素進(jìn)行分析構(gòu)建預(yù)測模型,例如:文獻(xiàn)[3]在多種不同空氣溫度和空氣濕度組合下連續(xù)測量多天糞便的NH3排放構(gòu)建雞糞NH3排放模型;文獻(xiàn)[4]通過對(duì)多種不同地板類型的豬舍的NH3排放進(jìn)行檢測和分析,得到不同地板類型豬舍的NH3排放系數(shù)。這些模型測量方法較復(fù)雜且受實(shí)際畜禽舍結(jié)構(gòu)影響大,實(shí)用性較差。另一種是基于人工智能算法的預(yù)測模型,例如,基于支持向量回歸機(jī)(Support vctor regression,SVR)、隨機(jī)森林(Random forest,RF)、極限隨機(jī)樹(Extra tree,ET)等傳統(tǒng)機(jī)器學(xué)習(xí)算法構(gòu)建的預(yù)測模型[5-7],基于卷積神經(jīng)網(wǎng)絡(luò)(Convolutional neural networks,CNN)、GRU等深度學(xué)習(xí)算法構(gòu)建的預(yù)測模型[2,8-13],將自回歸模型(Autoregressive model,AR)、ARIMA等統(tǒng)計(jì)學(xué)習(xí)算法與機(jī)器學(xué)習(xí)算法互相結(jié)合的混合預(yù)測模型[5,14-15]。與傳統(tǒng)理化統(tǒng)計(jì)模型相比,基于人工智能算法的模型具有預(yù)測精度高,可更好地提取氣體濃度變化特征的優(yōu)點(diǎn),已成為畜禽舍氣體環(huán)境控制領(lǐng)域研究的熱點(diǎn)。

畜禽舍內(nèi)有害氣體濃度隨時(shí)間不斷變化,具有典型的時(shí)序序列數(shù)據(jù)的特征[9]。在眾多的智能預(yù)測算法中,GRU模型不僅擅長提取時(shí)間序列中非線性長期依賴關(guān)系且訓(xùn)練速度較快。因此,許多學(xué)者關(guān)注GRU模型對(duì)時(shí)間序列數(shù)據(jù)的預(yù)測[16-17],同時(shí)探索將CNN、注意力機(jī)制與GRU模型相融合的混合網(wǎng)絡(luò)模型[9],以提高時(shí)序序列數(shù)據(jù)的關(guān)鍵特征及信息提取能力。在時(shí)序序列數(shù)據(jù)預(yù)測神經(jīng)網(wǎng)絡(luò)模型中,超參數(shù)是影響模型性能的重要因素。然而,依據(jù)傳統(tǒng)經(jīng)驗(yàn)選取超參數(shù)往往較難取得良好效果。因此,學(xué)者們開始嘗試使用麻雀搜索算法(Sparrow search algorithm,SSA)、粒子群優(yōu)化算法(Particle swarm optimization,PSO)、灰狼優(yōu)化算法(Grey wolf optimization,GWO)等優(yōu)化算法進(jìn)行神經(jīng)網(wǎng)絡(luò)超參數(shù)尋優(yōu)[2,18-24],來提升網(wǎng)絡(luò)模型預(yù)測性能。研究表明,通過模擬麻雀種群覓食和躲避天敵的行為而提出的麻雀搜索算法[21],其搜索精度和收斂速度均優(yōu)于PSO、GWO和遺傳算法(Genetic algorithm,GA),但是SSA仍然存在容易陷入局部最優(yōu)的問題[23],進(jìn)而導(dǎo)致模型預(yù)測精度低。

為了解決當(dāng)前有害氣體濃度預(yù)測模型存在的預(yù)測精度低的問題,本文構(gòu)建基于ISSA-GRU-ARIMA的豬舍有害氣體濃度預(yù)測模型,該模型使用改進(jìn)麻雀搜索算法對(duì)GRU模型超參數(shù)進(jìn)行尋優(yōu),以提升GRU模型非線性擬合能力,同時(shí)融合ARIMA模型對(duì)時(shí)間序列模型的殘差序列信息進(jìn)一步提取以提高模型預(yù)測精度。

1 材料與方法

1.1 數(shù)據(jù)來源

1.1.1豬舍結(jié)構(gòu)

本文數(shù)據(jù)來源于黑龍江省哈爾濱市尚志市亞布力驚哲養(yǎng)殖基地(44.78°N,128.47°E)的一個(gè)妊娠母豬舍,豬舍地面是水泥漏糞地板。豬舍尺寸為19.5 m×10.5 m×2.4 m,限位欄尺寸為2.2 m×0.65 m×1 m。豬舍采用機(jī)械通風(fēng),一側(cè)墻上安裝了3個(gè)風(fēng)機(jī),風(fēng)機(jī)1(fan1)和風(fēng)機(jī)2(fan2)尺寸均為0.8 m×0.8 m,額定功率為0.37 kW,風(fēng)機(jī)3(fan3)尺寸為1 m×1 m,額定功率為0.55 kW。舍內(nèi)飼養(yǎng)有37頭妊娠母豬,豬平均質(zhì)量約為165 kg,豬舍現(xiàn)場如圖1所示,平面結(jié)構(gòu)如圖2所示。

圖1 豬舍現(xiàn)場圖

圖2 豬舍平面結(jié)構(gòu)圖

1.1.2豬舍環(huán)境監(jiān)測

豬舍內(nèi)安裝了基于物聯(lián)網(wǎng)系統(tǒng)的自動(dòng)環(huán)境監(jiān)測及控制器,該控制器以STM32F103VET6芯片為核心。在豬舍內(nèi)設(shè)有3個(gè)環(huán)境監(jiān)測點(diǎn),每個(gè)監(jiān)測點(diǎn)分別安裝有溫濕度、NH3質(zhì)量濃度、CO2質(zhì)量濃度傳感器,傳感器距離地面高度為1.6 m,將3個(gè)環(huán)境監(jiān)測點(diǎn)監(jiān)測數(shù)據(jù)的平均值作為舍內(nèi)溫濕度、NH3質(zhì)量濃度、CO2質(zhì)量濃度取值。舍外布置1個(gè)溫濕度傳感器對(duì)舍外溫濕度進(jìn)行實(shí)時(shí)監(jiān)測,3個(gè)壓差傳感器分別布置在3個(gè)風(fēng)機(jī)正上方,對(duì)墻上風(fēng)機(jī)承受的壓力差進(jìn)行監(jiān)測,環(huán)境數(shù)據(jù)采樣時(shí)間間隔為1 min,采集到的環(huán)境數(shù)據(jù)通過Wi-Fi模塊上傳至物聯(lián)網(wǎng)云平臺(tái),傳感器輸出信號(hào)類型均為RS-485,通信協(xié)議為Modbus,傳感器具體參數(shù)如表1所示。通風(fēng)速率通過持續(xù)監(jiān)測風(fēng)機(jī)運(yùn)行時(shí)間、風(fēng)機(jī)承受的壓差變化以及現(xiàn)場實(shí)驗(yàn)數(shù)據(jù)擬合線性方程[25]計(jì)算得到,公式為

表1 傳感器技術(shù)參數(shù)

V=(0.534 2P1+112.9)e1+(0.521 7P2+106.5)e2+

(0.626 8P3+213.6)e3

(1)

式中Pi——風(fēng)機(jī)i所承受的壓差,Pa,i=1,2,3

ei——風(fēng)機(jī)i在1 min內(nèi)運(yùn)行的時(shí)間占比,%

V——豬舍風(fēng)機(jī)通風(fēng)速率,m3/min

豬舍投喂飼料和清理糞便均采用人工的方式,其中,喂食時(shí)間為每日08:00—08:30和14:00—14:30,清糞時(shí)間為每日08:30—09:00和14:30—15:00。

1.2 模型開發(fā)環(huán)境

本文模型開發(fā)硬件環(huán)境為Intel core i7 12700H處理器,Nvidia RTX3070顯卡,Windows 11 64位操作系統(tǒng)。模型算法開發(fā)采用Python 3.7編程語言結(jié)合Tensorflow框架和sklearn框架。

1.3 氣體濃度預(yù)測優(yōu)化模型構(gòu)建

本文構(gòu)建以GRU模型為核心的時(shí)序數(shù)據(jù)預(yù)測模型,提取豬舍環(huán)境數(shù)據(jù)的非線性特征,然后建立改進(jìn)的麻雀搜索算法對(duì)GRU模型的超參數(shù)進(jìn)行優(yōu)化,提升GRU模型的非線性特征擬合能力,再利用ARIMA模型的線性數(shù)據(jù)擬合能力來提取殘差數(shù)據(jù)的線性特征。建立的ISSA-GRU-ARIMA模型結(jié)合了ISSA的尋優(yōu)能力、GRU的非線性處理能力和ARIMA的線性處理能力,提升了模型的預(yù)測精度。將舍內(nèi)溫度、舍內(nèi)相對(duì)濕度、舍外溫度、舍外相對(duì)濕度、CO2質(zhì)量濃度、NH3質(zhì)量濃度、通風(fēng)速率的數(shù)據(jù)序列作為ISSA-GRU-ARIMA模型輸入,模型輸出為2個(gè)值,分別是1 h后的CO2質(zhì)量濃度和NH3質(zhì)量濃度。優(yōu)化模型構(gòu)建過程如下:

(1)將預(yù)處理后的環(huán)境數(shù)據(jù)劃分為訓(xùn)練集1、訓(xùn)練集2和測試集3部分,在訓(xùn)練集1上對(duì)GRU模型進(jìn)行訓(xùn)練,然后通過ISSA算法對(duì)模型超參數(shù)進(jìn)行尋優(yōu)得到最優(yōu)GRU模型,命名為ISSA-GRU模型。

(2)使用ISSA-GRU模型對(duì)訓(xùn)練集2上的氣體濃度進(jìn)行預(yù)測,分別得到NH3質(zhì)量濃度、CO2質(zhì)量濃度的觀測值與預(yù)測值之間的殘差序列DNH3和DCO2。

(3)以殘差序列DNH3和DCO2作為數(shù)據(jù)集構(gòu)建ARIMA模型,ARIMA模型可以提取DNH3和DCO2的變化規(guī)律。

(4)構(gòu)建ISSA-GRU-ARIMA模型對(duì)未來時(shí)刻的有害氣體濃度進(jìn)行預(yù)測,在測試集上對(duì)ISSA-GRU-ARIMA模型預(yù)測精度進(jìn)行測試。

1.3.1數(shù)據(jù)預(yù)處理及數(shù)據(jù)集劃分

選取從2022年8月10日到9月30日共52 d的環(huán)境數(shù)據(jù)(舍內(nèi)溫度、舍內(nèi)相對(duì)濕度、舍外溫度、舍外相對(duì)濕度、CO2質(zhì)量濃度、NH3質(zhì)量濃度、通風(fēng)速率)作為實(shí)驗(yàn)數(shù)據(jù)。對(duì)于傳感器采集的極個(gè)別時(shí)刻的缺失值采用線性插值法進(jìn)行補(bǔ)全;然后采取拉伊達(dá)準(zhǔn)則對(duì)數(shù)據(jù)進(jìn)行處理剔除其中的異常值,異常值被替換為相鄰數(shù)據(jù)的平均值。

1 h內(nèi)環(huán)境因素變化較小,因此將時(shí)間間隔為1 min的數(shù)據(jù)平均化處理為時(shí)間間隔為60 min的數(shù)據(jù),共得到1 248組數(shù)據(jù)。為了消除量綱的差異,將數(shù)據(jù)輸入模型前需要進(jìn)行歸一化,模型預(yù)測結(jié)果需要進(jìn)行反歸一化。

通過滑動(dòng)窗口構(gòu)建數(shù)據(jù)集,使用前3 h的數(shù)據(jù)序列對(duì)后1 h的有害氣體質(zhì)量濃度進(jìn)行預(yù)測。將實(shí)驗(yàn)數(shù)據(jù)劃分為訓(xùn)練集1、訓(xùn)練集2和測試集3部分,各部分占比分別是70%、20.4%和9.6%。訓(xùn)練集1用于ISSA-GRU模型訓(xùn)練,訓(xùn)練集2用于ISSA-GRU-ARIMA模型擬合,最后使用測試集對(duì)ISSA-GRU-ARIAM模型進(jìn)行測試。

1.3.2氣體質(zhì)量濃度預(yù)測GRU非線性模型

豬舍環(huán)境具有時(shí)序性、非線性、相互耦合的特性[8],尤其是舍內(nèi)NH3質(zhì)量濃度與CO2質(zhì)量濃度受多因素條件影響,因此需要構(gòu)建具有較強(qiáng)的非線性擬合能力和最大限度地提取時(shí)序數(shù)據(jù)特征的GRU神經(jīng)網(wǎng)絡(luò)模型來預(yù)測氣體質(zhì)量濃度變化。本文構(gòu)建的GRU神經(jīng)網(wǎng)絡(luò)模型由輸入層、隱藏層和全連接層3部分組成,其結(jié)構(gòu)如圖3所示。

圖3 GRU神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)

輸入層接收到舍內(nèi)溫度、舍內(nèi)相對(duì)濕度、舍外溫度、舍外相對(duì)濕度、CO2質(zhì)量濃度、NH3質(zhì)量濃度、通風(fēng)速率時(shí)序數(shù)據(jù),然后按照公式

(2)

式中zi——輸入的環(huán)境數(shù)據(jù)

z′i——?dú)w一化后的環(huán)境數(shù)據(jù)

zmax——輸入環(huán)境數(shù)據(jù)最大值

zmin——輸入環(huán)境數(shù)據(jù)最小值

進(jìn)行歸一化,歸一化后的數(shù)據(jù)直接傳入隱藏層。

隱藏層由2個(gè)GRU子層組成。兩個(gè)GRU子層的神經(jīng)元數(shù)量分別為486和315,GRU神經(jīng)網(wǎng)絡(luò)單元使用更新門和重置門兩個(gè)門來實(shí)現(xiàn)信息的記憶和傳遞[16],如圖4所示。更新門zt決定當(dāng)前節(jié)點(diǎn)的新增歷史信息比例,從而可以捕捉序列數(shù)據(jù)的短期依賴關(guān)系;重置門rt決定拋棄歷史信息的比例,可以更好地提取序列數(shù)據(jù)的長期依賴關(guān)系。時(shí)序序列多環(huán)境因子之間依賴關(guān)系的計(jì)算公式為

圖4 GRU神經(jīng)單元結(jié)構(gòu)

rt=σ(Wr[ht-1,at]+br)

(3)

zt=σ(Wz[ht-1,at]+bz)

(4)

(5)

(6)

式中at——當(dāng)前單元輸入值

ht-1——上一神經(jīng)單元狀態(tài)

ht——輸出狀態(tài)

Wr、Wz、Wh——重置門、更新門、當(dāng)前隱藏狀態(tài)的權(quán)值

br、bz、bh——重置門、更新門、當(dāng)前隱藏狀態(tài)的偏置值

σ——sigmod激活函數(shù)

全連接層每一個(gè)結(jié)點(diǎn)都與GRU子層的所有結(jié)點(diǎn)相連,通過GRU子層學(xué)習(xí)到的特征傳入全連接層。全連接層經(jīng)過非線性變換提取數(shù)據(jù)之間的依賴關(guān)系并輸出預(yù)測值,最后按照式(2)對(duì)輸出進(jìn)行反歸一化,全連接層的輸出為NH3質(zhì)量濃度預(yù)測值和CO2質(zhì)量濃度預(yù)測值。

GRU神經(jīng)網(wǎng)絡(luò)使用均方誤差(MSE)作為損失函數(shù),通過MSE計(jì)算誤差后反向傳播,然后通過Adam優(yōu)化器根據(jù)誤差不斷優(yōu)化權(quán)重,最終得到預(yù)測精度較高的模型。

1.3.3改進(jìn)麻雀搜索算法

由于GRU預(yù)測模型的超參數(shù)眾多,這些參數(shù)對(duì)網(wǎng)絡(luò)模型的性能影響較大,因此需要對(duì)GRU網(wǎng)絡(luò)模型超參數(shù)進(jìn)行優(yōu)化,以獲得最佳的模型性能。SSA是一種模擬麻雀覓食和躲避天敵行為而提出的群智能優(yōu)化算法[21],具有較好的全局探索能力和較快的收斂速度[22],有助于GRU預(yù)測模型對(duì)超參數(shù)進(jìn)行尋優(yōu),但是SSA存在容易陷入局部最小值的問題。

混沌序列具有規(guī)律性、遍歷性和隨機(jī)性,將其用于搜索算法可以有效提升尋優(yōu)效果,在混沌序列中,Tent序列相比常用的Logistic序列分布更加均勻且收斂速度更快;Tent混沌序列可以用于種群初始化以提高種群多樣性[26],還可以對(duì)局部極值進(jìn)行擾動(dòng)以擴(kuò)大搜索空間,從而有利于算法跳出局部最優(yōu)[24];高斯分布在期望附近概率最大,引入高斯變異可以對(duì)原始解周圍進(jìn)行重點(diǎn)搜索。因此,本文使用混沌序列、高斯變異和混沌擾動(dòng)對(duì)麻雀搜索算法進(jìn)行改進(jìn),來提升初始解分布的均勻性及局部尋優(yōu)效果,得到改進(jìn)的麻雀搜索算法。在ISSA中,一只麻雀代表一組解,包含時(shí)間窗口、學(xué)習(xí)率、批數(shù)目、第1層GRU子層神經(jīng)元個(gè)數(shù)和第2層GRU子層神經(jīng)元個(gè)數(shù)這些需要搜索的超參數(shù),麻雀個(gè)數(shù)代表解的個(gè)數(shù)。ISSA算法步驟如下:

(1)使用公式

(7)

式中O——粒子數(shù)量

rand(0,1)——[0,1]范圍內(nèi)的隨機(jī)數(shù)

產(chǎn)生Tent混沌序列S對(duì)種群進(jìn)行初始化,然后使用公式

(8)

umin——解序列第u維變量的最小值

Su——混沌序列的第u維變量

umax——解序列第u維變量的最大值

將其映射到原始解空間,產(chǎn)生的解維度為5,個(gè)數(shù)為10。

(2)計(jì)算種群中所有個(gè)體的適應(yīng)度fi,并對(duì)其進(jìn)行排序,適應(yīng)度是評(píng)判解優(yōu)劣的指標(biāo),本實(shí)驗(yàn)適應(yīng)度函數(shù)設(shè)置為均方誤差,適應(yīng)度越低表明優(yōu)化算法搜索到的解越優(yōu)秀。

(3)根據(jù)適應(yīng)度選出發(fā)現(xiàn)者和捕食者,然后分別按照公式

(9)

其中α∈(0,1]ST∈[0.5,1]R2∈[0,1]

式中t——迭代次數(shù)

itermax——最大迭代次數(shù)

α——隨機(jī)數(shù)

Xi,j——第i個(gè)解在第j維的位置信息

ST——預(yù)警值

R2——安全值,當(dāng)R2

Q——符合正態(tài)分布的隨機(jī)數(shù)

L——元素都為1的1×u的矩陣

(10)

其中

C+=CT(CCT)-1

式中XP——當(dāng)前發(fā)現(xiàn)者最優(yōu)位置

Xworst——全局最差位置

對(duì)其位置進(jìn)行更新。

(4)隨機(jī)選取種群中一定比例的個(gè)體作為警戒者,并按照公式

(11)

式中Xg——全局最優(yōu)位置

fg——全局最佳適應(yīng)度,fi>fg表示麻雀處于邊緣位置且容易被天敵發(fā)現(xiàn),fi=fg表示警戒者發(fā)現(xiàn)了敵人,需要前往安全位置

β——服從均值為0、方差為1的正態(tài)分布的隨機(jī)數(shù), 用于控制步長

K——隨機(jī)數(shù),表示麻雀移動(dòng)的方向,取(-1,1)

fworst——全局最差適應(yīng)度

ε——無限接近0的常數(shù),避免分母出現(xiàn)0

對(duì)其進(jìn)行更新。

(5)在進(jìn)行一輪迭代后,輸出當(dāng)前迭代的最小適應(yīng)度,以及對(duì)應(yīng)的GRU模型超參數(shù)尋優(yōu)值,然后重新計(jì)算所有解的適應(yīng)度fi和平均適應(yīng)度favg。

(6)當(dāng)fi

g=x(1+N(0,1))

(12)

式中x——需要進(jìn)行變異的解

N——隨機(jī)數(shù),均值為0,方差為1

g——高斯變異后的解

進(jìn)行高斯變異對(duì)此解周圍區(qū)域進(jìn)行重點(diǎn)搜索,然后將變異后的解與原始解進(jìn)行比較,兩者之中保留最優(yōu)的解。當(dāng)fi≥favg時(shí),使用公式

X′new=(X′+Xnew)/2

(13)

式中X′new——Tent混沌擾動(dòng)后產(chǎn)生的解

X′——需要進(jìn)行Tent混沌擾動(dòng)的解

Xnew——使用式(7)、(8)新產(chǎn)生的擾動(dòng)量

對(duì)此解進(jìn)行Tent混沌擾動(dòng),然后將擾動(dòng)后的解與原始解進(jìn)行比較,兩者之中保留最優(yōu)的解。

(7)對(duì)種群進(jìn)行更新后得到最優(yōu)適應(yīng)度以及最優(yōu)位置。

(8)當(dāng)算法達(dá)到最大迭代次數(shù)時(shí)停止算法,未達(dá)到時(shí)返回步驟(3)。

1.3.4ARIMA殘差預(yù)測模型

ISSA-GRU模型的有害氣體濃度預(yù)測精度雖然相比于單一GRU模型有所提升,但是通過本文實(shí)驗(yàn)研究發(fā)現(xiàn),ISSA-GRU模型預(yù)測得到的殘差序列仍然具備一定的規(guī)律,有必要對(duì)其進(jìn)行進(jìn)一步提取。ARIMA算法可以通過歷史數(shù)據(jù)對(duì)未來值進(jìn)行預(yù)測,不僅具有優(yōu)秀的時(shí)間序列提取能力,更能精確反映時(shí)間序列發(fā)展規(guī)律[14],因此使用ARIMA算法對(duì)歷史氣體質(zhì)量濃度殘差序列的有價(jià)值信息進(jìn)行提取和挖掘,然后對(duì)未來氣體質(zhì)量濃度殘差值進(jìn)行預(yù)測,以獲得更加精確的預(yù)測值。ARIMA(p,d,q)模型是由差分運(yùn)算和自回歸滑動(dòng)平均模型(Autoregressive moving average model,ARMA)模型組成,d指的是序列轉(zhuǎn)化成平穩(wěn)序列所需的差分次數(shù);ARMA(p,q)模型是由自回歸項(xiàng)和移動(dòng)平均項(xiàng)組成的,p指的是自回歸項(xiàng)階數(shù),q指的是移動(dòng)平均項(xiàng)階數(shù),ARMA(p,q)模型為

Nt=λ+α1Nt-1+α2Nt-2+…+αpNt-p+εt+

θ1εt-1+θ2εt-2+…+θqεt-q

(14)

式中Nt——?dú)埐钚蛄兄械趖個(gè)值

αi——自回歸項(xiàng)系數(shù)

θq——移動(dòng)平均項(xiàng)系數(shù)

εt——白噪聲序列中第t個(gè)值

λ——常數(shù)

ARIMA(p,d,q)模型建模步驟如下:

(1)使用ISSA-GRU模型對(duì)訓(xùn)練集2進(jìn)行預(yù)測,觀測值序列與預(yù)測值序列的差值為殘差序列,殘差序列為ARIMA模型的輸入。

(2)由于ARIMA模型的輸入序列必須為平穩(wěn)序列,使用單位根檢驗(yàn)對(duì)殘差時(shí)序數(shù)據(jù)的平穩(wěn)性進(jìn)行檢驗(yàn),當(dāng)P<0.05拒絕原假設(shè),該序列為平穩(wěn)序列,否則該序列為非平穩(wěn)序列,需要對(duì)序列進(jìn)行差分直到序列平穩(wěn),方能通過平穩(wěn)性檢驗(yàn)。

(3)對(duì)差分后的序列使用Ljung_Box檢驗(yàn)進(jìn)行白噪聲檢驗(yàn),當(dāng)P<0.05拒絕原假設(shè),該序列為非白噪聲序列,否則該序列為白噪聲序列,通過白噪聲檢驗(yàn)后進(jìn)行下一步。

(4)使用最小赤池信息準(zhǔn)則(Akaike information criterion,AIC)確定p和q的參數(shù),AIC取值最小的ARIMA模型為最優(yōu)模型。

(5)使用最優(yōu)的ARIMA(p,d,q)模型對(duì)殘差值進(jìn)行預(yù)測。

1.4 模型評(píng)價(jià)指標(biāo)

采用均方根誤差(RMSE)、平均絕對(duì)百分比誤差(MAPE)和決定系數(shù)R2對(duì)模型的性能進(jìn)行評(píng)估。

2 結(jié)果與分析

2.1 豬舍環(huán)境因子監(jiān)測結(jié)果及相關(guān)性分析

選取2022年8月10日到9月30日共52 d的預(yù)處理后環(huán)境數(shù)據(jù)作為實(shí)驗(yàn)數(shù)據(jù),豬舍環(huán)境因子變化曲線如圖5所示;其統(tǒng)計(jì)結(jié)果如表2所示。

表2 豬舍環(huán)境數(shù)據(jù)統(tǒng)計(jì)

圖5 環(huán)境因子變化曲線

豬舍內(nèi)有害氣體濃度受豬舍內(nèi)外溫度、濕度等多環(huán)境因素的影響,構(gòu)建有害氣體濃度預(yù)測模型時(shí)加入外部影響因素可以使模型更加充分學(xué)習(xí)到預(yù)測目標(biāo)的特征。使用Pearson’s系數(shù)分析法(P<0.05)對(duì)影響NH3質(zhì)量濃度和CO2質(zhì)量濃度的多環(huán)境因素進(jìn)行相關(guān)性分析,如表3所示。

表3 皮爾森系數(shù)相關(guān)性分析

豬舍內(nèi)NH3質(zhì)量濃度與CO2質(zhì)量濃度具有較強(qiáng)的正相關(guān)性,相關(guān)系數(shù)為0.72。豬舍內(nèi)NH3質(zhì)量濃度與舍內(nèi)、外溫度和通風(fēng)速率呈負(fù)相關(guān)性,相關(guān)系數(shù)分別為-0.49、-0.58和-0.36;NH3質(zhì)量濃度與舍內(nèi)、外濕度具有正相關(guān)性,相關(guān)系數(shù)分別為0.26和0.21。豬舍內(nèi)CO2質(zhì)量濃度與舍內(nèi)、外相對(duì)濕度具有正相關(guān)性,相關(guān)系數(shù)分別為0.41和0.30;CO2質(zhì)量濃度與舍內(nèi)、外溫度和通風(fēng)速率具有較強(qiáng)的負(fù)相關(guān)性,相關(guān)系數(shù)分別為-0.81、-0.62和-0.48。

相關(guān)性分析的結(jié)果證明了豬舍內(nèi)環(huán)境因子之間具有相關(guān)性,建模時(shí)加入外部影響因素可以減少模型預(yù)測誤差。因此,本文選取通風(fēng)速率、舍內(nèi)溫度、舍內(nèi)相對(duì)濕度、舍外溫度、舍外相對(duì)濕度、CO2質(zhì)量濃度作為NH3質(zhì)量濃度的外部影響因素;選取通風(fēng)速率、舍內(nèi)溫度、舍內(nèi)相對(duì)濕度、舍外溫度、舍外相對(duì)濕度、NH3質(zhì)量濃度作為CO2質(zhì)量濃度的外部影響因素。

2.2 ISSA-GRU模型構(gòu)建

在使用優(yōu)化算法對(duì)超參數(shù)進(jìn)行尋優(yōu)之前通過預(yù)實(shí)驗(yàn)確定GRU子層數(shù)目。實(shí)驗(yàn)結(jié)果如圖6所示,當(dāng)GRU子層數(shù)目為1時(shí),模型預(yù)測的RMSE分別為90.51 mg/m3和0.496 2 mg/m3;當(dāng)GRU子層數(shù)增加到2時(shí),模型預(yù)測的RMSE分別為83.15 mg/m3和0.381 mg/m3;當(dāng)GRU子層數(shù)繼續(xù)增加時(shí),預(yù)測精度卻下降,故本模型GRU子層數(shù)目設(shè)置為2。

圖6 不同GRU層數(shù)時(shí)模型預(yù)測的均方根誤差變化

使用ISSA算法對(duì)GRU神經(jīng)網(wǎng)絡(luò)的時(shí)間窗口、學(xué)習(xí)率、批數(shù)目、GRU子層1神經(jīng)元個(gè)數(shù)、GRU子層2神經(jīng)元個(gè)數(shù)這5個(gè)超參數(shù)進(jìn)行尋優(yōu),尋優(yōu)范圍如表4所示。ISSA對(duì)GRU模型超參數(shù)尋優(yōu)過程中優(yōu)化算法種群數(shù)目設(shè)置為10,迭代次數(shù)設(shè)置為50,預(yù)警值設(shè)置為0.8,發(fā)現(xiàn)者比例設(shè)置為80%,警戒者比例設(shè)置為20%。ISSA算法對(duì)GRU模型的超參數(shù)優(yōu)化結(jié)果如表4所示,當(dāng)時(shí)間窗口值為3,學(xué)習(xí)率為0.000 6,批數(shù)目為71,兩個(gè)GRU子層的神經(jīng)元數(shù)目分別為486和315時(shí)模型預(yù)測精度最優(yōu)。

表4 ISSA尋優(yōu)參數(shù)

為了更好評(píng)估ISSA算法對(duì)GRU模型預(yù)測精度的優(yōu)化,本文選取PSO算法、SSA算法和GWO算法進(jìn)行對(duì)比。如圖7所示。隨著迭代次數(shù)的增加4種算法都逐漸收斂,ISSA算法的適應(yīng)度最低(0.005 6),較PSO(0.006 2)、SSA(0.005 9)和GWO(0.006 0)3種算法的適應(yīng)度分別降低9.68%、5.08%和6.67%,因此說明ISSA算法具有更高的搜索精度。

圖7 4種不同算法適應(yīng)度變化曲線

2.3 ARIMA模型構(gòu)建

以訓(xùn)練集2為數(shù)據(jù)集,使用超參數(shù)優(yōu)化后的ISSA-GRU模型對(duì)舍內(nèi)NH3質(zhì)量濃度、CO2質(zhì)量濃度進(jìn)行預(yù)測,得到的監(jiān)測值與預(yù)測值的殘差序列如圖8所示,使用殘差序列構(gòu)建ARIMA模型。

圖8 殘差變化曲線

以NH3質(zhì)量濃度殘差序列為輸入,對(duì)ARIMA模型輸入序列的平穩(wěn)性進(jìn)行檢驗(yàn),P=8.85×10-14<0.05,該序列通過平穩(wěn)性檢驗(yàn),NH3質(zhì)量濃度殘差序列為平穩(wěn)序列,因此d=0。NH3質(zhì)量濃度殘差序列的白噪聲檢驗(yàn)P=5.44×10-4<0.05,NH3殘差序列為非白噪聲序列。以CO2質(zhì)量濃度殘差序列為輸入,對(duì)ARIMA模型輸入序列的平穩(wěn)性進(jìn)行檢驗(yàn),P=6.27×10-13<0.05,該序列通過平穩(wěn)性檢驗(yàn),CO2質(zhì)量濃度殘差序列為平穩(wěn)序列,因此d=0。CO2質(zhì)量濃度殘差序列的白噪聲檢驗(yàn)P=1.60×10-21<0.05,CO2質(zhì)量濃度殘差序列為非白噪聲序列。

使用網(wǎng)格搜索法確定p和q的取值,當(dāng)p取3且q取1時(shí)NH3質(zhì)量濃度殘差A(yù)RIMA預(yù)測模型AIC值最小,為447.232 2,此時(shí)ARIMA模型最優(yōu),NH3質(zhì)量濃度殘差最優(yōu)ARIMA模型為ARIMA(3,0,1),根據(jù)相同步驟可確定CO2質(zhì)量濃度殘差最優(yōu)ARIMA模型為ARIMA(4,0,1)。

2.4 ISSA-GRU-ARIMA模型預(yù)測結(jié)果驗(yàn)證

使用測試集對(duì)模型預(yù)測精度進(jìn)行測試,如圖9所示,ISSA-GRU-ARIMA模型對(duì)NH3質(zhì)量濃度和CO2質(zhì)量濃度的預(yù)測值變化趨勢與觀測值變化趨勢都非常接近,呈現(xiàn)出周期性正弦或余弦曲線波動(dòng)變化。NH3質(zhì)量濃度觀測值變化范圍為0.9～5.4 mg/m3,ISSA-GRU-ARIMA模型預(yù)測值變化范圍為 0.8～5.5 mg/m3,預(yù)測值與觀測值的最小誤差為0 mg/m3,最大誤差為1 mg/m3,平均誤差為0.2 mg/m3。CO2質(zhì)量濃度觀測值變化范圍為983.5～2 796.5 mg/m3,ISSA-GRU-ARIMA模型預(yù)測值變化范圍為978.5～2 823.5 mg/m3;預(yù)測值與觀測值的最小誤差為1 mg/m3,最大誤差為194 mg/m3,平均誤差為41 mg/m3。

圖9 不同模型預(yù)測值對(duì)比

2.4.1不同模型預(yù)測結(jié)果對(duì)比

為了進(jìn)一步驗(yàn)證模型預(yù)測效果,使用CNN、GRU、SVR、RF、PSO-GRU、GWO-GRU、SSA-GRU、ISSA-GRU、ISSA-CNN、ISSA-CNN-ARIMA模型進(jìn)行對(duì)比實(shí)驗(yàn),與本文所提模型ISSA-GRU-ARIMA使用相同數(shù)據(jù)集,實(shí)驗(yàn)結(jié)果如表5所示。CNN、GRU、SVR、RF是常見基礎(chǔ)模型,PSO-GRU、GWO-GRU、SSA-GRU是使用其他優(yōu)化算法對(duì)GRU模型進(jìn)行優(yōu)化,ISSA-CNN、ISSA-CNN-ARIMA是使用ISSA算法和ARIMA算法對(duì)CNN進(jìn)行優(yōu)化。

表5 模型預(yù)測精度對(duì)比

如表5所示,本文模型ISSA-GRU-ARIMA預(yù)測精度最優(yōu),NH3質(zhì)量濃度預(yù)測的RMSE、MAPE和R2分別為0.263 mg/m3、8.171%和0.928,CO2質(zhì)量濃度預(yù)測的RMSE、MAPE和R2分別為55.361 mg/m3、4.633%和0.985。傳統(tǒng)機(jī)器學(xué)習(xí)模型RF預(yù)測精度最差,NH3質(zhì)量濃度預(yù)測RMSE、MAPE和R2分別為0.722 mg/m3、20.592%和0.4978,CO2質(zhì)量濃度預(yù)測的RMSE和MAPE和R2分別為170.873 mg/m3、12.926%和0.859。本文模型與RF模型相比,NH3質(zhì)量濃度預(yù)測RMSE和MAPE分別下降63.51%和60.31%,R2提升86.41%,CO2質(zhì)量濃度預(yù)測的RMSE和MAPE分別下降67.58%和64.14%,R2提升14.65%。

ISSA-GRU模型與GRU模型相比,NH3和CO2質(zhì)量濃度預(yù)測的RMSE和MAPE分別降低21.69%和20.66%,17.25%和20.25%;R2分別提升了6.87%和1.28%,說明使用ISSA算法進(jìn)行尋優(yōu)顯著可提升GRU模型的預(yù)測精度。ISSA-GRU與PSO-GRU、GWO-GRU和SSA-GRU相比,NH3質(zhì)量濃度預(yù)測的RMSE分別下降18.97%、11.57%和7.27%,CO2質(zhì)量濃度預(yù)測的RMSE分別下降17.26%、17.24%和5.40%;NH3質(zhì)量濃度預(yù)測的MAPE分別下降26.48%、7.64%和14.91%,CO2質(zhì)量濃度預(yù)測的MAPE分別下降15.94%、14.77%和6.71%;NH3質(zhì)量濃度預(yù)測的R2分別上升5.63%、2.92%和1.69%,CO2質(zhì)量濃度預(yù)測的R2分別上升1.01%、1.00%和0.26%,說明ISSA相比其他優(yōu)化算法擁有更優(yōu)秀的尋優(yōu)能力。

ISSA-GRU-ARIMA模型與ISSA-GRU相比,NH3質(zhì)量濃度預(yù)測的RMSE和MAPE分別下降11.66%和12.13%,R2提升2.24%,CO2質(zhì)量濃度預(yù)測的RMSE和MAPE分別下降16.08和15.99%,R2提升0.64%,由此可見,融合線性統(tǒng)計(jì)ARIMA方法后的ISSA-GRU-ARIMA模型,可以充分挖掘ISSA-GRU預(yù)測殘差序列的有用信息,進(jìn)一步提升模型預(yù)測精度。

2.4.2結(jié)果分析

近年來,對(duì)豬舍氣體濃度進(jìn)行預(yù)測已成為豬舍環(huán)境控制研究的關(guān)注熱點(diǎn),國內(nèi)外學(xué)者針對(duì)畜禽舍內(nèi)不同氣體構(gòu)建了不同的預(yù)測模型。例如,文獻(xiàn)[7]提出的基于ANFIS的NH3質(zhì)量濃度預(yù)測模型R2為0.648; 文獻(xiàn)[9]構(gòu)建的融合注意力機(jī)制和CNN網(wǎng)絡(luò)的CNN-LSTM-AT模型,其CO2質(zhì)量濃度預(yù)測的R2為0.8170;文獻(xiàn)[27]使用不同時(shí)間間隔環(huán)境因子構(gòu)建基于長短期記憶網(wǎng)絡(luò)的CO2質(zhì)量濃度預(yù)測模型,其R2分別為0.859和0.656,文獻(xiàn)[28]使用多種機(jī)器學(xué)習(xí)模型對(duì)NH3質(zhì)量濃度進(jìn)行預(yù)測,梯度提升樹算法取得了最好的預(yù)測效果,R2為0.402。這些模型對(duì)NH3質(zhì)量濃度和CO2質(zhì)量濃度預(yù)測的R2范圍為0.4～0.9,本文所提出ISSA-GRU-ARIMA模型對(duì)NH3質(zhì)量濃度和CO2質(zhì)量濃度預(yù)測R2分別是0.928和0.985,優(yōu)于上述已有研究中的模型,可為豬舍氣體濃度預(yù)測及環(huán)境控制提供可行的方法。

綜上,本文ISSA-GRU-ARIMA氣體質(zhì)量濃度預(yù)測模型預(yù)測精度較優(yōu),可用于豬舍NH3質(zhì)量濃度和CO2質(zhì)量濃度預(yù)測。

3 結(jié)論

(1)通過ISSA算法優(yōu)化GRU模型的超參數(shù),可提升氣體質(zhì)量濃度預(yù)測模型的精度。ISSA-GRU模型NH3質(zhì)量濃度預(yù)測的RMSE、MAPE和R2分別為0.298 mg/m3、9.301%和0.907,CO2質(zhì)量濃度預(yù)測的分別為65.974 mg/m3、5.515%和0.979;與GRU(RMSE分別為0.381 mg/m3和83.151 mg/m3,MAPE分別為11.247%和6.916%,R2分別為0.849和0.967)相比,RMSE、MAPE分別降低21.69%和20.66%,17.25%和20.25%;R2分別提升6.87%和1.28%。

(2)融合統(tǒng)計(jì)學(xué)習(xí)ARIMA方法對(duì)ISSA-GRU模型進(jìn)行優(yōu)化,使得NH3質(zhì)量濃度和CO2質(zhì)量濃度預(yù)測的RMSE分別下降11.66%和16.08%,因此ARIMA模型對(duì)殘差進(jìn)行線性特征提取可以提升模型預(yù)測精度。

(3)提出的ISSA-GRU-ARIMA模型對(duì)NH3質(zhì)量度和CO2質(zhì)量濃度預(yù)測的R2分別為0.928和0.985,預(yù)測精度較高,預(yù)測結(jié)果可以為豬舍有害氣體濃度控制提供科學(xué)依據(jù)。