黃 超，廖玉芳，蔣元華，彭嘉棟

（1.湖南省氣候中心，長(zhǎng)沙 410008；2.湖南省氣象科學(xué)研究所，長(zhǎng)沙 410008；3.氣象防災(zāi)減災(zāi)湖南省重點(diǎn)實(shí)驗(yàn)室，長(zhǎng)沙 410008）

油茶（Camellia oleiferaAbel.）是中國(guó)特有的木本食用油料樹種。在中國(guó)，油茶、大豆、油菜和花生是主要的油料生產(chǎn)源，具有重要的經(jīng)濟(jì)效益。湖南省油茶種植面積位居全國(guó)第一，油茶是湖南省的重要經(jīng)濟(jì)作物之一，同時(shí)也是農(nóng)村脫貧的重要特色產(chǎn)業(yè)。同所有露天生產(chǎn)的農(nóng)作物一樣，氣象因素對(duì)油茶產(chǎn)量影響較大，開展基于氣象因子的油茶產(chǎn)量研究對(duì)提升油茶產(chǎn)業(yè)鏈效益極具現(xiàn)實(shí)意義。

目前基于氣象因子對(duì)作物產(chǎn)量的影響研究，主要采用統(tǒng)計(jì)方法［1］，早期的統(tǒng)計(jì)方法以線性模型為主，將氣象因子與產(chǎn)量的非線性關(guān)系簡(jiǎn)化成線性關(guān)系進(jìn)行產(chǎn)量預(yù)估，多元線性回歸模型和逐步回歸方法［2，3］被廣泛運(yùn)用。而氣象因子與產(chǎn)量之間存在著復(fù)雜的非線性關(guān)系，線性統(tǒng)計(jì)方法在處理簡(jiǎn)單因子時(shí)能夠滿足業(yè)務(wù)需求，但在面對(duì)大量數(shù)據(jù)的復(fù)雜關(guān)系時(shí)誤差相對(duì)較大，因此越來越多的研究開始采用非線性統(tǒng)計(jì)方法來提高預(yù)測(cè)準(zhǔn)確率，如聚類分析、神經(jīng)網(wǎng)絡(luò)、支持向量機(jī)等［4-6］方法都取得了一定的效果。氣象因子對(duì)油茶產(chǎn)量影響的研究目前主要采用以多元回歸方法為代表的線性統(tǒng)計(jì)方法，近年來主成分分析、神經(jīng)網(wǎng)絡(luò)等［7，8］方法也有應(yīng)用。

一般的統(tǒng)計(jì)方法需要研究者憑主觀認(rèn)識(shí)選取氣象因子進(jìn)行建模，即使使用逐步回歸方法也很難挑選到最優(yōu)的氣象因子對(duì)產(chǎn)量進(jìn)行分析［9，10］，而決策樹算法在一定程度上能解決這種問題。決策樹算法是數(shù)據(jù)挖掘中的一種分類算法，屬于非線性統(tǒng)計(jì)方法，能從大量數(shù)據(jù)中識(shí)別有用的規(guī)律，具備自動(dòng)挑選關(guān)鍵因子的優(yōu)勢(shì)，從而客觀地反映氣象因子與產(chǎn)量的相關(guān)關(guān)系。因此，本研究采用分類與回歸樹（Classification and regression tree，CART）和卡方自動(dòng)交叉檢驗(yàn)（Chi-Square automatic interaction detection，CHAID）兩種決策樹算法分別建立預(yù)測(cè)模型開展油茶產(chǎn)量預(yù)測(cè)。

1 數(shù)據(jù)與方法

1.1 資料來源與加工處理

1.1.1 資料來源 油茶鮮果產(chǎn)量數(shù)據(jù)來源于湖南省林業(yè)科學(xué)院的24個(gè)樣地實(shí)地測(cè)產(chǎn)（表1）。氣象資料來自湖南省97個(gè)地面氣象觀測(cè)站1961—2010年的觀測(cè)資料。

