王振龍,劉竹梅,呂海深,丁佳楠,陸云燕,王怡寧
1.安徽省(水利部淮委)水利科學研究院/水利水資源安徽省重點實驗室/五道溝水文水資源實驗站,安徽 蚌埠 233000;2.河海大學理學院,江蘇 南京 210098;3.南京水利科學研究院,江蘇 南京 210029
淮河流域是中國重要的農業生產基地,夏玉米是淮河流域主要的秋糧作物,對當地經濟發展及人民食品安全至關重要(張建軍等,2014;袁宏偉等,2018)。蒸散量是土壤蒸發量和作物蒸騰量之和,對地表能量循環和水循環有著非常重要的作用(Xu et al.,2005;趙靜等,2009;張淑杰等,2010;Lievens et al.,2017)。無水分脅迫時,作物系數是實際蒸散量與參考作物蒸散量的比值,是估算作物蒸散量的關鍵指標(汪順生等,2013;王維等,2015)。根據方法不同,可計算雙作物系數(Rosa et al.,2012;Ding et al.,2013)和單作物系數(李毅等,2018);根據時間尺度不同,可計算逐日作物系數(韓文霆等,2018)、旬作物系數(盧曉鵬等,2012)、月作物系數(云文麗等,2013)以及年平均作物系數(李波等,2020),進而可計算月平均蒸散量(涂安國等,2017)、各季節平均蒸散量(劉玉汐等,2019)或年平均蒸散量(吳文玉等,2013)。對于作物系數模型的建立,許多學者嘗試用優化算法進行研究:王怡寧等(2020)利用通徑分析篩選作物系數的影響因子,并結合BP神經網絡建立夏玉米作物系數模型,該模型對夏玉米蒸散量的預測誤差較小,準確率高于70%,但不足之處在于,沒有提出基于關鍵影響因子的蒸散量預測的具體概念性模型;金菊良等(2017)利用遺傳算法對無受旱的大豆作物系數及水分脅迫系數進行求解,建立了大豆作物系數模型,模型對大豆無受旱脅迫時的蒸散量估算誤差較小,均方根誤差均小于0.3 mm·d?1,但由于土壤水分脅迫系數也是由遺傳算法模擬得到,不是逐日實測值,故模型在估算受旱脅迫的大豆蒸散量時誤差增大;袁宏偉等(2018)利用遺傳算法優化太陽輻射參數,得到無受旱的玉米作物系數模型,該模型對無受旱玉米蒸散量的整體估算效果較好,優于FAO-56推薦值的估算結果,但在估算受旱玉米蒸散量時誤差明顯增大,且沒有提出對應的解決方法,使得模型總體上低估了玉米的蒸發蒸騰量。
遺傳算法是常用的尋優算法,可在搜索空間中全局搜索最優解(李華昌等,2005),但已有的研究一般只針對使用算法對各生育期作物系數進行尋優求解,未考慮到對算法輸出結果誤差的修正,簡單地尋優求解并不一定適合所有地區的數據,因此對算法輸出結果進行修正是有必要的。本文利用遺傳算法,并結合單作物系數法,建立五道溝地區夏玉米作物系數模型,并利用葉面積指數對模型進行修正,得到誤差更小的夏玉米作物系數及蒸散量估算模型。
實驗區位于安徽省蚌埠市固鎮縣五道溝水文實驗站,地處淮北平原南部,站內共有62套大型地中蒸滲儀群、10套大型稱重式蒸滲儀及自動高精度氣象站,設有蒸散發、氣象要素、潛水蒸發及土壤水分等觀測場地,地下水位埋深及年變幅為1—3 m,作物主要是冬小麥、大豆、夏玉米。土壤類型主要為砂姜黑土和黃潮土,分別占淮北平原土壤總面積的54%和33%,凋萎含水率為10%—13%,田間持水率為28%—30%(范月等,2020)。
本文研究土質為砂姜黑土,夏玉米生育期內逐日蒸散數據由稱重式蒸滲儀自動采集,蒸滲儀型號為FR101A,分辨率為0.025 mm,口徑面積4.0 m2,土柱高為4.0 m,蒸散數據每10分鐘記錄一次。蒸滲儀上方設有攝像機,可每時自動獲取作物生長圖像,便于劃分生長階段。蒸滲儀內10、30、50、100、200 cm埋深處各設有時域反射儀(TDR),通過土壤水分傳感器(型號為AV-EC5,分辨率為0.1%,精度為±3%,量程為 0—飽和,),每小時測量1次筒內土壤體積含水率。場地設有高精度氣象站,可每10分鐘獲取1次空氣濕度、風向、凈輻射等不同氣象數據,測量準確度為±5%,可用于計算參考作物蒸散量。
研究資料選取稱重式蒸滲儀地下水埋深為1 m、2018年6月22日—2018年10月8日的資料及同期氣象數據。由實測資料可知,2018年全年降雨量為1352.6 mm,最高日平均溫度為32.8 ℃。夏玉米全生育期內無明顯高溫干旱或洪澇等自然災害,作物僅接受自然降雨,無人工灌溉。蒸散量在8月8日、8月9日受人為活動干擾,數據異常,建模分析時應剔除并進行插補。
FAO推薦將生育期劃分為 4個生長階段:初期、發育期、中期和后期。根據實驗區夏玉米實際生長圖像,使用 PS軟件勾選作物像素以及土壤像素,像素之比便視為夏玉米的地面覆蓋率,將其各生長階段劃分見表1。

表1 夏玉米生長階段劃分Table 1 Summer maize growth stage division
用來計算 ET0的方法很多,最常用的是Penman-Monteith法,它將多種影響因素考慮在內,且在氣候差異較大的地區計算精度高,但該方法的缺點是需要的氣象資料較多,計算過程復雜,因此,聯合國糧農組織(FAO)修正了該公式,提出FAO Penman-Monteith公式來計算參考作物蒸散量(康燕霞等,2009)。公式如下:

式中,ET0為參考作物蒸散量(mm·d?1);Rn為地表凈輻射量[MJ·(m2·d)?1];G為土壤熱通量[MJ·(m2·d)?1];t為日平均氣溫(℃);u2為地面 2 m高處的平均風速(m·s?1);ps為飽和水汽壓(kPa);pa為實際水汽壓(kPa);Δ為飽和水汽壓與溫度曲線的斜率(kPa·℃?1);γ為干濕表常數(kPa·℃?1),各數據取自高精度氣象站(趙勇等,2005)。
土壤水分脅迫系數反映出土壤水分虧缺對作物蒸散發的抑制作用,本文用下式計算:

式中,Ksi為第i天的土壤水分脅迫系數;θi為第i天0—40 cm的平均土壤質量含水率;θ10i、θ30i為第i天10 cm、30 cm處的平均體積含水率;ρb為0—40 cm 的土壤容重(1.4 g·cm?3);θw為凋萎含水率,取10%;θf為田間持水率,取28%。
2.1.1 FAO單作物系數模型
根據有無土壤水分脅迫,分兩種情況具體討論:
在無土壤水分脅迫的條件時,作物的生長發育不會因為土壤中水分虧缺而受到抑制,其生長過程的影響因素主要為氣象條件、作物本身的生長發育特征以及土壤的物理特性,計算公式如下:

式中,ETc為實際蒸散量(mm·d?1);Kc為作物系數;ET0為參考作物蒸散量(mm·d?1)。
在有土壤水分脅迫的條件時,夏玉米的生長發育過程不僅會受到上述因素的影響,還會受到土壤水分虧缺的抑制作用,故添加水分脅迫系數Ks來調整實際蒸散量,計算公式如下:

2.1.2 遺傳算法模型
對于遺傳算法的原理,它是將實際問題的可能解集看作種群,每個個體都經過特定的基因編碼,并組成染色體。利用編碼產生初始種群后,采用適者生存、優勝劣汰的原理,根據適應度大小選擇存留的優秀個體,再組合交叉和變異,產生出新一代的種群,循環往復,逐代進化產生越來越好的近似解,它對結果沿多種路線平行搜索,會避免局部較優解,得到全局最優點,是一種全局最優化方法(喬均儉等,2007)。
算法流程為:設定初始參數;隨機產生初始種群;個體適應度評價;輪盤賭選擇存留個體;染色體交叉;染色體變異;對變異后的種群適應度值進行評估;終止條件的判斷,若終止,則輸出最優染色體及其適應度值,否則,返回輪盤賭再選擇新的存留個體。
本文遺傳算法模型的建立過程是:將夏玉米生育期分為4個階段,分別為初期、發育期、中期和后期。選擇初期、中期、后期訓練集中無水分脅迫的資料,優化變量為上述階段的作物系數,利用式(4)構造無水分脅迫時的目標函數,見式(6),再利用遺傳算法尋優求解最優作物系數:

約束條件為:

式中,n為各研究階段無土壤水分脅迫的天數;KCini、KCmin、KCend為初期、中期、后期的作物系數。
選用精度評價指標:平均絕對誤差Ema、均方根誤差Erms、相關系數r及準確率P;對于作物系數,定義絕對誤差在0.15、0.2、0.3以內的數據個數占總數據個數的比例為相應準確率,記為P0.15、P0.2和P0.3;對于蒸散量,定義絕對誤差在2、3、4 mm·d?1以內的數據個數占總數據個數的比例為相應準確率,記為P2、P3和P4。計算公式如下所示:

式中,yi為實測值;yi′為模型模擬值;為模擬值的平均值;為實測值的平均值;i為各研究階段的樣本個數,i=1, 2, …,m,m為各研究階段的樣本總數。
評價規則為:Ema、Erms越小,r越接近1,模型的擬合準確度越高,可信度越高,預測能力越強。
將夏玉米的各生長階段按7∶3的比例劃分訓練集與檢驗集,遺傳算法輸出初期、中期、后期作物系數最優解見表2,表中作物系數均為階段平均值。

表2 遺傳算法作物系數最優解Table 2 Optimal solution of crop coefficient based on genetic algorithm
由表2可知,初期作物系數最小,為0.280;中期作物系數最大,為0.611。整體上,從初期到中期,作物系數呈增加趨勢,從中期到后期,作物系數呈減小趨勢。繪制全生育期實測蒸散量與實際作物系數見圖 1,總體上,蒸散量與作物系數的變化過程是一致的。受人為活動的干擾,蒸散量在8月8日、8月9日顯示數據異常,蒸散量出現了極小值,9月3日、9月4日為大風天氣,董旭光等(2016)研究發現,對ET0影響較大的因子為風速,二者有顯著的正相關關系,這與王希等(2013)、張楊等(2018)的研究結果是一致的。故當風速較大時,蒸散量增加,ET0也會增大,此時作物系數會出現異常波動,所以在接下來建模分析時,應合理剔除異常數據,并用鄰近平均值進行插補。

圖1 全生育期蒸散量與作物系數Fig.1 Evapotranspiration and crop coefficient in whole growth period
3.1.1 作物系數的比較
由遺傳算法得到不同樣本模擬值與實測值的精度結果,見表 3。除發育期外,其他作物系數的計算誤差MAE、RMSE均小于0.2,說明利用遺傳算法得到的作物系數對實際值的擬合及預測精度較高,用該模型估算夏玉米蒸散量是可行的。

表3 不同樣本作物系數擬合精度值Table 3 Fitting accuracy of crop coefficients of different samples
FAO認為將初期與中期的作物系數進行線性插補,即可得到發育期作物系數的日變化值(宋迪等,2009)。本文針對發育期階段,若簡單地將遺傳算法輸出值進行線性插補,得到的發育期作物系數擬合值與真實值誤差較大,故需進一步對發育期作物系數線性插補值進行調整以提高精度。
3.1.2 發育期作物系數的修正
有學者研究表示作物系數與當地的植被指數存在聯系(Campos et al.,2017),受此啟發,繪制全生育期作物系數與葉面積指數實際值圖,見圖2。葉面積指數在全生育期內呈現先增后減的規律,其變化趨勢與作物系數變化基本一致。在夏玉米生長前期,即本文所劃分的初期與發育期階段,葉面積指數均呈增加的趨勢,尤其在發育期階段,葉面積指數增大明顯,玉米蒸騰作用變強,總蒸散量主要來自作物蒸騰,作物系數隨著作物的生長而變化,且隨著時間的推移,葉面積指數與作物系數呈現出較好的一致性,故引入葉面積指數和發育期天數來構建作物系數修正模型是合理的。

圖2 全生育期作物系數實際值與葉面積指數實際值Fig.2 Actual crop coefficient and actual LAI in whole growth period
考慮先將作物系數線性插補值與葉面積指數相乘,根據其變化特征,再將該乘積與某個指數函數、多項式函數、冪函數或有理函數等結合,選擇精度最高的函數,從而實現對發育期作物系數插補值的修正,模型式如下:

式中,Kc′為修正后的發育期作物系數;Kc″為發育期作物系數的線性插補值;LAI為發育期葉面積指數;f(x)為上述備選函數之一,x為發育期天數,x≤35 且x∈N*。
利用 MATLAB軟件,根據最小二乘法原理,對模型式中的f(x)進行擬合,得到各f(x)的殘差平方和SSE,擬合精度R2和均方根誤差RMSE見表4。

表4 不同f (x)的擬合精度值Table 4 Fitting accuracy of different f (x)
由表4可知,嘗試多種函數形式對f(x)進行擬合,R2均高于0.8,說明作物系數與葉面積指數、發育期天數之間確實存在顯著的函數關系。其中有理函數的擬合精度最高,R2達到0.92,故選擇有理函數加入模型式(11),最終修正模型式如下:

式中符號含義同式(11)。
利用式(12)修正作物系數模型,各精度值見表5。

表5 作物系數修正前后擬合精度值Table 5 Fitting accuracy of crop coefficient before and after adjustment
由表5可知,發育期的作物系數線性插補值經過模型修正后,其與實際值的誤差有明顯減小,且準確率P0.15、P0.2和P0.3有明顯提高,說明修正模型可用于發育期作物系數的調整。
3.1.3 修正后蒸散量的比較
圖3是訓練樣本中,遺傳算法預測蒸散量與實際值的對比。二者平均絕對誤差為0.786 mm·d?1,均方根誤差為 1.12 mm·d?1,P2為 91.9%,P3為97.3%,P4為100%,借助單作物系數法,修正后的遺傳算法模型得到的蒸散量與實際值擬合較好。

圖3 訓練集計算蒸散量與實際值比較Fig.3 Calculated ET compared with actual value in training set
圖4是檢驗樣本中,預測蒸散量與實際值的對比。二者的擬合程度較高,平均絕對誤差為 1.148 mm·d?1,均方根誤差為 1.61 mm·d?1,P2為 77.4%,P3為90.3%,P4為100%,說明本文模型對夏玉米蒸散量的預測效果較好。

圖4 檢驗集計算蒸散量與實際值比較Fig.4 Calculated ET compared with actual in testing set
表6是修正后遺傳算法模型模擬蒸散量的精度評價表。由表3、5、6可知,修正后的遺傳算法模型在擬合夏玉米作物系數、蒸散量方面精度較高,各誤差均保持在較小水平內。為進一步確定模型的適用性,本文將修正遺傳算法模型與FAO推薦作物系數模型的擬合精度進行比較。

表6 不同樣本蒸散量擬合精度值Table 6 Fitting accuracy of evapotranspiration of different samples
3.2.1 作物系數的對比
FAO推薦的夏玉米各生長階段作物系數分別為:初期0.3,中期1.2,后期0.6(劉鈺等,2000);計算全生育期逐日作物系數,兩種模型的擬合精度見表 7,可以看出,使用修正遺傳算法模型得到的作物系數擬合精度更高,可以滿足精度要求。

表7 全生育期不同模型作物系數擬合精度值Table 7 Fitting precision values of crop coefficients of different models during the whole growth period
3.2.2 蒸散量的對比
繪制上述兩種模型蒸散量模擬值與實際值的對比圖,見圖 5。結果表明,兩種模型擬合的蒸散量變化趨勢基本一致;FAO模型的誤差在初期、發育期較小,在中期、后期逐漸增大;修正后的遺傳算法模型,其模擬值更加接近實際值,究其原因,FAO推薦的夏玉米初期、中期、后期的作物系數分別為0.3、1.2、0.6,而本文根據實測數據擬合的初期、中期、后期作物系數分別為0.28、0.611、0.386,FAO推薦的作物系數除初期外明顯偏高,進而高估了夏玉米的實際蒸散量;相比之下,本文利用遺傳算法模型擬合的作物系數更加接近實際情況,且結合了復合函數對估算誤差偏高的生育階段進行優化,這使得本文模型對夏玉米蒸散量的估算誤差更小,可滿足一定的精度要求。全生育期內,修正遺傳算法模型的平均絕對誤差為0.893 mm·d?1,均方根誤差為 1.284 mm·d?1,P4為 100%,r為 0.863;FAO推薦作物系數模型的平均絕對誤差為 3.841 mm·d?1,均方根誤差為 5.28 mm·d?1,P4為 58.1%,r為0.654。本文模型在擬合夏玉米作物系數、蒸散量方面均有較高精度,可用于夏玉米蒸散量估算。

