魏慶偉，朱黎明，王福州

（1. 中國(guó)氣象局·河南省農(nóng)業(yè)氣象保障與應(yīng)用技術(shù)重點(diǎn)開(kāi)放實(shí)驗(yàn)室，河南鄭州 450003；2. 鶴壁市氣象局，河南鶴壁 458030；3. 揚(yáng)州大學(xué)水利科學(xué)與工程學(xué)院，江蘇揚(yáng)州 225009）

冬小麥?zhǔn)俏覈?guó)主要糧食作物之一，對(duì)穩(wěn)定我國(guó)糧食結(jié)構(gòu)以及維護(hù)糧食安全具有舉足輕重的地位[1?2]。作物發(fā)育時(shí)期反映了其生長(zhǎng)和發(fā)育規(guī)律[3?4]，及時(shí)準(zhǔn)確地獲取作物發(fā)育時(shí)期信息，可用于輔助監(jiān)測(cè)其生長(zhǎng)狀況，為科學(xué)管理作物提供有力依據(jù)[5?7]。對(duì)冬小麥發(fā)育時(shí)期的監(jiān)測(cè)是掌握冬小麥生長(zhǎng)發(fā)育狀態(tài)的必要手段，有利于對(duì)冬小麥的時(shí)空年際變化進(jìn)行合理分析，進(jìn)而對(duì)其估產(chǎn)模型進(jìn)行改進(jìn)[8?9]。因此，研究冬小麥關(guān)鍵發(fā)育時(shí)期的提取方法，可實(shí)現(xiàn)對(duì)冬小麥關(guān)鍵生長(zhǎng)階段的長(zhǎng)勢(shì)監(jiān)測(cè)，對(duì)估算冬小麥產(chǎn)量具有重要意義。

按照實(shí)現(xiàn)手段來(lái)劃分，當(dāng)前作物發(fā)育時(shí)期觀測(cè)的方法可分為人工觀測(cè)法和遙感監(jiān)測(cè)法。其中，傳統(tǒng)的人工觀測(cè)法通常是采用田間人工目視觀察，記錄作物生長(zhǎng)發(fā)育的季節(jié)變化和年際變化。人工觀測(cè)法簡(jiǎn)單易行、觀測(cè)準(zhǔn)確度較高，且觀測(cè)數(shù)據(jù)具有較好的連續(xù)性和穩(wěn)定性，但費(fèi)時(shí)費(fèi)力，受主觀影響較大[10?11]。遙感監(jiān)測(cè)法基于綠色植物光譜理論，認(rèn)為同一種作物在不同發(fā)育時(shí)期光譜反射特性不同，進(jìn)而反映作物生長(zhǎng)信息。遙感監(jiān)測(cè)法具有時(shí)效高、范圍寬、成本低等優(yōu)點(diǎn)，并能反映地面植被生長(zhǎng)發(fā)育的季節(jié)變化及年際變化等特點(diǎn)[12]。然而，衛(wèi)星遙感的光譜植被指數(shù)容易受云雨影響，導(dǎo)致監(jiān)測(cè)數(shù)據(jù)波動(dòng)。此外，受衛(wèi)星遙感時(shí)空分辨率影響，植被指數(shù)在時(shí)間上和空間上存在不連續(xù)現(xiàn)象[13]。近地遙感具有時(shí)間分辨率高、受大氣條件及云雨等天氣影響小、可野外長(zhǎng)時(shí)間固定觀測(cè)等優(yōu)勢(shì)，而且觀測(cè)數(shù)據(jù)的空間尺度適中，在很大程度上克服了衛(wèi)星遙感監(jiān)測(cè)面臨的上述問(wèn)題，是實(shí)現(xiàn)作物發(fā)育時(shí)期監(jiān)測(cè)的有力工具[14?16]。目前，已有的近地遙感研究多采用圖像傳感器的方式進(jìn)行植被物候期監(jiān)測(cè)[17?19]，鮮有采用光譜傳感器監(jiān)測(cè)農(nóng)作物發(fā)育時(shí)期的研究。

本研究以歸一化植被指數(shù)測(cè)量?jī)x（SRS-NDVI）為監(jiān)測(cè)手段，提出基于近地遙感的冬小麥關(guān)鍵發(fā)育時(shí)期監(jiān)測(cè)方法。首先，獲取SRS-NDVI 測(cè)量?jī)x監(jiān)測(cè)的NDVI 時(shí)序數(shù)據(jù)，并對(duì)數(shù)據(jù)進(jìn)行預(yù)處理。然后，運(yùn)用相關(guān)方法提取冬小麥關(guān)鍵發(fā)育時(shí)期，旨在提高冬小麥發(fā)育時(shí)期監(jiān)測(cè)的準(zhǔn)確性和全面性，為冬小麥發(fā)育時(shí)期的自動(dòng)化觀測(cè)提供新思路，為精準(zhǔn)農(nóng)業(yè)生產(chǎn)管理提供決策服務(wù)。

1 材料和方法

1.1 研究區(qū)域

研究區(qū)鶴壁市位于河南省北部，太行山東麓向華北平原過(guò)渡地帶，是黃淮海糧食生產(chǎn)核心區(qū)。研究區(qū)屬于暖溫帶半濕潤(rùn)型季風(fēng)氣候，年平均氣溫14.2～15.5 ℃，年平均降水量349.2～970.1 mm，年日照時(shí)數(shù)1 787.2～2 566.7 h，年平均無(wú)霜期220 d。所采用的SRS-NDVI測(cè)量?jī)x安裝于鶴壁市農(nóng)業(yè)氣象試驗(yàn)站試驗(yàn)田（圖1），站點(diǎn)所處農(nóng)田的農(nóng)業(yè)生產(chǎn)方式與該區(qū)農(nóng)業(yè)生產(chǎn)方式保持一致。

圖1 SRS-NDVI測(cè)量?jī)x布設(shè)環(huán)境Fig.1 Layout environment of SRS-NDVI measurement

1.2 數(shù)據(jù)獲取

研究區(qū)冬小麥常年10 月中旬播種，翌年6 月上旬成熟。為獲取冬小麥生長(zhǎng)季NDVI 曲線，選取2017 年10 月11 日—2018 年6 月10 日和2018 年10月11 日—2019 年6 月10 日的SRS-NDVI 測(cè)量?jī)x監(jiān)測(cè)數(shù)據(jù)。冬小麥發(fā)育時(shí)期人工觀測(cè)數(shù)據(jù)為鶴壁市農(nóng)業(yè)氣象試驗(yàn)站2018、2019年農(nóng)業(yè)氣象業(yè)務(wù)所觀測(cè)的冬小麥各發(fā)育日期。人工觀測(cè)發(fā)育期包括播種、出苗、三葉、分蘗、越冬始期、返青、起身、拔節(jié)、孕穗、抽穗、開(kāi)花、乳熟和成熟期，用于近地遙感NDVI監(jiān)測(cè)發(fā)育時(shí)期結(jié)果精度檢驗(yàn)。

1.3 觀測(cè)數(shù)據(jù)預(yù)處理

SRS-NDVI 測(cè)量?jī)x安裝在試驗(yàn)田植被冠層上方，包含向上和向下的光強(qiáng)檢測(cè)傳感器，該傳感器可同時(shí)檢測(cè)天空入射和冠層反射的近紅外（810 nm）和紅光（650 nm）的光照強(qiáng)度，通過(guò)入射和反射光強(qiáng)計(jì)算得到植被冠層近紅外與紅光的反射率，運(yùn)用NDVI 公式計(jì)算得到農(nóng)田尺度的作物冠層NDVI。研究所用SRS-NDVI 觀測(cè)儀時(shí)間分辨率為30 min，根據(jù)楊鈞森等[20]的研究結(jié)論，SRS-NDVI 觀測(cè)儀適宜測(cè)量條件的光照強(qiáng)度大于5.2 klx，在該條件下觀測(cè)數(shù)據(jù)具有較高的可靠性，且測(cè)量精度隨光照強(qiáng)度增加而提高。為確保觀測(cè)數(shù)據(jù)的可靠性，選取每天10:00—14:30 時(shí)段內(nèi)10 組觀測(cè)數(shù)據(jù)的平均值作為當(dāng)日冬小麥植被指數(shù)。

