許玉萍，劉鵬程，秦自成，伍曉陽(yáng)

（1.地理過(guò)程分析與模擬湖省重點(diǎn)實(shí)驗(yàn)室，湖武漢430079；2.華中師大學(xué) 城市與環(huán)境科學(xué)學(xué)院，湖武漢430079）

NDVI時(shí)序曲線形狀相似性模型的水稻提取方法

許玉萍1,2，劉鵬程1,2，秦自成1,2，伍曉陽(yáng)1,2

基于土地利用現(xiàn)狀圖與經(jīng)驗(yàn)觀測(cè)，提取標(biāo)準(zhǔn)水稻NDVI 時(shí)序曲線，利用傅立葉形狀描述子計(jì)算MoDIS-NDVI時(shí)序曲線與標(biāo)準(zhǔn)的水稻NDVI時(shí)序曲線形狀相似性距離，通過(guò)樣例數(shù)據(jù)探測(cè)未知像元與樣本的相似性距離閾值，從而判別雙季水稻種植區(qū)域。以江漢平原2010年的數(shù)據(jù)進(jìn)行實(shí)驗(yàn)，證明此方法識(shí)別的雙季水稻種植區(qū)域面積誤差為8.6%，總體精度為80%，較為理想。該方法將遙感光譜信息與幾何形狀的識(shí)別相結(jié)合，有效減少了個(gè)別時(shí)段光譜信息誤差引起的識(shí)別錯(cuò)誤，提高了識(shí)別水稻種植區(qū)域的有效性。

NDVI時(shí)序曲線；傅立葉形狀描述子；水稻種植區(qū)

水稻是世界上主要的糧食作物之一，在水稻主產(chǎn)區(qū)，水田灌溉往往要消耗超過(guò)80%的水資源[1]。稻田溫室氣體甲烷的排放超過(guò)大氣甲烷排放總量的10%[2]，是一種主要的溫室氣體排放源。及時(shí)、有效地獲取水稻種植空間分布信息，對(duì)我國(guó)農(nóng)業(yè)部門(mén)進(jìn)行科學(xué)指導(dǎo)及正確決策具有指導(dǎo)意義，也可以為水資源管理和溫室氣體排放估算提供依據(jù)[3-4]。

遙感技術(shù)的發(fā)展為大范圍水稻種植區(qū)的信息監(jiān)測(cè)提供了便捷、快速的應(yīng)用手段。鄭長(zhǎng)春等[5-7]利用MoDIS陸地表面反射率數(shù)據(jù)，結(jié)合水稻生長(zhǎng)期間物候特征，獲取研究區(qū)的NDVI、增強(qiáng)型植被指數(shù)、陸地表面水指數(shù)信息，構(gòu)建水稻種植提取模型。結(jié)果表明，此方法能取得精度較高的水稻種植數(shù)據(jù)。魏新彩等[8]利用HJ-1A/1B衛(wèi)星數(shù)據(jù)，采用類(lèi)似方法實(shí)現(xiàn)了水稻面積信息的提取，精度高于90%。張春桂等[9]利用MoDIS數(shù)據(jù)計(jì)算水稻葉面積指數(shù)LAI，實(shí)現(xiàn)了水稻面積的提取。顧曉鶴等[10]以SPoT5數(shù)據(jù)的水稻識(shí)別結(jié)果計(jì)算出圖像相似性指數(shù)，運(yùn)用SVM混合像元分解模型獲取較為精準(zhǔn)的樣本點(diǎn)，再對(duì)MoDIS時(shí)序數(shù)據(jù)進(jìn)行水稻面積測(cè)量，得到了較高的測(cè)量精度。

