基于NARX模型的參考作物蒸散發(fā)預(yù)測

武劍飛，康銀紅，宋鑫，梁友鵬

(四川農(nóng)業(yè)大學(xué)水利水電學(xué)院，四川 雅安 625014)

作物蒸散量是參與大氣水文循環(huán)的重要因子，同時也是農(nóng)業(yè)灌溉需水量制定的重要參數(shù).參考作物蒸散量(reference crop evapotranspiration，ET0)是國際上計算作物蒸散量的通用參數(shù)，準(zhǔn)確計算及預(yù)報ET0具有重要意義.目前關(guān)于ET0的確定主要分為3個方面：試驗測定、公式計算、數(shù)值模擬.運用蒸滲儀法測定的ET0精度好，但是儀器費用高且耗時長，只適用于一些特定研究，不適宜大范圍推廣[1-2].基于氣象參數(shù)建立的ET0計算公式是許多學(xué)者常用的方法，目前已發(fā)展出50多種計算方法.在1998年，聯(lián)合國糧農(nóng)組織(food and agriculture organization，F(xiàn)AO)把Penman-Monteith(PM)公式作為計算ET0的唯一標(biāo)準(zhǔn)方法，并且在全世界范圍內(nèi)都得到了很好應(yīng)用，但是該公式需要氣象參數(shù)多，計算復(fù)雜，對于一些氣象資料缺失的地區(qū)不具有適用性[3].因此，一些學(xué)者研究出了簡化的ET0計算模型，例如溫度法中的Hargreaves-Samani模型[4]，輻射法中的Makkink模型[5]、Priestley-Taylor模型[6]，經(jīng)驗法中的Irmak-Allen模型[7]等.這些簡化模型雖然能夠在氣象資料缺失的情況下計算ET0，但是這些模型都是在特定的地區(qū)發(fā)展而來的，對于不同的環(huán)境需要篩選不同的計算模型.

隨著計算機技術(shù)的發(fā)展以及仿真模擬的應(yīng)用，基于機器學(xué)習(xí)下的ET0模擬計算以及預(yù)測成為國內(nèi)外學(xué)者的研究熱點.ANTONOPOULOS等[8]、JAHANBANI等[9]運用人工神經(jīng)網(wǎng)絡(luò)(artificial neural network，ANN)模擬了ET0，結(jié)果表明ANN模擬的精度比簡化模型計算的精度高；ADAMALA[10]對印度不同地區(qū)建立廣義小波神經(jīng)網(wǎng)絡(luò)(wavelet-neural network，WNN)以模擬ET0，結(jié)果發(fā)現(xiàn)WNN模型比傳統(tǒng)的ANN模型、線性回歸模型和Hargreaves模型的模擬誤差小，具有良好的應(yīng)用能力.崔寧博等[11]基于思維進化算法(mind evolutionary algorithm)優(yōu)化的誤差反向傳波神經(jīng)網(wǎng)絡(luò)(back propagation neural network，BPNN)模型預(yù)測了西北干旱地區(qū)的ET0，發(fā)現(xiàn)該模型在相同氣象參數(shù)輸入的情況下，模擬精度遠(yuǎn)高于Hargreaves-Samani模型、Irmak模型和48-PM模型；馮禹等[12]基于只有溫度資料建立廣義回歸神經(jīng)網(wǎng)絡(luò)(generalized regression neural network，GRNN)，通過與Hargreaves模型、改進的Hargreaves模型比較，發(fā)現(xiàn)該模型在時空尺度上能夠很好地預(yù)測四川盆地的ET0；王升等[13]基于基因編程表達式算法(gene expression programming，GEP)以及徑向基神經(jīng)網(wǎng)絡(luò)算法(radial basis function-artificial neural network，RBFNN)建立了ET0預(yù)報模型，發(fā)現(xiàn)在氣象參數(shù)輸入相同的情況下GEP模型比RBFNN模型的計算精度更高.

