基于支持向量回歸的關中平原冬小麥估產研究

曾 妍，王 迪※，趙小娟

（1. 中國農業科學院農業資源與農業區劃研究所/農業農村部農業遙感重點實驗室，北京100081； 2. 青海省農牧業遙感中心，西寧810008）

0 引言

小麥在中國是僅次于水稻、玉米的第三大糧食作物［1］，小麥估產可為有關部門制定政策和經濟計劃提供依據，在糧食生產、流通和儲備等各環節的宏觀調控中發揮作用。中國是全球最大的小麥生產國和消費國，常年產量約占全球17%［2］，中國小麥產量變化對全球小麥產量的預測有重要影響，進而波及到國際糧價的穩定，因此，及時準確估測中國的小麥產量對保證全球糧食安全具有重要意義。

衛星遙感能及時高效地獲取大范圍地物信息，廣泛應用于農業監測領域，其中產量估算更是農情遙感監測的重要內容。目前，基于遙感信息的作物產量估測方法主要有3種：經驗模型或統計模型、機理模型和半經驗半機理模型［3-4］。機理模型能夠反映作物產量形成過程，但輸入參數多且校正困難，導致其應用受限。相較于機理模型，半經驗半機理模型雖有簡化，但模型輸入中需要的2個關鍵變量，光能利用效率和作物收獲指數，都難以在時空分布中進行準確地定量模擬［5］。統計模型方法盡管對作物產量形成的機理解釋性不強，但其結構簡單、可操作性強，契合作物估產的業務化應用需求，因此仍然是當前作物遙感估產的常用方法。

經驗模型以數學統計分析方法為基礎，包括線性模型和非線性模型。而糧食產量受眾多不確定因素的影響，是一個復雜的非線性系統，所以非線性估產模型的精度通常高于線性模型［6］。支持向量回歸（Support Vector Regression，SVR）方法是近年來廣泛應用的一種非線性回歸模型，其基本思想是，找到一個回歸平面，讓一個集合內所有數據到該平面的距離最近［7］。相較于神經網絡等機器學習方法，它的擬合和泛化推廣能力更優異［8-9］，在解決糧食估產等非線性、小樣本的實際問題上具有較大 優勢。

構建估產模型多采用與糧食產量具有高度相關性的參數，包括植被指數、太陽光合有效輻射、生物量等［10］，其中應用最廣泛的是能夠反映植被綠度狀態的歸一化植被指數（Normalized Difference Vegetation Index，NDVI）［11-13］。然而，除植被綠度狀態以外，影響農作物生長的土壤水分、氣溫等指標也與產量的形成緊密關聯。為此，文章綜合考慮作物生長過程中的水分脅迫指標、作物生長狀態指標，選取條件植被溫度指數（Vegetation Temperature Condition Index，VTCI）、葉面積指數（Leaf Area Index，LAI）作為模型輸入變量，采用支持向量回歸方法建立遙感估產模型，并對該方法的精度進行評價，為及時、準確的區域小麥估產提供新的解決途徑。

1 材料與方法

1.1 研究思路

回歸的基本思路是利用有限的樣本數據，建立輸入和輸出之間隱含的函數關系。該文在MATLAB環境下，通過LIBSVM3.23工具包［14］，采用支持向量回歸方法，基于所選的核函數和最優參數構建產量估測模型，找到各影響因素對產量最佳的映射關系，利用實際產量數據對模型模擬結果進行精度評價，技術路線如圖1所示。

1.2 研究區概況

研究區為陜西省中部關中平原（圖2a）的5個市，西安市、寶雞市、咸陽市、渭南市、銅川市，地處106°18′E～110°35′E，33°35′N～35°52′N（圖2b），總面積5.53萬km2。 關中平原三面環山，渭河橫跨平原，地勢西高東低，中部較為平坦寬闊，平均海拔約400 m。地屬大陸性季風氣候，多年平均氣溫13.9℃，多年平均降水量595.3 mm。該區域氣候溫潤、水系縱橫、土質疏松，是我國重要的糧食產區，也是重要的麥、棉產區，小麥占其耕地面積的50%左右。當地耕作制度為兩年三熟，作物種植模式主要是在灌溉區內進行冬小麥、夏玉米輪作，在旱作區域進行夏季休閑式的冬小麥 連作［15］。

圖2 研究區地理位置Fig.2 Location of study area

1.3 數據源

該文采用的遙感數據是研究區內2011—2016年冬小麥主要生育期（3月上旬至5月下旬）的Aqua-MODIS日地表反射率產品（MYD09GA）、日地表溫度產品（MYD11A1）、葉面積指數產品（MCD15A3）。所用的產量數據來源于2011—2016年的《陜西省統計年鑒》。將各市的小麥種植面積以及總產量匯總，表1中列出了計算后得到的對應年份各市冬小麥單產數據。

表1 2011—2016年研究區5市冬小麥單產Table 1 Yield per hectare of winter wheat of research area in five cities from 2011 to 2016 kg/hm2

1.4 研究方法

1.4.1 VTCI計算