目前這方面的研究主要是通過(guò)不同地物在生長(zhǎng)過(guò)程中光譜特征的差異來(lái)計(jì)算NDVI、增強(qiáng)型植被指數(shù)、陸地表面水指數(shù)等，從而識(shí)別水稻種植區(qū)。而有些影像產(chǎn)品波段數(shù)較少，無(wú)法進(jìn)行植被指數(shù)計(jì)算，如250 m空間分辨率的MoD09就無(wú)法計(jì)算陸地表面水指數(shù)等。本文將基于多時(shí)相植被指數(shù)值構(gòu)建光譜時(shí)相曲線，利用傅立葉形狀相似性識(shí)別模型來(lái)判別水稻種植區(qū)域。

1 形狀相似性距離模型提取方法

MoDIS-NDVI時(shí)間序列數(shù)據(jù)是一組能反映地物光譜隨時(shí)間變化情況的數(shù)值，相同地物在NDVI時(shí)間序列曲線變化上具有幾何形態(tài)的相似性。不同地物由于物候特征和光譜特性的差異性導(dǎo)致了時(shí)序曲線的形狀不同，如水體和植被的NDVI時(shí)序曲線。由于水體沒(méi)有生長(zhǎng)期，其N(xiāo)DVI值基本不隨時(shí)間變化而變化，而植被的NDVI值則會(huì)隨時(shí)間發(fā)生變化。依據(jù)這一特性，以江漢平原雙季稻種植區(qū)為例，利用NDVI時(shí)序數(shù)據(jù)和野外觀測(cè)經(jīng)驗(yàn)提取出標(biāo)準(zhǔn)的水稻NDVI生長(zhǎng)曲線，并建立NDVI時(shí)序數(shù)據(jù)的傅立葉形狀數(shù)學(xué)模型，計(jì)算與標(biāo)準(zhǔn)的水稻NDVI生長(zhǎng)曲線之間的形狀相似性距離，再根據(jù)土地利用類(lèi)型圖中水稻過(guò)渡區(qū)的形狀相似性距離閾值進(jìn)行識(shí)別，最后再用土地利用現(xiàn)狀圖進(jìn)行精度評(píng)價(jià)，具體技術(shù)流程如圖1所示。

圖1 技術(shù)流程圖

1.1 傅立葉形狀相似性數(shù)學(xué)模型

傅立葉形狀描述子是分析和識(shí)別地物形狀的重要方法之一，在圖像匹配中有較多應(yīng)用[11-12]，其基本思想是：如果物體的輪廓形狀是一條封閉的曲線，其輪廓形狀可以描述為以其周長(zhǎng)為周期的函數(shù)，這個(gè)周期函數(shù)可以用傅立葉級(jí)數(shù)形式表示，傅立葉級(jí)數(shù)中的一系列的級(jí)數(shù)被稱(chēng)為傅立葉形狀描述子。

由點(diǎn)列P0，P1，…，PN構(gòu)成多邊形（如圖2），以P0為起始點(diǎn)，其輪廓上任意一點(diǎn)P（s）都可以表達(dá)為曲線長(zhǎng)度s的函數(shù)：

iP0的曲線長(zhǎng)；Si≤s≤Si+1，0≤i≤N-1；Pi坐標(biāo)為（xi， yi）。函數(shù)的傅立葉級(jí)數(shù)表達(dá)式為：

圖2 多邊形上動(dòng)點(diǎn)的函數(shù)表達(dá)

其中，L為封閉曲線的周長(zhǎng)。由于傅立葉級(jí)數(shù)的子項(xiàng)具有很強(qiáng)的收斂性，當(dāng)系數(shù)項(xiàng)取得足夠多的階次時(shí)，可以將物體的形狀信息完全提取并恢復(fù)出來(lái)[13]。系數(shù)的表達(dá)式為：

研究表明，傅立葉形狀描述子會(huì)受形狀、位置、尺寸等因素的影響，為消除這些影響，需對(duì)傅立葉級(jí)數(shù)系數(shù)取模后進(jìn)行歸一化處理，得出歸一化傅立葉描述子。歸一化傅立葉描述子d（i）定義為：

歸一化后的傅立葉描述子可以計(jì)算任意兩個(gè)形狀的相似程度，比較m、n兩個(gè)形狀的相似性可以用歐氏距離來(lái)度量歸一化傅立葉描述子間的形狀差異：

