CN110363246A - 一种高时空分辨率植被指数ndvi的融合方法 - Google Patents

一种高时空分辨率植被指数ndvi的融合方法 Download PDF

Info

Publication number
CN110363246A
CN110363246A CN201910648630.2A CN201910648630A CN110363246A CN 110363246 A CN110363246 A CN 110363246A CN 201910648630 A CN201910648630 A CN 201910648630A CN 110363246 A CN110363246 A CN 110363246A
Authority
CN
China
Prior art keywords
ndvi
vegetation index
image
increment
index ndvi
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201910648630.2A
Other languages
English (en)
Other versions
CN110363246B (zh
Inventor
贾艳艳
邢学刚
李吉英
邢学军
赵昕
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Binzhou University
Original Assignee
Binzhou University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Binzhou University filed Critical Binzhou University
Priority to CN201910648630.2A priority Critical patent/CN110363246B/zh
Publication of CN110363246A publication Critical patent/CN110363246A/zh
Application granted granted Critical
Publication of CN110363246B publication Critical patent/CN110363246B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/23Clustering techniques
    • G06F18/232Non-hierarchical techniques
    • G06F18/2321Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/50Image enhancement or restoration using two or more images, e.g. averaging or subtraction
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing

Landscapes

  • Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Probability & Statistics with Applications (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Image Processing (AREA)

Abstract

本发明公开了一种高时空分辨率植被指数NDVI的融合方法,包括:A.利用像元降尺度算法对低分辨率植被指数NDVI影像进行降尺度处理;B.利用降尺度处理的植被指数NDVI影像和高分辨率的植被指数NDVI影像计算增量△NDVI,并统计增量△NDVI累积分布直方图进行聚类构建出增量△NDVI分类图;C.利用增量ΔNDVI分类图,构建高低分辨率增量ΔNDVI的线性混合模型,并利用岭回归方法分别计算出高分辨率增量D.计算tp时刻的高分辨率植被指数NDVI值。本发明利用增量△NDVI聚类生成分类图,充分利用了NDVI变化特征,相比于其他方法,本发明方法的融合效果更令人满意,为构建高时空分辨率的植被指数NDVI提供了一种行之有效的方法。

Description

一种高时空分辨率植被指数NDVI的融合方法
技术领域
本发明属于遥感影像数据处理领域,涉及一种高时空分辨率植被指数NDVI的融合方法, 尤其涉及一种在高分辨率地物分类数据获取费时困难且融合时期内地物发生变化的高时空分 辨率植被指数NDVI的融合方法。
背景技术
植被指数是卫星可见光和近红外波段多种波段组合的产物,其能够简单、有效的度量地 表植被状况,其中归一化植被指数(NDVI)能够反映土壤、潮湿地面、雪、枯叶等植被冠层 背景的影响,对地表植被的覆盖程度非常敏感,是检测和指示植被覆盖状况和动态的常用指 标之一。高时空分辨率的NDVI影像数据在监测植被覆盖变化、作物生长状况、识别地物类 型、作物估产、生物量估算、地表蒸散、土壤湿度监测、气候变化等方面具有重要的作用。
专利CN 105046648 B公开了一种构建高时空遥感数据的方法,其背景技术详细阐述了目 前卫星传感器数据存在的空间分辨率、时间分辨率和光谱分辨率3者之间的相互矛盾,就目 前技术水平而言,获得同时满足高空间分辨率、高时间分辨率和高光谱分辨率的卫星遥感数 据是不现实的。例如:MODIS(Moderate-resolution imagingspectroradiometer)影像可实现全 球范围的观测并能捕捉到地表日变化信息,但其较低的空间分辨率(250m-1000m)难以满足 景观破碎度或异质性较强地区的应用;与之相反,Landsat卫星提供空间分辨率为30m的多 光谱影像,而该卫星的重返周期为16天,且受天气的影响,很难获得连续有效的遥感影像, 难以满足植被动态监测的要求。
近些年,综合Landsat影像的高空间分辨率和MODIS影像的时间分辨率优势,国内外学 者提出了一系列获取高时空分辨率数据的遥感数据融合模型。现有的时空数据融合方法通常 遵循以下三个基本原理:线性混合模型、空间相似性和时间相似性。
线性混合模型最基本的假设是变量的线性可加性。低分辨率影像像元的像元值是所包含 的高分辨率影像像元值的平均值与偏差之和,进而可以通过线性混合模型把低分辨率影像的 像元和高分辨率影像的像元联系在一起。其中,偏差是两个传感器的系统差异,可由波段宽、 太阳几何参数和观测角度等因素引起。光学卫星遥感原始影像记录的是地表物体的辐射率, 使用者常常会把辐射率转化为反射率。尽管各个光谱波段的辐射率和反射率具有线性相加性, 而从原始数据获取而来的高级别产品并不具有线性相加性,如一些植被指数(例如NDVI,EVI)和陆地表面温度,但最近研究表明在不同的尺度上,用线性混合模型联系这些产品并 不会产生显著的错误。所以已有一些方法也用线性混合模型来融合这些高级别产品,如专利 CN 102831310 B公开的一种构建高时空分辨率数据的方法可以融合MODIS和Landsat的 NDVI产品,该方法包含这一基本原理。
空间相似性是大部分时空融合算法的一条原理。根据Tobler的地理学第一定律,任何事 物均是相关的,但越近的事物越相关。遥感卫星记录陆地表面地物的光谱信息,在临近区域 是相同地物的情况下,我们可以假设这些地物具有相似的光谱特征,所以光谱相似的临近像 元存在空间相似性。例如,森林区域临近的像元具有相同的树种和相似的组成结构(如树高 和茎干密度),导致这些像元具有高度相似的光谱特征和高度一致的时间变化特性。通常采用 计算像元之间相似性权值或利用土地覆盖类型数据来描述像元的空间相似性。
时间相似性是指后一时刻像元光谱值和前一时刻像元光谱值是相似的,换句话说,后一 时刻的像元值可以通过前一时刻像元值的函数来估算。就如一条曲线可以用若干条线段近似 的描述一样,虽然光谱特征及高级别光谱数据产品在长时间尺度上是非线性变化的,但可以 用若干个线性变化近似描述该地物光谱特征变化。
基于以上三个基本原理建立的模型在各自领域都取得了较好的效果,在植被指数方面, CN 102831310B公开的模型比经典STARFM模型和ESTARFM模型的融合精度和效率更好, 该模型升级版NDVI-LMGM的论述及评价见2015年发表于remote sensing期刊的文献Animproved method for producing high spatial-resolution NDVI time seriesdatasets with multi-temporal MODIS NDVI data and Landsat TM/ETM+images。但CN102831310B公开模型 实施过程中需要获得高分辨率地物分类影像,虽然有阈值分类法、有监督分类和非监督分类 可实现影像地物分类,但总体精度不高,而人机交互的目视影像地物分类方法受人的主观影 响较大,并且耗时费力;而且该模型构建过程中未考虑融合期间内地物覆盖度发生的变化, 如农作物收割造成t0、t1两个时刻区间内地物发生突变。
本发明通过降尺度处理的t0、tp、t1三个时刻植被指数NDVIFC数据和t0、t1两个时刻高分 辨率植被指数NDVIF数据设计三种植被指数NDVI增量ΔNDVI,并统计三种植被指数NDVI 增量ΔNDVI累积分布直方图进行聚类构建出植被指数NDVI增量ΔNDVI的分类图,该分类图 充分考虑了高/低分辨率植被指数变化特征;给出了混合像元宽度wmix和解窗口宽度wcal计算 方法,减少了用户需要设定的参数数量,提高了算法健壮性;相比于其他方法,本发明方法 的融合效果更令人满意,为构建高时空分辨率的植被指数NDVI提供了一种行之有效的方法。
发明内容
本发明针对NDVI-LMGM模型需要高分辨率地物分类影像及未考虑融合期间内地物覆 盖发生变化的问题,首先对低分率影像进行降尺度处理,通过降尺度处理的t0、tp、t1三个时 刻植被指数NDVIFC数据和t0、t1两个时刻高分辨率植被指数NDVIF数据设计三种ΔNDVI增量, 统计三种植被指数NDVI增量ΔNDVI累积分布直方图进行聚类构建出增量ΔNDVI的分类图, 然后利用植被指数NDVI增量ΔNDVI分类图,构建高低分辨率增量ΔNDVI的线性混合模型, 并利用岭回归方法分别计算出高分辨率植被指数NDVI增量最后 计算得到tp时刻的高分辨率植被指数值,实现数据融合。
具体地,本发明提出了一种高时空分辨率植被指数NDVI的融合方法,来融合通过Landsat 和MODIS获得的植被指数NDVI影像,产生每期MODIS所对应高分辨率植物指数NDVI影像。
本发明的技术方案说明如下。
一种高时空分辨率植被指数NDVI的融合方法,其特征在于:所述方法主要包括以下步 骤:
步骤A:利用像元降尺度算法对t0、tp、t1三个时刻低分辨率植被指数NDVIC影像进行降 尺度处理;
步骤B:利用降尺度处理的t0、tp、t1三个时刻植被指数NDVIFC影像和t0、t1两个时刻高 分辨率植被指数NDVIF影像计算t0至tp区间降尺度植被指数NDVIFC增量影像、tp至t1区间降尺度植被指数NDVIFC增量影像、t0至t1区间高分辨率植被指数NDVIF增量影像,并根据三种植被指数NDVI增量ΔNDVI影像累积分布直方图确定聚类中 心进行聚类构建出增量ΔNDVI分类图;
步骤C:利用增量ΔNDVI分类图,构建高低分辨率植被指数NDVI增量ΔNDVI的线性混合 模型,并利用岭回归方法分别计算出高分辨率植被指数NDVI增量
步骤D:利用高分辨率植被指数NDVI增量计算tp时刻的高分 辨率植被指数值;
所述步骤B、步骤C和步骤D的ti至tj区间植被指数NDVI增量的计算公式如 下:
式(1)中,为ti时刻植被指数NDVI值,为tj时刻植被指数NDVI值,ti与tj分别 为计算增量的开始时刻和结束时刻;
所述步骤A至步骤D中的t0、tp、t1三个时刻分别为基准时刻、预测时刻、基准时刻,t0、tp、 t1三个时刻均具有低分辨率植被指数NDVIC影像,只有在t0、t1两个时刻具有高分辨率植被指 数NDVIF影像;
所述步骤A至步骤D中的植被指数NDVI与植被指数NDVI增量ΔNDVI的上标F、FC、 C分别代表高分辨率值、低分辨率降尺度处理后值和低分辨率值。
所述步骤A中的像元降尺度算法,包括如下步骤:
步骤AA:以t0时刻的高分辨率植被指数NDVIF影像为母版,分别创建t0、tp、t1三个时刻低分辨率植被指数NDVIC的降尺度影像模板;
步骤AB:如果t0、tp、t1三个时刻降尺度影像模板的像元值为空值,则像元值 不做修改,如果t0、tp、t1三个时刻降尺度影像模板的像元值为非空值,则用降尺 度公式进行重采样,公式如下:
式(2)中,下标(x,y,t)表示t时刻(x,y)位置像元的植被指数NDVICF值,N为计算降尺度影 像像元与低分辨率影像地理区域重叠的低分辨率影像像元个数,S为降尺度影像像元地理面积, 为第i个与计算降尺度影像像元重叠的低分辨率影像像元地理重叠面积,为第i 个与计算降尺度影像像元重叠的低分辨率影像像元植被指数NDVIC值。
所述步骤B利用三种植被指数NDVI增量ΔNDVI影像累积分布直方图确定聚类中心进行 聚类构建出增量ΔNDVI分类图,包括如下步骤:
步骤BA:分别计算三种植被指数NDVI增量ΔNDVI影 像累积分布直方图;
步骤BB:分别将三种植被指数NDVI增量ΔNDVI累积分布平均划分为NUM等分,取累积分布每等分二分之一处的植被指数NDVI增量ΔNDVIi值为聚类特征分量中心,其中 i=1,2,…,NUM,NUM分类数由用户根据影像变化特征确定;依次取 的聚类特征分量中心构成特征空间矢量点,其中i,j=1,2,…,NUM,形成NUM*NUM个矢量点,并标注类别号,即为 NUM*NUM个类别的聚类中心点;
步骤BC:按两种植被指数NDVI增量ΔNDVI影像对应位置像元 ΔNDVI值构成1个2维矢量,依次计算该矢量与NUM*NUM个类别的聚类中心点 的欧式距离,其中i,j=1,2,…,NUM,取距离最小者的类别号为该像元类别号, 所有像元的类别号计算完成后,则构建出t0至tp区间植被指数NDVI增量ΔNDVI的分类图;
步骤BD:按两种植被指数NDVI增量ΔNDVI影像对应位置像元 ΔNDVI值构成1个2维矢量,依次计算该矢量与NUM*NUM个类别的聚类中心点 的欧式距离,其中i,j=1,2,…,NUM,取距离最小者的类别号为该像元类别号, 所有像元的类别号计算完成后,则构建出tp至t1区间植被指数NDVI增量ΔNDVI的分类图。
所述步骤C利用植被指数NDVI增量ΔNDVI分类图,构建高低分辨率植被指数NDVI增 量ΔNDVI的线性混合模型,并利用岭回归方法求高分辨率植被指数NDVI增量包括如下步骤:
步骤CA:根据NUM分类数计算混合像元宽度wmix,计算公式如下:
wmix=floor(sqrt(α*NUM))+1 (3)
式(3)中,floor()为下取整函数,返回不大于计算值的最大整数,sqrt()为开方根函数,α为 修正系数,取值范围1.5至4.0;;
步骤CB:根据NUM分类数计算求解窗口宽度wcal,计算公式如下:
wcal=floor(sqrt(β*NUM))+1 (4)
式(4)中,floor()为下取整函数,返回不大于计算值的最大整数,sqrt()为开方根函数,β为 修正系数,取值范围0.5至2.0;
步骤CC:根据两个区间植被指数NDVI增量ΔNDVI的分类图,提取每个混合像元内各 类别的丰度,得到类别丰度图,丰度计算公式如下:
式(5)中,为混合像元(i,j)在k区间c类别的丰度,k值为1或2,1代表t0至tp区间,2代表tp至t1区间,Qc为混合像元(i,j)中c类别的像元数,wmix为混合像元的宽度;
步骤CD:根据混合像元植被指数NDVI增量ΔNDVI等于该混合像元内各类别的植被指数 NDVI增量ΔNDVI与其丰度的线性组合原理,构建t0至tp区间高低分辨率植被指数NDVI增量 ΔNDVI的线性混合模型,模型如下:
式(6)中,为混合像元(i,j)在t0至tp区间降尺度影像像元植被指数NDVI增 量的平均值,fc 1(i,j)为t0至tp区间混合像元(i,j)中c类别的丰度,为 c类别在t0至tp区间高精度影像像元植被指数NDVI增量值,wcal为求解窗口宽度;
步骤CE:根据混合像元植被指数NDVI增量ΔNDVI等于该混合像元内各类别的植被指数 NDVI增量ΔNDVI与其丰度的线性组合原理,构建tp至t1区间高低分辨率植被指数NDVI增量 ΔNDVI的线性混合模型,模型如下:
式(7)中,为混合像元(i,j)在tp至t1区间降尺度影像像元植被指数NDVI增 量的平均值,fc 2(i,j)为tp至t1区间混合像元(i,j)中c类别的丰度,为 c类别在tp至t1区间高精度影像像元植被指数NDVI增量值,wcal为求解窗口宽度;
步骤CF:利用岭回归方法求解式(6)和式(7)得到高分辨率植被指数NDVI增量
所述步骤D利用高分辨率植被指数NDVI增量计算tp时刻的 高分辨率植被指数值的计算公式如下:
式(8)与(9)中,分别为t0和t1时刻高分辨率影像植被指数NDVIF值,w1为 权值系数,分别为t0至tp区间、tp至t1区间降尺度影像像元植被指数 NDVIFC增量。
本发明提供的一种高时空分辨率植被指数NDVI的融合方法可以有效解决NDVI-LMGM 模型需要高分辨率地物分类影像及未考虑融合期间内地物覆盖发生变化情况的问题。利用本 发明方法、NDVI-LMGM方法与ESTARFM方法融合Landsat8和MODIS影像数据构建高时 空分辨率植物指数NDVI影像。结果表明,本发明方法比NDVI-LMGM方法与ESTARFM方法具有更高的融合精度及更接近真实影像的目视效果,选择ESTARFM方法公开的测试数据进行测试,可得到本发明方法、NDVI-LMGM方法和ESTARFM方法的融合结果与真实测量 结果之间的相关系数R分别为0.90、0.87和0.86,均方根误差RMSE分别为0.0464、0.0537 和0.0561,绝对误差AAD分别为0.0296、0.036和0.038,本发明方法的融合效果更令人满 意,效果更佳。
附图说明
图1为本发明流程示意图;
图2为高低分辨率影像示意图;
图3为求解窗口与混合像元的组合示意图
图4为融合过程中输入的植被指数NDVI影像图
图5为本发明方法融合植被指数NDVI的结果图
图6为NDVI-LMGM方法融合植被指数NDVI的结果图
图7为ESTARFM方法融合植被指数NDVI的结果图
图8为真实观测植被指数NDVI图
图9为三种不同融合方法预测值与真实观测值的植被指数NDVI散点图
具体实施方式
为了对本发明的技术特征、应用目的和融合效果有更加清楚的理解,现对照附图说明本 发明的具体实施过程。但本领域的技术人员应该知道,以下实施例并不是对本发明技术方案 作的唯一限定,凡是在本发明技术方案精神实质下所做的任何等同变换或改动,均应视为属 于本发明的保护范围。
本发明提供的高时空分辨率植被指数NDVI的融合方法基于线性混合模型、空间相似性 和时间相似性三个基本原理来构建,本发明方法是对Landsat NDVI影像和MODISNDVI影 像进行融合的方法。本发明方法的主要创新点是以植被指数增量△NDVI累计分布直方图进 行聚类生成分类图,用于将Landsat植被指数NDVI影像(高空间分辨率,低时间分辨率) 与MODIS植被指数NDVI影像(低空间分辨率,高时间分辨率)进行数据融合,从t0、t1两个基准时刻的植被指数NDVI影像获得预测时刻tp的高空间分辨率的植被指数NDVI影像。
下面详细说明本发明方法的具体步骤。
步骤A、低分辨率植被指数NDVI影像进行降尺度处理。
融合的高低影像之间通常没有很好的对应关系,低分辨率影像像元对应多个高分辨率影 像像元。虽然一个低分辨率(500m*500m)影像像元与多个高分辨率(30m*30m)影像像元 对应,但是其中一些高分辨率(30m*30m)影像像元位于两个或多个低分辨率(500m*500m) 影像像元对应区域,如图2中的灰色背影影像像元。混合像元和求解窗口也与低空间分辨率 影像像元存在复杂的对应关系,如图3所示。为后期计算方便,本发明方法直接将低分辨率 (500m*500m)影像进行降尺度处理,降尺度处理后的影像与高分辨率影像的像元大小一致, 且像元的地理位置一一对应。降尺度处理过程具体如下:
步骤AA、获取t0时刻的高分辨率植被指数NDVI影像的地理坐标和分辨率数据,基于t0时 刻的高分辨率植被指数NDVI影像,分别创建t0、tp、t1三个时刻低分辨率植被指数NDVIC影像 的降尺度影像模板;
步骤AB、如果t0、tp、t1三个时刻降尺度影像模板的像元植被指数值为空值, 则像元值不做修改,如果t0、tp、t1三个时刻降尺度影像模板的像元植被指数值为 非空值,则用降尺度公式进行重采样,该公式如下:
式(1)中,下标(x,y,t)表示t时刻(x,y)位置像元的植被指数NDVICF值,N为计算降尺度影 像像元与低分辨率影像地理区域重叠的低分辨率影像像元个数,S为降尺度影像像元地理面积, 为第i个与计算降尺度影像像元重叠的低分辨率影像像元地理重叠面积,为第i 个与计算降尺度影像像元重叠的低分辨率影像像元植被指数NDVIC值。
其处理过程首先以高分辨影像为母版建立t0、tp、t1三个时刻低分辨率植被指数NDVI 影像降尺度影像模板,扫描三个时刻降尺度影像模板,利用式(1)对模板中非空像元值进行重 采样。
步骤B、构建植被指数NDVI增量ΔNDVI的分类图。
遥感影像地物识别和分类是遥感图像计算机信息提取的研究热点。在应用遥感技术解决 实际问题时,常常需要根据地物的特征进行归类,这一工作称为分类。目前主要采用人工目 视解释分类和计算机地物自动分类两种方式,人工目视解释分类依据影像色调和几何特征等 解释标志来识别地物类型,该分类方法受人的主观经验认识影响较大且费时费力效率低下; 计算机地物自动分类是通过模式识别理论,利用计算机将遥感影像自动分成若干地物类别的 方法,计算机分类效率高,但分类精度相对于人工目视分类较弱。
为了提高融合算法计算效率,NDVI-LMGM方法和ESTARFM方法推荐使用无监督聚类分类或阈值分类方法,这两种方法只针对地物分类信息构建,缺少对融合区间内地物变化特 征信息提取。本发明公开的高时空分辨率植被指数NDVI的融合方法以降尺度处理的t0、tp、 t1三个时刻植被指数NVDIFC数据和t0、t1两个时刻高分辨率NDVIF数据设计三种植被指数 NDVI增量ΔNDVI,这三种增量ΔNDVI为无监督聚类分类的特征分量。构建植被指数NDVI增量ΔNDVI分类图的过程如下:
步骤BA:分别计算三种植被指数NDVI增量ΔNDVI影 像累积分布直方图;
步骤BB:分别将三种植被指数NDVI增量ΔNDVI累积分布平均划分为NUM等分,取累积分布每等分二分之一处的植被指数NDVI增量ΔNDVIi值为聚类特征分量中心,其中 i=1,2,…,NUM,NUM分类数由用户根据影像变化特征确定;依次取 的聚类特征分量中心构成特征空间矢量点,其中i,j=1,2,…,NUM,形成NUM*NUM个矢量点,并标注类别号,即为 NUM*NUM个类别的聚类中心点;
步骤BC:按两种植被指数NDVI增量ΔNDVI影像对应位置像元 ΔNDVI值构成1个2维矢量,依次计算该矢量与NUM*NUM个类别的聚类中心点 的欧式距离,其中i,j=1,2,…,NUM,取距离最小者的类别号为该像元类别号, 所有像元的类别号计算完成后,则构建出t0至tp区间植被指数NDVI增量ΔNDVI的分类图;
步骤BD:按两种植被指数NDVI增量ΔNDVI影像对应位置像元 ΔNDVI值构成1个2维矢量,依次计算该矢量与NUM*NUM个类别的聚类中心点 的欧式距离,其中i,j=1,2,…,NUM,取距离最小者的类别号为该像元类别号, 所有像元的类别号计算完成后,则构建出tp至t1区间植被指数NDVI增量ΔNDVI的分类图。
步骤C、计算高分辨率植被指数NDVI增量
在遥感图像信息提取技术中,认为一个像元仅包含一类地物信息时称为纯像元,而对于 地物复杂或纹理区域的像元来说,一个像元往往包含多种地物信息,这类像元称为混合像元。 在已知高精度分类信息时,利用混合像元分解技术可获得相应高精度信息的光谱特征信息或 高级别产品信息。最常见的分解模型有线性混合模型、模糊模型、概率模型、几何光学模型 和随机几何模型等,其中线性混合模型原理简单、物理意义明确、效果高,目前应用最广。 线性混合模型是利用一个线性关系来表达影像中一个像元内的各地物类型、比例与地物的特 性。本发明方法、NDVI-LMGM方法和ESTARFM方法等高时空分辨率融合方法都基于线性 混合模型,并假设高分辨率影像像元为纯像元。本发明方法以线性混合模型计算高分辨率植 被指数NDVI增量的具体过程如下:
步骤CA:根据NUM分类数计算混合像元宽度wmix,计算公式如下:
wmix=floor(sqrt(α*NUM))+1 (2)
式(2)中,floor()为下取整函数,返回不大于计算值的最大整数,sqrt()为开方根函数,α为 修正系数,取值范围1.5至4.0,在技术效果分析案例中取值为2.0;
步骤CB:根据NUM分类数计算求解窗口宽度wcal,计算公式如下:
wcal=floor(sqrt(β*NUM))+1 (3)
式(3)中,floor()为下取整函数,返回不大于计算值的最大整数,sqrt()为开方根函数,β为 修正系数,取值范围0.5至2.0,在技术效果分析案例中取值为1.0;
步骤CC:根据两个区间植被指数NDVI增量ΔNDVI的分类图,提取每个混合像元内各 类别的丰度,得到类别丰度图,丰度计算公式如下:
式(4)中,为混合像元(i,j)在k区间c类别的丰度,k值为1或2,1代表t0至tp区间, 2代表tp至t1区间,Qc为混合像元(i,j)中c类别的像元数,wmix为混合像元的宽度;
步骤CD:根据混合像元植被指数NDVI增量ΔNDVI等于该混合像元内各类别的植被指数 NDVI增量ΔNDVI与其丰度的线性组合原理,构建t0至tp区间高低分辨率植被指数NDVI增量 ΔNDVI的线性混合模型,,模型如下:
式(5)中,为混合像元(i,j)在t0至tp区间降尺度影像像元植被指数NDVI增 量的平均值,fc 1(i,j)为t0至tp区间混合像元(i,j)中c类别的丰度,为 c类别在t0至tp区间高精度影像像元植被指数NDVI增量值,wcal为求解窗口宽度;
步骤CE:根据混合像元植被指数NDVI增量ΔNDVI等于该混合像元内各类别的植被指数 NDVI增量ΔNDVI与其丰度的线性组合原理,构建tp至t1区间高低分辨率植被指数NDVI增量 ΔNDVI的线性混合模型,模型如下:
式(6)中,为混合像元(i,j)在tp至t1区间降尺度影像像元植被指数NDVI增 量的平均值,fc 2(i,j)为tp至t1区间混合像元(i,j)中c类别的丰度,为 c类别在tp至t1区间高精度影像像元植被指数NDVI增量值,wcal为求解窗口宽度;
步骤CF:利用岭回归方法求解式(6)和式(7)得到高分辨率植被指数NDVI增量
步骤D、计算tp时刻的高分辨率植被指数值。
根据本发明方法步骤C计算获得的高分辨率植被指数NDVI增量结合t0和t1时刻高分辨率影像植被指数NDVI值能得到两个tp时刻 的高分辨率植被指数值,为提高计算精度,本发明采用加权平均方法计算最终植被指 数值,具体计算公式如下:
式(7)与(8)中,分别为t0和t1时刻高分辨率影像植被指数NDVIF值,w1为 权值系数,分别为t0至tp区间、tp至t1区间降尺度影像像元植被指数 NDVIFC增量。
如上所述,经过步骤A、B、C和D的计算后,可得到预测时刻tp的高分辨率植被指数NDVI影像,本发明方法比NDVI-LMGM方法与ESTARFM方法更充分使用高/低分辨率影像 像元特征值的变化信息,具有更高的融合精度。
为了更好地说明本发明方法的技术效果,将本发明方法与背景技术部分提到的NDVI-LMGM方法与ESTARFM方法进行对比,即针对本发明方法、NDVI-LMGM方法与 ESTARFM方法,分别进行植被指数NDVI影像数据融合测试。测试数据选择ESTARFM方 法公开的测试数据,该影像数据从网站:https://xiaolinzhu.weebly.com/open-source-code.html 下载,对测试数据进行植被指数NDVI计算后进行融合预测。图4a-e分别显示了融合过程中 输入的植被指数NDVI影像。
图5、图6和图7分别显示的是本发明方法、NDVI-LMGM方法与ESTARFM方法进行 融合植被指数NDVI的结果图,图8是7月11日真实观测植被指数NDVI影像。将三种方法 的融合结果与真实影像进行比较,从局部放大图上可以得到,NDVI-LMGM方法的融合结果 有明显的斑块问题;ESTARFM方法融合结果的局部放大图中标记区域位置影像比真实观测 植被指数NDVI影像暗色区域增大,预测NDVI偏低;本发明方法的融合结果与真实观测植 被指数NDVI影像最为接近,比以上两种融合方法的融合效果更佳。为客观定量比较三种融 合方法的效果,随机生成100000个点作为验证点,来绘制植被指数NDVI预测值与真值之间 的散点图;并通过统计分析发现,本发明方法、NDVI-LMGM方法和ESTARFM方法融合植 被指数NDVI预测值与真值之间的相关系数R分别为0.90、0.87和0.86,均方根误差RMSE 分别为0.0464、0.0537和0.0561,绝对误差AAD分别为0.0296、0.036和0.038,本发明方 法能够更加精确反映融合期间内地物变化情况,融合准确度最高。
综上所述,本发明的高时空分辨率植被指数NDVI的融合方法将融合期间内地物覆盖发 生变化情况的问题进行综合,设计植被指数NDVI增量为地物分类聚类方法的特征分量,有 效提高了数据融合的精度。本发明充分利用了融合区间内高/低分辨率影像像元之间的变化特 性,为生产高时空分辨率的植被指数NDVI提供了一种行之有效的方法。
以上所述仅是本发明的优选实施方式,应当指出,对于本技术领域的普通技术人员来说, 在不脱离本发明技术原理的前提下,还可以做出若干的改进和替换,这些改进和替换也应视 为本发明的保护范围。

Claims (5)

1.一种高时空分辨率植被指数NDVI的融合方法,其特征在于:所述方法主要包括以下步骤:
步骤A:利用像元降尺度算法对t0、tp、t1三个时刻低分辨率植被指数NDVIC影像进行降尺度处理;
步骤B:利用降尺度处理的t0、tp、t1三个时刻植被指数NDVIFC影像和t0、t1两个时刻高分辨率植被指数NDVIF影像计算t0至tp区间降尺度植被指数NDVIFC增量影像、tp至t1区间降尺度植被指数NDVIFC增量影像、t0至t1区间高分辨率植被指数NDVIF增量影像,并根据三种植被指数NDVI增量ΔNDVI影像累积分布直方图确定聚类中心进行聚类构建出增量ΔNDVI分类图;
步骤C:利用增量ΔNDVI分类图,构建高低分辨率植被指数NDVI增量ΔNDVI的线性混合模型,并利用岭回归方法分别计算出高分辨率植被指数NDVI增量
步骤D:利用高分辨率植被指数NDVI增量计算tp时刻的高分辨率植被指数值;
所述步骤B、步骤C和步骤D中的ti至tj区间植被指数NDVI增量的计算公式如下:
式(1)中,为ti时刻植被指数NDVI值,为tj时刻植被指数NDVI值,ti与tj分别为计算增量的开始时刻和结束时刻;
所述步骤A至步骤D中的t0、tp、t1三个时刻分别为基准时刻、预测时刻、基准时刻,t0、tp、t1三个时刻均具有低分辨率植被指数NDVIC影像,只有在t0、t1两个时刻具有高分辨率植被指数NDVIF影像;
所述步骤A至步骤D中的植被指数NDVI与植被指数NDVI增量ΔNDVI的上标F、FC、C分别代表高分辨率值、低分辨率降尺度处理后值和低分辨率值。
2.根据权利要求1所述的一种高时空分辨率植被指数NDVI的融合方法,其特征在于:所述步骤A中的像元降尺度算法,包括如下步骤:
步骤AA:以t0时刻的高分辨率植被指数NDVIF影像为母版,分别创建t0、tp、t1三个时刻低分辨率植被指数NDVIC的降尺度影像模板;
步骤AB:如果t0、tp、t1三个时刻降尺度影像模板的像元值为空值,则像元值不做修改,如果t0、tp、t1三个时刻降尺度影像模板的像元值为非空值,则用降尺度公式进行重采样,公式如下:
式(2)中,下标(x,y,t)表示t时刻(x,y)位置像元的植被指数NDVICF值,N为计算降尺度影像像元与低分辨率影像地理区域重叠的低分辨率影像像元个数,S为降尺度影像像元地理面积,为第i个与计算降尺度影像像元重叠的低分辨率影像像元地理重叠面积,为第i个与计算降尺度影像像元重叠的低分辨率影像像元植被指数NDVIC值。
3.根据权利要求1所述的一种高时空分辨率植被指数NDVI的融合方法,其特征在于:所述步骤B利用三种植被指数NDVI增量ΔNDVI影像累积分布直方图确定聚类中心进行聚类构建出增量ΔNDVI分类图,包括如下步骤:
步骤BA:分别计算三种植被指数NDVI增量ΔNDVI影像累积分布直方图;
步骤BB:分别将三种植被指数NDVI增量ΔNDVI累积分布平均划分为NUM等分,取累积分布每等分二分之一处的植被指数NDVI增量ΔNDVIi值为聚类特征分量中心,其中i=1,2,…,NUM,NUM分类数由用户根据影像变化特征确定;依次取 的聚类特征分量中心构成特征空间矢量点,其中i,j=1,2,…,NUM,形成NUM*NUM个矢量点,并标注类别号,即为NUM*NUM个类别的聚类中心点;
步骤BC:按两种植被指数NDVI增量ΔNDVI影像对应位置像元ΔNDVI值构成1个2维矢量,依次计算该矢量与NUM*NUM个类别的聚类中心点 的欧式距离,其中i,j=1,2,…,NUM,取距离最小者的类别号为该像元类别号,所有像元的类别号计算完成后,则构建出t0至tp区间植被指数NDVI增量ΔNDVI的分类图;
步骤BD:按两种植被指数NDVI增量ΔNDVI影像对应位置像元ΔNDVI值构成1个2维矢量,依次计算该矢量与NUM*NUM个类别的聚类中心点 的欧式距离,其中i,j=1,2,…,NUM,取距离最小者的类别号为该像元类别号,所有像元的类别号计算完成后,则构建出tp至t1区间植被指数NDVI增量ΔNDVI的分类图。
4.根据权利要求1所述的一种高时空分辨率植被指数NDVI的融合方法,其特征在于:所述步骤C利用植被指数NDVI增量ΔNDVI分类图,构建高低分辨率植被指数NDVI增量ΔNDVI的线性混合模型,并利用岭回归方法求高分辨率植被指数NDVI增量包括如下步骤:
步骤CA:根据NUM分类数计算混合像元宽度wmix,计算公式如下:
wmix=floor(sqrt(α*NUM))+1 (3)
式(3)中,floor()为下取整函数,返回不大于计算值的最大整数,sqrt()为开方根函数,α为修正系数,取值范围1.5至4.0;
步骤CB:根据NUM分类数计算求解窗口宽度wcal,计算公式如下:
wcal=floor(sqrt(β*NUM))+1 (4)
式(4)中,floor()为下取整函数,返回不大于计算值的最大整数,sqrt()为开方根函数,β为修正系数,取值范围0.5至2.0;
步骤CC:根据两个区间植被指数NDVI增量ΔNDVI的分类图,提取每个混合像元内各类别的丰度,得到类别丰度图,丰度计算公式如下:
式(5)中,为混合像元(i,j)在k区间c类别的丰度,k值为1或2,1代表t0至tp区间,2代表tp至t1区间,Qc为混合像元(i,j)中c类别的像元数,wmix为混合像元的宽度;
步骤CD:根据混合像元植被指数NDVI增量ΔNDVI等于该混合像元内各类别的植被指数NDVI增量ΔNDVI与其丰度的线性组合原理,构建t0至tp区间高低分辨率植被指数NDVI增量ΔNDVI的线性混合模型,模型如下:
式(6)中,为混合像元(i,j)在t0至tp区间降尺度影像像元植被指数NDVI增量的平均值,为t0至tp区间混合像元(i,j)中c类别的丰度,为c类别在t0至tp区间高精度影像像元植被指数NDVI增量值,wcal为求解窗口宽度;
步骤CE:根据混合像元植被指数NDVI增量ΔNDVI等于该混合像元内各类别的植被指数NDVI增量ΔNDVI与其丰度的线性组合原理,构建tp至t1区间高低分辨率植被指数NDVI增量ΔNDVI的线性混合模型,模型如下:
式(7)中,为混合像元(i,j)在tp至t1区间降尺度影像像元植被指数NDVI增量的平均值,为tp至t1区间混合像元(i,j)中c类别的丰度,为c类别在tp至t1区间高精度影像像元植被指数NDVI增量值,wcal为求解窗口宽度;
步骤CF:利用岭回归方法求解式(6)和式(7)得到高分辨率植被指数NDVI增量
5.根据权利要求1所述的一种高时空分辨率植被指数NDVI的融合方法,其特征在于:所述步骤D利用高分辨率植被指数NDVI增量计算tp时刻的高分辨率植被指数值的计算公式如下:
式(8)与(9)中,分别为t0和t1时刻高分辨率影像植被指数NDVIF值,w1为权值系数,分别为t0至tp区间、tp至t1区间降尺度影像像元植被指数NDVIFC增量。
CN201910648630.2A 2019-07-18 2019-07-18 一种高时空分辨率植被指数ndvi的融合方法 Active CN110363246B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910648630.2A CN110363246B (zh) 2019-07-18 2019-07-18 一种高时空分辨率植被指数ndvi的融合方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910648630.2A CN110363246B (zh) 2019-07-18 2019-07-18 一种高时空分辨率植被指数ndvi的融合方法

Publications (2)

Publication Number Publication Date
CN110363246A true CN110363246A (zh) 2019-10-22
CN110363246B CN110363246B (zh) 2023-05-09

Family

ID=68220342

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910648630.2A Active CN110363246B (zh) 2019-07-18 2019-07-18 一种高时空分辨率植被指数ndvi的融合方法

Country Status (1)

Country Link
CN (1) CN110363246B (zh)

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110909821A (zh) * 2019-12-03 2020-03-24 中国农业科学院农业资源与农业区划研究所 基于作物参考曲线进行高时空分辨率植被指数数据融合的方法
CN110991567A (zh) * 2019-12-27 2020-04-10 河南大学 一种由点到面的遥感瞬时地表温度数据检测方法
CN111523451A (zh) * 2020-04-22 2020-08-11 重庆邮电大学 一种构建高时空分辨率ndvi数据的方法
CN112052720A (zh) * 2020-07-16 2020-12-08 贵州师范学院 一种基于直方图聚类的高时空归一化植被指数ndvi的融合模型
CN112989286A (zh) * 2021-03-22 2021-06-18 自然资源部国土卫星遥感应用中心 一种时空信息融合的微波遥感土壤水分产品降尺度方法
CN113077384A (zh) * 2021-03-10 2021-07-06 中山大学 一种数据空间分辨率提高方法、装置、介质及终端设备
CN113327197A (zh) * 2021-05-10 2021-08-31 香港理工大学深圳研究院 遥感影像时空融合方法、智能终端及计算机可读存储介质
CN116258869A (zh) * 2023-01-10 2023-06-13 滁州学院 一种基于Sentinel-2遥感数据的毛竹林大小年边界线提取方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
AU2012101249A4 (en) * 2012-08-17 2012-09-20 Beijing Normal University Method for Generating High Spatial Resolution NDVI Time Series Data
CN105046648A (zh) * 2015-06-25 2015-11-11 北京师范大学 一种构建高时空遥感数据的方法
WO2016132161A1 (en) * 2015-02-16 2016-08-25 Kontoes Charalampos Method that detects areas of active fire hotspots in real-time, calculates the most probable ignition point and assesses fire probability indicators, using satellite images and fuel data.
CN107103584A (zh) * 2017-04-11 2017-08-29 北京师范大学 一种基于时空加权的生产高时空分辨率ndvi的方法
CN109753916A (zh) * 2018-12-28 2019-05-14 厦门理工学院 一种归一化差分植被指数尺度转换模型构建方法及装置

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
AU2012101249A4 (en) * 2012-08-17 2012-09-20 Beijing Normal University Method for Generating High Spatial Resolution NDVI Time Series Data
WO2016132161A1 (en) * 2015-02-16 2016-08-25 Kontoes Charalampos Method that detects areas of active fire hotspots in real-time, calculates the most probable ignition point and assesses fire probability indicators, using satellite images and fuel data.
CN105046648A (zh) * 2015-06-25 2015-11-11 北京师范大学 一种构建高时空遥感数据的方法
CN107103584A (zh) * 2017-04-11 2017-08-29 北京师范大学 一种基于时空加权的生产高时空分辨率ndvi的方法
CN109753916A (zh) * 2018-12-28 2019-05-14 厦门理工学院 一种归一化差分植被指数尺度转换模型构建方法及装置

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
苑惠丽;马荣华;李吉英;: "一种平原区园地遥感信息提取的新方法", 中国科学院大学学报 *
苑惠丽;马荣华;李吉英;余艳玲;: "基于Landsat-8 OLI的农作物信息提取研究――以安徽省蚌埠市为例", 金陵科技学院学报 *

Cited By (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110909821A (zh) * 2019-12-03 2020-03-24 中国农业科学院农业资源与农业区划研究所 基于作物参考曲线进行高时空分辨率植被指数数据融合的方法
CN110909821B (zh) * 2019-12-03 2020-07-28 中国农业科学院农业资源与农业区划研究所 基于作物参考曲线进行高时空分辨率植被指数数据融合的方法
CN110991567A (zh) * 2019-12-27 2020-04-10 河南大学 一种由点到面的遥感瞬时地表温度数据检测方法
CN111523451A (zh) * 2020-04-22 2020-08-11 重庆邮电大学 一种构建高时空分辨率ndvi数据的方法
CN111523451B (zh) * 2020-04-22 2023-05-30 重庆邮电大学 一种构建高时空分辨率ndvi数据的方法
CN112052720A (zh) * 2020-07-16 2020-12-08 贵州师范学院 一种基于直方图聚类的高时空归一化植被指数ndvi的融合模型
CN113077384A (zh) * 2021-03-10 2021-07-06 中山大学 一种数据空间分辨率提高方法、装置、介质及终端设备
CN113077384B (zh) * 2021-03-10 2022-04-29 中山大学 一种数据空间分辨率提高方法、装置、介质及终端设备
CN112989286A (zh) * 2021-03-22 2021-06-18 自然资源部国土卫星遥感应用中心 一种时空信息融合的微波遥感土壤水分产品降尺度方法
CN113327197A (zh) * 2021-05-10 2021-08-31 香港理工大学深圳研究院 遥感影像时空融合方法、智能终端及计算机可读存储介质
CN113327197B (zh) * 2021-05-10 2023-01-24 香港理工大学深圳研究院 遥感影像时空融合方法、智能终端及计算机可读存储介质
CN116258869A (zh) * 2023-01-10 2023-06-13 滁州学院 一种基于Sentinel-2遥感数据的毛竹林大小年边界线提取方法
CN116258869B (zh) * 2023-01-10 2023-08-18 滁州学院 一种基于Sentinel-2遥感数据的毛竹林大小年边界线提取方法

Also Published As

Publication number Publication date
CN110363246B (zh) 2023-05-09

Similar Documents

Publication Publication Date Title
CN110363246A (zh) 一种高时空分辨率植被指数ndvi的融合方法
Halme et al. Utility of hyperspectral compared to multispectral remote sensing data in estimating forest biomass and structure variables in Finnish boreal forest
CN104751478B (zh) 一种基于多特征融合的面向对象的建筑物变化检测方法
CN103364781B (zh) 基于遥感数据与地理信息系统的粮田地面参照点筛选方法
Liu et al. Estimating potato above-ground biomass by using integrated unmanned aerial system-based optical, structural, and textural canopy measurements
CN112381013B (zh) 基于高分辨率遥感影像的城市植被反演方法及系统
CN112348812B (zh) 林分年龄信息测量方法及装置
Guan et al. Modeling strawberry biomass and leaf area using object-based analysis of high-resolution images
WO2001033505A2 (en) Multi-variable model for identifying crop response zones in a field
CN106501186B (zh) 一种土壤含水量产品降尺度方法
CN110309780A (zh) 基于bfd-iga-svm模型的高分辨率影像房屋信息快速监督识别
Sun et al. Wheat head counting in the wild by an augmented feature pyramid networks-based convolutional neural network
CN114529097B (zh) 多尺度农作物物候期遥感降维预测方法
CN111582575A (zh) 一种多时空尺度下城市热环境形成发展主导因素识别方法
Zhou et al. Wheat phenology detection with the methodology of classification based on the time-series UAV images
CN112836725A (zh) 基于时序遥感数据的弱监督lstm循环神经网络稻田识别方法
CN111723711A (zh) 基于Pléiades和面向对象的地膜信息提取方法及系统
CN116403048B (zh) 一种基于多模态数据融合的农作物生长估计模型构建方法
Tian et al. Machine learning-based crop recognition from aerial remote sensing imagery
Zhang et al. Crop type mapping with temporal sample migration
AU2021101780A4 (en) Aboveground Biomass Estimation and Scale Conversion for Mean Regional Spectral Units
CN116863341B (zh) 基于时间序列卫星遥感影像的作物分类和识别方法和系统
Engstrom et al. Evaluating the Relationship between Contextual Features Derived from Very High Spatial Resolution Imagery and Urban Attributes: A Case Study in Sri Lanka
Huang et al. Information fusion approach for biomass estimation in a plateau mountainous forest using a synergistic system comprising UAS-based digital camera and LiDAR
Zhao et al. Improving object-oriented land use/cover classification from high resolution imagery by spectral similarity-based post-classification

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant