王鵬新 喬 琛 李 俐 周西嘉 許連香 胡亞京

(1.中國農(nóng)業(yè)大學(xué)信息與電氣工程學(xué)院， 北京 100083； 2.中國農(nóng)業(yè)大學(xué)土地科學(xué)與技術(shù)學(xué)院， 北京 100083)

0 引言

糧食的生產(chǎn)與安全對提升我國農(nóng)業(yè)的經(jīng)營管理水平、完善農(nóng)作物的種植結(jié)構(gòu)和確保糧食安全等具有重要意義[1-3]，因此及時準(zhǔn)確估測農(nóng)作物的產(chǎn)量有利于維持國家穩(wěn)定和促進經(jīng)濟發(fā)展。

干旱作為影響農(nóng)作物生長發(fā)育及產(chǎn)量的重要因素已受到廣泛的關(guān)注[4-5]。王鵬新等[6]在歸一化植被指數(shù)(Normalized difference vegetation index，NDVI)和地表溫度(Land surface temperature，LST)的散點圖呈三角形區(qū)域分布的基礎(chǔ)上，提出了條件植被溫度指數(shù)(Vegetation temperature condition index，VTCI)的干旱監(jiān)測方法，可用于反映農(nóng)作物生長過程中的水分虧缺信息。孫威等[7]對VTCI冷熱邊界的確定方法進行了完善，并驗證了利用 VTCI 進行干旱監(jiān)測的可行性。農(nóng)作物產(chǎn)量不僅受到水分脅迫的影響，還與作物的生長狀態(tài)密切相關(guān)[8]，葉面積指數(shù)(Leaf area index，LAI)通過對作物的光合速率和干物質(zhì)的積累量的反映能較好地表征農(nóng)作物的生長狀態(tài)[9]。

遙感技術(shù)對比傳統(tǒng)統(tǒng)計調(diào)查方法，憑借其覆蓋范圍廣、重訪周期短等獨特優(yōu)勢被廣泛應(yīng)用于作物長勢監(jiān)測及產(chǎn)量估測，同時對大規(guī)模的農(nóng)業(yè)生產(chǎn)調(diào)查、評價、監(jiān)測和管理具有獨特的作用[10]。隨著大數(shù)據(jù)技術(shù)的發(fā)展，遙感大數(shù)據(jù)將推動農(nóng)業(yè)遙感估產(chǎn)的發(fā)展[11]。目前，遙感估產(chǎn)方法中常采用的統(tǒng)計模型、機理模型和半機理模型等均能夠較好地對作物進行估產(chǎn)[12]。但由于在實際應(yīng)用中機理和半機理模型存在需要輸入較多的參數(shù)問題，因此機理和半機理模型存在一定的局限性。而統(tǒng)計模型其估產(chǎn)的精確度依賴于選取遙感影像的時相，對作物的生長和產(chǎn)量形成的機理解釋性不強[12]，因此在實際應(yīng)用中同樣具有一定局限性。在作物生長過程中，經(jīng)常受到各種因素影響，同時，這些因素在作物不同生育時期產(chǎn)生不同的影響，即使采用相同的估算方式也會得到不同的結(jié)果與估測精度。雖然其估產(chǎn)精度不同，但是這些模型可以提供不同的有用信息，如將一些方法丟棄，就會失去有用信息，從而影響估產(chǎn)的精度[13]。因此，可以將單個估測模型提供的有用信息進行組合，以提高估產(chǎn)的精度。組合預(yù)測(Combination forecasting，CF)將不同預(yù)測模型進行有效組合，可視為對無限逼近真實數(shù)據(jù)生成過程的有效補充[14]。

BATES等[15]首次提出了組合預(yù)測理論，這一預(yù)測理論在提高精度的同時更充分利用了預(yù)測樣本所表達的信息，受到了國內(nèi)外廣大學(xué)者重視。組合預(yù)測根據(jù)各單一模型的信息貢獻度，進而基于計算得到的單一模型權(quán)重構(gòu)建組合預(yù)測模型，從而實現(xiàn)減少預(yù)測誤差、提高預(yù)測精度的目標(biāo)[16]。特別是在時間序列中分析真實數(shù)據(jù)生成過程中，通常具有區(qū)制轉(zhuǎn)換(Regime shift)或參數(shù)漂移等特性。組合預(yù)測方法的引入，可減少由參數(shù)錯誤或模型錯誤帶來的預(yù)測誤差[17]，甚至在單一預(yù)測結(jié)果存在有偏性的情況下，通過組合能產(chǎn)生具有無偏性的預(yù)測結(jié)果[18]。因此，組合模型具有普遍性的優(yōu)點，最終的預(yù)測結(jié)果更接近實際數(shù)據(jù)。根據(jù)組合預(yù)測方法綜合手段的不同，可分為權(quán)重綜合法和區(qū)域綜合法兩種。與區(qū)域綜合法相比，權(quán)重綜合法得到了廣泛的應(yīng)用。目前，組合預(yù)測現(xiàn)在已經(jīng)在多個領(lǐng)域內(nèi)得到了應(yīng)用[19-20]，但是在農(nóng)業(yè)領(lǐng)域還少有報道。

本文以河北中部平原為研究區(qū)域，選取與玉米長勢和產(chǎn)量密切相關(guān)的VTCI和LAI為特征變量，采用極限梯度提升樹(Extreme gradient boosting，XGBoost)和隨機森林(Random forest，RF)兩種機器學(xué)習(xí)算法模型分別估測研究區(qū)域的玉米單產(chǎn)，并借鑒經(jīng)濟學(xué)中合作博弈論Shapley值利益分配理論，通過組合預(yù)測模型中的權(quán)重合成思想，以單一預(yù)測模型的均方誤差為基礎(chǔ)確定單一模型權(quán)重得到組合預(yù)測模型，以期為玉米長勢監(jiān)測和產(chǎn)量估測提供新方法。

1 材料與方法

1.1 研究區(qū)概況

河北中部平原屬黃淮海平原組成部分，是我國重要的玉米種植和生產(chǎn)基地[21]，位于114°32′～117°36′E與36°57′～39°50′N之間，包括石家莊等5個市的53個縣(區(qū))，其土地面積約為5.30×104km2。河北中部平原地處典型溫帶大陸性季風(fēng)氣候區(qū)，四季分明、雨熱同期。河北省中部平原年降水量350～800 mm。降水時空分布不均，南方的降水量比北方高，夏季降水量高于冬季，豐水年降水量與枯水年降水量相差較大。該地區(qū)的耕種管理制度為一年兩熟，結(jié)合該地區(qū)玉米實際生長狀況，本文將河北中部平原玉米生長劃分為4個生育時期：7月上旬至7月中旬的出苗-拔節(jié)期、7月下旬至8月上旬的拔節(jié)-抽穗期、8月中旬至9月上旬的抽穗-乳熟期和9月中旬至9月下旬的乳熟-成熟期。根據(jù)王鵬新等[22]提出的作物分類方法，進而獲得研究區(qū)域2010—2018年玉米種植區(qū)分布圖，其中河北中部平原2017年玉米種植區(qū)分布圖如圖1所示。

1.2 數(shù)據(jù)獲取與處理

1.2.1VTCI時間序列的生成

選取2010—2018年每年7—9月空間分辨率和時間分辨率分別為1 000 m、1 d的MODIS日地表溫度產(chǎn)品(MYD11A1)及日地表反射率產(chǎn)品(MYD09GA)。利用MRT對日地表溫度和日地表反射率產(chǎn)品進行預(yù)處理之后，得到研究區(qū)域的日LST和日NDVI產(chǎn)品，運用最大值合成法生成旬LST與NDVI的最大值合成產(chǎn)品；基于生成的多年某一旬的NDVI和LST最大值合成產(chǎn)品再次使用最大值合成技術(shù)生成多年的旬NDVI和LST的最大值合成產(chǎn)品；基于生成的多年某一旬的LST最大值合成產(chǎn)品利用最小值合成技術(shù)逐像素提取最小值，計算得到多年旬LST的最大—最小值合成產(chǎn)品，并以此通過計算生成河北中部平原旬VTCI時間序列數(shù)據(jù)。VTCI的計算公式為[6,23]

(1)

(2)

LSTmax(NDVIi)=a+bNDVIi

(3)

LSTmin(NDVIi)=a′+b′NDVIi

(4)

式中NDVI——歸一化植被指數(shù)

ρNIR、ρred——近紅外、紅光波段的反射率

LST——地表溫度

LSTmax(NDVIi)、LSTmin(NDVIi)——當(dāng)NDVIi為某一特定值時所有像素地表溫度最大值和最小值，即熱邊界和冷邊界

LST(NDVIi)——研究區(qū)域內(nèi)某一像素的NDVI值為NDVIi時的地表溫度

a、b、a′、b′——LST和NDVI散點圖近似得到的待定系數(shù)

1.2.2LAI時間序列的生成

選取2010—2018年7—9月的MODIS葉面積指數(shù)產(chǎn)品(MCD15A3H)，其空間分辨率和時間分辨率分別為500 m、4 d。利用上包絡(luò)線Savitzky-Golay濾波對經(jīng)MRT處理后得到的原始LAI產(chǎn)品進行平滑處理以消除云層覆蓋、大氣溶膠等因素引起的數(shù)據(jù)驟降的現(xiàn)象[23]，濾波處理后的葉面積指數(shù)變化趨于平穩(wěn)且更加符合玉米的生長物侯特征，解決了原始數(shù)據(jù)存在的質(zhì)量問題。由于LAI與VTCI時間分辨率和取值范圍不同，因此對濾波后的葉面積指數(shù)進行歸一化處理，并通過觀察多時相MODIS的MCD15A3H原始數(shù)據(jù)相元統(tǒng)計直方圖將取值范圍設(shè)置為0～7，進而計算得到研究區(qū)域2010—2018年玉米主要生育時期的LAI時間序列數(shù)據(jù)。

1.2.3玉米生育時期VTCI和LAI的計算及異常點數(shù)據(jù)處理

基于研究區(qū)玉米生育時期的劃分結(jié)果，將生育時期內(nèi)所包含的多旬VTCI和LAI平均值作為研究區(qū)域玉米該生育時期逐像素的VTCI和LAI值；通過疊加河北中部平原各縣(區(qū))行政邊界圖，將各縣(區(qū))所包含的玉米生育時期逐像素VTCI和LAI值的平均值作為生育時期縣(區(qū))尺度的VTCI和LAI值。同時在構(gòu)建回歸模型時，剔除加權(quán)VTCI與LAI殘差置信區(qū)間在[-4 000,4 000] kg/km2以外的異常點數(shù)據(jù)。

1.3 研究方法

1.3.1極限梯度提升算法

極限梯度提升算法(XGBoost)是一種基于梯度提升的決策樹(Gradient boosting decision tree)集成算法[24]，其基分類器主要包含分類和回歸樹(Classification and regression tree，CART)。本文玉米單產(chǎn)估測是回歸問題，所以其基模型選擇為回歸樹。含有K棵決策樹XGBoost的樹集成模型定義為

(5)

xi——樣本所對應(yīng)的特征變量VTCI與LAI

fk——第k個決策樹的預(yù)測函數(shù)

樹集成優(yōu)化模型可以定義為

(6)

(7)

l(yi,i)=(yi-i)2

式中l(wèi)(yi,i)——損失函數(shù)，即均方誤差

Ω(f)——正則化項

γ——復(fù)雜度參數(shù)

T——樹中葉子節(jié)點個數(shù)

λ——固定系數(shù)

ω——葉子節(jié)點量化權(quán)重向量

求解式(6)的優(yōu)化問題，即求解CART樹的結(jié)構(gòu)。通過保留訓(xùn)練好的前t-1輪樹模型不變，在第t輪時添加一個新的預(yù)測函數(shù)，迭代計算得到最終預(yù)測結(jié)果[25]

(8)

ft(xi)——第t輪加入的新的預(yù)測函數(shù)

XGBoost模型可以通過輸出玉米4個生育時期VTCI或LAI的特征變量重要性來評估不同特征變量對玉米產(chǎn)量的影響程度。在XGBoost中常用基于增益、覆蓋度、頻率的特征重要性指標(biāo)進行特征重要性評價，其中基于增益(gain)表示各特征變量對XGBoost模型中每棵樹采取每個特征變量的貢獻而計算的模型相對貢獻度為

(9)

式中G——VTCI與LAI基于增益獲取的特征重要性得分值

GL、GR——左、右子樹中的梯度

HL、HR——左、右子樹中的二階梯度

c——給定的臨界值

1.3.2隨機森林算法

隨機森林(RF)是一種基于袋裝法(bagging)理論實現(xiàn)的集成學(xué)習(xí)算法，其基模型是決策樹模型。隨機森林具有極好的準(zhǔn)確率，同時可以評估各特征在回歸問題上的重要性。隨機森林將若干沒有聯(lián)系的回歸決策樹{y(x,θk),k=1,2，…,K}構(gòu)成K棵集成決策樹[26]，可以定義為

(10)

式中y(x)——估測研究區(qū)域玉米的單產(chǎn)

x——特征變量的輸入，即VTCI與LAI

文獻[26]基于RF算法確定了河北中部平原玉米各個生育時期的權(quán)重，進而構(gòu)建了加權(quán)VTCI和LAI與玉米單產(chǎn)之間的雙參數(shù)估產(chǎn)模型，結(jié)果顯示，雙參數(shù)估產(chǎn)模型精度相對較高，達到極顯著的水平(P<0.001)，且基于隨機森林算法確定的玉米各生育時期權(quán)重均較為合理。因此，本文引用其基于隨機森林算法確定的河北中部平原玉米各生育時期的VTCI和LAI權(quán)重進行估產(chǎn)模型的構(gòu)建。

1.3.3組合預(yù)測模型

XGBoost算法與RF算法均為統(tǒng)計理論的機器學(xué)習(xí)方法，但二者各有優(yōu)勢。XGBoost算法通過引入正則化項在很大程度上避免了過擬合問題的出現(xiàn)，同時采用優(yōu)化算法降低了問題計算復(fù)雜度，在處理大規(guī)模數(shù)據(jù)集方面具有明顯的優(yōu)勢，已被廣泛應(yīng)用于不同領(lǐng)域的研究[27-29]。相對于XGBoost，RF更適合處理高維數(shù)據(jù)，對噪聲數(shù)據(jù)的容忍性較高。通過結(jié)合XGBoost算法和RF算法的優(yōu)點，參考經(jīng)濟學(xué)中合作博弈論的Shapley值理論，通過計算確定組合預(yù)測模型總誤差在單一預(yù)測模型之間的分布，并在此基礎(chǔ)上確定各單一預(yù)測模型的權(quán)重，具體步驟如下：

(1)求取Shapley值

(11)

(12)

式中i——預(yù)測模型序號，即XGBoost與RF預(yù)測模型序號

E′i——第i個預(yù)測模型的Shapley值，即XGBoost與RF預(yù)測模型計算得到的均方誤差

s——包含模型XGBoost與RF的所有集合

|s|——預(yù)測模型個數(shù)

E(s)——集合s的收益

E(s-{i})——集合s中去除成員i后的收益

(2)計算權(quán)重

組合預(yù)測模型中若某單一預(yù)測模型所獲得的誤差分配值越大，表示預(yù)測精度越低，在組合預(yù)測模型中的權(quán)重就越小。基于此原則，預(yù)測模型i的權(quán)重λi定義為

(13)

最終的組合預(yù)測模型可描述為

(14)

式中Y——組合預(yù)測模型估測的玉米單產(chǎn)

Yi——模型XGBoost與模型RF估測的玉米單產(chǎn)

以通過計算得到的2010—2017年玉米4個生育時期的VTCI和LAI數(shù)據(jù)作為特征變量，相對應(yīng)的玉米單產(chǎn)數(shù)據(jù)作為目標(biāo)變量(每年53組數(shù)據(jù)，共424組數(shù)據(jù))，以2010—2017年(除2012年)數(shù)據(jù)作為訓(xùn)練集合，2012年數(shù)據(jù)作為測試集合。由于特征變量與目標(biāo)變量之間差值較大，特征變量值的范圍在0～1之間，而目標(biāo)變量值在3 000～9 000 kg/km2之間。由于特征變量和目標(biāo)變量處于不同區(qū)間的值域差異可能對結(jié)果造成不同的影響，取對數(shù)后不會改變數(shù)據(jù)的性質(zhì)和相對關(guān)系，對訓(xùn)練數(shù)據(jù)集和測試數(shù)據(jù)集中的目標(biāo)變量進行取對數(shù)處理，并通過對估產(chǎn)模型輸出的估測數(shù)據(jù)進行取對數(shù)還原，得到最終的估測單產(chǎn)數(shù)據(jù)。采用均方根誤差(RMSE)、平均相對誤差(MRE)和決定系數(shù)(R2)等指標(biāo)對模型的精度進行評價。