當(dāng)Dis=0時(shí)，認(rèn)為兩個(gè)形狀完全相似。Dis越大，物體的形狀差異越大。

1.2 研究區(qū)概況

江漢平原位于長(zhǎng)江流域中下游地區(qū)，該區(qū)域耕地面積為70%，其中水田面積近60%。水系眾多，雨熱同期，光熱充足，潮土和水稻土的比例占耕地面積的90%以上，適合水稻種植，且以雙季稻種植為主。

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

本文使用數(shù)據(jù)包括江漢平原矢量圖、中國(guó)1∶100 000土地利用類(lèi)型圖、MoD13Q1影像數(shù)據(jù)，其中中國(guó)1∶100 000比例尺土地利用現(xiàn)狀監(jiān)測(cè)數(shù)據(jù)來(lái)源于中國(guó)科學(xué)院資源環(huán)境科學(xué)數(shù)據(jù)中心（http：//www. resdc.cn）。將這些數(shù)據(jù)統(tǒng)一轉(zhuǎn)換成Albers等面積投影和WGS84坐標(biāo)系，利用江漢平原矢量圖進(jìn)行裁剪。

水稻在不同的生長(zhǎng)期，相應(yīng)的光譜特征也隨之變化，根據(jù)水稻不同生長(zhǎng)時(shí)期的光譜特征及物候特征來(lái)選擇一定時(shí)期的影像[14]。根據(jù)研究區(qū)物候歷，湖北省早稻生長(zhǎng)季時(shí)段為3月中旬～7月中旬，晚稻為6月中旬～10月下旬；考慮到水稻抽穗揚(yáng)花至成熟是影響水稻產(chǎn)量與品質(zhì)的關(guān)鍵時(shí)期，所以選取了水稻關(guān)鍵生長(zhǎng)期5月～9月的影像數(shù)據(jù)[15]。5月上旬為雙季早稻的移栽期，7月上旬到中旬為早稻成熟期，7月中旬晚稻開(kāi)始移栽，到9月中旬晚稻進(jìn)入成熟期。

本文采用2010年MoDIS的16 d合成植被指數(shù)（MoD13Q1），數(shù)據(jù)產(chǎn)品從NASA數(shù)據(jù)站https：// ladsweb.nascom.nasa.gov/data/下載獲得，覆蓋江漢平原的影像軌道號(hào)為h27v05、h27v06，空間分辨率為250 m，根據(jù)研究區(qū)水稻物候特征選取了4月下旬～9月中旬的10期影像數(shù)據(jù)。MoDIS數(shù)據(jù)為HDF格式，投影為正弦曲線投影。使用MRT軟件進(jìn)行影像的拼接，提取NDVI圖層并轉(zhuǎn)換為tiff數(shù)據(jù)格式，然后轉(zhuǎn)換為Albers等面積投影和WGS84坐標(biāo)系，最后進(jìn)行裁剪。

1.4 標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線提取

根據(jù)土地利用現(xiàn)狀圖與野外觀察經(jīng)驗(yàn)，在研究區(qū)內(nèi)選取20個(gè)典型的水稻種植點(diǎn)作為水稻純像元點(diǎn)，分布如圖3。從NDVI時(shí)序數(shù)據(jù)中提取這20個(gè)點(diǎn)各時(shí)期的NDVI值，求取各期NDVI值的平均值，得到NDVI均值連成曲線作為標(biāo)準(zhǔn)水稻生長(zhǎng)曲線，如圖4所示。

圖3 水稻純像元樣本點(diǎn)分布

圖4 MODIS-NDVI時(shí)間序列的標(biāo)準(zhǔn)水稻生長(zhǎng)曲線

