基于雙參數(shù)和Morlet多時間尺度特性的冬小麥單產(chǎn)估測

張 悅 王鵬新 張樹譽 梅樹立 李紅梅 陳 弛

(1.中國農(nóng)業(yè)大學信息與電氣工程學院， 北京 100083； 2.陜西省氣象局， 西安 710014)

0 引言

與傳統(tǒng)的監(jiān)測方式相比，遙感技術具有獲取信息量大、快速、覆蓋面積大的優(yōu)勢，是及時掌握農(nóng)業(yè)資源、作物長勢、農(nóng)業(yè)災害等信息的最佳手段[1]，廣泛應用于作物估產(chǎn)研究。利用遙感進行作物產(chǎn)量估算的方法主要包括：遙感統(tǒng)計估產(chǎn)模型、干物質-產(chǎn)量模型和作物模型模擬等[2]。遙感統(tǒng)計模型主要指利用產(chǎn)量與遙感參數(shù)進行相關統(tǒng)計分析從而實現(xiàn)作物估產(chǎn)，常用的遙感參數(shù)包括歸一化植被指數(shù)(Normalized difference vegetation index，NDVI)[3]、葉面積指數(shù)(Leaf area index，LAI)[4]、植被狀態(tài)指數(shù)(Vegetation condition index，VCI)[5]等。干物質-產(chǎn)量模型通過遙感信息數(shù)據(jù)估測作物地上干物質質量，再依據(jù)作物干物質質量與果實部分間的關系得到作物產(chǎn)量[6]。遙感作物模型模擬是將遙感數(shù)據(jù)作為模型校正的數(shù)據(jù)源之一，對作物模型進行參數(shù)本地化標定后，在氣象、土壤、種植信息等數(shù)據(jù)的驅動下進行作物生長模擬和產(chǎn)量估算[7-10]。相較于遙感統(tǒng)計估產(chǎn)模型，由于上述兩種模型方法所需的部分驅動參數(shù)在大范圍內定量獲取目前存在一定困難，導致上述兩種方法在大區(qū)域應用中仍然受到一定限制。因此，本文選用遙感統(tǒng)計估產(chǎn)模型進行單產(chǎn)估測。

根據(jù)遙感統(tǒng)計估產(chǎn)模型中輸入遙感參數(shù)的個數(shù)可以分為基于單參數(shù)的統(tǒng)計估產(chǎn)模型和基于多參數(shù)(兩個及兩個以上)的統(tǒng)計估產(chǎn)模型。對基于單參數(shù)的遙感統(tǒng)計估產(chǎn)模型，已廣泛應用于作物長勢監(jiān)測、田間環(huán)境監(jiān)測及產(chǎn)量估測[11-16]。對于多參數(shù)遙感統(tǒng)計模型，研究人員利用其對水稻[17]、油菜[18]、玉米[19]等作物產(chǎn)量進行估測和預測。作物生長是一個由多種因素共同作用的復雜過程，相較于單參數(shù)的統(tǒng)計模型，基于多參數(shù)的統(tǒng)計模型更能全面考慮不同因素對作物生長及最終產(chǎn)量的影響，因此，本文選擇與作物生長過程密切相關的VTCI和LAI 作為遙感參數(shù)進行估產(chǎn)模型的構建與分析。

目前，我國對作物估產(chǎn)多立足于單時間尺度進行研究，而多時間尺度有助于分析不同時間尺度對最終結果的影響程度，揭示隱含在序列中的潛在規(guī)律。小波分析作為一種有效的時間序列分析方法，具有時頻多分辨率特征，可以對時間序列的內部結構給以更加詳細的說明[20-24]。因此，將小波分析用于作物估產(chǎn)具有明顯的優(yōu)勢。

本文以陜西省關中平原為研究區(qū)域，以冬小麥主要生育期VTCI、LAI和單產(chǎn)時間序列作為研究對象，開展基于Morlet的連續(xù)小波變換和交叉小波變換，通過計算小波互相關度，進一步分析在不同時間尺度下各生育時期VTCI、LAI和單產(chǎn)間的相關性，確定VTCI、LAI在不同生育時期的權重，進而分別構建基于加權VTCI、加權LAI的單參數(shù)和雙參數(shù)冬小麥估產(chǎn)模型，以期獲得更好的能用于作物產(chǎn)量估測的模型。

1 材料與方法

1.1 研究區(qū)域

關中平原位于陜西省中部的渭河流域，西起寶雞，東至潼關，北到陜北黃土高原，南達秦嶺，位于106°22′～110°24′E，33°57′～35°39′N之間，其行政區(qū)域包括西安市、銅川市、寶雞市、咸陽市、渭南市5個省轄市和楊凌示范區(qū)(圖1)。該地區(qū)地勢西高東低，土壤肥沃，光照條件較好，種植模式主要為冬小麥與夏玉米輪作，是陜西省主要的糧食生產(chǎn)基地和全國重要的商品糧產(chǎn)區(qū)。關中平原為典型大陸性季風氣候區(qū)，年平均溫度6～13℃，年均降水量550～700 mm，受季風氣候影響，降水主要集中在7—9月，且在空間分布不均勻，年際變化與年內變化均較大，因此冬小麥在全生育時期內存在不同程度的水分不足現(xiàn)象。

1.2 數(shù)據(jù)來源與預處理

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

采用的遙感數(shù)據(jù)為MODIS的地表覆蓋類型產(chǎn)品(MCD12Q1)、日地表溫度產(chǎn)品(MYD11A1)、日地表反射率產(chǎn)品(MYD09GA)和葉面積指數(shù)產(chǎn)品(MCD15A3H)；所采用的2011—2017年冬小麥單產(chǎn)數(shù)據(jù)來自陜西省各市統(tǒng)計局發(fā)布的統(tǒng)計年鑒；所使用的降水量數(shù)據(jù)從陜西省氣象局獲得。由于銅川市位于關中平原向陜北黃土高原的過渡地帶，冬小麥面積相對較小，且主要分布在其南部的渭北旱塬，因此，選用關中平原2011—2018年其余4市共24縣的主要生育期旬時間尺度的VTCI、LAI和冬小麥單產(chǎn)數(shù)據(jù)進行相關研究。

1.2.2冬小麥種植區(qū)域提取

利用MODIS Terra+Aqua三級地表覆蓋類型年度全球500 m產(chǎn)品MCD12Q1提取冬小麥種植面積。根據(jù)國際地圈生物圈計劃(IGBP)的分類方案，將MCD12Q1地表覆蓋類型產(chǎn)品與研究區(qū)行政邊界矢量圖進行疊加，得到關中平原縣(區(qū))的小麥種植區(qū)分布圖。

1.2.3生育時期VTCI計算

基于2011—2018年每年3—5月Aqua-MODIS的日地表溫度產(chǎn)品(MYD11A1)和日地表反射率產(chǎn)品(MYD09GA)，得到每年關中平原日LST和日NDVI；應用最大值合成法，以旬為時間尺度，生成每年3—5月的NDVI和LST最大合成產(chǎn)品；利用多年某一旬的NDVI和LST的最大值合成產(chǎn)品，采用最大值合成技術分別生成多年的旬NDVI和LST的最大值合成產(chǎn)品；基于每年3—5月以旬為單位的LST最大值合成產(chǎn)品，逐像素取最小值，生成多年旬LST最大-最小值合成產(chǎn)品，并以此計算旬時間尺度的VTCI[15,25]為

(1)

其中

Lmax(Ni)=a+bNi

(2)

(3)

式中L(Ni)——研究區(qū)域內某一像素的NDVI值為Ni時的地表溫度

Lmax(Ni)、Lmin(Ni)——研究區(qū)域內當NDVI值等于Ni時所有像素地表溫度的最大值和最小值，稱作熱邊界和冷邊界

a、b、a′、b′——待定系數(shù)，由研究區(qū)域NDVI和LST的散點圖近似獲得

結合關中平原冬小麥的生長情況，將冬小麥越冬后的生育時期劃分為返青期(3月上旬—中旬)、拔節(jié)期(3月下旬—4月中旬)、抽穗-灌漿期(4月下旬—5月上旬)和乳熟期(5月中旬—下旬)，并將這4個生育時期稱為冬小麥主要生育期[26]。取某一生育時期包含的多旬VTCI的平均值作為該生育時期的VTCI值；依據(jù)關中平原各縣的冬小麥種植區(qū)分布圖，取各縣域內種植區(qū)所包含像素的VTCI平均值作為該區(qū)域該年該生育時期的VTCI值。

1.2.4生育時期LAI計算

選取研究區(qū)域2011—2018年每年3—5月MODIS MCD15A3H產(chǎn)品進行葉面積指數(shù)提取，該產(chǎn)品是基于Terra和Aqua衛(wèi)星的MODIS傳感器獲得的，與MOD15A2和MYD15A2產(chǎn)品相比，MCD15A3H產(chǎn)品既有較高的時間分辨率(4 d)又有較高的空間分辨率(500 m)，有利于作物長勢監(jiān)測及產(chǎn)量估測。該產(chǎn)品由于云和大氣等因素的影響存在數(shù)據(jù)驟降現(xiàn)象，因此通過上包絡線Savitzky-Golay(S-G)濾波對原始葉面積指數(shù)產(chǎn)品進行平滑處理[27]，經(jīng)上包絡線S-G濾波平滑處理后的葉面積指數(shù)更加符合冬小麥實際生長情況。

為使LAI與VTCI具有相同的取值范圍，將S-G濾波后的LAI進行歸一化處理；將冬小麥各旬所包含的多時相LAI的最大值作為各旬的LAI值，取某一生育時期包含的多旬LAI的平均值作為該生育時期的LAI值；依據(jù)關中平原各縣的冬小麥種植區(qū)分布圖，取各縣域內種植區(qū)所包含像素的LAI平均值作為該區(qū)域該年該生育時期的LAI值。