動態(tài)神經(jīng)網(wǎng)絡(luò)是一種具有動態(tài)反饋性的神經(jīng)網(wǎng)絡(luò)，能夠有效實現(xiàn)復(fù)雜的動態(tài)系統(tǒng)建模.常見的動態(tài)神經(jīng)網(wǎng)絡(luò)有帶外部輸入非線性自回歸(nonlinear autoregressive models with exogenous inputs，NARX)神經(jīng)網(wǎng)絡(luò)、非線性自回歸(nonlinear autoregressive，NAR)神經(jīng)網(wǎng)絡(luò)等，目前廣泛應(yīng)用于電力負(fù)荷預(yù)測[14]、鍋爐溫度預(yù)測[15]等方面，對于ET0的預(yù)測還未有報道.

1 材料與方法

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

選取四川西南山地的攀枝花站點作為研究站點，氣象數(shù)據(jù)包括日最高溫(Tmax)、日最低溫(Tmin)、日照時數(shù)(n)、相對濕度(RH)和10 m高度處的風(fēng)速(采用FAO風(fēng)速輪廓線關(guān)系[3]換算為2 m高度處的風(fēng)速u2)的2010—2018年的逐日氣象資料.其中，2010—2017年的數(shù)據(jù)作為訓(xùn)練集，2018年的數(shù)據(jù)作為測試集；數(shù)據(jù)來自國家氣象信息中心，且經(jīng)過較好的質(zhì)量控制.

1.2 參考作物計算模型

選取Penman-Monteith模型計算的ET0作為標(biāo)準(zhǔn)值，選用Hargreaves-Samani 模型、Irmak-Allen模型、Makkink模型以及Priestley-Taylor模型進行對比分析，具體計算公式見表1.

表1中ET0為參考作物蒸發(fā)蒸騰量，mm/d；Δ為飽和水氣壓曲線斜率，kPa/℃；Rn為作物表面凈輻射，MJ/(m2·d)；Ra為地外輻射，MJ/(m2·d)；Rs為太陽輻射，MJ/(m2·d)；G為土壤熱通量，MJ/(m2·d)，以天計算的土壤熱通量為0；γ為濕度常數(shù)，kPa/℃；Tmean為日平均溫度，℃；es為飽和水氣壓，kPa；ea為實際水氣壓，kPa；λ為氣化潛熱，2.45 MJ/kg.

表1 參考作物蒸散量計算公式

1.3 NARX模型

NARX是一種非線性動態(tài)循環(huán)網(wǎng)絡(luò)模型，由輸入層、隱含層、輸出層及輸入和輸出延遲構(gòu)成[16]，模型結(jié)構(gòu)見圖1.

圖1 NARX模型結(jié)構(gòu)圖

模型帶有外部輸入變量，能夠接受來自輸出神經(jīng)元的反饋，反映過去輸入與輸出的歷史信息，是一種具有記憶功能的動態(tài)神經(jīng)網(wǎng)絡(luò)[17-18].NARX同時引進外部輸入和輸出兩種時間序列，運用輸出y(t)和另一個時間序列x(t)的歷史值以模擬預(yù)測未來時間內(nèi)y(t)，模型數(shù)學(xué)表達形式為

y=f[x(t-1),…,x(t-m),y(t-1),y(t-2),…,(t-n)]，

(1)

式中：m，n分別為時間序列x(t)和y(t)的延遲階數(shù)；t為時間.

從實現(xiàn)方式來說，耗材信息化建設(shè)的本質(zhì)是信息技術(shù)與醫(yī)療行業(yè)的有機融合，其主導(dǎo)因素是信息技術(shù)的運用。但隨著醫(yī)院管理需求的變化，其對應(yīng)的應(yīng)用需求也要進行適時調(diào)整，這就決定了在信息化過程中，信息技術(shù)與信息管理是既獨立又密切相連的特殊關(guān)系，充分保障信息化技術(shù)手段平穩(wěn)發(fā)揮作用是醫(yī)院管理部門的重要職責(zé)，通過技術(shù)手段實現(xiàn)醫(yī)院運營管理一體化，能更好地發(fā)揮信息資源的高效利用［3］。

NARX模型的輸入及輸出延時階數(shù)設(shè)為1，隱含層的個數(shù)為10，訓(xùn)練函數(shù)為“trainlm”，利用Matlab軟件進行模擬預(yù)測.