從標(biāo)準(zhǔn)水稻生長(zhǎng)曲線圖中可以看出2個(gè)時(shí)期的NDVI高峰值——6月中旬和8月初，體現(xiàn)出研究區(qū)水稻種植的主要特征，即研究區(qū)內(nèi)雙季稻早稻從4月下旬開(kāi)始移栽到水田，因此NDVI值較低，隨著水稻的返青，NDVI值逐漸增加。稻田時(shí)間序列光譜反射率顯示，NDVI最大值出現(xiàn)在水稻抽穗期，原因是水稻開(kāi)始從營(yíng)養(yǎng)生長(zhǎng)階段轉(zhuǎn)向生殖生長(zhǎng)階段[16]，即對(duì)應(yīng)早稻生長(zhǎng)期的6月初到中旬；之后水稻葉片開(kāi)始衰老枯黃，生物量下降，NDVI值隨之下降，直到早稻7月初成熟收割，NDVI值降到最低；隨著晚稻的移栽和逐漸生長(zhǎng)，NDVI值相應(yīng)地再次經(jīng)歷上升，一直到8月初抽穗期的高峰后再開(kāi)始下降。

1.5 NDVI時(shí)序曲線形狀描述及相似性距離計(jì)算

NDVI時(shí)序曲線是一條未封閉的曲線，為了能夠使用傅立葉形狀描述子模型進(jìn)行相似性分析，對(duì)其以首尾點(diǎn)連線并作鏡像處理，構(gòu)成封閉曲線，鏡像部分與原要素形狀完全相同，可視為原要素的對(duì)偶形狀[17]。圖5為標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線經(jīng)過(guò)鏡像處理后的多邊形以及對(duì)不同NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線相似性比較示意圖，圖中藍(lán)色多邊形為NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻生長(zhǎng)曲線形狀相似距離為0.044 008的鏡像處理，黃色多邊形相似距離為0.033 606，灰色多邊形相似距離為0.026 318，紅色多邊形相似距離為0.015 237。圖中可以看出，兩個(gè)多邊形形狀越接近，相似性距離Dis值越小。根據(jù)模型計(jì)算后得到的NDVI時(shí)序曲線形狀相似性距離如圖6所示。

圖5 不同NDVI時(shí)序曲線鏡像處理后形狀相似性比較示意圖

圖6 NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻生長(zhǎng)曲線形狀相似性距離

1.6 形狀相似性距離閾值確定及地物重分類(lèi)

研究區(qū)內(nèi)水稻與非水稻的過(guò)渡區(qū)域?qū)嶋H上是其他容易與水稻混淆的地物，這類(lèi)地物的NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線形狀相似程度介于相似距離上限值與下限值之間。在這些過(guò)渡區(qū)域選取樣本點(diǎn)，提取NDVI時(shí)序曲線，計(jì)算其與標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線的相似性距離，并選取中位數(shù)作為水稻分類(lèi)閾值。

在土地利用現(xiàn)狀圖上選取26個(gè)水稻與非水稻過(guò)渡區(qū)域的樣本點(diǎn)，提取NDVI時(shí)序曲線值。計(jì)算這26條NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻生長(zhǎng)曲線的相似性距離，為了不受極大極小值對(duì)數(shù)據(jù)的影響，取26個(gè)形狀相似性距離的中位數(shù)0.023 92作為閾值。利用閾值對(duì)地物進(jìn)行重分類(lèi)后的結(jié)果如圖7所示。

圖7 MODIS數(shù)據(jù)提取水稻信息空間分布圖

2 精度評(píng)價(jià)

經(jīng)統(tǒng)計(jì)，利用形狀相似性距離模型提取出的水稻種植面積是10 427.8 km2，土地利用類(lèi)型圖上水稻種植面積為11 412.8 km2，二者誤差為8.6%。從土地利用類(lèi)型圖中隨機(jī)抽取部分樣點(diǎn)，對(duì)結(jié)果作精度分析（表1），總體精度為80%，Kappa系數(shù)為94.4%。

表 1 分類(lèi)精度評(píng)價(jià)/%