1.2.5時間序列生成

考慮到時間序列過長的問題，依據(jù)關中平原各縣分布的方位，分為3次，每次均勻不重復，選擇8縣為1組，分別依次構建各縣2011—2017年各生育時期VTCI、LAI和單產(chǎn)時間序列。最終形成均包含4個生育時期的3組VTCI、LAI時間序列及其與每組相對應的單產(chǎn)時間序列，共計27個時間序列，各時間序列長度均為56 a，時間間隔為1 a。

1.3 研究方法

1.3.1Morlet連續(xù)小波變換

1.3.1.1Morlet小波

利用小波對時間序列進行分析時，為使處理后得到平滑連續(xù)的小波振幅，通常使用非正交小波對時間序列進行變換分析。同時在利用小波分析時希望得到更多信息，而復數(shù)小波具有實部和虛部，可以得到時間序列的振幅和相位信息，因此選用復數(shù)且具有非正交性的小波。Morlet小波既具有非正交性，同時是由Gaussian函數(shù)調制的指數(shù)復數(shù)小波，能夠滿足本文對時間序列的分析要求。Morlet小波表示為[28]

(4)

對應的頻率域小波函數(shù)為

(5)

其中

式中t——時間

ω0——無量綱角頻率

ω——角頻率s——伸縮尺度

H(ω)——Heavyside函數(shù)

當ω0=6時滿足容許性條件，在該條件下，Morlet小波的伸縮尺度s與傅里葉分析中周期T的關系為：T=1.03 s，兩者差距較小，因此在本文中假設兩者相等。

1.3.1.2小波變換

利用Morlet小波函數(shù)對各生育時期VTCI、LAI和單產(chǎn)時間序列進行連續(xù)小波變換[28]，計算式為

(6)

其中

s=s02jdj(j=0,1,…,J)

(7)

(8)

式中Wn(s)——小波變換函數(shù)

dt——時間序列的時間間隔，取dt=1

xn′——某一特定時間指數(shù)下的時間序列

n′——某一特定的時間指數(shù)

n——時間指數(shù)

M——時間序列數(shù)據(jù)個數(shù)

s0——最小尺度，取s0=2dt

dj——離散尺度間距，取dj=0.125

J——時間尺度個數(shù)

1.3.1.3小波全譜

通過伸縮小波尺度s并沿著時間指數(shù)n進行局部化，獲得在不同時間尺度下對應的能量密度在時間域中的分布，即通過小波功率譜|Wn(s)|2顯示時間序列在時頻域的振動能量[28]。將其在周期上進行時間平均得到小波全譜(Global wavelet spectrum)

(9)

小波全譜可以用來反映時間序列波動能量隨尺度的分布情況。通過各生育時期VTCI、LAI與單產(chǎn)時間序列的小波全譜圖，可以識別每個時間序列的周期波動特征及波動強度，其振蕩能量峰值對應的周期為該序列的主振蕩周期，表明在該周期下小波系數(shù)的振蕩最強烈，最能代表原序列的周期性變化規(guī)律。

1.3.2交叉小波變換及顯著性檢驗

1.3.2.1交叉小波功率譜

交叉小波變換是將小波變換和交叉譜分析相結合而產(chǎn)生的一種信號分析技術，是表征兩個時間序列在不同時間尺度上關聯(lián)程度與位相關系的重要指標[29-30]，交叉小波功率譜重點突出序列之間在時頻域中高能量區(qū)的相互關系，反映時間序列間的共振周期。各生育時期VTCI時間序列或LAI時間序列與單產(chǎn)時間序列y之間的小波交叉譜為[31]

Wxny(s)=Wxn(s)Wy(s)*

(10)

式中Wxn(s)——VTCI序列或LAI序列的小波變換系數(shù)

Wy(s)*——單產(chǎn)序列的小波變換系數(shù)的復共軛

交叉小波功率譜(|Wxny(s)|)越大，說明兩序列在該時間尺度的相關性越強，用以確定冬小麥不同生育時期VTCI或LAI與單產(chǎn)序列之間的共振周期。

1.3.2.2顯著性檢驗

交叉小波的置信度水平來源于兩個χ2分布的小波功率譜乘積的平方根[28]，概率計算式為

(11)

式中P——概率v——自由度

σxn、σy——生育時期時間序列xn和單產(chǎn)時間序列y的標準差

Zv(p)——P的置信水平

選用的Morlet是復數(shù)小波，v=2，則Z2(95%)=3.999。如果|Wxny(s)|>P，說明對應的周期是顯著共振周期。

1.3.3小波互相關分析方法

各生育時期VTCI或LAI時間序列xn與單產(chǎn)時間序列y之間的小波互相關系數(shù)[32]可表示為

(12)

式中CWRn(s)——小波互相關系數(shù)