條件植被溫度指數（VTCI）是基于歸一化植被指數（NDVI）和地表溫度（LST）特征空間呈三角形區域分布的特點提出的，能體現植物生長過程中受到的水分脅迫，是一種近實時的干旱監測方法，計算方法為［16］：

式（1）～（3）中，LSTmax（NDVIi）、LSTmin（NDVIi）分別被稱作熱、冷邊界，表示在研究區內，當NDVIi值為某一特定值時，所有像素地表溫度的最大值和最小值；LST（NDVIi）表示某一像素的NDVI值為NDVIi時的地表溫度；a，b，a′，b′是待定系數，由研究區內的NDVI和LST散點圖近似得到。

計算VTCI采用的遙感數據是研究區內2011—2016年冬小麥主要生育期（3月上旬至5月下旬）的Aqua-MODIS日地表反射率產品（MYD09GA）、日地表溫度產品（MYD11A1）。利用MODIS重投影工具MRT（MODIS Reprojection Tools）對原始數據進行拼接、重采樣、投影轉換等預處理操作后，得到研究區日LST和日地表反射率數據。利用日地表反射率數據計算日NDVI，用最大值合成方法分別生成逐像素的旬LST和旬NDVI最大值合成數據，并按照VTCI定義計算出旬VTCI。

將冬小麥的主要生育期劃分為4個生育時期，結合關中平原冬小麥越冬后的生長情況，將返青期定為3月上旬至3月中旬，拔節期為3月下旬至4月中旬，抽穗—灌漿期為4月下旬至5月上旬，乳熟期為5月中旬至5月下旬［17］。計算每個生育期內各旬VTCI的平均值作為該生育期的VTCI值，計算研究區2011—2016年5市各生育期的VTCI。

1.4.2 LAI計算

采用2011—2016年冬小麥主要生育期內的葉面積指數產品（MCD15A3）計算LAI，這一產品每4 d合成1次，時間分辨率較高，空間分辨率為500 m，有利于農作物長勢監測。同樣對影像進行投影轉換、裁剪等預處理，利用S-G濾波器消除云和大氣等噪聲，改善數據質量。然后采用取最大值的方法逐像素地生成旬LAI，取各生育期內所包含的多旬LAI的均值分別作為該生育期LAI，實現LAI、VTCI時間分辨率的 統一。

1.4.3 SVR模型構建

①支持向量回歸

支持向量回歸的基本思想是采用一個非線性映射把數據映射到高維特征空間，然后在高維特征空間構造回歸估計函數，再映回到原空間。通過定義適當的核函數k（x，xi）來實現非線性變換，無需得到確切的映射函數Φ（x），給計算帶來了極大方便［18］。回歸方程為：

式（4）中，ai-a*

i是各支持向量的系數。

②格式轉換及歸一化

通過運行FormatDataLibsvm.xls，將數據格式轉化為LIBSVM規定的輸入格式。為避免數據之間的量級差別，對VTCI、LAI、單產數據分別采用歸一化處理，歸一化公式 為：

式（5）中，Xmax和Xmin分別為數據中的最大值和最小值，經過歸一化消除了模型輸入由于量綱和單位不同造成的影響，使樣本數據更加適應回歸建模與分析。

③核函數選擇

支持向量回歸包括線性和非線性，在線性回歸的基礎上引入核函數，得到非線性回歸。支持向量回歸將高維空間內積運算化簡成了從低維空間進行輸入的核函數計 算［19-20］，選擇的核函數的類型決定了支持向量回歸的很多特性，從某種程度上直接影響著回歸模型的擬合能力和預測精度。因此，核函數的選擇是用支持向量回歸預測糧食產量的一個關鍵步驟。

常用的核函數有多項式核函數、徑向基（RBF）核函數（即高斯核函數）、Sigmoid核函數等，選用不同的核函數可以構造出不同的支持向量回歸模型。

許多實驗研究表明［21-22］，當缺少先驗知識時，選用徑向基核函數訓練建模的效果較好，所得模型的總體性能較高，而且徑向基核函數只含一個未知參數σ，易于進行優化，所以該研究選擇的核函數為徑向基核函數。公式為：

④參數尋優

在用支持向量回歸解決實際問題時，參數的選擇對回歸模型的性能影響很大。應用RBF核函數時，懲罰因子C和參數σ與回歸模型的學習性能息息相關。懲罰因子C用于控制模型對誤差范圍以外的樣本的懲罰程度，也可以說是對離群點的重視程度，懲罰因子過高會造成過擬合。徑向基核函數的參數σ反映了訓練樣本數據的范圍或分布，確定局部鄰域的寬度。

該研究選用LIBSVM中的ε-SVR方法，用網格搜索和n-fold交叉驗證的方法進行參數尋優，確定C和σ值。當有不同的C和σ都對應最高的精度時，把參數C最小的那組C和σ作為最佳參數，避免懲罰參數C太高造成過學習狀態。

1.4.4 模型精度評價