為了驗(yàn)證模型的實(shí)用性，選取了一塊典型水稻種植區(qū)（如圖8），用疊置算法比較利用形狀相似性距離模型提取出的水稻種植區(qū)域與土地利用類(lèi)型圖中水稻種植區(qū)域，圖9是通過(guò)空間疊加得到的同為水稻的像元分布圖。統(tǒng)計(jì)得到，兩者同為水稻種植的面積為636 km2，分別占形狀相似性距離模型提取水稻面積和土地利用圖水稻面積的57%和61%。

圖8 典型區(qū)基于形狀相似性距離模型水稻信息提取與土地利用圖水稻空間分布對(duì)比

圖9 典型區(qū)同為水稻信息空間分布圖

3 結(jié) 語(yǔ)

本研究通過(guò)計(jì)算NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻生長(zhǎng)曲線形狀相似性距離，探測(cè)水稻與非水稻時(shí)序曲線形狀相似性距離閾值，以提取水稻種植區(qū)域。這種方法能夠充分利用MoDIS時(shí)間序列優(yōu)勢(shì)，在某些遙感影像產(chǎn)品波段數(shù)較少，無(wú)法進(jìn)行其他植被指數(shù)運(yùn)算的情況下，提取水稻的種植信息；而且通過(guò)傅立葉形狀分析能夠減少個(gè)別時(shí)相的植被指數(shù)的粗差對(duì)水稻種植區(qū)域提取的影響。如果水稻重要生長(zhǎng)階段的NDVI值出現(xiàn)異常，光譜信息的粗差會(huì)使基于光譜特征提取的水稻種植區(qū)域信息產(chǎn)生錯(cuò)誤。圖10中水稻NDVI時(shí)序曲線在晚稻移栽時(shí)期（7月17日）的NDVI值達(dá)到了0.8以上，與正常水稻移栽期的NDVI值不符，但其計(jì)算的與標(biāo)準(zhǔn)水稻生長(zhǎng)曲線形狀相似性距離為0.015 237，能夠從概率上識(shí)別為水稻地物?？梢?jiàn)，在計(jì)算形狀相似性距離值時(shí)不會(huì)因個(gè)別時(shí)段的形狀相似性距離值偏大而偏大，而是在整體的計(jì)算值上進(jìn)行比較，此方法能夠減少水稻生長(zhǎng)期個(gè)別時(shí)段光譜信息粗差引起的識(shí)別錯(cuò)誤。

圖10 一個(gè)NDVI時(shí)序曲線鏡像處理后的形狀相似性示意

從提取后的水稻種植面積精度分析來(lái)看，研究區(qū)內(nèi)湖泊水系分布較密集的區(qū)域水稻面積信息提取精度較低，這與湖泊水系分割水稻田，迫使水稻田斑塊變小而變得較難識(shí)別有關(guān)。

利用形狀相似性距離模型提取的水稻種植面積與土地利用類(lèi)型圖中的水稻面積信息結(jié)果較為一致，二者誤差在8.6%，提取的總體精度達(dá)到80%，較為理想，可以作為以后水稻種植面積信息提取的新方法。

本研究尚有一些不足之處，如MoDIS時(shí)序影像雖然利用算法減弱了云霧等噪聲對(duì)地面真實(shí)信息的影響，但是因?yàn)樘幚矸椒ㄝ^為粗糙，難以達(dá)到理想效果；選取的土地利用類(lèi)型圖的數(shù)據(jù)本身也存在誤差，今后可選取更高分辨率的影像圖作為地面真實(shí)信息。

[1] FAoSTAT.Statistical Database of the Food and Agricultural organization of the United Nations[R]. 2001

[2] PRATHER M, EHHALT D. Atmospheric Chemistry and Greenhouse Gases[M]. Cambridge, U.K.： Cambridge University, 2001

[3] JANG M W, CHoI J Y, LEE J J. A Spatial Reasoning Approach to Estimating Paddy rice Water Demand in Hwanghaenam-do,North Korea[J].Agricultural Water Management,2007, 89(3)：185-198

