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

許玉萍，劉鵬程，秦自成，伍曉陽

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

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

許玉萍1,2，劉鵬程1,2，秦自成1,2，伍曉陽1,2

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

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

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

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

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

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

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

圖1 技術(shù)流程圖

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

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

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

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

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

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

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

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

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

1.2 研究區(qū)概況

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

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

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

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

本文采用2010年MoDIS的16 d合成植被指數(shù)（MoD13Q1），數(shù)據(jù)產(chǎn)品從NASA數(shù)據(jù)站https：// ladsweb.nascom.nasa.gov/data/下載獲得，覆蓋江漢平原的影像軌道號為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)水稻生長曲線，如圖4所示。

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

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

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

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

NDVI時(shí)序曲線是一條未封閉的曲線，為了能夠使用傅立葉形狀描述子模型進(jìn)行相似性分析，對其以首尾點(diǎn)連線并作鏡像處理，構(gòu)成封閉曲線，鏡像部分與原要素形狀完全相同，可視為原要素的對偶形狀[17]。圖5為標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線經(jīng)過鏡像處理后的多邊形以及對不同NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻NDVI時(shí)序曲線相似性比較示意圖，圖中藍(lán)色多邊形為NDVI時(shí)序曲線與標(biāo)準(zhǔn)水稻生長曲線形狀相似距離為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)水稻生長曲線形狀相似性距離

1.6 形狀相似性距離閾值確定及地物重分類

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

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

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

2 精度評價(jià)

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

表 1 分類精度評價(jià)/%

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

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

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

3 結(jié) 語

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

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

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

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

本研究尚有一些不足之處，如MoDIS時(shí)序影像雖然利用算法減弱了云霧等噪聲對地面真實(shí)信息的影響，但是因?yàn)樘幚矸椒ㄝ^為粗糙，難以達(dá)到理想效果；選取的土地利用類型圖的數(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] 鄭長春,王秀珍,黃敬峰.多時(shí)相MoDIS影像的浙江省水稻種植面積信息提取方法研究[J].浙江大學(xué)學(xué)報(bào)(農(nóng)業(yè)與生命科學(xué)版),2008,35(1)：98-104

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

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

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

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

[10] 顧曉鶴,韓立建,張錦水,等.基于相似性分析的中低分辨率復(fù)合水稻種植面積測量法[J].中國農(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] 艾廷華,帥赟,李精忠.基于形狀相似性識別的空間查詢[J].測繪學(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] 鄧愛娟,劉敏,萬素琴,等.湖北省雙季稻生長季降水及洪澇變化特征[J].長江流域資源與環(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] 劉鵬程,羅靜,艾廷華,等.基于線要素綜合的形狀相似性評價(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)目來源：國家自然科學(xué)基金資助項(xiàng)目（41371183、41531180）；中央高校基本科研資助項(xiàng)目（CCNU15A02004）。