1.4 NAR模型

NAR是一種基于時間序列且具有反饋和記憶功能的動態(tài)神經(jīng)網(wǎng)絡(luò)模型，模型的輸入和輸出并不是靜態(tài)的映射方式，而是每個輸出基于當(dāng)前系統(tǒng)的動態(tài)結(jié)果綜合得到[15].模型由時間序列y(t)的歷史值預(yù)測未來時段內(nèi)的y(t)值，數(shù)學(xué)表達式為

y(t)=f[y(t-1),…,y(t-n)].

(2)

模型由輸入層、隱含層、輸出層和滯后層組成.NAR模型的隱含層個數(shù)默認(rèn)取10，取滯后階數(shù)默認(rèn)值1∶2，訓(xùn)練函數(shù)為“trainlm”.

1.5 模型評價方法

選取均方根誤差(RMSE)、平均絕對偏差(MAE)、平均偏差(MBE)、納什效率系數(shù)(NSE)、綜合評價指標(biāo)(GPI)評價各模型的精度，計算公式為

(3)

(4)

(5)

(6)

(7)

RMSE，MAE和MBE越接近0，表明模型的精度越好；NSE越接近1，表明模型效率高.GPI采用排名形式呈現(xiàn)，排名越靠前，模型的性能越好.

2 結(jié)果分析

2.1 不同氣象因子組合下的NARX模型預(yù)測

基于Tmax，Tmin，RH，u2，n這5種與ET0相關(guān)的氣象因子進行隨機組合，建立11組NARX模型，分析不同氣象因子組合下的模型模擬效果，具體結(jié)果見表2.

表2 不同氣象因子組合下NARX模型的評價指標(biāo)

當(dāng)輸入5個氣象因子時，模型NARX-1的RMSE，MAE，MBE分別為0.425 mm/d，0.320 mm/d，0.069 mm/d，NES為0.920，GPI排名第11，模型模擬效果在所有組合中效果最差，這與馮禹等[19]、魏俊等[20]對ELM模型研究結(jié)果(輸入5個氣象因子的模型模擬效果最優(yōu))有所差異，這可能是由于某個因子對ET0的影響已包含在另一因子中(例如風(fēng)速可通過溫度反映對ET0的影響)，而該因子能夠充分反映對ET0的影響，NARX模型包含了過多的冗余信息.過多的輸入維度使得NARX模型趨于復(fù)雜，泛化能力下降，模擬精度較低.

當(dāng)輸入4個氣象因子時，模型NARX精度相比模型NARX-1有所提高，但是NARX-2和NARX-3之間差異比較大.NARX-2(缺失RH)的RMSE，MBE分別為0.363 mm/d和0.035 mm/d，比NARX-3(缺失n)的RMSE，MBE分別減少了0.082 mm/d和0.032 mm/d，而NSE增幅為3.2%，兩者GPI排名分別為第8、第9.說明日照時數(shù)對ET0的影響比相對濕度大，這與馮禹等[21]對四川地區(qū)蒸散發(fā)研究中得出日照時數(shù)是四川地區(qū)的主要影響因子的結(jié)論一致.

當(dāng)輸入3個氣象因子時，模型NARX精度差異明顯.NARX-5(缺失RH，n)精度最高，它與NARX-4(缺失RH，u2)的RMSE均低于0.370 mm/d，NSE大于0.940，而NARX-6(缺失n，u2)的RMSE高于0.400 mm/d，NSE低于0.930，這說明n和u2對ET0有著重要的影響.除了NARX-6對ET0低估外(MBE<0),其他模型均高估了ET0.對比模型NARX-5(輸入Tmax，Tmin，u2)和NARX-3(輸入Tmax，Tmin，RH，u2)可以發(fā)現(xiàn)，在缺失RH時，模型NARX-5的RMSE下降幅度為25%，NSE增加了4.1%，GPI排名從第9升至第2，說明相對濕度對該地區(qū)ET0影響較小.