該文選取決定系數（Coefficient of Determination，R2）、均方根誤差（Root Mean Square Error，RMSE）、平均絕對百分比誤差（Mean Absolute Percentage Error，MAPE）作為模型擬合程度的評價指標，具體公式為：

式（7）～（9）中，i表示第i個樣本點數據，Yi為第i個實際產量，單位為kg/hm2；Ei為根據模型擬合的冬小麥產量估算值，單位為kg/hm2；為實際產量的均值，單位為kg/hm2；為模型估測產量的均值，單位為kg/hm2。

2 結果與分析

2.1 回歸結果及誤差分析

研究以2011—2015年研究區的數據作為訓練集，采用SVR算法構建冬小麥估產模型，圖3為2011—2015年研究區各市的實際產量和回歸產量，可以看出預測產量和對應的實際產量波動的趨勢大致相符。以實際產量和估測產量分別為橫縱坐標繪制散點圖（圖4），添加趨勢線和45°參考線，趨勢線對比參考線偏離較小，可見模型對測試集的擬合效果較好。計算模型對訓練集樣本的預測精度，回歸模型具有較高的決定系數R2=0.88，估測產量與實際產量間的RMSE為336.39 kg/hm2，MAPE為6.12%。

對訓練集的擬合學習說明了SVR構建的模型能夠以較小的誤差擬合歷史數據。采用構建好的模型預測2016年5市冬小麥單產，表2為預測結果及誤差，可見5市數據共同構建的模型進行單產估測時有一定誤差，其中面積較小的銅川市的相對誤差較小，僅有1.28%。分別對2011—2015年5市冬小麥產量數據進行建模預測（表3），預測結果的RMSE平均值為232.63 kg/hm2，其中寶雞市、渭南市的MAPE低至2.98%、3.37%，5市MAPE的平均值為5.74%，可見模型精度較高。

圖3 2011—2015年5市冬小麥實際產量與回歸產量對比Fig.3 Comparison of actual yield and regression yield from 2011 to 2015

圖4 訓練集擬合情況Fig.4 Fitting effect of training set

表2 2016年5市冬小麥回歸預測結果及誤差Table 2 Regression prediction results and errors

表3 5市分別構建回歸模型的預測結果及精度Table 3 Prediction results and accuracy of regression models constructed by different municipalities

2.2 參數敏感性分析

通過觀察不同參數對預測精度的影響，判斷對預測精度提高幫助最大的參數。圖5、圖6給出了C和σ值變化時3個誤差評價指標的相應變化。結果顯示，發生同等程度的變化時，核參數σ引起的MAPE、RMSE、R2的改變都比懲罰因子C引起的改變更大。因此，核參數σ比懲罰因子C更重要。對2個參數進行耦合分析，發現該研究中當懲罰因子C為2、核參數σ為0.354時模型精度最高。

圖5 不同參數C條件下的回歸模型精度變化Fig.5 The accuracy change of regression model under different parameter C

3 結論與討論

3.1 結論

該研究計算了2種與植被生長狀況息息相關的參數——VTCI、LAI，將冬小麥各個生育期的VTCI和LAI作為輸入變量，劃分訓練集和測試集，用支持向量回歸算法對這些數據進行學習和訓練，從而找到其內部隱含的函數關系。尋找最優的模型參數，建立冬小麥各生育期的VTCI、LAI與冬小麥產量之間的支持向量回歸估產模型。

圖6 不同參數σ條件下的回歸模型精度變化Fig.6 The accuracy change of regression model under different parameter σ

支持向量回歸模型較好地反映了輸入的影響因素VTCI、LAI同冬小麥單產之間的復雜的非線性映射關系，建立的回歸模型決定系數達到0.88，模型擬合較為理想，此外，模型參數的選取對模型精度有明顯的影響，徑向基核參數σ比懲罰因子C對預測結果影響更大。

該文采用支持向量機方法對關中平原的冬小麥產量進行預測，在樣本較少的情況下，建立的回歸模型精度較高，能夠推廣到其他地區，也可以為其他作物的產量預測研究提供參考，該方法在糧食產量預測領域有良好的應用前景。

3.2 討論

該文僅選用了VTCI和LAI 2個指標，只能在一定程度上描述作物本身的生長狀態。實際上，影響糧食產量的外部因素眾多，還包括很多非遙感信息，如自然災害、勞動主體的積極性、農業科技水平、農業的基礎設施等。該研究缺少足夠的先驗知識，模型設置過于理想化，并未充分考慮多層次的因素，若要得到更為科學合理的估產結果，仍需繼續探索與糧食產量密切相關的指標。

該研究核函數的選定依照前人的經驗，而無充分的理論支撐，因此核函數的選用是未來需要進一步研究的內容。隨著新技術的不斷涌現，支持向量回歸和其他理論方法的結合值得深入發掘，例如參數尋優可采用粒子群優化算法［23］、遺傳算法，核函數可選用混合核函數［24］等。

此外，該研究采用統計年鑒的產量數據進行驗證，但產量數據不可避免地存在著差錯或漏報［25］, 因此如何有效驗證模型精度也是有待解決的問題之一。