許玉萍,劉鵬程,秦自成,伍曉陽
(1.地理過程分析與模擬湖省重點實驗室,湖武漢430079;2.華中師大學 城市與環境科學學院,湖武漢430079)
NDVI時序曲線形狀相似性模型的水稻提取方法
許玉萍1,2,劉鵬程1,2,秦自成1,2,伍曉陽1,2

基于土地利用現狀圖與經驗觀測,提取標準水稻NDVI 時序曲線,利用傅立葉形狀描述子計算MoDIS-NDVI時序曲線與標準的水稻NDVI時序曲線形狀相似性距離,通過樣例數據探測未知像元與樣本的相似性距離閾值,從而判別雙季水稻種植區域。以江漢平原2010年的數據進行實驗,證明此方法識別的雙季水稻種植區域面積誤差為8.6%,總體精度為80%,較為理想。該方法將遙感光譜信息與幾何形狀的識別相結合,有效減少了個別時段光譜信息誤差引起的識別錯誤,提高了識別水稻種植區域的有效性。
NDVI時序曲線;傅立葉形狀描述子;水稻種植區
水稻是世界上主要的糧食作物之一,在水稻主產區,水田灌溉往往要消耗超過80%的水資源[1]。稻田溫室氣體甲烷的排放超過大氣甲烷排放總量的10%[2],是一種主要的溫室氣體排放源。及時、有效地獲取水稻種植空間分布信息,對我國農業部門進行科學指導及正確決策具有指導意義,也可以為水資源管理和溫室氣體排放估算提供依據[3-4]。
遙感技術的發展為大范圍水稻種植區的信息監測提供了便捷、快速的應用手段。鄭長春等[5-7]利用MoDIS陸地表面反射率數據,結合水稻生長期間物候特征,獲取研究區的NDVI、增強型植被指數、陸地表面水指數信息,構建水稻種植提取模型。結果表明,此方法能取得精度較高的水稻種植數據。魏新彩等[8]利用HJ-1A/1B衛星數據,采用類似方法實現了水稻面積信息的提取,精度高于90%。張春桂等[9]利用MoDIS數據計算水稻葉面積指數LAI,實現了水稻面積的提取。顧曉鶴等[10]以SPoT5數據的水稻識別結果計算出圖像相似性指數,運用SVM混合像元分解模型獲取較為精準的樣本點,再對MoDIS時序數據進行水稻面積測量,得到了較高的測量精度。
目前這方面的研究主要是通過不同地物在生長過程中光譜特征的差異來計算NDVI、增強型植被指數、陸地表面水指數等,從而識別水稻種植區。而有些影像產品波段數較少,無法進行植被指數計算,如250 m空間分辨率的MoD09就無法計算陸地表面水指數等。本文將基于多時相植被指數值構建光譜時相曲線,利用傅立葉形狀相似性識別模型來判別水稻種植區域。
MoDIS-NDVI時間序列數據是一組能反映地物光譜隨時間變化情況的數值,相同地物在NDVI時間序列曲線變化上具有幾何形態的相似性。不同地物由于物候特征和光譜特性的差異性導致了時序曲線的形狀不同,如水體和植被的NDVI時序曲線。由于水體沒有生長期,其NDVI值基本不隨時間變化而變化,而植被的NDVI值則會隨時間發生變化。依據這一特性,以江漢平原雙季稻種植區為例,利用NDVI時序數據和野外觀測經驗提取出標準的水稻NDVI生長曲線,并建立NDVI時序數據的傅立葉形狀數學模型,計算與標準的水稻NDVI生長曲線之間的形狀相似性距離,再根據土地利用類型圖中水稻過渡區的形狀相似性距離閾值進行識別,最后再用土地利用現狀圖進行精度評價,具體技術流程如圖1所示。

圖1 技術流程圖
1.1 傅立葉形狀相似性數學模型
傅立葉形狀描述子是分析和識別地物形狀的重要方法之一,在圖像匹配中有較多應用[11-12],其基本思想是:如果物體的輪廓形狀是一條封閉的曲線,其輪廓形狀可以描述為以其周長為周期的函數,這個周期函數可以用傅立葉級數形式表示,傅立葉級數中的一系列的級數被稱為傅立葉形狀描述子。
由點列P0,P1,…,PN構成多邊形(如圖2),以P0為起始點,其輪廓上任意一點P(s)都可以表達為曲線長度s的函數:

iP0的曲線長;Si≤s≤Si+1,0≤i≤N-1;Pi坐標為(xi, yi)。函數的傅立葉級數表達式為:


圖2 多邊形上動點的函數表達
其中,L為封閉曲線的周長。由于傅立葉級數的子項具有很強的收斂性,當系數項取得足夠多的階次時,可以將物體的形狀信息完全提取并恢復出來[13]。系數的表達式為:

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

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