冬小麥生長(zhǎng)季內(nèi)農(nóng)田的NDVI 時(shí)序曲線應(yīng)平穩(wěn)變化，受冬季或早春時(shí)積雪、云與其他不利大氣條件的影響，NDVI 觀測(cè)值降低，使得時(shí)間序列NDVI曲線在某一時(shí)間段內(nèi)呈現(xiàn)不規(guī)則波動(dòng)，從而影響了植被參數(shù)反演和信息提取的精度。NDVI 時(shí)序曲線中的突變點(diǎn)與植被自身的生長(zhǎng)規(guī)律不相符，屬噪聲點(diǎn)，需要降低噪音，重構(gòu)NDVI時(shí)序。

目前噪音去除方法主要有鄰域差值分析法、濾波法、最大值合成法和曲線擬合法等算法模型[21?22]。結(jié)合冬小麥生長(zhǎng)過(guò)程N(yùn)DVI 變化規(guī)律，本研究采用鄰域差值分析法重構(gòu)NDVI 時(shí)序，并運(yùn)用S-G 濾波法進(jìn)行平滑處理。

將i與i+1 連續(xù)2 d 觀測(cè)的植被指數(shù)值按照冬小麥生長(zhǎng)的時(shí)間順序（從播種期到收獲期）作差值處理，兩者差值≥0.1 時(shí)，取較大值作為第(i+1)天的NDVI 值。重復(fù)上述步驟，直至所有相鄰差值均<0.1，處理停止。表達(dá)式為

式（1）中，Ni為第i天NDVI 觀測(cè)值，Ni+1為第（i+1）天NDVI觀測(cè)值。

S-G 濾波法由Savitzky 和Golay 于1964 年提出，之后被廣泛地運(yùn)用于數(shù)據(jù)流平滑除噪，是一種在時(shí)域內(nèi)基于局域多項(xiàng)式最小二乘法擬合的濾波方法[23?24]。這種濾波法最大的特點(diǎn)是在濾除噪聲的同時(shí)可以確保數(shù)據(jù)的形狀、寬度不變。S-G 濾波器由輸入的濾波窗口寬度和多項(xiàng)式次數(shù)生成，通過(guò)對(duì)原始序列進(jìn)行卷積計(jì)算來(lái)平滑數(shù)據(jù)，其表達(dá)式為

式（2）中，Yj*是時(shí)間j處濾波后的NDVI值，Yj+i是NDVI 原始序列中第j處第i個(gè)值，Ci是第i個(gè)NDVI值濾波時(shí)的系數(shù)，N是卷積的數(shù)目，也等于濾波窗口的寬度（2m+1）。

1.4 NDVI時(shí)間序列函數(shù)擬合

多項(xiàng)式擬合是用一個(gè)多項(xiàng)式展開(kāi)去擬合包含數(shù)個(gè)分析格點(diǎn)的一小塊分析區(qū)域中的所有觀測(cè)點(diǎn)，得到觀測(cè)數(shù)據(jù)的客觀分析場(chǎng)，展開(kāi)系數(shù)用最小二乘法擬合確定。結(jié)合冬小麥生長(zhǎng)季NDVI 時(shí)序數(shù)據(jù)特點(diǎn)，采用四次多項(xiàng)式函數(shù)擬合冬小麥的生長(zhǎng)曲線。擬合多項(xiàng)式為

式（3）中，t為冬小麥生育累積天數(shù)（d），y（t）為冬小麥生育累積天數(shù)為t時(shí)的NDVI 值，a、b、c、d、e為擬合多項(xiàng)式的系數(shù)。

1.5 特征點(diǎn)提取方法