σ2(Wxn(s))——尺度s對應的VTCI與單產(chǎn)時間序列小波系數(shù)的方差

σ2(Wy(s))——尺度s對應的LAI與單產(chǎn)時間序列小波系數(shù)的方差

Cov(Wxn,Wy)s——尺度s對應的VTCI或LAI與單產(chǎn)序列小波系數(shù)的協(xié)方差

基于CWRn(s)，定義各生育時期VTCI或LAI與單產(chǎn)的小波互相關度(Wavelet cross-correlation degree，WCCD)[32]為

(13)

式中f(CWRn(s))——尺度s下各生育時期VTCI或LAI與單產(chǎn)的CWRn(s)的權重

1.4 數(shù)據(jù)分析指標

選用決定系數(shù)(R2)和歸一化均方根誤差(Normalized root mean square error，NRMSE)作為精度評價指標對估產(chǎn)模型精度進行評價。

2 結果與分析

2.1 主振蕩周期分析

運用Morlet小波對2011—2017年關中平原3個組別的各生育時期VTCI、LAI和單產(chǎn)時間序列進行多尺度特征分析獲得小波功率譜，并將小波功率譜在各時間尺度上取平均，得到小波全譜。依據(jù)不同生育時期VTCI、LAI及冬小麥單產(chǎn)在不同時間尺度下全譜的變化，以第1組為例進行分析說明(圖2)，峰值對應的時間尺度為該時間序列的主振蕩周期。

2.1.1生育時期VTCI與冬小麥單產(chǎn)的主振蕩周期分析

利用小波變換對第1組各生育時期VTCI和單產(chǎn)進行多尺度分析(圖2)，確定主振蕩周期。在返青期，VTCI全譜在2～3 a、3～4 a和6～8 a處達到峰值，即存在2～3 a、3～4 a和6～8 a的主振蕩周期；在拔節(jié)期，VTCI有2～3 a、3～4 a及6～8 a的主振蕩周期，其中，全譜在3～4 a處達到最高點，全譜振蕩最為強烈；在抽穗-灌漿期，VTCI呈2～3 a、3～4 a和10～12 a的主振蕩周期，且隨著時間尺度的增大，VTCI全譜在14 a后呈上升趨勢；乳熟期VTCI全譜在時間尺度2～3 a、3～4 a和6～8 a下達到峰值，即存在2～3 a、3～4 a和6～8 a的主振蕩周期；單產(chǎn)存在2～3 a、4～5 a和10～12 a的主振蕩周期。通過以上分析可知，在不同的生育時期，VTCI在整體上表現(xiàn)出主振蕩周期的一致性，各生育時期的VTCI和單產(chǎn)均存在2～3 a的主振蕩周期。

2.1.2生育時期LAI與冬小麥單產(chǎn)的主振蕩周期分析

通過對第1組各生育時期LAI和單產(chǎn)全譜隨時間尺度變化(圖2)的特征進行分析可知，返青期，LAI在整個時間尺度上存在2～3 a、6～8 a和10～11 a的主振蕩周期；拔節(jié)期，LAI分別在 2～3 a、3～5 a、6～8 a和11～12 a處達到峰值，其中在10～12 a處振蕩最為強烈；抽穗-灌漿期，LAI存在3～4 a、6～8 a和11～12 a的主振蕩周期；乳熟期LAI全譜存在3～4 a和6～8 a的主振蕩周期，單產(chǎn)存在2～3 a、4～5 a和10～12 a的主振蕩周期。在各生育時期LAI均有3 a和6～8 a的主振蕩周期，且在14 a后隨時間尺度的增大，各生育時期LAI的小波功率不斷增大，說明在不同的生育時期，LAI全譜隨時間尺度的變化具有相似性，即在各生育時期LAI有相似的主振蕩周期。

通過對比不同生育時期的VTCI和LAI全譜與單產(chǎn)時間序列全譜隨時間尺度的變化發(fā)現(xiàn)，相較于VTCI，LAI與單產(chǎn)時間序列在不同時間尺度上的變化更為相似，因此，在不同的生育時期LAI與單產(chǎn)隨時間尺度的變化更具有一致性。

由以上結果可知，各生育時期VTCI和LAI與單產(chǎn)序列間均有2～3 a的主振蕩周期，但小波功率譜的方法主要針對單個時間序列進行分析，對于兩序列間的相關性無法直觀表達，因此，為進一步驗證上述結果，利用交叉功率譜確定各生育時期VTCI、LAI與冬小麥單產(chǎn)間的共振周期，分析VTCI和LAI與單產(chǎn)間的相互關系。

2.2 共振周期分析

2.2.1生育時期VTCI與冬小麥單產(chǎn)的共振周期分析