當(dāng)輸入2個氣象因子時，NARX-7的RMSE比NARX-8低0.050 mm/d，MAE減少0.063 mm/d，MBE降低58%，NSE增加1.7%，GPI排名第4，精度明顯優(yōu)于NARX-8.對比NARX-8(缺失Tmax，Tmin，RH)和NARX-4(缺失RH，u2)和NARX-6(缺失n，u2)，減少輸入?yún)?shù)后模型維度降低，精度并不一定提高.這主要與氣象因子對ET0的貢獻程度有關(guān).

當(dāng)輸入1個氣象因子時，模型NARX-11(輸入RH)精度最低.NARX-9(輸入u2)與NARX-10(輸入n)在所有模型中精度均較高.其中，模型NARX-9的GPI排名第1，NSE為0.964，而RMSE和MAE均低于0.300 mm/d.進一步說明日照時數(shù)在該地區(qū)對ET0影響較大.而基于較少氣象因子建立的BPNN模型[22]和WNN模型[23]的模擬精度較差，這說明NARX模型泛化能力更強，模擬效果更好.

2.2 模型NARX與NAR精度比較

為了比較模型NAR與NARX的精度，選取文中精度較差的模型NARX-1，分別繪制NARX-1和NAR的預(yù)測值與PM計算值的線性擬合圖，如圖2所示.由圖可知，模型NARX-1的預(yù)測值與PM計算值擬合效果最好，趨勢線斜率趨近于1，R2>0.900，與PM計算值有很好的相關(guān)性.模型NAR的預(yù)測值分布較為散亂，斜率為0.706 3，R2遠(yuǎn)小于0.700，模擬誤差比較大.因此，模型NARX比NAR的模擬精度更高.

2.3 NARX模型與其他模型的精度比較

以PM模型計算的ET0作為標(biāo)準(zhǔn)值，比較基于Tmax和Tmin的NARX-7和Hargreaves-Samani模型以及基于Tmax，Tmin和n的NARX-4，Irmak-Allen，Makkink，Priestley-Taylor模型的精度.模擬精度如圖3所示.各模型與PM模型的擬合效果之間差異具有統(tǒng)計學(xué)意義.模型Irmak-Allen和Makkink的ET0整體在1∶1線以下，其中Irmak-Allen擬合效果最差(擬合斜率為0.5168，R2為0.771)；模型Hargreaves-Samani和Priestley-Taylor的ET0在1∶1線附近分布較為散亂，R2在0.800左右；模型NARX-7和NARX-4的ET0均勻分布在1∶1線附近，且R2在0.950以上，擬合效果最優(yōu).

圖2 NARX-1，NAR模型與PM模型的線性擬合

圖3 PM模型與其他模型的線性擬合

進一步對各模型的精度進行統(tǒng)計分析.圖4為NARX模型與其他計算模型精度比較，圖中σ為精度比較誤差.由圖可以直觀看出各模型誤差情況.模型Irmak-Allen在所有模型中誤差較高，RMSE大于1.000 mm/d，NSE小于0.200，精度最差；在輸入氣象因子相同的情況下，NARX模型的誤差值均小于簡化模型的誤差值，模型Irmak-Allen和Makkink的MBE均小于-0.500 mm/d，表明計算的ET0值偏低較多.對比 NARX-7(輸入Tmax，Tmin)和Hargreaves-Samani(輸入Tmax，Tmin)模型，NARX-7的RMSE，MAE和MBE均低于0.500 mm/d，NSE接近于1，而Hargreaves-Samani的各項誤差遠(yuǎn)大于0.600 mm/d，NSE低于0.700，比NARX-7的精度低很多.相比較輸入相同氣象因子的模型NARX-4，Irmak-Allen，Makkink和Priestley-Taylor，除了模型NARX-4和Priestley-Taylor的MBE比較接近，其他情況下模型NARX-4的MAE和RMSE比其他3個模型誤差小得多，NARX-4的NSE在0.900以上，其他模型的NSE均低于0.700，說明NARX-4模擬精度遠(yuǎn)比其他3個模型精度高.

以上分析表明，在日照時數(shù)、濕度、風(fēng)速氣象資料缺失的情況下，NARX模型的模擬精度比簡化模型更具有優(yōu)勢，可作為該地區(qū)預(yù)測ET0的優(yōu)良模型.