提取作物發(fā)育時(shí)期常采用的方法主要有動(dòng)態(tài)閾值法、曲線速率法和極值法等。在農(nóng)作物生長(zhǎng)周期中，當(dāng)從一個(gè)生長(zhǎng)階段向另一個(gè)生長(zhǎng)階段轉(zhuǎn)換時(shí)，函數(shù)曲線從一個(gè)近線性趨勢(shì)轉(zhuǎn)向另一個(gè)近線性趨勢(shì)，通過(guò)計(jì)算植被指數(shù)曲線速率極值可以確定冬小麥的發(fā)育時(shí)期，曲線速率公式為：

式（4）中，dx表示自變量x的單位變化量，dy表示函數(shù)y單位變化量，y′表示自變量x增加單位dx時(shí)，函數(shù)y的變化速率。

動(dòng)態(tài)閾值法使用植被指數(shù)定義了作物營(yíng)養(yǎng)生長(zhǎng)和生殖生長(zhǎng)的開(kāi)始期，即在作物生長(zhǎng)曲線上升階段，將距離植被指數(shù)最小值為曲線增幅10%的位置定義成作物營(yíng)養(yǎng)生長(zhǎng)的開(kāi)始期；在作物生長(zhǎng)曲線下降階段，將距離植被指數(shù)最小值為曲線增幅90%的位置定義成生殖生長(zhǎng)的開(kāi)始期[25]。

不同區(qū)域冬小麥生長(zhǎng)規(guī)律不同，動(dòng)態(tài)閾值法不能有效全面提取發(fā)育時(shí)期特征點(diǎn)。因此，本研究提出廣義動(dòng)態(tài)閾值法，即植被指數(shù)最小值與最大值差值的百分比與冬小麥發(fā)育時(shí)期進(jìn)行匹配，確定為某一發(fā)育時(shí)期的特征點(diǎn)閾值，并以此為基礎(chǔ)進(jìn)行發(fā)育時(shí)期監(jiān)測(cè)及對(duì)比分析。

2 結(jié)果與分析

2.1 NDVI時(shí)間序列曲線重構(gòu)

經(jīng)鄰域差值分析法重構(gòu)NDVI 時(shí)序數(shù)據(jù)、運(yùn)用S-G 濾波法對(duì)重構(gòu)NDVI 進(jìn)行平滑處理后，2017—2018 年度生長(zhǎng)季冬小麥NDVI 時(shí)序數(shù)據(jù)如圖2 所示，2018—2019 年度生長(zhǎng)季冬小麥NDVI 時(shí)序數(shù)據(jù)如圖3 所示。為了更直觀體現(xiàn)NDVI 時(shí)間序列曲線預(yù)處理效果，選取2019 年冬小麥生育累積天數(shù)為120～140 d的數(shù)據(jù)處理結(jié)果進(jìn)行展示，如圖4所示。

圖2 2017—2018年度冬小麥生長(zhǎng)季NDVI時(shí)間序列Fig.2 NDVI time series of winter wheat in the whole growth stage in 2017—2018

圖3 2018—2019年度冬小麥生長(zhǎng)季NDVI 時(shí)間序列Fig.3 NDVI time series of winter wheat in the whole growth stage in 2018—2019

圖4 2019年冬小麥NDVI時(shí)間序列（生育累積天數(shù)120～140 d）Fig.4 NDVI time series of winter wheat in 2019（120—140 d of growing period）

2.2 冬小麥關(guān)鍵發(fā)育時(shí)期特征點(diǎn)提取

根據(jù)常年衛(wèi)星遙感監(jiān)測(cè)的冬小麥生長(zhǎng)季NDVI變化規(guī)律：出苗后NDVI 開(kāi)始逐漸增長(zhǎng)，三葉期至越冬期前為快速增長(zhǎng)階段，越冬期NDVI 基本保持不變，返青期后NDVI 再次快速上升，抽穗期達(dá)到最大，之后開(kāi)始逐步下降，乳熟期后NDVI 快速下降，成熟收獲后達(dá)到最小。為準(zhǔn)確提取冬小麥關(guān)鍵發(fā)育時(shí)期特征點(diǎn)，以冬小麥生育累積天數(shù)第83 天為界，將預(yù)處理后冬小麥生長(zhǎng)季NDVI 時(shí)序數(shù)據(jù)分為兩部分，分別進(jìn)行關(guān)鍵發(fā)育時(shí)期特征點(diǎn)提取。