采用交叉小波變換分析各生育時期VTCI與單產(chǎn)序列在不同時間尺度的關聯(lián)性及位相關系，突出時間序列之間在時頻域中高能量區(qū)的相互關系。通過Morlet小波函數(shù)對每組各生育時期VTCI和單產(chǎn)時間序列進行交叉小波能量譜計算，揭示兩時間序列的位相關系，得到序列之間的共振周期并以第1組為例進行分析說明(圖3)。圖中粗的黑實線圈閉區(qū)域為功率譜值通過置信水平為95%的標準背景譜檢驗；細的黑實線為小波邊界效應影響錐線，其包絡區(qū)域為有效值。圖中箭頭表示位相關系，→表示單產(chǎn)與VTCI同位相，說明兩者為正相關關系；←表示單產(chǎn)與VTCI反位相，說明兩者為負相關關系；↗、↙表示VTCI滯后單產(chǎn)變化且分別為正相關關系、負相關關系；↘、↖表示VTCI超前單產(chǎn)變化且分別為正相關關系、負相關關系。在返青期，VTCI與單產(chǎn)間存在2～3 a和8～10 a的共振周期，其中在8～10 a的相關性更顯著，在此時間尺度下，返青期VTCI與單產(chǎn)間的關系表現(xiàn)為VTCI滯后于單產(chǎn)變化的正相關關系，說明在返青期，水分脅迫對于單產(chǎn)的影響存在滯后效應，對產(chǎn)量的影響可延期至拔節(jié)期；在拔節(jié)期，VTCI與單產(chǎn)間存在3 a左右的共振周期，同時主要表現(xiàn)為VTCI與單產(chǎn)同位相的正相關關系，并通過了95%的顯著性檢驗，說明在拔節(jié)期，VTCI對單產(chǎn)的影響具有實時性且影響較大；在抽穗-灌漿期，VTCI與單產(chǎn)存在3 a和8～10 a的共振周期，主要表現(xiàn)為VTCI滯后于單產(chǎn)變化的正相關關系并通過了95%的顯著性檢驗，說明VTCI變化晚于產(chǎn)量的變化，對產(chǎn)量的影響可延期至乳熟期，且在抽穗-灌漿期VTCI對單產(chǎn)的影響較為顯著；在乳熟期，VTCI與單產(chǎn)存在3～4 a的共振周期，VTCI與單產(chǎn)間的相互關系主要表現(xiàn)為同位相正相關關系，說明VTCI與單產(chǎn)兩者間變化具有同步性。

通過對各生育時期VTCI與單產(chǎn)間的共振周期分析可知，在不同的時間尺度下，各生育時期的VTCI與單產(chǎn)間的共振周期不同，但在拔節(jié)期均表現(xiàn)出VTCI對單產(chǎn)的影響較為顯著的特性。

2.2.2生育時期LAI與冬小麥單產(chǎn)的共振周期分析

對3組各生育時期LAI與冬小麥單產(chǎn)進行交叉小波能量譜計算，并以第1組為例進行分析說明(圖4)。圖中粗的黑實線圈閉區(qū)域為功率譜值通過置信水平為95%的標準背景譜檢驗；細的黑實線為小波邊界效應影響錐線，其包絡區(qū)域為有效值。圖中箭頭表示位相關系，→表示單產(chǎn)與LAI同位相，說明兩者為正相關關系；←表示單產(chǎn)與LAI反位相，說明兩者為負相關關系；↗、↙表示LAI滯后單產(chǎn)變化且分別為正相關關系、負相關關系；↘、↖表示LAI超前單產(chǎn)變化且分別為正相關關系、負相關關系。返青期LAI與單產(chǎn)明顯存在11～12 a的共振周期，且表現(xiàn)為LAI超前于單產(chǎn)變化的正相關關系并通過了95%的顯著性檢驗，說明在返青期，LAI的變化先于單產(chǎn)的變化，表明冬小麥早期的長勢對后期的長勢有顯著的影響；在拔節(jié)期，LAI與單產(chǎn)主要存在10～11 a的共振周期，且LAI先于單產(chǎn)變化，同樣表明拔節(jié)期的長勢對后期的長勢有較大的影響；在抽穗-灌漿期，LAI與單產(chǎn)存在9～10 a共振周期，兩者之間是同位相正相關關系，說明在抽穗-灌漿期，LAI與單產(chǎn)同時變化且為正相關；在乳熟期，LAI與單產(chǎn)存在3 a和10 a的共振周期，均表現(xiàn)為兩者同位相正相關關系。

通過對各生育時期LAI與單產(chǎn)間共振周期的分析可知，在不同的時間尺度下，各生育時期的LAI與單產(chǎn)的相關關系不同，但基本上表現(xiàn)出在返青期和拔節(jié)期，LAI超前于單產(chǎn)變化，在抽穗-灌漿期和乳熟期LAI與單產(chǎn)同時變化，表明冬小麥返青后，返青期和拔節(jié)期的長勢對后期的長勢有較大的影響，抽穗-灌漿期和乳熟期的長勢對單產(chǎn)有較大的影響。

2.3 小波互相關分析

2.3.1生育時期VTCI與單產(chǎn)特征尺度和相關度分析