圖5 全生育期不同方法ET計算值與實際值比較Fig.5 Calculated ET compared with actual value based on different ways during the whole growth period
在無土壤水分脅迫條件下,作物系數是作物實際蒸散量與參考作物蒸散量的比值,其數值可反映出因作物種類差異、土壤水肥條件、作物生長發育狀況等差異對農田蒸散量的影響(Alberto et al.,2014)。本文采用單作物系數法計算夏玉米實際作物系數,相比于雙作物系數法,單作物系數法可使用較少的數據而得到較高精度的作物系數估計(曹永強等,2019)。本文將夏玉米全生育期劃分為4個生長階段,并利用無受旱脅迫的實測蒸散資料和同期氣象數據,結合遺傳算法得到無受旱脅迫下的玉米初期、中期和后期作物系數,根據葉面積指數和發育期天數構造修正模型,對誤差較大的發育期作物系數進行調整,得到夏玉米各生長階段的作物系數估算值,進一步可估算夏玉米的蒸散量。
由本文模型得到的夏玉米作物系數,全生育期呈現出先增后降的變動趨勢,這與肖然等(2019)的研究結果相一致,且估算蒸散量的平均絕對誤差小于1 mm·d?1,模擬效果較好。進一步將該模型與FAO推薦作物系數模型進行對比,在估算蒸散量方面,FAO推薦作物系數模型在初期的估算誤差較小,中期、后期估算誤差較大;而本文模型的估算值更加符合實際蒸散量。一方面,在生長初期,玉米植株矮小,葉片不發達,農田蒸散量主要來自于土壤蒸散,植株蒸騰占比很少,且總蒸散量很小,所以在生長初期,無論哪種模型的估算誤差都較小;到發育期,隨著光熱條件逐漸充足,玉米生長發育迅速,葉片快速伸展,由于在中期、后期,FAO推薦的作物系數值偏大,故從該階段開始,FAO推薦作物系數模型的估算誤差開始增大,而本文模型借助復合函數對算法輸出的作物系數進行了修正,縮小了誤差,優勢得以體現;到中期和后期,FAO推薦的作物系數明顯偏高,使得蒸散量估算誤差進一步增大,而本文模型立足于實測資料,對蒸散量估算的均方根誤差小于1.5 mm·d?1,可滿足一定的精度要求,可作為估算夏玉米蒸散量的一種方法。另一方面,FAO推薦的作物系數值是針對于作物生長階段均無水分脅迫、耕作和田間水分管理條件良好,作物能獲得最大收成時的值,而作物實際生長狀況會受多種因素干擾,這就使得FAO推薦的值會高于實際值,從而增大了誤差(袁宏偉等,2018)。
本文修正遺傳算法模型對夏玉米蒸散量的估算誤差可滿足一定的精度要求,不失為一種估算夏玉米蒸散量的新思路,但本文計算是基于地下水埋深為1 m的砂姜黑土實測數據,由于不同土質、不同地下水埋深會對作物系數產生不同的影響,不同作物的生長發育狀況也不盡相同,故基于不同土質、不同地下水埋深、不同作物生長情況下的蒸散量估算需進一步研究。
(1)本文構建了基于遺傳算法的夏玉米蒸散估算模型,并利用葉面積指數構造修正模型,對發育期作物系數進行數值修正,模型對夏玉米全生育期作物系數的擬合、預測精度均較高。根據修正遺傳算法模型,訓練集中作物系數模擬值與實際值的MAE為0.083,RMSE為0.112,P0.3為97.3%;檢驗集中兩者MAE為0.101,RMSE為0.805,P0.3為93.5%。本文模型可以提高作物系數的計算精度,可用于夏玉米作物系數的模擬。
(2)基于修正模型得到的夏玉米蒸散量計算精度較高。訓練集中蒸散量模擬值與實際值的 MAE為 0.786 mm·d?1,RMSE 為 1.120 mm·d?1,r為 0.898;檢驗集中兩者MAE為1.148 mm·d?1,RMSE為1.610 mm·d?1,r為0.813,蒸散量的擬合精度均高于 FAO推薦作物系數模型。
(3)與FAO推薦作物系數模型相比,全生育期內,修正遺傳算法模型作物系數和蒸散量的擬合誤差、準確率均有明顯提高;對全生育期內逐日作物系數與蒸散量的計算,本文提出的修正遺傳算法模型模擬精度符合精度要求。