[4] GAo J F, PAN G X, JIANG X S, et al. Land-use Induced Changes in Topsoil organic Carbon Stock of Paddy Fields Using MoDIS and TM/ETM Analysis： A Case Study of Wujiang County,China[J]. Journal of Environmental Sciences, 2008, 20(7)： 852-858

[5] 鄭長(zhǎng)春,王秀珍,黃敬峰.多時(shí)相MoDIS影像的浙江省水稻種植面積信息提取方法研究[J].浙江大學(xué)學(xué)報(bào)(農(nóng)業(yè)與生命科學(xué)版),2008,35(1)：98-104

[6] 馮銳,張玉書(shū),錢(qián)永蘭,等.基于多時(shí)相MoDIS數(shù)據(jù)的東北地區(qū)一季稻面積提取[J].生態(tài)學(xué)雜志,2011,30(11)：2 570-2 576

[7] 張莉,吳文斌,左麗君,等.基于EoS/MoDIS數(shù)據(jù)的南方水稻面積提取技術(shù)[J].中國(guó)農(nóng)業(yè)資源與區(qū)劃,2011,32(4)：39-44

[8] 魏新彩,王新生,劉海,等.HJ衛(wèi)星圖像水稻種植面積的識(shí)別分析[J].地球信息科學(xué)學(xué)報(bào),2012,14(3)：382-388

[9] 張春桂,林晶,吳振海,等.基于MoDIS數(shù)據(jù)的水稻種植面積監(jiān)測(cè)方法研究[J].自然資源學(xué)報(bào),2007,22(1)：1-8

[10] 顧曉鶴,韓立建,張錦水,等.基于相似性分析的中低分辨率復(fù)合水稻種植面積測(cè)量法[J].中國(guó)農(nóng)業(yè)科學(xué),2008,41(4)：978-985

[11] ZAHN C T, RoSKIES R Z. Fourier Descriptors for Plane Closed Curves[J].IEEE Transactions on Computers,1972,21(3)：269-281

[12] LIU P, LI X, LIU W, et al. Fourier-based Multi-Scale Representation and Progressive Transmission of Cartographic Curves on the Internet[J]. Cartography and Geographic Information Science, 2015：1-15

[13] 艾廷華,帥赟,李精忠.基于形狀相似性識(shí)別的空間查詢(xún)[J].測(cè)繪學(xué)報(bào),2009,38(4)：356-362

[14] 于文穎,馮銳,紀(jì)瑞鵬,等.基于MoDIS數(shù)據(jù)的水稻種植面積提取研究進(jìn)展[J].氣象與環(huán)境學(xué)報(bào),2011,27(2)：56-61

[15] 鄧愛(ài)娟,劉敏,萬(wàn)素琴,等.湖北省雙季稻生長(zhǎng)季降水及洪澇變化特征[J].長(zhǎng)江流域資源與環(huán)境,2012(增刊)：173-178

[16] SHIBAYAMA M,AKIYAMA T. Seasonal Visible, Near-infrared and Mid-infrared Spectra of Rice Canopies in Relation to LAI and Aboveground Dry Phytomass[J]. Remote Sensing of Environment, 1989, 27(2)： 119-127

[17] 劉鵬程,羅靜,艾廷華,等.基于線要素綜合的形狀相似性評(píng)價(jià)模型[J].武漢大學(xué)學(xué)報(bào)（信息科學(xué)版）,2012,37(1)：114-117

P237.3

B

1672-4623（2016）08-0056-05

10.3969/j.issn.1672-4623.2016.08.019

許玉萍，碩士研究生，研究方向?yàn)榭臻g數(shù)據(jù)挖掘與可視化表達(dá)。

2015-12-04。

項(xiàng)目來(lái)源：國(guó)家自然科學(xué)基金資助項(xiàng)目（41371183、41531180）；中央高?；究蒲匈Y助項(xiàng)目（CCNU15A02004）。