小波互相關系數(shù)可以反映兩個時間序列在整個時間域中不同尺度的相關程度。利用Morlet小波獲取不同組別各生育時期VTCI和冬小麥單產(chǎn)之間的小波互相關系數(shù)并以第1組為例進行說明(圖5)。通過分析可知，隨著時間尺度的變化，各生育時期的VTCI與單產(chǎn)間的小波互相關系數(shù)發(fā)生改變且存在負相關系數(shù)。依據(jù)VTCI值越小，旱情越嚴重，產(chǎn)量越低的客觀規(guī)律，選擇各生育時期VTCI與單產(chǎn)呈正相關對應的主振蕩周期和共振周期作為分析各生育時期VTCI對單產(chǎn)的相對重要程度的特征時間尺度，通過主振蕩周期確定的特征時間尺度為：返青期2～3 a、3～4 a和7～8 a，拔節(jié)期為2～3 a、6～7 a，抽穗-灌漿期2～3 a和10～11 a，乳熟期為3 a；通過共振周期確定的特征時間尺度為：返青期8～10 a，拔節(jié)期為3 a，抽穗-灌漿期3 a和8～9 a，乳熟期為3～4 a。

分別計算在不同特征時間尺度下的加權期望，得出各生育時期VTCI與單產(chǎn)間的小波互相關度并進行歸一化，得到各生育時期VTCI權重并對不同組別中得到的各生育時期的VTCI權重取平均值作為最終各生育時期VTCI的權重(表1)。可知，經(jīng)過小波變換和交叉小波變換得到的各生育時期VTCI與單產(chǎn)的相關性大小相一致，即拔節(jié)期VTCI與單產(chǎn)間的小波互相關度最大，其次為抽穗-灌漿期，返青期和乳熟期較小。拔節(jié)期是冬小麥根、莖、葉生長的主要階段，對土壤中水分的吸收利用最為迫切，土壤中水分虧缺將顯著影響根莖葉干物質及植株干物質質量的累積速率，從而影響冬小麥最終的長勢及產(chǎn)量。冬小麥在抽穗-灌漿期主要進行生長方式的轉變，由營養(yǎng)生長轉向生殖生長，水分的虧缺會顯著影響其進行光合作用的速率，減少淀粉、蛋白質和有機質的合成，造成冬小麥粒重明顯降低，因此抽穗-灌漿期對最終產(chǎn)量的影響位于第2位。在乳熟期階段，穗粒結構已經(jīng)形成，對一定的水分虧缺表現(xiàn)出較強的忍受力；返青期冬小麥的葉、莖、根等器官增長較為緩慢且干物質量積累不大，其水分虧缺對株高、最終的分蘗、葉面積及干物質累積量的影響

表1 冬小麥各生育時期VTCI與單產(chǎn)小波互相關度及權重Tab.1 Wavelet cross-correlation degrees of yield and VTCIs at four growth stages of winter wheat and weights of VTCIs

2.3.2生育時期LAI與單產(chǎn)的特征尺度和相關度分析

利用Morlet小波計算不同組別各生育時期LAI與單產(chǎn)之間的小波互相關系數(shù)并以第1組為例進行分析說明(圖5)。可以看出，各生育時期的LAI與單產(chǎn)間的小波互相關系數(shù)存在負數(shù)。依據(jù)在一定范圍內，冬小麥葉面積指數(shù)與單產(chǎn)呈正相關的先驗知識，選擇各生育時期LAI與單產(chǎn)呈正相關對應的主振蕩周期和共振周期作為分析各生育時期LAI對單產(chǎn)的相對重要程度的特征時間尺度，通過主振蕩周期確定的特征時間尺度為：返青期2～3 a、6～7 a和10～11 a，拔節(jié)期為2～3 a、4～5 a、6～7 a和11～12 a，抽穗-灌漿期為3～4 a、6～7 a和11～12 a，乳熟期為3 a和6～7 a；通過共振周期確定的特征時間尺度為：返青期為11～12 a，拔節(jié)期為10～11 a，抽穗-灌漿期為9～10 a，乳熟期為3 a和10 a。

基于Morlet確定的時間尺度，分別對不同特征時間尺度的小波互相關系數(shù)求解其加權期望值，得到4個生育時期LAI的小波互相關度，進行歸一化，得到各生育時期LAI權重，并對不同組別中得到的各生育時期的LAI權重取平均值作為最終各生育時期LAI的權重(表2)。可以看出，基于主振蕩周期和共振周期得出的小波互相關度均表現(xiàn)為LAI在抽穗-灌漿期和乳熟期與冬小麥單產(chǎn)的小波互相關度大于其在返青期和拔節(jié)期與單產(chǎn)的小波互相關度，說明在冬小麥的生長過程中LAI在抽穗-灌漿期、乳熟期對于單產(chǎn)的影響較大，而在返青期和拔節(jié)期，LAI對于冬小麥單產(chǎn)的影響較小。在抽穗-灌漿期和乳熟期，冬小麥主要進行生殖生長，主要決定了籽粒質量，相較于其他階段，與單產(chǎn)間的相關性更大。