1.1.2 資料處理 對(duì)于基礎(chǔ)氣象數(shù)據(jù)，將其分為47項(xiàng)氣象指標(biāo)（表2），并采用1961—2010年的均值和標(biāo)準(zhǔn)差進(jìn)行標(biāo)準(zhǔn)化處理。對(duì)于油茶測(cè)產(chǎn)數(shù)據(jù)，因無明顯趨勢(shì)變化，所以直接進(jìn)行標(biāo)準(zhǔn)化處理。

表1 測(cè)產(chǎn)點(diǎn)數(shù)據(jù)信息

油茶物候期主要包括春梢萌動(dòng)期、春梢生長(zhǎng)期、花芽分化前期、花芽現(xiàn)形期、夏梢生長(zhǎng)期、花芽成熟期、開花期、果實(shí)第一次膨大期、果實(shí)膨大高峰期、油脂轉(zhuǎn)化和積累高峰期、果實(shí)成熟期共11個(gè)時(shí)期，此外根據(jù)油茶生育期跨年的特點(diǎn)，將整個(gè)生育期和結(jié)果年數(shù)據(jù)單獨(dú)進(jìn)行分析，并綜合所有物候期數(shù)據(jù)，將物候期劃分為13個(gè)（表3）。

表2 氣象指標(biāo)名稱

表3 物候期劃分

1.2 決策樹基本原理

決策樹算法是數(shù)據(jù)挖掘算法中的一種白箱分類方法，擅長(zhǎng)處理非線性問題，主要包含C4.5、分類回歸樹（CART）、卡方自動(dòng)交叉檢驗(yàn)（CHAID）等經(jīng)典算法［11，12］。C4.5算 法 只 能 處 理 離 散 數(shù) 據(jù)，CART、CHAID算法可以處理連續(xù)數(shù)據(jù)，因此本研究選取CART和CHAID兩種算法對(duì)油茶產(chǎn)量和氣象數(shù)據(jù)進(jìn)行分析。

決策樹算法以遞歸的形式不斷對(duì)節(jié)點(diǎn)進(jìn)行分割，并且通過預(yù)先定義的分離規(guī)則和分類優(yōu)度來確定分割的閾值，直至達(dá)到終止條件并形成決策樹［13，14］。CART算法是通過計(jì)算分割過程中節(jié)點(diǎn)N包含樣本的預(yù)測(cè)變量y的冗余平方和來實(shí)現(xiàn)的，具體公式如下：

式中，μ=，n為節(jié)點(diǎn)N所包含的樣本數(shù)。

CHAID算法是通過卡方值來確定最佳分割點(diǎn)的［15］，并且屬于多變量分析。CHAID根據(jù)卡方值的大小順序進(jìn)行分類。它以因變量為根結(jié)點(diǎn)，計(jì)算分類的卡方值χ2，對(duì)每個(gè)自變量進(jìn)行分類，具體公式如下：

式中，A i為i水平的觀察頻數(shù)，E i為i水平的期望頻數(shù)，n為總頻數(shù)，p i為i水平的期望頻率，E i=n p i，k為單元格數(shù)。當(dāng)n較大時(shí)，χ2統(tǒng)計(jì)量近似服從k-1個(gè)自由度的卡方分布。

1.3 油茶預(yù)測(cè)模型質(zhì)量評(píng)價(jià)方法

對(duì)于建立的油茶產(chǎn)量預(yù)測(cè)模型，采用平均相對(duì)誤差、趨勢(shì)準(zhǔn)確率以及產(chǎn)量偏多（少）三成和五成準(zhǔn)確率指標(biāo)評(píng)價(jià)其質(zhì)量?jī)?yōu)劣。具體定義如下。

1）平均相對(duì)誤差：

式中，r表示平均相對(duì)誤差，n為樣本數(shù)，x i為模擬產(chǎn)量，x'i為實(shí)際產(chǎn)量。

2）趨勢(shì)準(zhǔn)確率：

式中，n為樣本數(shù)，m i=x i為模擬產(chǎn)量，x'i為實(shí)際產(chǎn)量為平均產(chǎn)量。

3）產(chǎn)量偏多（少）三成、五成準(zhǔn)確率：

式中，i為樣本數(shù)，x i為模擬產(chǎn)量，x'i為實(shí)際產(chǎn)量，為平均產(chǎn)量。當(dāng)計(jì)算產(chǎn)量偏多（少）三成和五成準(zhǔn)確率時(shí)，C分別取0.3和0.5。

2 結(jié)果與分析

首先對(duì)產(chǎn)量數(shù)據(jù)和氣象因子進(jìn)行標(biāo)準(zhǔn)化處理，對(duì)標(biāo)準(zhǔn)化處理后的氣象因子和產(chǎn)量數(shù)據(jù)進(jìn)行相關(guān)分析，篩選出通過0.05置信水平檢驗(yàn)的氣象因子；然后將篩選后的氣象因子數(shù)據(jù)代入模型中進(jìn)行計(jì)算，得到總的預(yù)測(cè)結(jié)果，最后將數(shù)據(jù)劃分為13個(gè)物候期，并分別使用模型進(jìn)行計(jì)算，得到各個(gè)物候期的擬合產(chǎn)量。建模流程如圖1所示。

圖1 建模流程

2.1 區(qū)域產(chǎn)量模型

將2010—2016年的24個(gè)測(cè)產(chǎn)點(diǎn)數(shù)據(jù)合并為一個(gè)數(shù)據(jù)序列，基于所有的氣象指標(biāo)進(jìn)行油茶產(chǎn)量模擬。為了防止過擬合，在參數(shù)設(shè)置時(shí)需要保證每個(gè)葉節(jié)點(diǎn)的樣本總量不小于總樣本數(shù)的5%，同時(shí)對(duì)決策樹采取后剪枝策略來減少樹的分支［16］。

基于24個(gè)測(cè)產(chǎn)點(diǎn)標(biāo)準(zhǔn)化后數(shù)據(jù)建立的最優(yōu)CART決策樹模型見圖2，以節(jié)點(diǎn)1為例來解釋，N表示節(jié)點(diǎn)中的樣本數(shù)，cumt051 and 0_9＞0.050表示判斷條件（開花期5℃以上活動(dòng)積溫的標(biāo)準(zhǔn)化后數(shù)據(jù)大于0.050），如果滿足該條件就進(jìn)入節(jié)點(diǎn)3繼續(xù)判斷，不滿足則進(jìn)入節(jié)點(diǎn)2，以此類推，達(dá)到?jīng)Q策樹終止條件后停止分類，將樣本平均值作為模擬值輸出，形成判斷油茶產(chǎn)量的規(guī)則。同時(shí)，決策樹算法是從眾多氣象因子中選取關(guān)鍵因子組成決策樹，且越排在樹的上層重要性越高。由圖2可見，開花期5℃以上活動(dòng)積溫（cumt051 and 0_9）、春梢萌動(dòng)期25 mm以上降水日數(shù)（rda0251 and 0_3）、果實(shí)成熟期無日照天數(shù)（sunnod2 and 0_13）、開花期平均最高氣溫（tmmean2 and 0_9）、果實(shí)第一次膨大期小于等于0℃的低溫日數(shù)（tnd0002 and 0_10）、春梢萌動(dòng)期累積降水日數(shù)（rdaccu2 and 0_3）、果實(shí)第一次膨大期氣溫日較差（ranget2 and 0_10）是影響油茶產(chǎn)量的關(guān)鍵氣象因子。

圖2 基于區(qū)域產(chǎn)量標(biāo)準(zhǔn)化數(shù)據(jù)建立的CART最優(yōu)決策樹模型

分別使用CART和CHAID算法進(jìn)行建模，最優(yōu)模型的相對(duì)誤差分別為36.00%、38.30%，趨勢(shì)準(zhǔn)確率分別為81.20%、85.10%，結(jié)果表明，直接對(duì)所有站點(diǎn)數(shù)據(jù)建模相對(duì)誤差較大，準(zhǔn)確率不高，這是由于湖南省各地區(qū)地形、氣候條件、油茶品種差異較大，具有區(qū)域種植的特點(diǎn)，單一數(shù)學(xué)模型無法準(zhǔn)確對(duì)各地區(qū)產(chǎn)量與氣象因子關(guān)系進(jìn)行識(shí)別，因此需要分站點(diǎn)、分區(qū)域進(jìn)行建模，從而提高準(zhǔn)確率。

2.2 測(cè)產(chǎn)點(diǎn)產(chǎn)量模型

由于測(cè)產(chǎn)數(shù)據(jù)樣本偏少，為保證樣本數(shù)據(jù)足夠多，采取合并相似站點(diǎn)數(shù)據(jù)的方法來合并區(qū)域產(chǎn)量數(shù)據(jù)集，即對(duì)某一單站與其他23站的油茶產(chǎn)量序列進(jìn)行相關(guān)性分析，取相關(guān)系數(shù)大小前6位的測(cè)站點(diǎn)數(shù)據(jù)合并成新的數(shù)據(jù)集。將2010—2015年的24個(gè)測(cè)站點(diǎn)的油茶產(chǎn)量數(shù)據(jù)采取上述操作進(jìn)行區(qū)域合并，作為模型的訓(xùn)練集，將2016年的油茶產(chǎn)量數(shù)據(jù)作為驗(yàn)證集對(duì)模型進(jìn)行驗(yàn)證。

分別對(duì)24個(gè)測(cè)產(chǎn)點(diǎn)數(shù)據(jù)進(jìn)行建模，根據(jù)相對(duì)誤差從15個(gè)模型中選取最優(yōu)模型后計(jì)算各項(xiàng)評(píng)價(jià)指標(biāo)，得到兩種決策樹方法的最優(yōu)模型質(zhì)量結(jié)果（表4）。由表4可見，CART和CHAID的相對(duì)誤差分別為8.80%、14.30%，趨勢(shì)準(zhǔn)確率分別為97.40%、92.20%，均優(yōu)于基于逐步回歸方法的26.00%的相對(duì)誤差和87.30%的趨勢(shì)準(zhǔn)確率。

表4 兩種決策樹方法的最優(yōu)模型質(zhì)量結(jié)果（單位：%）

使用兩種決策樹算法對(duì)不同物候期數(shù)據(jù)建模，得到不同時(shí)段模型的相對(duì)誤差和準(zhǔn)確率。由圖3和圖4可見，兩種算法在各個(gè)物候期的準(zhǔn)確率分布較為一致，基于11個(gè)物候期產(chǎn)量數(shù)據(jù)建模的模型平均相對(duì)誤差較小，趨勢(shì)準(zhǔn)確率較高，說明在數(shù)據(jù)充足的條件下，決策樹算法能夠更好地識(shí)別氣象指標(biāo)與產(chǎn)量的關(guān)系。開花期、春梢萌動(dòng)期相對(duì)誤差同樣較低，表明該物候期氣象條件對(duì)油茶產(chǎn)量有較大影響，而春梢生長(zhǎng)期、花芽分化前期、花芽現(xiàn)形期、夏梢生長(zhǎng)期相對(duì)誤差較大，趨勢(shì)準(zhǔn)確率較低，說明這4個(gè)物候期的氣象條件對(duì)油茶產(chǎn)量影響較小。

圖3 各物候期模擬產(chǎn)量相對(duì)誤差

圖4 各物候期模擬產(chǎn)量趨勢(shì)準(zhǔn)確率

由圖5可知，24個(gè)測(cè)產(chǎn)點(diǎn)中，兩種算法中最優(yōu)模型的最大平均相對(duì)誤差為17.70%，最小平均相對(duì)誤差為0.50%，其中CART算法有16個(gè)站點(diǎn)相對(duì)誤差小于10.00%，CHAID算法模擬性能相對(duì)較差，僅5個(gè)站點(diǎn)相對(duì)誤差小于10.00%。為了分析兩種決策樹算法相對(duì)誤差在各個(gè)區(qū)域的分布情況，圖5給出了兩種算法的相對(duì)誤差的分布，可以看出兩種算法均在湘中東部地區(qū)有相對(duì)誤差的低值區(qū)，而湘南地區(qū)相對(duì)誤差較大。

