NDVI时序曲线形状相似性模型的水稻提取方法

许玉萍，刘鹏程，秦自成，伍晓阳

（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就无法计算陆地表面水指数等。本文将基于多时相植被指数值构建光谱时相曲线，利用傅立叶形状相似性识别模型来判别水稻种植区域。

1 形状相似性距离模型提取方法

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数据提取水稻信息空间分布图

2 精度评价

经统计，利用形状相似性距离模型提取出的水稻种植面积是10 427.8 km2，土地利用类型图上水稻种植面积为11 412.8 km2，二者误差为8.6%。从土地利用类型图中随机抽取部分样点，对结果作精度分析（表1），总体精度为80%，Kappa系数为94.4%。

表 1 分类精度评价/%

为了验证模型的实用性，选取了一块典型水稻种植区（如图8），用叠置算法比较利用形状相似性距离模型提取出的水稻种植区域与土地利用类型图中水稻种植区域，图9是通过空间叠加得到的同为水稻的像元分布图。统计得到，两者同为水稻种植的面积为636 km2，分别占形状相似性距离模型提取水稻面积和土地利用图水稻面积的57%和61%。

图8 典型区基于形状相似性距离模型水稻信息提取与土地利用图水稻空间分布对比

图9 典型区同为水稻信息空间分布图

3 结 语

本研究通过计算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）；中央高校基本科研资助项目（CCNU15A02004）。