2.4 單產(chǎn)估測模型的構建與應用

2.4.1單產(chǎn)估測模型的構建與精度分析

根據(jù)主振蕩周期和共振周期確定的冬小麥與各生育時期VTCI的權重，分別建立基于加權VTCI與小麥單產(chǎn)間的一元線性回歸模型(表3)。分析可知，基于主振蕩周期和共振周期建立的單產(chǎn)估測模型決定系數(shù)(R2)分別達到了0.259和0.263，并通過顯著性檢驗，達極顯著水平。基于主振蕩周期和共振周期建立的加權VTCI與單產(chǎn)的估測模型歸一化均方根誤差(NRMSE)分別為16.88%和16.83%，表明基于共振周期確定的產(chǎn)量估測模型精度高于基于主振蕩周期所確定的產(chǎn)量估測模型，說明利用交叉小波變換能夠更好地發(fā)現(xiàn)各生育時期VTCI與單產(chǎn)間的相互關系。

表3 單產(chǎn)估測模型與精度評價結果Tab.3 Yield estimation models and accuracy evaluation results

利用主振蕩周期和共振周期在特征時間尺度下確定的冬小麥各生育時期LAI的權重，計算關中平原的冬小麥加權LAI，分別建立加權LAI與小麥單產(chǎn)間的估產(chǎn)模型(表3)。基于主振蕩周期和共振周期建立的單產(chǎn)估測模型R2分別為0.520和0.522，并通過顯著性檢驗。同時，基于主振蕩周期和共振周期建立的加權LAI與單產(chǎn)的估測模型NRMSE分別為13.58%和13.56%，表明基于共振周期確定的產(chǎn)量估測模型精度與基于主振蕩周期所確定的產(chǎn)量估測模型的精度相近。通過與上述利用加權VTCI構建的估產(chǎn)模型精度對比分析可知，利用加權LAI構建的估產(chǎn)模型精度更高。

根據(jù)Morlet確定的主振蕩周期和共振周期分別建立基于雙參數(shù)的加權VTCI和加權LAI與小麥單產(chǎn)間的線性回歸模型(表3)。結果表明，基于主振蕩周期和共振周期確定的加權VTCI和加權LAI與小麥單產(chǎn)的相關性達到顯著性水平(P<0.05)，建立的單產(chǎn)估測模型R2分別達到0.531和0.533，均高于基于單產(chǎn)數(shù)構建的估產(chǎn)模型精度。基于主振蕩周期和共振周期建立的加權VTCI和LAI與單產(chǎn)的估測模型歸一化均方根誤差分別為13.52%和13.40%，表明基于雙參數(shù)構建的估產(chǎn)模型誤差均小于單參數(shù)構建的估產(chǎn)模型且基于共振周期確定的產(chǎn)量估測模型精度比基于主振蕩周期所確定的產(chǎn)量估測模型精度略高。兩者雖估產(chǎn)結果基本相近，但研究VTCI和LAI與冬小麥單產(chǎn)間相互關系的思路不同，基于主振蕩周期確定的估產(chǎn)模型重點研究單個時間序列的振蕩變化情況，而共振周期重點考慮各生育時期VTCI或LAI與單產(chǎn)間的相關關系，因此，使用共振周期構建的VTCI和LAI雙參數(shù)的估產(chǎn)模型更為合理。

2.4.2冬小麥單產(chǎn)估測

利用2011—2017年數(shù)據(jù)確定共振周期構建的加權VTCI和加權LAI雙參數(shù)的估產(chǎn)模型對關中平原2011—2018年單產(chǎn)進行逐像素估測(圖6)。通過分析可知，關中平原產(chǎn)量分布具有明顯的空間分布特征，即西部地區(qū)產(chǎn)量最高，中部地區(qū)次之，東部地區(qū)產(chǎn)量最小。其中西部和中部地區(qū)是關中平原種植作物的主要區(qū)域，估產(chǎn)空間分布特征與實際情況相符。在西部地區(qū)，2011—2018年冬小麥單產(chǎn)估測范圍為3 800～6 200 kg/hm2，其中，最高產(chǎn)量是2015年鳳翔縣，為6 171.88 kg/hm2，最低產(chǎn)量是2013年眉縣，為3 839.01 kg/hm2，平均4 901.80 kg/hm2；在中部地區(qū)，2011—2018年單產(chǎn)估測范圍為3 500～5 900 kg/hm2，其中單產(chǎn)最低是2013年永壽縣，為3 516.16 kg/hm2，最高產(chǎn)量是2017年閻良區(qū)，為5 821.86 kg/hm2，平均單產(chǎn)4 411.49 kg/hm2；在東部地區(qū)，2011—2018年冬小麥單產(chǎn)平均估測范圍為3 000～4 400 kg/hm2，其中，最高產(chǎn)量是2015年蒲城縣，為4 394.17 kg/hm2，最低產(chǎn)量是2013年澄城縣，為3 047.18 kg/hm2，平均單產(chǎn)3 666.65 kg/hm2。在對2011—2018年估產(chǎn)時發(fā)現(xiàn)，2013年產(chǎn)量最低。2013年冬小麥在生育期內研究區(qū)域平均降水量為243.78 mm，明顯低于2011—2018年生育時期內平均降水量(284.88 mm)，進一步驗證了本文所得出的2013年產(chǎn)量最低的結論。