根據(jù)觀測(cè)數(shù)據(jù)所呈現(xiàn)的特點(diǎn)，使用廣義動(dòng)態(tài)閾值法提取冬小麥關(guān)鍵發(fā)育時(shí)期的特征點(diǎn)。首先，將NDVI 曲線上升階段的曲線增幅10%的位置歸類為出苗期和返青期，將越冬前曲線增幅65%的位置歸類為分蘗期，將越冬后曲線增幅90%的位置歸類為拔節(jié)期，再將NDVI 曲線下降階段的曲線增幅90%的位置歸類為乳熟期。然后，使用變化速率法，在NDVI 擬合曲線增長(zhǎng)速率最大位置處，提取冬小麥三葉期、起身期；使用極值法，提取越冬始期、抽穗期和成熟期。最終，將1 月1 日前NDVI 擬合函數(shù)曲線最大值歸類為越冬始期，1 月1 日后NDVI 時(shí)序曲線最大值歸類為抽穗期，1 月1 日后冬小麥抽穗前的最小值歸類為成熟期。2018、2019 年冬小麥各發(fā)育時(shí)期特征點(diǎn)提取結(jié)果如圖5—8所示。

圖5 2017年10月11日—12月31日冬小麥NDVI時(shí)序特征曲線Fig.5 NDVI time series characteristic curve of winter wheat from October 11 to December 31，2017

圖6 2018年1月1日—6月10日冬小麥NDVI 時(shí)序特征曲線Fig.6 NDVI time series characteristic curve of winter wheat from January 1 to June 10，2018

圖7 2018年10月11日—12月31日冬小麥NDVI 時(shí)序特征曲線Fig.7 NDVI time series characteristic curve of winter wheat from October 11 to December 31，2018

圖8 2019年1月1日—6月10日冬小麥NDVI 時(shí)序特征曲線Fig.8 NDVI time series characteristic curve of winter wheat from January 1 to June 10，2019

2.3 冬小麥發(fā)育時(shí)期提取結(jié)果與誤差分析

根據(jù)NDVI 曲線特征點(diǎn)提取冬小麥發(fā)育時(shí)期監(jiān)測(cè)日期，并與人工觀測(cè)日期進(jìn)行誤差對(duì)比分析，如表1 所示。由表1 可以看出，2017—2018 年度冬小麥分蘗期（2017 年11 月8 日）NDVI 監(jiān)測(cè)與人工觀測(cè)日期最吻合，絕對(duì)誤差為0；越冬始期（2017 年12 月23 日）的誤差最大，較人工觀測(cè)日期晚8 d；其他發(fā)育時(shí)期絕對(duì)誤差均≤4 d，平均誤差2.7 d。2018—2019 年度冬小麥三葉期（2018 年11 月10 日）、乳熟期（2019 年5 月22 日）NDVI 監(jiān)測(cè)與人工觀測(cè)日期最吻合，絕對(duì)誤差為0；出苗期（2018年10月30日）、抽穗期（2019 年4 月14 日）誤差較大，NDVI 監(jiān)測(cè)日期較人工觀測(cè)日期誤差分別為6、8 d；其他發(fā)育時(shí)期絕對(duì)誤差均≤2 d，平均誤差2.3 d。綜上，NDVI 時(shí)序數(shù)據(jù)監(jiān)測(cè)發(fā)育時(shí)期平均誤差為2.5 d。

表1 冬小麥發(fā)育時(shí)期人工與NDVI監(jiān)測(cè)結(jié)果對(duì)比Tab.1 Comparison of winter wheat growth stage between artificial and NDVI monitoring datec

續(xù)表1 冬小麥發(fā)育時(shí)期人工與NDVI監(jiān)測(cè)結(jié)果對(duì)比Tab.1（Continued） Comparison of winter wheat growth stage between artificial and NDVI monitoring datec

3 結(jié)論與討論