2 結(jié)果與分析

2.1 估產(chǎn)模型構(gòu)建

構(gòu)建XGBoost估產(chǎn)模型與RF估產(chǎn)模型，借鑒Shapley值理論，計算得到XGBoost與RF估產(chǎn)模型的Shapley值。根據(jù)Shapley值確定XGBoost估產(chǎn)模型與RF估產(chǎn)模型權(quán)重，進而完成對組合預(yù)測模型的構(gòu)建。

2.1.1XGBoost與RF的估產(chǎn)模型構(gòu)建

基于XGBoost與RF算法輸出的特征重要性值，通過歸一化計算得到研究區(qū)域玉米各生育時期VTCI與LAI的權(quán)重(表1)?？梢钥闯?，玉米生育后期的LAI權(quán)重高于玉米生育前期LAI的權(quán)重，表明玉米生育前期的LAI對玉米產(chǎn)量影響程度較弱，玉米生育后期的LAI對玉米產(chǎn)量影響程度較強。原因可能是玉米LAI變化規(guī)律呈前期增長緩慢、中期增長快速、后期下降緩慢的偏峰曲線，其中玉米在出苗-拔節(jié)期和拔節(jié)-抽雄期主要以分化莖葉的營養(yǎng)生長為主，此階段玉米生長迅速，葉片迅速增多增大；以抽雄-乳熟期為界，玉米進入以生殖生長為主的生育后期，光合產(chǎn)物的分配模式主要以果穗為中心，是玉米產(chǎn)量形成的重要時期。因此玉米生育后期的LAI對玉米產(chǎn)量的影響程度強于玉米生育前期。拔節(jié)-乳熟期的VTCI權(quán)重高于出苗-拔節(jié)期和乳熟-成熟期的權(quán)重，表明拔節(jié)-乳熟期的VTCI對玉米產(chǎn)量的影響程度較強，出苗-拔節(jié)期和乳熟-成熟期的VTCI對玉米產(chǎn)量的影響程度較弱。原因可能是在拔節(jié)-抽雄期、抽雄-乳熟期玉米進入營養(yǎng)生長的高峰期，此階段玉米生長迅速，對土壤中的水分吸收也最為急迫，若發(fā)生水分脅迫將減緩玉米根莖葉的生長發(fā)育，降低光合作用對玉米干物質(zhì)積累速率，減少蛋白質(zhì)和有機質(zhì)的合成，造成玉米粒重明顯降低，最終影響玉米產(chǎn)量，因此拔節(jié)-抽雄期、抽雄-乳熟期的VTCI對玉米產(chǎn)量的影響程度較強。而在出苗-拔節(jié)期，玉米由于植株矮小，對水分的需求量較少，在乳熟-成熟期，玉米處于生育后期，生長變緩，對一定的水分脅迫表現(xiàn)出一定的忍受力，因此出苗-拔節(jié)期、乳熟-成熟期VTCI對玉米產(chǎn)量的影響程度較弱。綜上所述，基于XGBoost算法與RF算法確定的研究區(qū)域玉米各生育時期的權(quán)重均較為合理。

表1 玉米各生育時期的權(quán)重Tab.1 Weight results of each growth stage of maize

分別將XGBoost估產(chǎn)模型與RF估產(chǎn)模型得到的玉米各生育時期的特征權(quán)重進行加權(quán)VTCI和LAI計算，進而構(gòu)建基于加權(quán)VTCI和LAI與玉米單產(chǎn)之間的回歸估測模型，將2010—2017年(除2012年)數(shù)據(jù)代入XGBoost估產(chǎn)模型與RF估產(chǎn)模型，進行可視化分析(圖2)。可以看出，XGBoost估產(chǎn)模型R2為0.31，均方根誤差為940.91 kg/km2，RF估產(chǎn)模型R2為0.30，均方根誤差為947.50 kg/km2，XGBoost估產(chǎn)模型與RF估產(chǎn)模型均通過顯著性檢驗(P<0.001)。結(jié)果表明，XGBoost估產(chǎn)模型與RF估產(chǎn)模型玉米單產(chǎn)估測精度均較高，可用于研究區(qū)域各縣(區(qū))玉米的單產(chǎn)估測。

2.1.2組合預(yù)測模型的構(gòu)建