圖4 NARX模型與其他計算模型精度比較

3 討 論

基于氣象因子建立的Hargreaves-Samani和Makkink等各種ET0簡化模型雖然能夠很好地解決某些氣象資料缺失下的ET0獲取問題，但是這些簡化模型普遍存在精度較差且模型具有區(qū)域性問題，需要進行篩選且修正后才能應(yīng)用到該地區(qū).隨著機器模擬技術(shù)的發(fā)展，基于各種機器學(xué)習(xí)條件下的ET0模型不需要具體的公式，僅建立氣象因子與ET0之間的關(guān)系進行模擬，極大地提高了ET0的計算精度.NARX模型與靜態(tài)神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)相比，輸入層增加了歷史值y(t)，每次y(t)反饋到輸入層，從而達到長期記憶的能力，而且對y(t)影響因素較大的輸入變量能夠很好反映變量間的非線性關(guān)系.

研究表明，在氣象因子缺失條件相同的情況下，NARX模型比簡化的ET0計算模型的精度更高.通過建立不同氣象因子下的NARX模型發(fā)現(xiàn)，基于單個氣象因子的NARX模型的精度(RMSE低于0.400 mm/d，NSE在0.940以上)高于其他一些機器學(xué)習(xí)模型，例如張皓杰等[24]對西北地區(qū)建立極限學(xué)習(xí)機(extreme learning machine，ELM)模型時，基于單個氣象因子的ELM模型精度最低(RMSE高于1.300 mm/d，NSE在0.500以下)；KISI等[25]對Ankara和Kirikkale地區(qū)比較了模型ELM，WANN和ANN在不同氣象因子組合下的精度，結(jié)果發(fā)現(xiàn)基于少量氣象因子建立的這些模型精度均較低(RMSE均高于1.200 mm/d，NSE低于0.600).這主要是由于模型結(jié)構(gòu)不同導(dǎo)致，模型NARX比模型ELM和WANN在輸入層增加了一個外部輸入反饋，使得模型可以根據(jù)歷史值調(diào)整輸出，從而逼近實際值.

文中應(yīng)用NARX模型僅進行了單站點短期ET0的預(yù)測，但對于多站點長序列的ET0預(yù)測以及和其他神經(jīng)網(wǎng)絡(luò)模型對比，還有待后續(xù)進一步地研究.

4 結(jié) 論

以攀枝花為研究站點，基于5個氣象因子作為輸入?yún)?shù)，PM模型計算的ET0作為輸出，構(gòu)建了動態(tài)神經(jīng)網(wǎng)絡(luò)中NARX和NAR模型，并且與其他模型進行對比分析，得到以下主要結(jié)論.

1) 不同氣象因子組合下的NARX模型對比發(fā)現(xiàn)，輸入單個氣象因子的模型精度普遍高于其他組合，基于u2的NARX-9模擬精度最高，其RMSE為0.285 mm/d，MAE和MBE分別為0.237 mm/d和0.019 mm/d，NSE為0.964，GPI排名第1，能夠很好地實現(xiàn)該地區(qū)缺失氣象資料下的ET0預(yù)測.

2) 將模型NARX-1與NAR對比發(fā)現(xiàn)，帶有外部輸入的NARX-1模型精度高于無外部輸入的NAR模型.

3) 通過模型NARX-4，NARX-7和模型Hargreaves-Samani，Irmak-Allen，Makkink以及Priestley-Taylor對比分析發(fā)現(xiàn)，基于Tmax，Tmin，n時，模型NARX-4的RMSE，MAE，MBE均在0.500 mm/d以下，NSE高于0.950，而模型Irmak-Allen，Makkink，Priestley-Taylor的RMSE，MAE均在0.600 mm/d以上，NSE均在0.700以下，精度相對較差；基于Tmax，Tmin時，模型NARX-7(RMSE<0.500 mm/d，NSE>0.940)精度高于模型Hargreaves-Samani(RMSE>0.600 mm/d，NSE<0.700).因此，在氣象資料缺失條件下，NARX模型相較于其他簡化計算模型的精度高，更適合用于模擬計算ET0.