選取2016年油茶產(chǎn)量數(shù)據(jù)對(duì)模型進(jìn)行驗(yàn)證，由圖6可知，CART和CHAID兩種算法的最優(yōu)模型的最小相對(duì)誤差分別為0.40%、0.03%。24個(gè)區(qū)域站點(diǎn)中，CART算法有15個(gè)站點(diǎn)相對(duì)誤差在10.00%以內(nèi)，2個(gè)測(cè)產(chǎn)點(diǎn)相對(duì)誤差較大。CHAID算法有11個(gè)站點(diǎn)相對(duì)誤差在10.00%以內(nèi)，3個(gè)測(cè)產(chǎn)點(diǎn)相對(duì)誤差偏大。此外，從各站點(diǎn)的相對(duì)誤差分布來看，CART算法在湘中一帶有較高的準(zhǔn)確率，CHAID在湘南有較高的準(zhǔn)確率。

在建模過程中，決策樹算法選取對(duì)產(chǎn)量產(chǎn)生主要影響的氣象因子，對(duì)模型選取的氣象指標(biāo)頻率進(jìn)行排序，得到影響油茶產(chǎn)量的關(guān)鍵氣象指標(biāo)。綜合兩種方法各物候期的模擬準(zhǔn)確率并結(jié)合蔣元華等［17］的研究，開花期、果實(shí)第一次膨大期、油脂轉(zhuǎn)化和積累高峰期是油茶生長(zhǎng)的3個(gè)關(guān)鍵物候期。關(guān)鍵物候期中決策樹算法挑選頻率排名前五的氣象因子見圖7，由圖7可知，溫度類指標(biāo)在油茶生長(zhǎng)關(guān)鍵期起著最重要的作用。開花期0℃以上積溫和平均最高氣溫的入選頻率最高；在果實(shí)第一次膨大期、油脂轉(zhuǎn)化和積累高峰期，兩種決策樹算法挑選的因子比較接近，氣溫日較差、平均最低氣溫和高溫日數(shù)分別排在各自物候期的前列。

圖5 最優(yōu)模型相對(duì)誤差分布

圖6 2016年驗(yàn)證產(chǎn)量數(shù)據(jù)最小相對(duì)誤差分布

圖7 關(guān)鍵物候期中頻率前五位的氣象因子

3 小結(jié)與討論

CART和CHAID兩種決策樹算法對(duì)歷史產(chǎn)量數(shù)據(jù)模擬的平均相對(duì)誤差分別為8.80%、14.30%，趨勢(shì)準(zhǔn)確率分別為97.40%、92.20%，均優(yōu)于逐步回歸方法。對(duì)2016年油茶產(chǎn)量進(jìn)行驗(yàn)證，24個(gè)區(qū)域站點(diǎn)中CART算法有15個(gè)站點(diǎn)相對(duì)誤差在10.00%以內(nèi)，CHAID算法有11個(gè)站點(diǎn)相對(duì)誤差在10.00%以內(nèi)。湘中東部地區(qū)平原地帶氣象臺(tái)站分布密集，氣象數(shù)據(jù)能真實(shí)反映測(cè)產(chǎn)點(diǎn)情況，整體準(zhǔn)確率較高，高海拔地區(qū)氣象站點(diǎn)往往離測(cè)產(chǎn)點(diǎn)距離較遠(yuǎn)，誤差偏大。

決策樹挑選的重要?dú)庀笠蜃又校瑴囟阮愔笜?biāo)在油茶生長(zhǎng)關(guān)鍵期起著最重要的作用，開花期0℃以上積溫和平均最高氣溫對(duì)產(chǎn)量影響程度最高，主要原因是低溫天氣花粉開裂受到抑制，溫度會(huì)影響昆蟲進(jìn)行授粉，果實(shí)第一次膨大期、油脂轉(zhuǎn)化和積累高峰期的重要?dú)庀笠蜃臃謩e為氣溫日較差、平均最低氣溫和高溫日數(shù)，主要與高溫不利于果實(shí)增長(zhǎng)和油脂積累有關(guān)。