為進一步提高玉米單產(chǎn)估測精度，進行組合預(yù)測模型的構(gòu)建。需要確定XGBoost與RF估產(chǎn)模型的權(quán)重，基于2010—2017年(除2012年)由單一估產(chǎn)模型輸出的各縣(區(qū))的估測單產(chǎn)數(shù)據(jù)(共371組)，分別計算XGBoost估產(chǎn)模型、RF估產(chǎn)模型和組合預(yù)測模型的均方誤差，將各模型的均方誤差代入式(11)，將得到的XGBoost與RF估產(chǎn)模型的Shapley值代入式(13)，則可以得到XGBoost與RF估產(chǎn)模型的權(quán)重分別為0.48與0.52，則組合預(yù)測模型為

Y=0.48YXGBoost+0.52YRF

(15)

通過計算，組合預(yù)測模型的R2為0.32，模型通過顯著性檢驗且達顯著水平(P<0.001)。同時，組合預(yù)測模型均方根誤差為939.81 kg/km2，與單一XGBoost估產(chǎn)模型與RF估產(chǎn)模型相比，組合預(yù)測模型決定系數(shù)與均方根誤差均得到提升。結(jié)果表明，組合模型的玉米單產(chǎn)估測精度優(yōu)于單一估產(chǎn)模型，可用于河北中部平原的玉米單產(chǎn)估測。

2.2 玉米估測單產(chǎn)的精度評價

將2012年數(shù)據(jù)分別代入XGBoost估產(chǎn)模型、RF估產(chǎn)模型與組合預(yù)測模型，將玉米的的估測單產(chǎn)與實際單產(chǎn)數(shù)據(jù)進行數(shù)據(jù)可視化分析(圖3)。對比分析，可以看出XGBoost估產(chǎn)模型、RF估產(chǎn)模型和組合預(yù)測模型的估測單產(chǎn)與實際單產(chǎn)均呈顯著的正相關(guān)關(guān)系。同時，R2均不小于0.50，可以反映實際單產(chǎn)的波動均有50%以上能被估測單產(chǎn)的波動所描述，即玉米的實際單產(chǎn)與估測單產(chǎn)之間的誤差較小。其中，組合預(yù)測模型R2達到0.52，均方根誤差為831.14 kg/km2，平均相對誤差為9.86%，對比單一XGBoost估產(chǎn)模型與RF估產(chǎn)模型，組合預(yù)測模型的決定系數(shù)、均方根誤差與平均相對誤差均得到提升。進一步表明，基于Shapley值的組合預(yù)測模型的估產(chǎn)精度優(yōu)于單一估產(chǎn)模型。

2.3 區(qū)域玉米單產(chǎn)的估測與時空變化分析

將組合估產(chǎn)模型應(yīng)用于河北中部平原2010—2018年逐像素玉米單產(chǎn)估測(圖4)。結(jié)果表明，研究區(qū)域玉米的估測單產(chǎn)隨年際變化呈現(xiàn)先減少后增加的波動變化趨勢。其中，2010—2014年玉米估測單產(chǎn)整體呈現(xiàn)逐年下降的趨勢，并在2014年玉米平均單產(chǎn)達到最低，平均單產(chǎn)在6 000 kg/km2左右，原因可能是2014年河北中部平原降水較少且發(fā)生階段性干旱，導(dǎo)致玉米單產(chǎn)減少嚴(yán)重；2010、2011、2013年平均估測單產(chǎn)相差不大，均在6 350 kg/km2左右。對比2010—2014年玉米估測單產(chǎn)年際變化趨勢，2015—2018年的玉米單產(chǎn)整體呈逐年上升的趨勢，其中2018年玉米平均估測單產(chǎn)最高，約為6 500 kg/km2，2016年和2017年玉米平均估測單產(chǎn)次之，約為6 400 kg/km2，2015年玉米平均估測單產(chǎn)達到最低，約為6 200 kg/km2。