準(zhǔn)確、高效地獲取冬小麥發(fā)育時(shí)期信息，對(duì)于提高田間精細(xì)管理以及預(yù)測(cè)產(chǎn)量等具有重要意義。本研究利用SRS-NDVI測(cè)量?jī)x觀測(cè)了鶴壁農(nóng)業(yè)氣象試驗(yàn)站2017—2018、2018—2019年度冬小麥生長(zhǎng)季的NDVI 時(shí)序數(shù)據(jù)，采用鄰域差值分析法重構(gòu)NDVI時(shí)序，并運(yùn)用S-G 濾波法進(jìn)行平滑處理。以冬小麥生育累積天數(shù)第83天為界，分別使用廣義動(dòng)態(tài)閾值法、曲線速率法和極值法提取冬小麥關(guān)鍵發(fā)育時(shí)期的特征點(diǎn)。研究結(jié)果表明，鄰域差值分析法可有效去除觀測(cè)數(shù)據(jù)的明顯異常值，NDVI時(shí)序數(shù)據(jù)經(jīng)S-G濾波法處理后，與冬小麥的實(shí)際生長(zhǎng)過(guò)程更為吻合，為作物發(fā)育時(shí)期的準(zhǔn)確提取提供了更精準(zhǔn)的數(shù)據(jù)支持。

本研究監(jiān)測(cè)了冬小麥從出苗到成熟共計(jì)10 個(gè)發(fā)育時(shí)期類別，較衛(wèi)星遙感監(jiān)測(cè)能夠更為準(zhǔn)確、全面監(jiān)測(cè)冬小麥生長(zhǎng)季內(nèi)的關(guān)鍵發(fā)育時(shí)期，是由于近地遙感SRS-NDVI 觀測(cè)儀時(shí)間分辨率（30 min）明顯高于衛(wèi)星遙感（8～30 d），受大氣條件及云雨等天氣影響小，能夠監(jiān)測(cè)冬小麥NDVI 的細(xì)微變化規(guī)律，進(jìn)一步提高冬小麥NDVI 變化精度所致。與衛(wèi)星遙感監(jiān)測(cè)冬小麥在返青期植被指數(shù)增長(zhǎng)速率最大的結(jié)論不同[7?8]，本研究發(fā)現(xiàn)近地遙感監(jiān)測(cè)冬小麥起身期NDVI 增長(zhǎng)速率最大。同時(shí)，研究發(fā)現(xiàn)，越冬期NDVI 最低值與成熟期NDVI 值吻合度非常高，可以作為判斷冬小麥成熟期的依據(jù)。2018、2019 年冬小麥生長(zhǎng)季NDVI 時(shí)序數(shù)據(jù)監(jiān)測(cè)發(fā)育時(shí)期平均誤差分別為2.7、2.3 d。在已有的研究中，如楊琳等[5]利用MODIS NDVI 數(shù)據(jù)提取了江蘇省冬小麥物候期，識(shí)別結(jié)果平均誤差5.3～8.6 d；孔令寅等[7]基于關(guān)鍵發(fā)育期的冬小麥長(zhǎng)勢(shì)遙感監(jiān)測(cè)中，提取的全生育期平均絕對(duì)誤差為11.2 d；黃健熙等[8]基于遙感和積溫?cái)?shù)據(jù)提取了冬小麥生育期，平均誤差為3.8～7.4 d。因此，與已有基于衛(wèi)星遙感監(jiān)測(cè)冬小麥發(fā)育時(shí)期的研究相比，本研究采用近地遙感監(jiān)測(cè)冬小麥發(fā)育時(shí)期精度明顯提升。本研究方法能夠滿足冬小麥發(fā)育時(shí)期觀測(cè)要求，可為冬小麥發(fā)育時(shí)期的自動(dòng)化觀測(cè)方法研究提供新思路。由于數(shù)據(jù)限制，本研究?jī)H對(duì)單一站點(diǎn)的冬小麥關(guān)鍵生育時(shí)期監(jiān)測(cè)方法進(jìn)行了研究；而且不同方法對(duì)不同區(qū)域冬小麥生育期監(jiān)測(cè)的影響敏感程度不同，監(jiān)測(cè)方法的普適性未得到充分且有效驗(yàn)證，今后仍需進(jìn)一步驗(yàn)證與分析。