當Dis=0時,認為兩個形狀完全相似。Dis越大,物體的形狀差異越大。
1.2 研究區概況
江漢平原位于長江流域中下游地區,該區域耕地面積為70%,其中水田面積近60%。水系眾多,雨熱同期,光熱充足,潮土和水稻土的比例占耕地面積的90%以上,適合水稻種植,且以雙季稻種植為主。
1.3 數據獲取與處理
本文使用數據包括江漢平原矢量圖、中國1∶100 000土地利用類型圖、MoD13Q1影像數據,其中中國1∶100 000比例尺土地利用現狀監測數據來源于中國科學院資源環境科學數據中心(http://www. resdc.cn)。將這些數據統一轉換成Albers等面積投影和WGS84坐標系,利用江漢平原矢量圖進行裁剪。
水稻在不同的生長期,相應的光譜特征也隨之變化,根據水稻不同生長時期的光譜特征及物候特征來選擇一定時期的影像[14]。根據研究區物候歷,湖北省早稻生長季時段為3月中旬~7月中旬,晚稻為6月中旬~10月下旬;考慮到水稻抽穗揚花至成熟是影響水稻產量與品質的關鍵時期,所以選取了水稻關鍵生長期5月~9月的影像數據[15]。5月上旬為雙季早稻的移栽期,7月上旬到中旬為早稻成熟期,7月中旬晚稻開始移栽,到9月中旬晚稻進入成熟期。
本文采用2010年MoDIS的16 d合成植被指數(MoD13Q1),數據產品從NASA數據站https:// ladsweb.nascom.nasa.gov/data/下載獲得,覆蓋江漢平原的影像軌道號為h27v05、h27v06,空間分辨率為250 m,根據研究區水稻物候特征選取了4月下旬~9月中旬的10期影像數據。MoDIS數據為HDF格式,投影為正弦曲線投影。使用MRT軟件進行影像的拼接,提取NDVI圖層并轉換為tiff數據格式,然后轉換為Albers等面積投影和WGS84坐標系,最后進行裁剪。
1.4 標準水稻NDVI時序曲線提取
根據土地利用現狀圖與野外觀察經驗,在研究區內選取20個典型的水稻種植點作為水稻純像元點,分布如圖3。從NDVI時序數據中提取這20個點各時期的NDVI值,求取各期NDVI值的平均值,得到NDVI均值連成曲線作為標準水稻生長曲線,如圖4所示。

圖3 水稻純像元樣本點分布

圖4 MODIS-NDVI時間序列的標準水稻生長曲線
從標準水稻生長曲線圖中可以看出2個時期的NDVI高峰值——6月中旬和8月初,體現出研究區水稻種植的主要特征,即研究區內雙季稻早稻從4月下旬開始移栽到水田,因此NDVI值較低,隨著水稻的返青,NDVI值逐漸增加。稻田時間序列光譜反射率顯示,NDVI最大值出現在水稻抽穗期,原因是水稻開始從營養生長階段轉向生殖生長階段[16],即對應早稻生長期的6月初到中旬;之后水稻葉片開始衰老枯黃,生物量下降,NDVI值隨之下降,直到早稻7月初成熟收割,NDVI值降到最低;隨著晚稻的移栽和逐漸生長,NDVI值相應地再次經歷上升,一直到8月初抽穗期的高峰后再開始下降。
1.5 NDVI時序曲線形狀描述及相似性距離計算
NDVI時序曲線是一條未封閉的曲線,為了能夠使用傅立葉形狀描述子模型進行相似性分析,對其以首尾點連線并作鏡像處理,構成封閉曲線,鏡像部分與原要素形狀完全相同,可視為原要素的對偶形狀[17]。圖5為標準水稻NDVI時序曲線經過鏡像處理后的多邊形以及對不同NDVI時序曲線與標準水稻NDVI時序曲線相似性比較示意圖,圖中藍色多邊形為NDVI時序曲線與標準水稻生長曲線形狀相似距離為0.044 008的鏡像處理,黃色多邊形相似距離為0.033 606,灰色多邊形相似距離為0.026 318,紅色多邊形相似距離為0.015 237。圖中可以看出,兩個多邊形形狀越接近,相似性距離Dis值越小。根據模型計算后得到的NDVI時序曲線形狀相似性距離如圖6所示。

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

圖6 NDVI時序曲線與標準水稻生長曲線形狀相似性距離
1.6 形狀相似性距離閾值確定及地物重分類
研究區內水稻與非水稻的過渡區域實際上是其他容易與水稻混淆的地物,這類地物的NDVI時序曲線與標準水稻NDVI時序曲線形狀相似程度介于相似距離上限值與下限值之間。在這些過渡區域選取樣本點,提取NDVI時序曲線,計算其與標準水稻NDVI時序曲線的相似性距離,并選取中位數作為水稻分類閾值。
在土地利用現狀圖上選取26個水稻與非水稻過渡區域的樣本點,提取NDVI時序曲線值。計算這26條NDVI時序曲線與標準水稻生長曲線的相似性距離,為了不受極大極小值對數據的影響,取26個形狀相似性距離的中位數0.023 92作為閾值。利用閾值對地物進行重分類后的結果如圖7所示。

圖7 MODIS數據提取水稻信息空間分布圖
經統計,利用形狀相似性距離模型提取出的水稻種植面積是10 427.8 km2,土地利用類型圖上水稻種植面積為11 412.8 km2,二者誤差為8.6%。從土地利用類型圖中隨機抽取部分樣點,對結果作精度分析(表1),總體精度為80%,Kappa系數為94.4%。

表 1 分類精度評價/%
為了驗證模型的實用性,選取了一塊典型水稻種植區(如圖8),用疊置算法比較利用形狀相似性距離模型提取出的水稻種植區域與土地利用類型圖中水稻種植區域,圖9是通過空間疊加得到的同為水稻的像元分布圖。統計得到,兩者同為水稻種植的面積為636 km2,分別占形狀相似性距離模型提取水稻面積和土地利用圖水稻面積的57%和61%。

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

圖9 典型區同為水稻信息空間分布圖
本研究通過計算NDVI時序曲線與標準水稻生長曲線形狀相似性距離,探測水稻與非水稻時序曲線形狀相似性距離閾值,以提取水稻種植區域。這種方法能夠充分利用MoDIS時間序列優勢,在某些遙感影像產品波段數較少,無法進行其他植被指數運算的情況下,提取水稻的種植信息;而且通過傅立葉形狀分析能夠減少個別時相的植被指數的粗差對水稻種植區域提取的影響。如果水稻重要生長階段的NDVI值出現異常,光譜信息的粗差會使基于光譜特征提取的水稻種植區域信息產生錯誤。圖10中水稻NDVI時序曲線在晚稻移栽時期(7月17日)的NDVI值達到了0.8以上,與正常水稻移栽期的NDVI值不符,但其計算的與標準水稻生長曲線形狀相似性距離為0.015 237,能夠從概率上識別為水稻地物??梢?,在計算形狀相似性距離值時不會因個別時段的形狀相似性距離值偏大而偏大,而是在整體的計算值上進行比較,此方法能夠減少水稻生長期個別時段光譜信息粗差引起的識別錯誤。

圖10 一個NDVI時序曲線鏡像處理后的形狀相似性示意
從提取后的水稻種植面積精度分析來看,研究區內湖泊水系分布較密集的區域水稻面積信息提取精度較低,這與湖泊水系分割水稻田,迫使水稻田斑塊變小而變得較難識別有關。
利用形狀相似性距離模型提取的水稻種植面積與土地利用類型圖中的水稻面積信息結果較為一致,二者誤差在8.6%,提取的總體精度達到80%,較為理想,可以作為以后水稻種植面積信息提取的新方法。
本研究尚有一些不足之處,如MoDIS時序影像雖然利用算法減弱了云霧等噪聲對地面真實信息的影響,但是因為處理方法較為粗糙,難以達到理想效果;選取的土地利用類型圖的數據本身也存在誤差,今后可選取更高分辨率的影像圖作為地面真實信息。
[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] 鄭長春,王秀珍,黃敬峰.多時相MoDIS影像的浙江省水稻種植面積信息提取方法研究[J].浙江大學學報(農業與生命科學版),2008,35(1):98-104
[6] 馮銳,張玉書,錢永蘭,等.基于多時相MoDIS數據的東北地區一季稻面積提取[J].生態學雜志,2011,30(11):2 570-2 576
[7] 張莉,吳文斌,左麗君,等.基于EoS/MoDIS數據的南方水稻面積提取技術[J].中國農業資源與區劃,2011,32(4):39-44
[8] 魏新彩,王新生,劉海,等.HJ衛星圖像水稻種植面積的識別分析[J].地球信息科學學報,2012,14(3):382-388
[9] 張春桂,林晶,吳振海,等.基于MoDIS數據的水稻種植面積監測方法研究[J].自然資源學報,2007,22(1):1-8
[10] 顧曉鶴,韓立建,張錦水,等.基于相似性分析的中低分辨率復合水稻種植面積測量法[J].中國農業科學,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].測繪學報,2009,38(4):356-362
[14] 于文穎,馮銳,紀瑞鵬,等.基于MoDIS數據的水稻種植面積提取研究進展[J].氣象與環境學報,2011,27(2):56-61
[15] 鄧愛娟,劉敏,萬素琴,等.湖北省雙季稻生長季降水及洪澇變化特征[J].長江流域資源與環境,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] 劉鵬程,羅靜,艾廷華,等.基于線要素綜合的形狀相似性評價模型[J].武漢大學學報(信息科學版),2012,37(1):114-117
P237.3
B
1672-4623(2016)08-0056-05
10.3969/j.issn.1672-4623.2016.08.019
許玉萍,碩士研究生,研究方向為空間數據挖掘與可視化表達。
2015-12-04。
項目來源:國家自然科學基金資助項目(41371183、41531180);中央高?;究蒲匈Y助項目(CCNU15A02004)。