3 討論

選擇與作物長勢密切相關的條件植被溫度指數(shù)(VTCI)和葉面積指數(shù)(LAI)作為研究指數(shù)，Morlet小波作為函數(shù)，利用小波變換和交叉小波變換分別分析不同時間尺度下冬小麥各生育時期VTCI和LAI與單產(chǎn)時間序列間的主振蕩周期和共振周期，從而分別構建基于加權VTCI、加權LAI的單參數(shù)和雙參數(shù)估產(chǎn)模型。經(jīng)過比較得出，基于雙參數(shù)的估產(chǎn)模型R2均大于0.53，明顯高于王鵬新等[24]在市域尺度利用VTCI構建的估產(chǎn)模型(R2=0.437)和本文利用VTCI、LAI構建的單參數(shù)估產(chǎn)模型的決定系數(shù)，原因是基于雙參數(shù)的估產(chǎn)模型能更充分反映冬小麥在生長過程中水分脅迫和葉面積指數(shù)等因素對于產(chǎn)量的影響，但在作物生長過程中受到多種因素的影響，在今后的研究中還需進一步充分考慮其他因素的影響。

通過將本文獲得的VTCI各生育時期權重與王鵬新等[24]利用小波分析進行市域尺度VTCI多尺度分析得出的權重比較發(fā)現(xiàn)，其所確定的VTCI各生育時期權重間的差異性高于本文基于縣域尺度得出的VTCI權重，原因可能為，在利用小波對縣域尺度VTCI進行多尺度分析時，VTCI在空間維度中復雜性明顯增加，進而增加了小波進行多尺度分析時的誤差，最終對VTCI權重的計算造成了影響。

在利用小波變換對VTCI和LAI進行多尺度分析實現(xiàn)特征提取時，發(fā)現(xiàn)基于LAI構建的產(chǎn)量估測模型精度高于基于VTCI構建的估產(chǎn)模型，原因可能是LAI在整個生育期內有明顯的先增大后減小的變化趨勢和主要特征，而VTCI在整個生育期間變化規(guī)律不明顯且存在大量的細節(jié)特征，因此，在利用Morlet小波對VTCI和LAI進行多尺度分析時，對有明顯主要特征的LAI的表達效果更好，相比之下對VTCI進行特征提取較難。因此，在今后的研究過程中將針對時間序列的變化特性選擇更為適宜的小波函數(shù)進行多尺度分析以提高估產(chǎn)精度。此外，本研究采用MCD12Q1產(chǎn)品中國際地圈生物圈計劃(IGBP)的分類方案對研究區(qū)冬小麥種植區(qū)域進行提取。然而，由于該方案僅對“農(nóng)用地或農(nóng)用地/自然植被”類型進行分類，因此，在對冬小麥分布信息提取時可能產(chǎn)生誤差，進而對冬小麥的估產(chǎn)結果造成一定影響，今后需要進一步對作物分布提取信息支持下的冬小麥單產(chǎn)估測進行精細化研究。

4 結論

(1)通過小波變換和交叉小波變換分別確定冬小麥各生育時期VTCI、LAI與單產(chǎn)的主振蕩周期和共振周期基本具有一致性，不同生育時期VTCI、LAI與單產(chǎn)間存在不同的主振蕩周期和共振周期，但基本上均在2～3 a有明顯的振蕩周期。利用交叉小波分析各生育時期LAI與單產(chǎn)之間的關系發(fā)現(xiàn)，在返青期和拔節(jié)期，LAI先于單產(chǎn)變化，在抽穗-灌漿期和乳熟期，LAI與單產(chǎn)同步變化，因此，LAI在后兩個生育時期對單產(chǎn)的影響更為顯著。

(2)利用小波變換和交叉小波變換獲得的權重，建立加權VTCI、加權LAI單參數(shù)和雙參數(shù)3種單產(chǎn)估測模型，結果表明，基于加權VTCI、加權LAI雙參數(shù)構建的單產(chǎn)估測模型精度均高于單參數(shù)估產(chǎn)模型，基于共振周期確定的雙參數(shù)估產(chǎn)模型精度略高于基于主振蕩周期構建的雙參數(shù)估產(chǎn)模型，說明基于雙參數(shù)的共振周期建立的估產(chǎn)模型能夠更充分反映作物在生長過程中受到的各種因素的影響，能更好地反映作物生長的實際情況。