研究區(qū)域玉米的單產(chǎn)估測空間上呈現(xiàn)西部地區(qū)玉米估測單產(chǎn)最高，南部和北部地區(qū)玉米單產(chǎn)次之，東部地區(qū)玉米單產(chǎn)最低的分布特征。河北中部平原的北部地區(qū)中，2012年玉米的估測單產(chǎn)最高，2014年玉米估測單產(chǎn)最低，分別約為6 000 kg/km2和4 500 kg/km2，2015、2016、2017年玉米估測單產(chǎn)相差不大，均約為5 000 kg/km2，2011年和2018年玉米的估測單產(chǎn)略低于2012年；東部地區(qū)中，2016年玉米的估測單產(chǎn)最高，2014年玉米的估測單產(chǎn)最低，分別約5 000 kg/km2和4 500 kg/km2，2016—2018年東部地區(qū)的玉米估測單產(chǎn)相差不大，均約為4 700 kg/km2。南部地區(qū)玉米的估測單產(chǎn)分別在2012年和2014年達到最高和最低，分別約6 500、4 500 kg/km2，其余年份玉米估測單產(chǎn)相差不大。西部地區(qū)中，2012、2015年玉米的估測單產(chǎn)較高且相差不大，均約為6 700 kg/km2，2010、2014年玉米估測單產(chǎn)較低，均約為4 500 kg/km2。

3 討論

經(jīng)過對不同估產(chǎn)模型特點的研究分析發(fā)現(xiàn)，在玉米單產(chǎn)估測過程中，不同估產(chǎn)模型對特征因素的信息提取方式不完全一樣。同時以往單一估產(chǎn)模型均注重特征因素對產(chǎn)量的影響，而缺少考慮估產(chǎn)模型對特征因素的信息提取以及對估測產(chǎn)量的影響。本文在考慮VTCI和LAI作物長勢監(jiān)測指標(biāo)的基礎(chǔ)上，也充分考慮了估產(chǎn)模型對特征因素的信息提取以及對估測產(chǎn)量的影響，利用RF對噪聲數(shù)據(jù)的容忍性較高以及XGBoost較低的計算復(fù)雜度等優(yōu)點，借鑒經(jīng)濟學(xué)中Shapley值理論計算得到XGBoost和RF估產(chǎn)模型的權(quán)重，進而構(gòu)建組合估產(chǎn)模型，實現(xiàn)了單一估產(chǎn)模型的綜合利用。在今后的研究中，可以在Shapley值理論的基礎(chǔ)上嘗試加入其他估產(chǎn)模型，以期進一步全面地提取特征因素的信息，避免單一估產(chǎn)模型在估產(chǎn)過程中造成有用信息的浪費，從而實現(xiàn)玉米單產(chǎn)估測精度的進一步提升。

本文雖然選取了與玉米長勢和單產(chǎn)密切相關(guān)的VTCI和LAI作為特征變量，但玉米的生長發(fā)育除了受到水分脅迫和生長狀態(tài)的影響之外，還受到其他自然因素和人為因素的影響，例如溫度、土壤肥力、田間管理水平等因素。因此在基于VTCI與LAI作物生長指標(biāo)的基礎(chǔ)上，未來研究應(yīng)進一步綜合考慮與玉米單產(chǎn)相關(guān)性較大的其他因素。此外，基于Shapley值構(gòu)建的組合預(yù)測模型對河北中部平原各縣(區(qū))玉米單產(chǎn)估測精度雖較高，但缺少對農(nóng)學(xué)先驗知識的考慮，今后研究中可以通過主觀賦權(quán)法進一步修正XGBoost估產(chǎn)模型和RF估產(chǎn)模型權(quán)重，使模型權(quán)重更優(yōu)。

4 結(jié)論

(1)通過借鑒組合預(yù)測模型思想，利用經(jīng)濟學(xué)合作博弈論Shapley值利益分配理論確定了XGBoost估產(chǎn)模型與RF估產(chǎn)模型的權(quán)重，進而構(gòu)建組合預(yù)測模型。結(jié)果表明，基于Shapley值的組合預(yù)測模型估測精度較高，且優(yōu)于單一估產(chǎn)模型的精度。

(2)將組合預(yù)測模型應(yīng)用于研究區(qū)域2010—2018年逐像素的玉米單產(chǎn)估測。結(jié)果表明，從研究區(qū)域玉米單產(chǎn)估測的年際變化看，河北中部平原2010—2018年各縣(區(qū))的玉米估測單產(chǎn)波動變化，但總體呈先減少后增加的趨勢；從玉米單產(chǎn)估測的空間分布看，整體呈西部地區(qū)最高，南部和北部地區(qū)次之，東部地區(qū)最低的特征。